Merge branch 'master' into 943-xBBCSD-read-of-unitialised-value

This commit is contained in:
langou
2026-03-25 13:21:59 -06:00
committed by GitHub
6319 changed files with 101044 additions and 54638 deletions
+1 -1
View File
@@ -33,4 +33,4 @@ build_script:
- cmake --build .
test_script:
- ctest -j2
- ctest -j2 --output-on-failure
+65
View File
@@ -0,0 +1,65 @@
using BinaryBuilder, Pkg
haskey(ENV, "BLAS_LAPACK_RELEASE") || error("The environment variable BLAS_LAPACK_RELEASE is not defined.")
haskey(ENV, "BLAS_LAPACK_COMMIT") || error("The environment variable BLAS_LAPACK_COMMIT is not defined.")
haskey(ENV, "BLAS_LAPACK_URL") || error("The environment variable BLAS_LAPACK_URL is not defined.")
name = "blas_lapack"
version = VersionNumber(ENV["BLAS_LAPACK_RELEASE"])
# Collection of sources required to complete build
sources = [
GitSource(ENV["BLAS_LAPACK_URL"], ENV["BLAS_LAPACK_COMMIT"])
]
# Bash recipe for building across all platforms
script = raw"""
cd ${WORKSPACE}/srcdir/lapack
# FortranCInterface_VERIFY fails on macOS, but it's not actually needed for the current build
sed -i 's/FortranCInterface_VERIFY/# FortranCInterface_VERIFY/g' ./CBLAS/CMakeLists.txt
sed -i 's/FortranCInterface_VERIFY/# FortranCInterface_VERIFY/g' ./LAPACKE/include/CMakeLists.txt
mkdir build && cd build
cmake .. \
-DCBLAS=ON \
-DLAPACKE=ON \
-DCMAKE_INSTALL_PREFIX="$prefix" \
-DCMAKE_FIND_ROOT_PATH="$prefix" \
-DCMAKE_TOOLCHAIN_FILE="${CMAKE_TARGET_TOOLCHAIN}" \
-DCMAKE_BUILD_TYPE=Release \
-DBUILD_SHARED_LIBS=OFF \
-DBUILD_INDEX64_EXT_API=OFF \
-DTEST_FORTRAN_COMPILER=OFF \
-DLAPACKE_WITH_TMG=OFF
make -j${nproc}
make install
install_license $WORKSPACE/srcdir/lapack/LICENSE
"""
# These are the platforms we will build for by default, unless further
# platforms are passed in on the command line
platforms = supported_platforms()
platforms = expand_gfortran_versions(platforms)
# The products that we will ensure are always built
products = [
FileProduct("lib/libblas.a", :libblas_a),
FileProduct("lib/libcblas.a", :libcblas_a),
FileProduct("lib/liblapack.a", :liblapack_a),
FileProduct("lib/liblapacke.a", :liblapacke_a),
# LibraryProduct("libblas", :libblas),
# LibraryProduct("libcblas", :libcblas),
# LibraryProduct("liblapack", :liblapack),
# LibraryProduct("liblapacke", :liblapacke),
]
# Dependencies that must be installed before this package can be built
dependencies = [
Dependency(PackageSpec(name="CompilerSupportLibraries_jll", uuid="e66e0078-7015-5450-92f7-15fbd957f2ae")),
]
# Build the tarballs, and possibly a `build.jl` as well.
build_tarballs(ARGS, name, version, sources, script, platforms, products, dependencies; julia_compat="1.6")
+90
View File
@@ -0,0 +1,90 @@
# Version
haskey(ENV, "BLAS_LAPACK_RELEASE") || error("The environment variable BLAS_LAPACK_RELEASE is not defined.")
version = VersionNumber(ENV["BLAS_LAPACK_RELEASE"])
version2 = ENV["BLAS_LAPACK_RELEASE"]
package = "blas_lapack"
platforms = [
("aarch64-apple-darwin-libgfortran5" , "lib", "dylib"),
# ("aarch64-linux-gnu-libgfortran3" , "lib", "so" ),
# ("aarch64-linux-gnu-libgfortran4" , "lib", "so" ),
("aarch64-linux-gnu-libgfortran5" , "lib", "so" ),
# ("aarch64-linux-musl-libgfortran3" , "lib", "so" ),
# ("aarch64-linux-musl-libgfortran4" , "lib", "so" ),
# ("aarch64-linux-musl-libgfortran5" , "lib", "so" ),
# ("powerpc64le-linux-gnu-libgfortran3" , "lib", "so" ),
# ("powerpc64le-linux-gnu-libgfortran4" , "lib", "so" ),
# ("powerpc64le-linux-gnu-libgfortran5" , "lib", "so" ),
# ("x86_64-apple-darwin-libgfortran3" , "lib", "dylib"),
# ("x86_64-apple-darwin-libgfortran4" , "lib", "dylib"),
("x86_64-apple-darwin-libgfortran5" , "lib", "dylib"),
# ("x86_64-linux-gnu-libgfortran3" , "lib", "so" ),
# ("x86_64-linux-gnu-libgfortran4" , "lib", "so" ),
("x86_64-linux-gnu-libgfortran5" , "lib", "so" ),
# ("x86_64-linux-musl-libgfortran3" , "lib", "so" ),
# ("x86_64-linux-musl-libgfortran4" , "lib", "so" ),
# ("x86_64-linux-musl-libgfortran5" , "lib", "so" ),
# ("x86_64-unknown-freebsd-libgfortran3", "lib", "so" ),
# ("x86_64-unknown-freebsd-libgfortran4", "lib", "so" ),
# ("x86_64-unknown-freebsd-libgfortran5", "lib", "so" ),
# ("x86_64-w64-mingw32-libgfortran3" , "bin", "dll" ),
# ("x86_64-w64-mingw32-libgfortran4" , "bin", "dll" ),
("x86_64-w64-mingw32-libgfortran5" , "bin", "dll" ),
]
for (platform, libdir, ext) in platforms
tarball_name = "$package.v$version.$platform.tar.gz"
if isfile("products/$(tarball_name)")
# Unzip the tarball generated by BinaryBuilder.jl
isdir("products/$platform") && rm("products/$platform", recursive=true)
mkdir("products/$platform")
run(`tar -xzf products/$(tarball_name) -C products/$platform`)
if isfile("products/$platform/deps.tar.gz")
# Unzip the tarball of the dependencies
run(`tar -xzf products/$platform/deps.tar.gz -C products/$platform`)
# Copy the license of each dependency
for folder in readdir("products/$platform/deps/licenses")
cp("products/$platform/deps/licenses/$folder", "products/$platform/share/licenses/$folder")
end
rm("products/$platform/deps/licenses", recursive=true)
# Copy the shared library of each dependency
for file in readdir("products/$platform/deps")
cp("products/$platform/deps/$file", "products/$platform/$libdir/$file")
end
# Remove the folder used to unzip the tarball of the dependencies
rm("products/$platform/deps", recursive=true)
rm("products/$platform/deps.tar.gz", recursive=true)
end
# Create the archives *_binaries
isfile("$(package)_binaries.$version2.$platform.tar.gz") && rm("$(package)_binaries.$version2.$platform.tar.gz")
isfile("$(package)_binaries.$version2.$platform.zip") && rm("$(package)_binaries.$version2.$platform.zip")
cd("products/$platform")
# Create a folder with the version number of the package
mkdir("$(package)_binaries.$version2")
for folder in ("include", "share", "lib")
cp(folder, "$(package)_binaries.$version2/$folder")
end
cd("$(package)_binaries.$version2")
if ext == "dll"
run(`zip -r --symlinks ../../../$(package)_binaries.$version2.$platform.zip include share lib`)
else
run(`tar -czf ../../../$(package)_binaries.$version2.$platform.tar.gz include share lib`)
end
cd("../../..")
# Remove the folder used to unzip the tarball generated by BinaryBuilder.jl
rm("products/$platform", recursive=true)
else
@warn("The tarball for the platform $platform was not generated!")
end
end
+5 -5
View File
@@ -75,12 +75,12 @@ jobs:
- name: Install ninja-build tool
uses: seanmiddleditch/gha-setup-ninja@16b940825621068d98711680b6c3ff92201f8fc0 # v3
- name: Use GCC-11 on MacOS
- name: Use GCC-14 on MacOS
if: ${{ matrix.os == 'macos-latest' }}
run: >
cmake -B build -G Ninja
-D CMAKE_C_COMPILER="gcc-11"
-D CMAKE_Fortran_COMPILER="gfortran-11"
-D CMAKE_C_COMPILER="gcc-14"
-D CMAKE_Fortran_COMPILER="gfortran-14"
-D USE_FLAT_NAMESPACE:BOOL=ON
- name: Special flags for Windows
@@ -221,7 +221,7 @@ jobs:
- name: Test
working-directory: ${{github.workspace}}/build
run: |
ctest -C ${{env.BUILD_TYPE}} --schedule-random -j2 -T memcheck > memcheck.out
ctest -C ${{env.BUILD_TYPE}} --output-on-failure --schedule-random -j2 -T memcheck > memcheck.out
cat memcheck.out
if tail -n 1 memcheck.out | grep -q "Memory checking results:"; then
exit 0
@@ -237,4 +237,4 @@ jobs:
fi
done
exit 0
fi
fi
+2 -2
View File
@@ -90,8 +90,8 @@ jobs:
echo "DOCSDIR = ${{github.workspace}}/DOCS" >> make.inc
- name: Alias for GCC compilers
run: |
sudo ln -s $(which gcc-11) /usr/local/bin/gcc
sudo ln -s $(which gfortran-11) /usr/local/bin/gfortran
sudo ln -s $(which gcc-14) /usr/local/bin/gcc
sudo ln -s $(which gfortran-14) /usr/local/bin/gfortran
- name: Install
run: |
make -s -j2 all
+239
View File
@@ -0,0 +1,239 @@
name: Release
on:
push:
# Sequence of patterns matched against refs/tags
tags:
- 'v*' # Push events to matching v*, i.e. v1.0, v2023.11.15
jobs:
build-linux-x64:
name: blas / lapack -- Linux (x86_64) -- Release ${{ github.ref_name }}
runs-on: ubuntu-latest
steps:
- name: Checkout lapack
uses: actions/checkout@v4
- name: Install Julia
uses: julia-actions/setup-julia@v2
with:
version: "1.7"
arch: x64
- name: Set the environment variables BINARYBUILDER_AUTOMATIC_APPLE, BLAS_LAPACK_RELEASE, BLAS_LAPACK_COMMIT
shell: bash
run: |
echo "BINARYBUILDER_AUTOMATIC_APPLE=true" >> $GITHUB_ENV
echo "BLAS_LAPACK_RELEASE=${{ github.ref_name }}" >> $GITHUB_ENV
echo "BLAS_LAPACK_COMMIT=${{ github.sha }}" >> $GITHUB_ENV
echo "BLAS_LAPACK_URL=https://github.com/${{ github.repository }}.git" >> $GITHUB_ENV
- name: Cross-compilation of blas / lapack -- x86_64-linux-gnu-libgfortran5
run: |
julia --color=no -e 'using Pkg; Pkg.add("BinaryBuilder")'
julia --color=no .github/julia/build_tarballs.jl x86_64-linux-gnu-libgfortran5 --verbose
- name: Archive artifact
run: julia --color=no .github/julia/generate_binaries.jl
- name: Upload artifact
uses: actions/upload-artifact@v4
with:
name: blas_lapack_binaries.${{ github.ref_name }}.x86_64-linux-gnu-libgfortran5.tar.gz
path: ./blas_lapack_binaries.${{ github.ref_name }}.x86_64-linux-gnu-libgfortran5.tar.gz
build-linux-aarch64:
name: blas / lapack -- Linux (aarch64) -- Release ${{ github.ref_name }}
runs-on: ubuntu-latest
steps:
- name: Checkout lapack
uses: actions/checkout@v4
- name: Install Julia
uses: julia-actions/setup-julia@v2
with:
version: "1.7"
arch: x64
- name: Set the environment variables BINARYBUILDER_AUTOMATIC_APPLE, BLAS_LAPACK_RELEASE, BLAS_LAPACK_COMMIT
shell: bash
run: |
echo "BINARYBUILDER_AUTOMATIC_APPLE=true" >> $GITHUB_ENV
echo "BLAS_LAPACK_RELEASE=${{ github.ref_name }}" >> $GITHUB_ENV
echo "BLAS_LAPACK_COMMIT=${{ github.sha }}" >> $GITHUB_ENV
echo "BLAS_LAPACK_URL=https://github.com/${{ github.repository }}.git" >> $GITHUB_ENV
- name: Cross-compilation of blas / lapack -- aarch64-linux-gnu-libgfortran5
run: |
julia --color=no -e 'using Pkg; Pkg.add("BinaryBuilder")'
julia --color=no .github/julia/build_tarballs.jl aarch64-linux-gnu-libgfortran5 --verbose
- name: Archive artifact
run: julia --color=no .github/julia/generate_binaries.jl
- name: Upload artifact
uses: actions/upload-artifact@v4
with:
name: blas_lapack_binaries.${{ github.ref_name }}.aarch64-linux-gnu-libgfortran5.tar.gz
path: ./blas_lapack_binaries.${{ github.ref_name }}.aarch64-linux-gnu-libgfortran5.tar.gz
build-windows-x64:
name: blas / lapack -- Windows (x86_64) -- Release ${{ github.ref_name }}
runs-on: ubuntu-latest
steps:
- name: Checkout lapack
uses: actions/checkout@v4
- name: Install Julia
uses: julia-actions/setup-julia@v2
with:
version: "1.7"
arch: x64
- name: Set the environment variables BINARYBUILDER_AUTOMATIC_APPLE, BLAS_LAPACK_RELEASE, BLAS_LAPACK_COMMIT
shell: bash
run: |
echo "BINARYBUILDER_AUTOMATIC_APPLE=true" >> $GITHUB_ENV
echo "BLAS_LAPACK_RELEASE=${{ github.ref_name }}" >> $GITHUB_ENV
echo "BLAS_LAPACK_COMMIT=${{ github.sha }}" >> $GITHUB_ENV
echo "BLAS_LAPACK_URL=https://github.com/${{ github.repository }}.git" >> $GITHUB_ENV
- name: Cross-compilation of blas / lapack -- x86_64-w64-mingw32-libgfortran5
run: |
julia --color=no -e 'using Pkg; Pkg.add("BinaryBuilder")'
julia --color=no .github/julia/build_tarballs.jl x86_64-w64-mingw32-libgfortran5 --verbose
- name: Archive artifact
run: julia --color=no .github/julia/generate_binaries.jl
- name: Upload artifact
uses: actions/upload-artifact@v4
with:
name: blas_lapack_binaries.${{ github.ref_name }}.x86_64-w64-mingw32-libgfortran5.zip
path: ./blas_lapack_binaries.${{ github.ref_name }}.x86_64-w64-mingw32-libgfortran5.zip
build-mac-x64:
name: blas / lapack -- macOS (x86_64) -- Release ${{ github.ref_name }}
runs-on: ubuntu-latest
steps:
- name: Checkout lapack
uses: actions/checkout@v4
- name: Install Julia
uses: julia-actions/setup-julia@v2
with:
version: "1.7"
arch: x64
- name: Set the environment variables BINARYBUILDER_AUTOMATIC_APPLE, BLAS_LAPACK_RELEASE, BLAS_LAPACK_COMMIT
shell: bash
run: |
echo "BINARYBUILDER_AUTOMATIC_APPLE=true" >> $GITHUB_ENV
echo "BLAS_LAPACK_RELEASE=${{ github.ref_name }}" >> $GITHUB_ENV
echo "BLAS_LAPACK_COMMIT=${{ github.sha }}" >> $GITHUB_ENV
echo "BLAS_LAPACK_URL=https://github.com/${{ github.repository }}.git" >> $GITHUB_ENV
- name: Cross-compilation of blas / lapack -- x86_64-apple-darwin-libgfortran5
run: |
julia --color=no -e 'using Pkg; Pkg.add("BinaryBuilder")'
julia --color=no .github/julia/build_tarballs.jl x86_64-apple-darwin-libgfortran5 --verbose
- name: Archive artifact
run: julia --color=no .github/julia/generate_binaries.jl
- name: Upload artifact
uses: actions/upload-artifact@v4
with:
name: blas_lapack_binaries.${{ github.ref_name }}.x86_64-apple-darwin-libgfortran5.tar.gz
path: ./blas_lapack_binaries.${{ github.ref_name }}.x86_64-apple-darwin-libgfortran5.tar.gz
build-mac-aarch64:
name: blas / lapack -- macOS (aarch64) -- Release ${{ github.ref_name }}
runs-on: ubuntu-latest
steps:
- name: Checkout lapack
uses: actions/checkout@v4
- name: Install Julia
uses: julia-actions/setup-julia@v2
with:
version: "1.7"
arch: x64
- name: Set the environment variables BINARYBUILDER_AUTOMATIC_APPLE, BLAS_LAPACK_RELEASE, BLAS_LAPACK_COMMIT
shell: bash
run: |
echo "BINARYBUILDER_AUTOMATIC_APPLE=true" >> $GITHUB_ENV
echo "BLAS_LAPACK_RELEASE=${{ github.ref_name }}" >> $GITHUB_ENV
echo "BLAS_LAPACK_COMMIT=${{ github.sha }}" >> $GITHUB_ENV
echo "BLAS_LAPACK_URL=https://github.com/${{ github.repository }}.git" >> $GITHUB_ENV
- name: Cross-compilation of blas / lapack -- aarch64-apple-darwin-libgfortran5
run: |
julia --color=no -e 'using Pkg; Pkg.add("BinaryBuilder")'
julia --color=no .github/julia/build_tarballs.jl aarch64-apple-darwin-libgfortran5 --verbose
- name: Archive artifact
run: julia --color=no .github/julia/generate_binaries.jl
- name: Upload artifact
uses: actions/upload-artifact@v4
with:
name: blas_lapack_binaries.${{ github.ref_name }}.aarch64-apple-darwin-libgfortran5.tar.gz
path: ./blas_lapack_binaries.${{ github.ref_name }}.aarch64-apple-darwin-libgfortran5.tar.gz
release:
name: Create Release and Upload Binaries
needs: [build-windows-x64, build-linux-x64, build-linux-aarch64, build-mac-x64, build-mac-aarch64]
runs-on: ubuntu-latest
steps:
- name: Checkout lapack
uses: actions/checkout@v4
- name: Download artifacts
uses: actions/download-artifact@v4
with:
path: .
- name: Create GitHub Release
run: |
gh release create ${{ github.ref_name }} \
--title "${{ github.ref_name }}" \
--notes "" \
--verify-tag
env:
GH_TOKEN: ${{ secrets.GITHUB_TOKEN }}
- name: Upload Linux (x86_64) artifact
run: |
gh release upload ${{ github.ref_name }} \
blas_lapack_binaries.${{ github.ref_name }}.x86_64-linux-gnu-libgfortran5.tar.gz/blas_lapack_binaries.${{ github.ref_name }}.x86_64-linux-gnu-libgfortran5.tar.gz#blas_lapack.${{ github.ref_name }}.linux.x86_64.tar.gz
env:
GH_TOKEN: ${{ secrets.GITHUB_TOKEN }}
- name: Upload Linux (aarch64) artifact
run: |
gh release upload ${{ github.ref_name }} \
blas_lapack_binaries.${{ github.ref_name }}.aarch64-linux-gnu-libgfortran5.tar.gz/blas_lapack_binaries.${{ github.ref_name }}.aarch64-linux-gnu-libgfortran5.tar.gz#blas_lapack.${{ github.ref_name }}.linux.aarch64.tar.gz
env:
GH_TOKEN: ${{ secrets.GITHUB_TOKEN }}
- name: Upload Mac (x86_64) artifact
run: |
gh release upload ${{ github.ref_name }} \
blas_lapack_binaries.${{ github.ref_name }}.x86_64-apple-darwin-libgfortran5.tar.gz/blas_lapack_binaries.${{ github.ref_name }}.x86_64-apple-darwin-libgfortran5.tar.gz#blas_lapack.${{ github.ref_name }}.mac.x86_64.tar.gz
env:
GH_TOKEN: ${{ secrets.GITHUB_TOKEN }}
- name: Upload Mac (aarch64) artifact
run: |
gh release upload ${{ github.ref_name }} \
blas_lapack_binaries.${{ github.ref_name }}.aarch64-apple-darwin-libgfortran5.tar.gz/blas_lapack_binaries.${{ github.ref_name }}.aarch64-apple-darwin-libgfortran5.tar.gz#blas_lapack.${{ github.ref_name }}.mac.aarch64.tar.gz
env:
GH_TOKEN: ${{ secrets.GITHUB_TOKEN }}
- name: Upload Windows (x86_64) artifact
run: |
gh release upload ${{ github.ref_name }} \
blas_lapack_binaries.${{ github.ref_name }}.x86_64-w64-mingw32-libgfortran5.zip/blas_lapack_binaries.${{ github.ref_name }}.x86_64-w64-mingw32-libgfortran5.zip#blas_lapack.${{ github.ref_name }}.windows.x86_64.zip
env:
GH_TOKEN: ${{ secrets.GITHUB_TOKEN }}
+4 -4
View File
@@ -32,12 +32,12 @@ jobs:
steps:
- name: "Checkout code"
uses: actions/checkout@c85c95e3d7251135ab7dc9ce3241c5835cc595a9 # v3.5.3
uses: actions/checkout@d632683dd7b4114ad314bca15554477dd762a938 # tag=v4.2.0
with:
persist-credentials: false
- name: "Run analysis"
uses: ossf/scorecard-action@08b4669551908b1024bb425080c797723083c031 # v2.2.0
uses: ossf/scorecard-action@62b2cac7ed8198b15735ed49ab1e5cf35480ba46 # v2.4.0
with:
results_file: results.sarif
results_format: sarif
@@ -59,7 +59,7 @@ jobs:
# Upload the results as artifacts (optional). Commenting out will disable uploads of run results in SARIF
# format to the repository Actions tab.
- name: "Upload artifact"
uses: actions/upload-artifact@0b7f8abb1508181956e8e162db84b466c27e18ce # v3.1.2
uses: actions/upload-artifact@b4b15b8c7c6ac21ea08fcf65892d2ee8f75cf882 # v4.4.3
with:
name: SARIF file
path: results.sarif
@@ -67,6 +67,6 @@ jobs:
# Upload the results to GitHub's code scanning dashboard.
- name: "Upload to code-scanning"
uses: github/codeql-action/upload-sarif@f9a7c6738f28efb36e31d49c53a201a9c5d6a476 # v2.14.2
uses: github/codeql-action/upload-sarif@662472033e021d55d94146f66f6058822b0b39fd # v3.27.0
with:
sarif_file: results.sarif
+4
View File
@@ -43,3 +43,7 @@ build*
DOCS/man
DOCS/explore-html
output_err
# Mod files from compilation in SRC
SRC/la_constants.mod
SRC/la_xisnan.mod
+25 -15
View File
@@ -29,18 +29,18 @@
# Level 1 BLAS
#---------------------------------------------------------
set(SBLAS1 isamax.f sasum.f saxpy.f scopy.f sdot.f snrm2.f90
set(SBLAS1 isamax.f sasum.f saxpy.f saxpby.f scopy.f sdot.f snrm2.f90
srot.f srotg.f90 sscal.f sswap.f sdsdot.f srotmg.f srotm.f)
set(CBLAS1 scabs1.f scasum.f scnrm2.f90 icamax.f caxpy.f ccopy.f
set(CBLAS1 scabs1.f scasum.f scnrm2.f90 icamax.f90 caxpy.f caxpby.f ccopy.f
cdotc.f cdotu.f csscal.f crotg.f90 cscal.f cswap.f csrot.f)
set(DBLAS1 idamax.f dasum.f daxpy.f dcopy.f ddot.f dnrm2.f90
set(DBLAS1 idamax.f dasum.f daxpy.f daxpby.f dcopy.f ddot.f dnrm2.f90
drot.f drotg.f90 dscal.f dsdot.f dswap.f drotmg.f drotm.f)
set(DB1AUX sscal.f isamax.f)
set(ZBLAS1 dcabs1.f dzasum.f dznrm2.f90 izamax.f zaxpy.f zcopy.f
set(ZBLAS1 dcabs1.f dzasum.f dznrm2.f90 izamax.f90 zaxpy.f zaxpby.f zcopy.f
zdotc.f zdotu.f zdscal.f zrotg.f90 zscal.f zswap.f zdrot.f)
set(CB1AUX
@@ -49,7 +49,7 @@ set(CB1AUX
sswap.f)
set(ZB1AUX
icamax.f idamax.f
icamax.f90 idamax.f
cgemm.f cherk.f cscal.f ctrsm.f
dasum.f daxpy.f dcopy.f ddot.f dgemm.f dgemv.f dnrm2.f90 drot.f dscal.f
dswap.f
@@ -65,7 +65,8 @@ set(ALLBLAS lsame.f xerbla.f xerbla_array.f)
#---------------------------------------------------------
set(SBLAS2 sgemv.f sgbmv.f ssymv.f ssbmv.f sspmv.f
strmv.f stbmv.f stpmv.f strsv.f stbsv.f stpsv.f
sger.f ssyr.f sspr.f ssyr2.f sspr2.f)
sger.f ssyr.f sspr.f ssyr2.f sspr2.f
sskewsymv.f sskewsyr2.f)
set(CBLAS2 cgemv.f cgbmv.f chemv.f chbmv.f chpmv.f
ctrmv.f ctbmv.f ctpmv.f ctrsv.f ctbsv.f ctpsv.f
@@ -73,7 +74,8 @@ set(CBLAS2 cgemv.f cgbmv.f chemv.f chbmv.f chpmv.f
set(DBLAS2 dgemv.f dgbmv.f dsymv.f dsbmv.f dspmv.f
dtrmv.f dtbmv.f dtpmv.f dtrsv.f dtbsv.f dtpsv.f
dger.f dsyr.f dspr.f dsyr2.f dspr2.f)
dger.f dsyr.f dspr.f dsyr2.f dspr2.f
dskewsymv.f dskewsyr2.f)
set(ZBLAS2 zgemv.f zgbmv.f zhemv.f zhbmv.f zhpmv.f
ztrmv.f ztbmv.f ztpmv.f ztrsv.f ztbsv.f ztpsv.f
@@ -82,15 +84,17 @@ set(ZBLAS2 zgemv.f zgbmv.f zhemv.f zhbmv.f zhpmv.f
#---------------------------------------------------------
# Level 3 BLAS
#---------------------------------------------------------
set(SBLAS3 sgemm.f ssymm.f ssyrk.f ssyr2k.f strmm.f strsm.f)
set(SBLAS3 sgemm.f ssymm.f ssyrk.f ssyr2k.f strmm.f strsm.f sgemmtr.f
sskewsymm.f sskewsyr2k.f)
set(CBLAS3 cgemm.f csymm.f csyrk.f csyr2k.f ctrmm.f ctrsm.f
chemm.f cherk.f cher2k.f)
chemm.f cherk.f cher2k.f cgemmtr.f)
set(DBLAS3 dgemm.f dsymm.f dsyrk.f dsyr2k.f dtrmm.f dtrsm.f)
set(DBLAS3 dgemm.f dsymm.f dsyrk.f dsyr2k.f dtrmm.f dtrsm.f dgemmtr.f
dskewsymm.f dskewsyr2k.f)
set(ZBLAS3 zgemm.f zsymm.f zsyrk.f zsyr2k.f ztrmm.f ztrsm.f
zhemm.f zherk.f zher2k.f)
zhemm.f zherk.f zher2k.f zgemmtr.f)
set(SOURCES)
@@ -124,9 +128,9 @@ if(BUILD_INDEX64_EXT_API)
#Add _64 suffix to all Fortran functions via macros
foreach(F IN LISTS SOURCES_64_F)
if(CMAKE_Fortran_COMPILER_ID STREQUAL "NAG")
set(COPT_64_F -fpp)
set_source_files_properties(${F} PROPERTIES COMPILE_FLAGS "-fpp")
else()
set(COPT_64_F -cpp)
set_source_files_properties(${F} PROPERTIES COMPILE_FLAGS "-cpp")
endif()
file(STRINGS ${F} ${F}.lst)
list(FILTER ${F}.lst INCLUDE REGEX "subroutine|SUBROUTINE|external|EXTERNAL|function|FUNCTION")
@@ -137,10 +141,10 @@ if(BUILD_INDEX64_EXT_API)
string(REGEX REPLACE "^[a-zA-Z0-9_ *]*(subroutine|SUBROUTINE|external|EXTERNAL|function|FUNCTION)[ ]*[*]?" "" FUNC ${FUNC})
string(REGEX REPLACE "[(][a-zA-Z0-9_, )]*$" "" FUNC ${FUNC})
string(STRIP ${FUNC} FUNC)
list(APPEND COPT_64_F "-D${FUNC}=${FUNC}_64")
list(APPEND COPT_64_F "${FUNC}=${FUNC}_64")
endforeach()
list(REMOVE_DUPLICATES COPT_64_F)
set_source_files_properties(${F} PROPERTIES COMPILE_OPTIONS "${COPT_64_F}")
set_source_files_properties(${F} PROPERTIES COMPILE_DEFINITIONS "${COPT_64_F}")
endforeach()
endif()
@@ -148,6 +152,12 @@ add_library(${BLASLIB}
$<TARGET_OBJECTS:${BLASLIB}_obj>
$<$<BOOL:${BUILD_INDEX64_EXT_API}>: $<TARGET_OBJECTS:${BLASLIB}_64_obj>>)
# For flang, use C linker instead of Fortran linker to avoid macOS-specific flags
# that CMake adds (tested CMake 4.2).
if (CMAKE_Fortran_COMPILER_ID STREQUAL "LLVMFlang")
set_target_properties (${BLASLIB} PROPERTIES LINKER_LANGUAGE C)
endif()
set_target_properties(
${BLASLIB} PROPERTIES
VERSION ${LAPACK_VERSION}
+14 -10
View File
@@ -69,19 +69,19 @@ all: $(BLASLIB)
# Comment out the next 6 definitions if you already have
# the Level 1 BLAS.
#---------------------------------------------------------
SBLAS1 = isamax.o sasum.o saxpy.o scopy.o sdot.o snrm2.o \
SBLAS1 = isamax.o sasum.o saxpy.o saxpby.o scopy.o sdot.o snrm2.o \
srot.o srotg.o sscal.o sswap.o sdsdot.o srotmg.o srotm.o
$(SBLAS1): $(FRC)
CBLAS1 = scabs1.o scasum.o scnrm2.o icamax.o caxpy.o ccopy.o \
CBLAS1 = scabs1.o scasum.o scnrm2.o icamax.o caxpy.o caxpby.o ccopy.o \
cdotc.o cdotu.o csscal.o crotg.o cscal.o cswap.o csrot.o
$(CBLAS1): $(FRC)
DBLAS1 = idamax.o dasum.o daxpy.o dcopy.o ddot.o dnrm2.o \
DBLAS1 = idamax.o dasum.o daxpy.o daxpby.o dcopy.o ddot.o dnrm2.o \
drot.o drotg.o dscal.o dsdot.o dswap.o drotmg.o drotm.o
$(DBLAS1): $(FRC)
ZBLAS1 = dcabs1.o dzasum.o dznrm2.o izamax.o zaxpy.o zcopy.o \
ZBLAS1 = dcabs1.o dzasum.o dznrm2.o izamax.o zaxpy.o zaxpby.o zcopy.o \
zdotc.o zdotu.o zdscal.o zrotg.o zscal.o zswap.o zdrot.o
$(ZBLAS1): $(FRC)
@@ -105,7 +105,8 @@ $(ALLBLAS): $(FRC)
#---------------------------------------------------------
SBLAS2 = sgemv.o sgbmv.o ssymv.o ssbmv.o sspmv.o \
strmv.o stbmv.o stpmv.o strsv.o stbsv.o stpsv.o \
sger.o ssyr.o sspr.o ssyr2.o sspr2.o
sger.o ssyr.o sspr.o ssyr2.o sspr2.o \
sskewsymv.o sskewsyr2.o
$(SBLAS2): $(FRC)
CBLAS2 = cgemv.o cgbmv.o chemv.o chbmv.o chpmv.o \
@@ -115,7 +116,8 @@ $(CBLAS2): $(FRC)
DBLAS2 = dgemv.o dgbmv.o dsymv.o dsbmv.o dspmv.o \
dtrmv.o dtbmv.o dtpmv.o dtrsv.o dtbsv.o dtpsv.o \
dger.o dsyr.o dspr.o dsyr2.o dspr2.o
dger.o dsyr.o dspr.o dsyr2.o dspr2.o \
dskewsymv.o dskewsyr2.o
$(DBLAS2): $(FRC)
ZBLAS2 = zgemv.o zgbmv.o zhemv.o zhbmv.o zhpmv.o \
@@ -127,18 +129,20 @@ $(ZBLAS2): $(FRC)
# Comment out the next 4 definitions if you already have
# the Level 3 BLAS.
#---------------------------------------------------------
SBLAS3 = sgemm.o ssymm.o ssyrk.o ssyr2k.o strmm.o strsm.o
SBLAS3 = sgemm.o ssymm.o ssyrk.o ssyr2k.o strmm.o strsm.o sgemmtr.o \
sskewsymm.o sskewsyr2k.o
$(SBLAS3): $(FRC)
CBLAS3 = cgemm.o csymm.o csyrk.o csyr2k.o ctrmm.o ctrsm.o \
chemm.o cherk.o cher2k.o
chemm.o cherk.o cher2k.o cgemmtr.o
$(CBLAS3): $(FRC)
DBLAS3 = dgemm.o dsymm.o dsyrk.o dsyr2k.o dtrmm.o dtrsm.o
DBLAS3 = dgemm.o dsymm.o dsyrk.o dsyr2k.o dtrmm.o dtrsm.o dgemmtr.o \
dskewsymm.o dskewsyr2k.o
$(DBLAS3): $(FRC)
ZBLAS3 = zgemm.o zsymm.o zsyrk.o zsyr2k.o ztrmm.o ztrsm.o \
zhemm.o zherk.o zher2k.o
zhemm.o zherk.o zher2k.o zgemmtr.o
$(ZBLAS3): $(FRC)
ALLOBJ = $(SBLAS1) $(SBLAS2) $(SBLAS3) $(DBLAS1) $(DBLAS2) $(DBLAS3) \
+144
View File
@@ -0,0 +1,144 @@
*> \brief \b CAXPBY
*
* =========== DOCUMENTATION ===========
*
* Online html documentation available at
* http://www.netlib.org/lapack/explore-html/
*
* Definition:
* ===========
*
* SUBROUTINE CAXPBY(N,CA,CX,INCX,CB,CY,INCY)
*
* .. Scalar Arguments ..
* COMPLEX CA,CB
* INTEGER INCX,INCY,N
* ..
* .. Array Arguments ..
* COMPLEX CX(*),CY(*)
* ..
*
*
*> \par Purpose:
* =============
*>
*> \verbatim
*>
*> CAXPBY constant times a vector plus constant times a vector.
*>
*> Y = ALPHA * X + BETA * Y
*>
*> \endverbatim
*
* Arguments:
* ==========
*
*> \param[in] N
*> \verbatim
*> N is INTEGER
*> number of elements in input vector(s)
*> \endverbatim
*>
*> \param[in] CA
*> \verbatim
*> CA is COMPLEX
*> On entry, CA specifies the scalar alpha.
*> \endverbatim
*>
*> \param[in] CX
*> \verbatim
*> CX is COMPLEX array, dimension ( 1 + ( N - 1 )*abs( INCX ) )
*> \endverbatim
*>
*> \param[in] INCX
*> \verbatim
*> INCX is INTEGER
*> storage spacing between elements of CX
*> \endverbatim
*>
*> \param[in] CB
*> \verbatim
*> CB is COMPLEX
*> On entry, CB specifies the scalar beta.
*> \endverbatim
*>
*> \param[in,out] CY
*> \verbatim
*> CY is COMPLEX array, dimension ( 1 + ( N - 1 )*abs( INCY ) )
*> \endverbatim
*>
*> \param[in] INCY
*> \verbatim
*> INCY is INTEGER
*> storage spacing between elements of CY
*> \endverbatim
*
* Authors:
* ========
*
*> \author Univ. of Tennessee
*> \author Univ. of California Berkeley
*> \author Univ. of Colorado Denver
*> \author NAG Ltd.
*> \author Martin Koehler, MPI Magdeburg
*
*> \ingroup axpby
*
* =====================================================================
SUBROUTINE CAXPBY(N,CA,CX,INCX,CB,CY,INCY)
IMPLICIT NONE
*
* -- Reference BLAS level1 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
*
* .. Scalar Arguments ..
COMPLEX CA, CB
INTEGER INCX,INCY,N
* ..
* .. Array Arguments ..
COMPLEX CX(*),CY(*)
* ..
* .. External Subroutines ..
EXTERNAL CSCAL
*
* =====================================================================
*
* .. Local Scalars ..
INTEGER I,IX,IY
* ..
IF (N.LE.0) RETURN
IF (CA .EQ. (0.0,0.0) .AND. CB.NE.(0.0,0.0)) THEN
CALL CSCAL(N,CB, CY, INCY)
RETURN
END IF
IF (INCX.EQ.1 .AND. INCY.EQ.1) THEN
*
* code for both increments equal to 1
*
DO I = 1,N
CY(I) = CB*CY(I) + CA*CX(I)
END DO
ELSE
*
* code for unequal increments or equal increments
* not equal to 1
*
IX = 1
IY = 1
IF (INCX.LT.0) IX = (-N+1)*INCX + 1
IF (INCY.LT.0) IY = (-N+1)*INCY + 1
DO I = 1,N
CY(IY) = CB*CY(IY) + CA*CX(IX)
IX = IX + INCX
IY = IY + INCY
END DO
END IF
*
RETURN
*
* End of CAXBPY
*
END
+1
View File
@@ -85,6 +85,7 @@
*>
* =====================================================================
SUBROUTINE CAXPY(N,CA,CX,INCX,CY,INCY)
IMPLICIT NONE
*
* -- Reference BLAS level1 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -78,6 +78,7 @@
*>
* =====================================================================
SUBROUTINE CCOPY(N,CX,INCX,CY,INCY)
IMPLICIT NONE
*
* -- Reference BLAS level1 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -80,6 +80,7 @@
*>
* =====================================================================
COMPLEX FUNCTION CDOTC(N,CX,INCX,CY,INCY)
IMPLICIT NONE
*
* -- Reference BLAS level1 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -80,6 +80,7 @@
*>
* =====================================================================
COMPLEX FUNCTION CDOTU(N,CX,INCX,CY,INCY)
IMPLICIT NONE
*
* -- Reference BLAS level1 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -187,6 +187,7 @@
* =====================================================================
SUBROUTINE CGBMV(TRANS,M,N,KL,KU,ALPHA,A,LDA,X,INCX,
+ BETA,Y,INCY)
IMPLICIT NONE
*
* -- Reference BLAS level2 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+29 -7
View File
@@ -35,6 +35,16 @@
*>
*> alpha and beta are scalars, and A, B and C are matrices, with op( A )
*> an m by k matrix, op( B ) a k by n matrix and C an m by n matrix.
*>
*> Note: if alpha and/or beta is zero, some parts of the matrix-matrix
*> operations are not performed. This results in the following NaN/Inf
*> propagation quirks:
*>
*> 1. If alpha is zero, NaNs or Infs in A or B do not affect the result.
*> 2. If both alpha and beta are zero, then a zero matrix is returned in C,
*> irrespective of any NaNs or Infs in A, B or C.
*> 3. If only beta is zero, alpha*op( A )*op( B ) is returned, irrespective
*> of any NaNs or Infs in C.
*> \endverbatim
*
* Arguments:
@@ -92,7 +102,9 @@
*> \param[in] ALPHA
*> \verbatim
*> ALPHA is COMPLEX
*> On entry, ALPHA specifies the scalar alpha.
*> On entry, ALPHA specifies the scalar alpha. If ALPHA is zero the
*> values in A and B do not affect the result. This also means that
*> NaN/Inf propagation from A and B is inhibited if ALPHA is zero.
*> \endverbatim
*>
*> \param[in] A
@@ -102,7 +114,10 @@
*> Before entry with TRANSA = 'N' or 'n', the leading m by k
*> part of the array A must contain the matrix A, otherwise
*> the leading k by m part of the array A must contain the
*> matrix A.
*> matrix A, except if ALPHA is zero.
*> If ALPHA is zero, none of the values in A affect the result, even
*> if they are NaN/Inf. This also implies that if ALPHA is zero,
*> the matrix elements of A need not be initialized by the caller.
*> \endverbatim
*>
*> \param[in] LDA
@@ -121,7 +136,10 @@
*> Before entry with TRANSB = 'N' or 'n', the leading k by n
*> part of the array B must contain the matrix B, otherwise
*> the leading n by k part of the array B must contain the
*> matrix B.
*> matrix B, except if ALPHA is zero.
*> If ALPHA is zero, none of the values in B affect the result, even
*> if they are NaN/Inf. This also implies that if ALPHA is zero,
*> the matrix elements of B need not be initialized by the caller.
*> \endverbatim
*>
*> \param[in] LDB
@@ -136,16 +154,19 @@
*> \param[in] BETA
*> \verbatim
*> BETA is COMPLEX
*> On entry, BETA specifies the scalar beta. When BETA is
*> supplied as zero then C need not be set on input.
*> On entry, BETA specifies the scalar beta. If BETA is zero the
*> values in C do not affect the result. This also means that
*> NaN/Inf propagation from C is inhibited if BETA is zero.
*> \endverbatim
*>
*> \param[in,out] C
*> \verbatim
*> C is COMPLEX array, dimension ( LDC, N )
*> Before entry, the leading m by n part of the array C must
*> contain the matrix C, except when beta is zero, in which
*> case C need not be set on entry.
*> contain the matrix C, except if beta is zero.
*> If beta is zero, none of the values in C affect the result, even
*> if they are NaN/Inf. This also implies that if beta is zero,
*> the matrix elements of C need not be initialized by the caller.
*> On exit, the array C is overwritten by the m by n matrix
*> ( alpha*op( A )*op( B ) + beta*C ).
*> \endverbatim
@@ -185,6 +206,7 @@
* =====================================================================
SUBROUTINE CGEMM(TRANSA,TRANSB,M,N,K,ALPHA,A,LDA,B,LDB,
+ BETA,C,LDC)
IMPLICIT NONE
*
* -- Reference BLAS level3 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+569
View File
@@ -0,0 +1,569 @@
*> \brief \b CGEMMTR
*
* =========== DOCUMENTATION ===========
*
* Online html documentation available at
* http://www.netlib.org/lapack/explore-html/
*
* Definition:
* ===========
*
* SUBROUTINE CGEMMTR(UPLO,TRANSA,TRANSB,N,K,ALPHA,A,LDA,B,LDB,BETA,
* C,LDC)
*
* .. Scalar Arguments ..
* COMPLEX ALPHA,BETA
* INTEGER K,LDA,LDB,LDC,N
* CHARACTER TRANSA,TRANSB, UPLO
* ..
* .. Array Arguments ..
* COMPLEX A(LDA,*),B(LDB,*),C(LDC,*)
* ..
*
*
*> \par Purpose:
* =============
*>
*> \verbatim
*>
*> CGEMMTR performs one of the matrix-matrix operations
*>
*> C := alpha*op( A )*op( B ) + beta*C,
*>
*> where op( X ) is one of
*>
*> op( X ) = X or op( X ) = X**T,
*>
*> alpha and beta are scalars, and A, B and C are matrices, with op( A )
*> an n by k matrix, op( B ) a k by n matrix and C an n by n matrix.
*> Thereby, the routine only accesses and updates the upper or lower
*> triangular part of the result matrix C. This behaviour can be used if
*> the resulting matrix C is known to be Hermitian or symmetric.
*> \endverbatim
*
* Arguments:
* ==========
*
*> \param[in] UPLO
*> \verbatim
*> UPLO is CHARACTER*1
*> On entry, UPLO specifies whether the lower or the upper
*> triangular part of C is access and updated.
*>
*> UPLO = 'L' or 'l', the lower triangular part of C is used.
*>
*> UPLO = 'U' or 'u', the upper triangular part of C is used.
*> \endverbatim
*
*> \param[in] TRANSA
*> \verbatim
*> TRANSA is CHARACTER*1
*> On entry, TRANSA specifies the form of op( A ) to be used in
*> the matrix multiplication as follows:
*>
*> TRANSA = 'N' or 'n', op( A ) = A.
*>
*> TRANSA = 'T' or 't', op( A ) = A**T.
*>
*> TRANSA = 'C' or 'c', op( A ) = A**H.
*> \endverbatim
*>
*> \param[in] TRANSB
*> \verbatim
*> TRANSB is CHARACTER*1
*> On entry, TRANSB specifies the form of op( B ) to be used in
*> the matrix multiplication as follows:
*>
*> TRANSB = 'N' or 'n', op( B ) = B.
*>
*> TRANSB = 'T' or 't', op( B ) = B**T.
*>
*> TRANSB = 'C' or 'c', op( B ) = B**H.
*> \endverbatim
*>
*> \param[in] N
*> \verbatim
*> N is INTEGER
*> On entry, N specifies the number of rows and columns of
*> the matrix C, the number of columns of op(B) and the number
*> of rows of op(A). N must be at least zero.
*> \endverbatim
*>
*> \param[in] K
*> \verbatim
*> K is INTEGER
*> On entry, K specifies the number of columns of the matrix
*> op( A ) and the number of rows of the matrix op( B ). K must
*> be at least zero.
*> \endverbatim
*>
*> \param[in] ALPHA
*> \verbatim
*> ALPHA is COMPLEX.
*> On entry, ALPHA specifies the scalar alpha.
*> \endverbatim
*>
*> \param[in] A
*> \verbatim
*> A is COMPLEX array, dimension ( LDA, ka ), where ka is
*> k when TRANSA = 'N' or 'n', and is n otherwise.
*> Before entry with TRANSA = 'N' or 'n', the leading n by k
*> part of the array A must contain the matrix A, otherwise
*> the leading k by m part of the array A must contain the
*> matrix A.
*> \endverbatim
*>
*> \param[in] LDA
*> \verbatim
*> LDA is INTEGER
*> On entry, LDA specifies the first dimension of A as declared
*> in the calling (sub) program. When TRANSA = 'N' or 'n' then
*> LDA must be at least max( 1, n ), otherwise LDA must be at
*> least max( 1, k ).
*> \endverbatim
*>
*> \param[in] B
*> \verbatim
*> B is COMPLEX array, dimension ( LDB, kb ), where kb is
*> n when TRANSB = 'N' or 'n', and is k otherwise.
*> Before entry with TRANSB = 'N' or 'n', the leading k by n
*> part of the array B must contain the matrix B, otherwise
*> the leading n by k part of the array B must contain the
*> matrix B.
*> \endverbatim
*>
*> \param[in] LDB
*> \verbatim
*> LDB is INTEGER
*> On entry, LDB specifies the first dimension of B as declared
*> in the calling (sub) program. When TRANSB = 'N' or 'n' then
*> LDB must be at least max( 1, k ), otherwise LDB must be at
*> least max( 1, n ).
*> \endverbatim
*>
*> \param[in] BETA
*> \verbatim
*> BETA is COMPLEX.
*> On entry, BETA specifies the scalar beta. When BETA is
*> supplied as zero then C need not be set on input.
*> \endverbatim
*>
*> \param[in,out] C
*> \verbatim
*> C is COMPLEX array, dimension ( LDC, N )
*> Before entry, the leading n by n part of the array C must
*> contain the matrix C, except when beta is zero, in which
*> case C need not be set on entry.
*> On exit, the upper or lower triangular part of the matrix
*> C is overwritten by the n by n matrix
*> ( alpha*op( A )*op( B ) + beta*C ).
*> \endverbatim
*>
*> \param[in] LDC
*> \verbatim
*> LDC is INTEGER
*> On entry, LDC specifies the first dimension of C as declared
*> in the calling (sub) program. LDC must be at least
*> max( 1, n ).
*> \endverbatim
*
* Authors:
* ========
*
*> \author Martin Koehler
*
*> \ingroup gemmtr
*
*> \par Further Details:
* =====================
*>
*> \verbatim
*>
*> Level 3 Blas routine.
*>
*> -- Written on 19-July-2023.
*> Martin Koehler, MPI Magdeburg
*> \endverbatim
*>
* =====================================================================
SUBROUTINE CGEMMTR(UPLO,TRANSA,TRANSB,N,K,ALPHA,A,LDA,B,LDB,
+ BETA,C,LDC)
IMPLICIT NONE
*
* -- Reference BLAS level3 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
*
* .. Scalar Arguments ..
COMPLEX ALPHA,BETA
INTEGER K,LDA,LDB,LDC,N
CHARACTER TRANSA,TRANSB,UPLO
* ..
* .. Array Arguments ..
COMPLEX A(LDA,*),B(LDB,*),C(LDC,*)
* ..
*
* =====================================================================
*
* .. External Functions ..
LOGICAL LSAME
EXTERNAL LSAME
* ..
* .. External Subroutines ..
EXTERNAL XERBLA
* ..
* .. Intrinsic Functions ..
INTRINSIC CONJG,MAX
* ..
* .. Local Scalars ..
COMPLEX TEMP
INTEGER I,INFO,J,L,NROWA,NROWB,ISTART, ISTOP
LOGICAL CONJA,CONJB,NOTA,NOTB,UPPER
* ..
* .. Parameters ..
COMPLEX ONE
PARAMETER (ONE= (1.0E+0,0.0E+0))
COMPLEX ZERO
PARAMETER (ZERO= (0.0E+0,0.0E+0))
* ..
*
* Set NOTA and NOTB as true if A and B respectively are not
* conjugated or transposed, set CONJA and CONJB as true if A and
* B respectively are to be transposed but not conjugated and set
* NROWA and NROWB as the number of rows of A and B respectively.
*
NOTA = LSAME(TRANSA,'N')
NOTB = LSAME(TRANSB,'N')
CONJA = LSAME(TRANSA,'C')
CONJB = LSAME(TRANSB,'C')
IF (NOTA) THEN
NROWA = N
ELSE
NROWA = K
END IF
IF (NOTB) THEN
NROWB = K
ELSE
NROWB = N
END IF
UPPER = LSAME(UPLO, 'U')
*
* Test the input parameters.
*
INFO = 0
IF ((.NOT. UPPER) .AND. (.NOT. LSAME(UPLO, 'L'))) THEN
INFO = 1
ELSE IF ((.NOT.NOTA) .AND. (.NOT.CONJA) .AND.
+ (.NOT.LSAME(TRANSA,'T'))) THEN
INFO = 2
ELSE IF ((.NOT.NOTB) .AND. (.NOT.CONJB) .AND.
+ (.NOT.LSAME(TRANSB,'T'))) THEN
INFO = 3
ELSE IF (N.LT.0) THEN
INFO = 4
ELSE IF (K.LT.0) THEN
INFO = 5
ELSE IF (LDA.LT.MAX(1,NROWA)) THEN
INFO = 8
ELSE IF (LDB.LT.MAX(1,NROWB)) THEN
INFO = 10
ELSE IF (LDC.LT.MAX(1,N)) THEN
INFO = 13
END IF
IF (INFO.NE.0) THEN
CALL XERBLA('CGEMMTR',INFO)
RETURN
END IF
*
* Quick return if possible.
*
IF (N.EQ.0) RETURN
*
* And when alpha.eq.zero.
*
IF (ALPHA.EQ.ZERO) THEN
IF (BETA.EQ.ZERO) THEN
DO 20 J = 1,N
IF (UPPER) THEN
ISTART = 1
ISTOP = J
ELSE
ISTART = J
ISTOP = N
END IF
DO 10 I = ISTART, ISTOP
C(I,J) = ZERO
10 CONTINUE
20 CONTINUE
ELSE
DO 40 J = 1,N
IF (UPPER) THEN
ISTART = 1
ISTOP = J
ELSE
ISTART = J
ISTOP = N
END IF
DO 30 I = ISTART, ISTOP
C(I,J) = BETA*C(I,J)
30 CONTINUE
40 CONTINUE
END IF
RETURN
END IF
*
* Start the operations.
*
IF (NOTB) THEN
IF (NOTA) THEN
*
* Form C := alpha*A*B + beta*C.
*
DO 90 J = 1,N
IF (UPPER) THEN
ISTART = 1
ISTOP = J
ELSE
ISTART = J
ISTOP = N
END IF
IF (BETA.EQ.ZERO) THEN
DO 50 I = ISTART, ISTOP
C(I,J) = ZERO
50 CONTINUE
ELSE IF (BETA.NE.ONE) THEN
DO 60 I = ISTART, ISTOP
C(I,J) = BETA*C(I,J)
60 CONTINUE
END IF
DO 80 L = 1,K
TEMP = ALPHA*B(L,J)
DO 70 I = ISTART, ISTOP
C(I,J) = C(I,J) + TEMP*A(I,L)
70 CONTINUE
80 CONTINUE
90 CONTINUE
ELSE IF (CONJA) THEN
*
* Form C := alpha*A**H*B + beta*C.
*
DO 120 J = 1,N
IF (UPPER) THEN
ISTART = 1
ISTOP = J
ELSE
ISTART = J
ISTOP = N
END IF
DO 110 I = ISTART, ISTOP
TEMP = ZERO
DO 100 L = 1,K
TEMP = TEMP + CONJG(A(L,I))*B(L,J)
100 CONTINUE
IF (BETA.EQ.ZERO) THEN
C(I,J) = ALPHA*TEMP
ELSE
C(I,J) = ALPHA*TEMP + BETA*C(I,J)
END IF
110 CONTINUE
120 CONTINUE
ELSE
*
* Form C := alpha*A**T*B + beta*C
*
DO 150 J = 1,N
IF (UPPER) THEN
ISTART = 1
ISTOP = J
ELSE
ISTART = J
ISTOP = N
END IF
DO 140 I = ISTART, ISTOP
TEMP = ZERO
DO 130 L = 1,K
TEMP = TEMP + A(L,I)*B(L,J)
130 CONTINUE
IF (BETA.EQ.ZERO) THEN
C(I,J) = ALPHA*TEMP
ELSE
C(I,J) = ALPHA*TEMP + BETA*C(I,J)
END IF
140 CONTINUE
150 CONTINUE
END IF
ELSE IF (NOTA) THEN
IF (CONJB) THEN
*
* Form C := alpha*A*B**H + beta*C.
*
DO 200 J = 1,N
IF (UPPER) THEN
ISTART = 1
ISTOP = J
ELSE
ISTART = J
ISTOP = N
END IF
IF (BETA.EQ.ZERO) THEN
DO 160 I = ISTART,ISTOP
C(I,J) = ZERO
160 CONTINUE
ELSE IF (BETA.NE.ONE) THEN
DO 170 I = ISTART, ISTOP
C(I,J) = BETA*C(I,J)
170 CONTINUE
END IF
DO 190 L = 1,K
TEMP = ALPHA*CONJG(B(J,L))
DO 180 I = ISTART, ISTOP
C(I,J) = C(I,J) + TEMP*A(I,L)
180 CONTINUE
190 CONTINUE
200 CONTINUE
ELSE
*
* Form C := alpha*A*B**T + beta*C
*
DO 250 J = 1,N
IF (UPPER) THEN
ISTART = 1
ISTOP = J
ELSE
ISTART = J
ISTOP = N
END IF
IF (BETA.EQ.ZERO) THEN
DO 210 I = ISTART, ISTOP
C(I,J) = ZERO
210 CONTINUE
ELSE IF (BETA.NE.ONE) THEN
DO 220 I = ISTART, ISTOP
C(I,J) = BETA*C(I,J)
220 CONTINUE
END IF
DO 240 L = 1,K
TEMP = ALPHA*B(J,L)
DO 230 I = ISTART, ISTOP
C(I,J) = C(I,J) + TEMP*A(I,L)
230 CONTINUE
240 CONTINUE
250 CONTINUE
END IF
ELSE IF (CONJA) THEN
IF (CONJB) THEN
*
* Form C := alpha*A**H*B**H + beta*C.
*
DO 280 J = 1,N
IF (UPPER) THEN
ISTART = 1
ISTOP = J
ELSE
ISTART = J
ISTOP = N
END IF
DO 270 I = ISTART, ISTOP
TEMP = ZERO
DO 260 L = 1,K
TEMP = TEMP + CONJG(A(L,I))*CONJG(B(J,L))
260 CONTINUE
IF (BETA.EQ.ZERO) THEN
C(I,J) = ALPHA*TEMP
ELSE
C(I,J) = ALPHA*TEMP + BETA*C(I,J)
END IF
270 CONTINUE
280 CONTINUE
ELSE
*
* Form C := alpha*A**H*B**T + beta*C
*
DO 310 J = 1,N
IF (UPPER) THEN
ISTART = 1
ISTOP = J
ELSE
ISTART = J
ISTOP = N
END IF
DO 300 I = ISTART, ISTOP
TEMP = ZERO
DO 290 L = 1,K
TEMP = TEMP + CONJG(A(L,I))*B(J,L)
290 CONTINUE
IF (BETA.EQ.ZERO) THEN
C(I,J) = ALPHA*TEMP
ELSE
C(I,J) = ALPHA*TEMP + BETA*C(I,J)
END IF
300 CONTINUE
310 CONTINUE
END IF
ELSE
IF (CONJB) THEN
*
* Form C := alpha*A**T*B**H + beta*C
*
DO 340 J = 1,N
IF (UPPER) THEN
ISTART = 1
ISTOP = J
ELSE
ISTART = J
ISTOP = N
END IF
DO 330 I = ISTART, ISTOP
TEMP = ZERO
DO 320 L = 1,K
TEMP = TEMP + A(L,I)*CONJG(B(J,L))
320 CONTINUE
IF (BETA.EQ.ZERO) THEN
C(I,J) = ALPHA*TEMP
ELSE
C(I,J) = ALPHA*TEMP + BETA*C(I,J)
END IF
330 CONTINUE
340 CONTINUE
ELSE
*
* Form C := alpha*A**T*B**T + beta*C
*
DO 370 J = 1,N
IF (UPPER) THEN
ISTART = 1
ISTOP = J
ELSE
ISTART = J
ISTOP = N
END IF
DO 360 I = ISTART, ISTOP
TEMP = ZERO
DO 350 L = 1,K
TEMP = TEMP + A(L,I)*B(J,L)
350 CONTINUE
IF (BETA.EQ.ZERO) THEN
C(I,J) = ALPHA*TEMP
ELSE
C(I,J) = ALPHA*TEMP + BETA*C(I,J)
END IF
360 CONTINUE
370 CONTINUE
END IF
END IF
*
RETURN
*
* End of CGEMMTR
*
END
+1
View File
@@ -157,6 +157,7 @@
*>
* =====================================================================
SUBROUTINE CGEMV(TRANS,M,N,ALPHA,A,LDA,X,INCX,BETA,Y,INCY)
IMPLICIT NONE
*
* -- Reference BLAS level2 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -127,6 +127,7 @@
*>
* =====================================================================
SUBROUTINE CGERC(M,N,ALPHA,X,INCX,Y,INCY,A,LDA)
IMPLICIT NONE
*
* -- Reference BLAS level2 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -127,6 +127,7 @@
*>
* =====================================================================
SUBROUTINE CGERU(M,N,ALPHA,X,INCX,Y,INCY,A,LDA)
IMPLICIT NONE
*
* -- Reference BLAS level2 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -184,6 +184,7 @@
*>
* =====================================================================
SUBROUTINE CHBMV(UPLO,N,K,ALPHA,A,LDA,X,INCX,BETA,Y,INCY)
IMPLICIT NONE
*
* -- Reference BLAS level2 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -188,6 +188,7 @@
*>
* =====================================================================
SUBROUTINE CHEMM(SIDE,UPLO,M,N,ALPHA,A,LDA,B,LDB,BETA,C,LDC)
IMPLICIT NONE
*
* -- Reference BLAS level3 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -151,6 +151,7 @@
*>
* =====================================================================
SUBROUTINE CHEMV(UPLO,N,ALPHA,A,LDA,X,INCX,BETA,Y,INCY)
IMPLICIT NONE
*
* -- Reference BLAS level2 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -132,6 +132,7 @@
*>
* =====================================================================
SUBROUTINE CHER(UPLO,N,ALPHA,X,INCX,A,LDA)
IMPLICIT NONE
*
* -- Reference BLAS level2 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -147,6 +147,7 @@
*>
* =====================================================================
SUBROUTINE CHER2(UPLO,N,ALPHA,X,INCX,Y,INCY,A,LDA)
IMPLICIT NONE
*
* -- Reference BLAS level2 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -194,6 +194,7 @@
*>
* =====================================================================
SUBROUTINE CHER2K(UPLO,TRANS,N,K,ALPHA,A,LDA,B,LDB,BETA,C,LDC)
IMPLICIT NONE
*
* -- Reference BLAS level3 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -170,6 +170,7 @@
*>
* =====================================================================
SUBROUTINE CHERK(UPLO,TRANS,N,K,ALPHA,A,LDA,BETA,C,LDC)
IMPLICIT NONE
*
* -- Reference BLAS level3 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -146,6 +146,7 @@
*>
* =====================================================================
SUBROUTINE CHPMV(UPLO,N,ALPHA,AP,X,INCX,BETA,Y,INCY)
IMPLICIT NONE
*
* -- Reference BLAS level2 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -127,6 +127,7 @@
*>
* =====================================================================
SUBROUTINE CHPR(UPLO,N,ALPHA,X,INCX,AP)
IMPLICIT NONE
*
* -- Reference BLAS level2 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -142,6 +142,7 @@
*>
* =====================================================================
SUBROUTINE CHPR2(UPLO,N,ALPHA,X,INCX,Y,INCY,AP)
IMPLICIT NONE
*
* -- Reference BLAS level2 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -86,6 +86,7 @@
!
! =====================================================================
subroutine CROTG( a, b, c, s )
implicit none
integer, parameter :: wp = kind(1.e0)
!
! -- Reference BLAS level1 routine --
+1
View File
@@ -75,6 +75,7 @@
*>
* =====================================================================
SUBROUTINE CSCAL(N,CA,CX,INCX)
IMPLICIT NONE
*
* -- Reference BLAS level1 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -95,6 +95,7 @@
*
* =====================================================================
SUBROUTINE CSROT( N, CX, INCX, CY, INCY, C, S )
IMPLICIT NONE
*
* -- Reference BLAS level1 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -75,6 +75,7 @@
*>
* =====================================================================
SUBROUTINE CSSCAL(N,SA,CX,INCX)
IMPLICIT NONE
*
* -- Reference BLAS level1 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -78,6 +78,7 @@
*>
* =====================================================================
SUBROUTINE CSWAP(N,CX,INCX,CY,INCY)
IMPLICIT NONE
*
* -- Reference BLAS level1 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -186,6 +186,7 @@
*>
* =====================================================================
SUBROUTINE CSYMM(SIDE,UPLO,M,N,ALPHA,A,LDA,B,LDB,BETA,C,LDC)
IMPLICIT NONE
*
* -- Reference BLAS level3 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -185,6 +185,7 @@
*>
* =====================================================================
SUBROUTINE CSYR2K(UPLO,TRANS,N,K,ALPHA,A,LDA,B,LDB,BETA,C,LDC)
IMPLICIT NONE
*
* -- Reference BLAS level3 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -164,6 +164,7 @@
*>
* =====================================================================
SUBROUTINE CSYRK(UPLO,TRANS,N,K,ALPHA,A,LDA,BETA,C,LDC)
IMPLICIT NONE
*
* -- Reference BLAS level3 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -183,6 +183,7 @@
*>
* =====================================================================
SUBROUTINE CTBMV(UPLO,TRANS,DIAG,N,K,A,LDA,X,INCX)
IMPLICIT NONE
*
* -- Reference BLAS level2 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -186,6 +186,7 @@
*>
* =====================================================================
SUBROUTINE CTBSV(UPLO,TRANS,DIAG,N,K,A,LDA,X,INCX)
IMPLICIT NONE
*
* -- Reference BLAS level2 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -139,6 +139,7 @@
*>
* =====================================================================
SUBROUTINE CTPMV(UPLO,TRANS,DIAG,N,AP,X,INCX)
IMPLICIT NONE
*
* -- Reference BLAS level2 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -141,6 +141,7 @@
*>
* =====================================================================
SUBROUTINE CTPSV(UPLO,TRANS,DIAG,N,AP,X,INCX)
IMPLICIT NONE
*
* -- Reference BLAS level2 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -174,6 +174,7 @@
*>
* =====================================================================
SUBROUTINE CTRMM(SIDE,UPLO,TRANSA,DIAG,M,N,ALPHA,A,LDA,B,LDB)
IMPLICIT NONE
*
* -- Reference BLAS level3 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -144,6 +144,7 @@
*>
* =====================================================================
SUBROUTINE CTRMV(UPLO,TRANS,DIAG,N,A,LDA,X,INCX)
IMPLICIT NONE
*
* -- Reference BLAS level2 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -177,6 +177,7 @@
*>
* =====================================================================
SUBROUTINE CTRSM(SIDE,UPLO,TRANSA,DIAG,M,N,ALPHA,A,LDA,B,LDB)
IMPLICIT NONE
*
* -- Reference BLAS level3 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -146,6 +146,7 @@
*>
* =====================================================================
SUBROUTINE CTRSV(UPLO,TRANS,DIAG,N,A,LDA,X,INCX)
IMPLICIT NONE
*
* -- Reference BLAS level2 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -68,6 +68,7 @@
*>
* =====================================================================
DOUBLE PRECISION FUNCTION DASUM(N,DX,INCX)
IMPLICIT NONE
*
* -- Reference BLAS level1 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+149
View File
@@ -0,0 +1,149 @@
*> \brief \b DAXPBY
*
* =========== DOCUMENTATION ===========
*
* Online html documentation available at
* http://www.netlib.org/lapack/explore-html/
*
* Definition:
* ===========
*
* SUBROUTINE DAXPBY(N,DA,DX,INCX,DB,DY,INCY)
*
* .. Scalar Arguments ..
* DOUBLE PRECISION DA,DB
* INTEGER INCX,INCY,N
* ..
* .. Array Arguments ..
* DOUBLE PRECISION DX(*),DY(*)
* ..
*
*
*> \par Purpose:
* =============
*>
*> \verbatim
*>
*> DAXPBY constant times a vector plus constant times a vector.
*>
*> Y = ALPHA * X + BETA * Y
*>
*> \endverbatim
*
* Arguments:
* ==========
*
*> \param[in] N
*> \verbatim
*> N is INTEGER
*> number of elements in input vector(s)
*> \endverbatim
*>
*> \param[in] DA
*> \verbatim
*> DA is DOUBLE PRECISION
*> On entry, DA specifies the scalar alpha.
*> \endverbatim
*>
*> \param[in] DX
*> \verbatim
*> DX is DOUBLE PRECISION array, dimension ( 1 + ( N - 1 )*abs( INCX ) )
*> \endverbatim
*>
*> \param[in] INCX
*> \verbatim
*> INCX is INTEGER
*> storage spacing between elements of DX
*> \endverbatim
*>
*> \param[in] DB
*> \verbatim
*> DB is DOUBLE PRECISION
*> On entry, DB specifies the scalar beta.
*> \endverbatim
*>
*> \param[in,out] DY
*> \verbatim
*> DY is DOUBLE PRECISION array, dimension ( 1 + ( N - 1 )*abs( INCY ) )
*> \endverbatim
*>
*> \param[in] INCY
*> \verbatim
*> INCY is INTEGER
*> storage spacing between elements of DY
*> \endverbatim
*
* Authors:
* ========
*
*> \author Univ. of Tennessee
*> \author Univ. of California Berkeley
*> \author Univ. of Colorado Denver
*> \author NAG Ltd.
*> \author Martin Koehler, MPI Magdeburg
*
*> \ingroup axpby
*
* =====================================================================
SUBROUTINE DAXPBY(N,DA,DX,INCX,DB,DY,INCY)
IMPLICIT NONE
*
* -- Reference BLAS level1 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
*
* .. Scalar Arguments ..
DOUBLE PRECISION DA,DB
INTEGER INCX,INCY,N
* ..
* .. Array Arguments ..
DOUBLE PRECISION DX(*),DY(*)
* ..
* .. External Subroutines
EXTERNAL DSCAL
*
* =====================================================================
*
* .. Local Scalars ..
INTEGER I,IX,IY,M,MP1
* ..
* .. Intrinsic Functions ..
INTRINSIC MOD
* ..
IF (N.LE.0) RETURN
* Scale if DA.EQ.0
IF (DA.EQ.0.0D0 .AND. DB.NE.0.0D0) THEN
CALL DSCAL(N, DB, DY, INCY)
RETURN
END IF
IF (INCX.EQ.1 .AND. INCY.EQ.1) THEN
*
* code for both increments equal to 1
*
*
*
DO I = 1,N
DY(I) = DB*DY(I) + DA*DX(I)
END DO
ELSE
*
* code for unequal increments or equal increments
* not equal to 1
*
IX = 1
IY = 1
IF (INCX.LT.0) IX = (-N+1)*INCX + 1
IF (INCY.LT.0) IY = (-N+1)*INCY + 1
DO I = 1,N
DY(IY) = DB*DY(IY) + DA*DX(IX)
IX = IX + INCX
IY = IY + INCY
END DO
END IF
RETURN
*
* End of DAXPBY
*
END
+1
View File
@@ -86,6 +86,7 @@
*>
* =====================================================================
SUBROUTINE DAXPY(N,DA,DX,INCX,DY,INCY)
IMPLICIT NONE
*
* -- Reference BLAS level1 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -44,6 +44,7 @@
*
* =====================================================================
DOUBLE PRECISION FUNCTION DCABS1(Z)
IMPLICIT NONE
*
* -- Reference BLAS level1 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -79,6 +79,7 @@
*>
* =====================================================================
SUBROUTINE DCOPY(N,DX,INCX,DY,INCY)
IMPLICIT NONE
*
* -- Reference BLAS level1 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -79,6 +79,7 @@
*>
* =====================================================================
DOUBLE PRECISION FUNCTION DDOT(N,DX,INCX,DY,INCY)
IMPLICIT NONE
*
* -- Reference BLAS level1 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -185,6 +185,7 @@
* =====================================================================
SUBROUTINE DGBMV(TRANS,M,N,KL,KU,ALPHA,A,LDA,X,INCX,
+ BETA,Y,INCY)
IMPLICIT NONE
*
* -- Reference BLAS level2 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+35 -7
View File
@@ -35,6 +35,16 @@
*>
*> alpha and beta are scalars, and A, B and C are matrices, with op( A )
*> an m by k matrix, op( B ) a k by n matrix and C an m by n matrix.
*>
*> Note: if alpha and/or beta is zero, some parts of the matrix-matrix
*> operations are not performed. This results in the following NaN/Inf
*> propagation quirks:
*>
*> 1. If alpha is zero, NaNs or Infs in A or B do not affect the result.
*> 2. If both alpha and beta are zero, then a zero matrix is returned in C,
*> irrespective of any NaNs or Infs in A, B or C.
*> 3. If only beta is zero, alpha*op( A )*op( B ) is returned, irrespective
*> of any NaNs or Infs in C.
*> \endverbatim
*
* Arguments:
@@ -51,6 +61,9 @@
*> TRANSA = 'T' or 't', op( A ) = A**T.
*>
*> TRANSA = 'C' or 'c', op( A ) = A**T.
*>
*> Note: TRANSA = 'C' is supported for the sake of API consistency
*> between all ?GEMM variants.
*> \endverbatim
*>
*> \param[in] TRANSB
@@ -64,6 +77,9 @@
*> TRANSB = 'T' or 't', op( B ) = B**T.
*>
*> TRANSB = 'C' or 'c', op( B ) = B**T.
*>
*> Note: TRANSB = 'C' is supported for the sake of API consistency
*> between all ?GEMM variants.
*> \endverbatim
*>
*> \param[in] M
@@ -92,7 +108,9 @@
*> \param[in] ALPHA
*> \verbatim
*> ALPHA is DOUBLE PRECISION.
*> On entry, ALPHA specifies the scalar alpha.
*> On entry, ALPHA specifies the scalar alpha. If ALPHA is zero the
*> values in A and B do not affect the result. This also means that
*> NaN/Inf propagation from A and B is inhibited if ALPHA is zero.
*> \endverbatim
*>
*> \param[in] A
@@ -102,7 +120,10 @@
*> Before entry with TRANSA = 'N' or 'n', the leading m by k
*> part of the array A must contain the matrix A, otherwise
*> the leading k by m part of the array A must contain the
*> matrix A.
*> matrix A, except if ALPHA is zero.
*> If ALPHA is zero, none of the values in A affect the result, even
*> if they are NaN/Inf. This also implies that if ALPHA is zero,
*> the matrix elements of A need not be initialized by the caller.
*> \endverbatim
*>
*> \param[in] LDA
@@ -121,7 +142,10 @@
*> Before entry with TRANSB = 'N' or 'n', the leading k by n
*> part of the array B must contain the matrix B, otherwise
*> the leading n by k part of the array B must contain the
*> matrix B.
*> matrix B, except if ALPHA is zero.
*> If ALPHA is zero, none of the values in B affect the result, even
*> if they are NaN/Inf. This also implies that if ALPHA is zero,
*> the matrix elements of B need not be initialized by the caller.
*> \endverbatim
*>
*> \param[in] LDB
@@ -136,16 +160,19 @@
*> \param[in] BETA
*> \verbatim
*> BETA is DOUBLE PRECISION.
*> On entry, BETA specifies the scalar beta. When BETA is
*> supplied as zero then C need not be set on input.
*> On entry, BETA specifies the scalar beta. If BETA is zero the
*> values in C do not affect the result. This also means that
*> NaN/Inf propagation from C is inhibited if BETA is zero.
*> \endverbatim
*>
*> \param[in,out] C
*> \verbatim
*> C is DOUBLE PRECISION array, dimension ( LDC, N )
*> Before entry, the leading m by n part of the array C must
*> contain the matrix C, except when beta is zero, in which
*> case C need not be set on entry.
*> contain the matrix C, except if beta is zero.
*> If beta is zero, none of the values in C affect the result, even
*> if they are NaN/Inf. This also implies that if beta is zero,
*> the matrix elements of C need not be initialized by the caller.
*> On exit, the array C is overwritten by the m by n matrix
*> ( alpha*op( A )*op( B ) + beta*C ).
*> \endverbatim
@@ -185,6 +212,7 @@
* =====================================================================
SUBROUTINE DGEMM(TRANSA,TRANSB,M,N,K,ALPHA,A,LDA,B,LDB,
+ BETA,C,LDC)
IMPLICIT NONE
*
* -- Reference BLAS level3 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+431
View File
@@ -0,0 +1,431 @@
*> \brief \b DGEMMTR
*
* =========== DOCUMENTATION ===========
*
* Online html documentation available at
* http://www.netlib.org/lapack/explore-html/
*
* Definition:
* ===========
*
* SUBROUTINE DGEMMTR(UPLO,TRANSA,TRANSB,N,K,ALPHA,A,LDA,B,LDB,BETA,
* C,LDC)
*
* .. Scalar Arguments ..
* DOUBLE PRECISION ALPHA,BETA
* INTEGER K,LDA,LDB,LDC,N
* CHARACTER TRANSA,TRANSB, UPLO
* ..
* .. Array Arguments ..
* DOUBLE PRECISION A(LDA,*),B(LDB,*),C(LDC,*)
* ..
*
*
*> \par Purpose:
* =============
*>
*> \verbatim
*>
*> DGEMMTR performs one of the matrix-matrix operations
*>
*> C := alpha*op( A )*op( B ) + beta*C,
*>
*> where op( X ) is one of
*>
*> op( X ) = X or op( X ) = X**T,
*>
*> alpha and beta are scalars, and A, B and C are matrices, with op( A )
*> an n by k matrix, op( B ) a k by n matrix and C an n by n matrix.
*> Thereby, the routine only accesses and updates the upper or lower
*> triangular part of the result matrix C. This behaviour can be used if
*> the resulting matrix C is known to be symmetric.
*> \endverbatim
*
* Arguments:
* ==========
*
*> \param[in] UPLO
*> \verbatim
*> UPLO is CHARACTER*1
*> On entry, UPLO specifies whether the lower or the upper
*> triangular part of C is access and updated.
*>
*> UPLO = 'L' or 'l', the lower triangular part of C is used.
*>
*> UPLO = 'U' or 'u', the upper triangular part of C is used.
*> \endverbatim
*
*> \param[in] TRANSA
*> \verbatim
*> TRANSA is CHARACTER*1
*> On entry, TRANSA specifies the form of op( A ) to be used in
*> the matrix multiplication as follows:
*>
*> TRANSA = 'N' or 'n', op( A ) = A.
*>
*> TRANSA = 'T' or 't', op( A ) = A**T.
*>
*> TRANSA = 'C' or 'c', op( A ) = A**T.
*> \endverbatim
*>
*> \param[in] TRANSB
*> \verbatim
*> TRANSB is CHARACTER*1
*> On entry, TRANSB specifies the form of op( B ) to be used in
*> the matrix multiplication as follows:
*>
*> TRANSB = 'N' or 'n', op( B ) = B.
*>
*> TRANSB = 'T' or 't', op( B ) = B**T.
*>
*> TRANSB = 'C' or 'c', op( B ) = B**T.
*> \endverbatim
*>
*> \param[in] N
*> \verbatim
*> N is INTEGER
*> On entry, N specifies the number of rows and columns of
*> the matrix C, the number of columns of op(B) and the number
*> of rows of op(A). N must be at least zero.
*> \endverbatim
*>
*> \param[in] K
*> \verbatim
*> K is INTEGER
*> On entry, K specifies the number of columns of the matrix
*> op( A ) and the number of rows of the matrix op( B ). K must
*> be at least zero.
*> \endverbatim
*>
*> \param[in] ALPHA
*> \verbatim
*> ALPHA is DOUBLE PRECISION.
*> On entry, ALPHA specifies the scalar alpha.
*> \endverbatim
*>
*> \param[in] A
*> \verbatim
*> A is DOUBLE PRECISION array, dimension ( LDA, ka ), where ka is
*> k when TRANSA = 'N' or 'n', and is n otherwise.
*> Before entry with TRANSA = 'N' or 'n', the leading n by k
*> part of the array A must contain the matrix A, otherwise
*> the leading k by m part of the array A must contain the
*> matrix A.
*> \endverbatim
*>
*> \param[in] LDA
*> \verbatim
*> LDA is INTEGER
*> On entry, LDA specifies the first dimension of A as declared
*> in the calling (sub) program. When TRANSA = 'N' or 'n' then
*> LDA must be at least max( 1, n ), otherwise LDA must be at
*> least max( 1, k ).
*> \endverbatim
*>
*> \param[in] B
*> \verbatim
*> B is DOUBLE PRECISION array, dimension ( LDB, kb ), where kb is
*> n when TRANSB = 'N' or 'n', and is k otherwise.
*> Before entry with TRANSB = 'N' or 'n', the leading k by n
*> part of the array B must contain the matrix B, otherwise
*> the leading n by k part of the array B must contain the
*> matrix B.
*> \endverbatim
*>
*> \param[in] LDB
*> \verbatim
*> LDB is INTEGER
*> On entry, LDB specifies the first dimension of B as declared
*> in the calling (sub) program. When TRANSB = 'N' or 'n' then
*> LDB must be at least max( 1, k ), otherwise LDB must be at
*> least max( 1, n ).
*> \endverbatim
*>
*> \param[in] BETA
*> \verbatim
*> BETA is DOUBLE PRECISION.
*> On entry, BETA specifies the scalar beta. When BETA is
*> supplied as zero then C need not be set on input.
*> \endverbatim
*>
*> \param[in,out] C
*> \verbatim
*> C is DOUBLE PRECISION array, dimension ( LDC, N )
*> Before entry, the leading n by n part of the array C must
*> contain the matrix C, except when beta is zero, in which
*> case C need not be set on entry.
*> On exit, the upper or lower triangular part of the matrix
*> C is overwritten by the n by n matrix
*> ( alpha*op( A )*op( B ) + beta*C ).
*> \endverbatim
*>
*> \param[in] LDC
*> \verbatim
*> LDC is INTEGER
*> On entry, LDC specifies the first dimension of C as declared
*> in the calling (sub) program. LDC must be at least
*> max( 1, n ).
*> \endverbatim
*
* Authors:
* ========
*
*> \author Martin Koehler
*
*> \ingroup gemmtr
*
*> \par Further Details:
* =====================
*>
*> \verbatim
*>
*> Level 3 Blas routine.
*>
*> -- Written on 19-July-2023.
*> Martin Koehler, MPI Magdeburg
*> \endverbatim
*>
* =====================================================================
SUBROUTINE DGEMMTR(UPLO,TRANSA,TRANSB,N,K,ALPHA,A,LDA,B,LDB,
+ BETA,C,LDC)
IMPLICIT NONE
*
* -- Reference BLAS level3 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
*
* .. Scalar Arguments ..
DOUBLE PRECISION ALPHA,BETA
INTEGER K,LDA,LDB,LDC,N
CHARACTER TRANSA,TRANSB,UPLO
* ..
* .. Array Arguments ..
DOUBLE PRECISION A(LDA,*),B(LDB,*),C(LDC,*)
* ..
*
* =====================================================================
*
* .. External Functions ..
LOGICAL LSAME
EXTERNAL LSAME
* ..
* .. External Subroutines ..
EXTERNAL XERBLA
* ..
* .. Intrinsic Functions ..
INTRINSIC MAX
* ..
* .. Local Scalars ..
DOUBLE PRECISION TEMP
INTEGER I,INFO,J,L,NROWA,NROWB, ISTART, ISTOP
LOGICAL NOTA,NOTB, UPPER
* ..
* .. Parameters ..
DOUBLE PRECISION ONE,ZERO
PARAMETER (ONE=1.0D+0,ZERO=0.0D+0)
* ..
*
* Set NOTA and NOTB as true if A and B respectively are not
* transposed and set NROWA and NROWB as the number of rows of A
* and B respectively.
*
NOTA = LSAME(TRANSA,'N')
NOTB = LSAME(TRANSB,'N')
IF (NOTA) THEN
NROWA = N
ELSE
NROWA = K
END IF
IF (NOTB) THEN
NROWB = K
ELSE
NROWB = N
END IF
UPPER = LSAME(UPLO, 'U')
*
* Test the input parameters.
*
INFO = 0
IF ((.NOT. UPPER) .AND. (.NOT. LSAME(UPLO, 'L'))) THEN
INFO = 1
ELSE IF ((.NOT.NOTA) .AND. (.NOT.LSAME(TRANSA,'C')) .AND.
+ (.NOT.LSAME(TRANSA,'T'))) THEN
INFO = 2
ELSE IF ((.NOT.NOTB) .AND. (.NOT.LSAME(TRANSB,'C')) .AND.
+ (.NOT.LSAME(TRANSB,'T'))) THEN
INFO = 3
ELSE IF (N.LT.0) THEN
INFO = 4
ELSE IF (K.LT.0) THEN
INFO = 5
ELSE IF (LDA.LT.MAX(1,NROWA)) THEN
INFO = 8
ELSE IF (LDB.LT.MAX(1,NROWB)) THEN
INFO = 10
ELSE IF (LDC.LT.MAX(1,N)) THEN
INFO = 13
END IF
IF (INFO.NE.0) THEN
CALL XERBLA('DGEMMTR',INFO)
RETURN
END IF
*
* Quick return if possible.
*
IF (N.EQ.0) RETURN
*
* And if alpha.eq.zero.
*
IF (ALPHA.EQ.ZERO) THEN
IF (BETA.EQ.ZERO) THEN
DO 20 J = 1,N
IF (UPPER) THEN
ISTART = 1
ISTOP = J
ELSE
ISTART = J
ISTOP = N
END IF
DO 10 I = ISTART, ISTOP
C(I,J) = ZERO
10 CONTINUE
20 CONTINUE
ELSE
DO 40 J = 1,N
IF (UPPER) THEN
ISTART = 1
ISTOP = J
ELSE
ISTART = J
ISTOP = N
END IF
DO 30 I = ISTART, ISTOP
C(I,J) = BETA*C(I,J)
30 CONTINUE
40 CONTINUE
END IF
RETURN
END IF
*
* Start the operations.
*
IF (NOTB) THEN
IF (NOTA) THEN
*
* Form C := alpha*A*B + beta*C.
*
DO 90 J = 1,N
IF (UPPER) THEN
ISTART = 1
ISTOP = J
ELSE
ISTART = J
ISTOP = N
END IF
IF (BETA.EQ.ZERO) THEN
DO 50 I = ISTART, ISTOP
C(I,J) = ZERO
50 CONTINUE
ELSE IF (BETA.NE.ONE) THEN
DO 60 I = ISTART, ISTOP
C(I,J) = BETA*C(I,J)
60 CONTINUE
END IF
DO 80 L = 1,K
TEMP = ALPHA*B(L,J)
DO 70 I = ISTART, ISTOP
C(I,J) = C(I,J) + TEMP*A(I,L)
70 CONTINUE
80 CONTINUE
90 CONTINUE
ELSE
*
* Form C := alpha*A**T*B + beta*C
*
DO 120 J = 1,N
IF (UPPER) THEN
ISTART = 1
ISTOP = J
ELSE
ISTART = J
ISTOP = N
END IF
DO 110 I = ISTART, ISTOP
TEMP = ZERO
DO 100 L = 1,K
TEMP = TEMP + A(L,I)*B(L,J)
100 CONTINUE
IF (BETA.EQ.ZERO) THEN
C(I,J) = ALPHA*TEMP
ELSE
C(I,J) = ALPHA*TEMP + BETA*C(I,J)
END IF
110 CONTINUE
120 CONTINUE
END IF
ELSE
IF (NOTA) THEN
*
* Form C := alpha*A*B**T + beta*C
*
DO 170 J = 1,N
IF (UPPER) THEN
ISTART = 1
ISTOP = J
ELSE
ISTART = J
ISTOP = N
END IF
IF (BETA.EQ.ZERO) THEN
DO 130 I = ISTART,ISTOP
C(I,J) = ZERO
130 CONTINUE
ELSE IF (BETA.NE.ONE) THEN
DO 140 I = ISTART,ISTOP
C(I,J) = BETA*C(I,J)
140 CONTINUE
END IF
DO 160 L = 1,K
TEMP = ALPHA*B(J,L)
DO 150 I = ISTART,ISTOP
C(I,J) = C(I,J) + TEMP*A(I,L)
150 CONTINUE
160 CONTINUE
170 CONTINUE
ELSE
*
* Form C := alpha*A**T*B**T + beta*C
*
DO 200 J = 1,N
IF (UPPER) THEN
ISTART = 1
ISTOP = J
ELSE
ISTART = J
ISTOP = N
END IF
DO 190 I = ISTART, ISTOP
TEMP = ZERO
DO 180 L = 1,K
TEMP = TEMP + A(L,I)*B(J,L)
180 CONTINUE
IF (BETA.EQ.ZERO) THEN
C(I,J) = ALPHA*TEMP
ELSE
C(I,J) = ALPHA*TEMP + BETA*C(I,J)
END IF
190 CONTINUE
200 CONTINUE
END IF
END IF
*
RETURN
*
* End of DGEMMTR
*
END
+1
View File
@@ -155,6 +155,7 @@
*>
* =====================================================================
SUBROUTINE DGEMV(TRANS,M,N,ALPHA,A,LDA,X,INCX,BETA,Y,INCY)
IMPLICIT NONE
*
* -- Reference BLAS level2 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -127,6 +127,7 @@
*>
* =====================================================================
SUBROUTINE DGER(M,N,ALPHA,X,INCX,Y,INCY,A,LDA)
IMPLICIT NONE
*
* -- Reference BLAS level2 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+3 -2
View File
@@ -85,11 +85,12 @@
!> \endverbatim
!>
! =====================================================================
function DNRM2( n, x, incx )
function DNRM2( n, x, incx )
implicit none
integer, parameter :: wp = kind(1.d0)
real(wp) :: DNRM2
!
! -- Reference BLAS level1 routine (version 3.9.1) --
! -- Reference BLAS level1 routine --
! -- Reference BLAS is a software package provided by Univ. of Tennessee, --
! -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
! March 2021
+1
View File
@@ -89,6 +89,7 @@
*>
* =====================================================================
SUBROUTINE DROT(N,DX,INCX,DY,INCY,C,S)
IMPLICIT NONE
*
* -- Reference BLAS level1 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -89,6 +89,7 @@
!
! =====================================================================
subroutine DROTG( a, b, c, s )
implicit none
integer, parameter :: wp = kind(1.d0)
!
! -- Reference BLAS level1 routine --
+5
View File
@@ -38,6 +38,10 @@
*> H=( ) ( ) ( ) ( )
*> (DH21 DH22), (DH21 1.D0), (-1.D0 DH22), (0.D0 1.D0).
*> SEE DROTMG FOR A DESCRIPTION OF DATA STORAGE IN DPARAM.
*>
*> IF DFLAG IS NOT ONE OF THE LISTED ABOVE, THE BEHAVIOR IS UNDEFINED.
*> NANS IN DFLAG MAY NOT PROPAGATE TO THE OUTPUT.
*>
*> \endverbatim
*
* Arguments:
@@ -93,6 +97,7 @@
*
* =====================================================================
SUBROUTINE DROTM(N,DX,INCX,DY,INCY,DPARAM)
IMPLICIT NONE
*
* -- Reference BLAS level1 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -87,6 +87,7 @@
*
* =====================================================================
SUBROUTINE DROTMG(DD1,DD2,DX1,DY1,DPARAM)
IMPLICIT NONE
*
* -- Reference BLAS level1 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -181,6 +181,7 @@
*>
* =====================================================================
SUBROUTINE DSBMV(UPLO,N,K,ALPHA,A,LDA,X,INCX,BETA,Y,INCY)
IMPLICIT NONE
*
* -- Reference BLAS level2 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -76,6 +76,7 @@
*>
* =====================================================================
SUBROUTINE DSCAL(N,DA,DX,INCX)
IMPLICIT NONE
*
* -- Reference BLAS level1 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -116,6 +116,7 @@
*>
* =====================================================================
DOUBLE PRECISION FUNCTION DSDOT(N,SX,INCX,SY,INCY)
IMPLICIT NONE
*
* -- Reference BLAS level1 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+365
View File
@@ -0,0 +1,365 @@
*> \brief \b DSKEWSYMM
*
* =========== DOCUMENTATION ===========
*
* Online html documentation available at
* http://www.netlib.org/lapack/explore-html/
*
* Definition:
* ===========
*
* SUBROUTINE DSKEWSYMM(SIDE,UPLO,M,N,ALPHA,A,LDA,B,LDB,BETA,C,LDC)
*
* .. Scalar Arguments ..
* DOUBLE PRECISION ALPHA,BETA
* INTEGER LDA,LDB,LDC,M,N
* CHARACTER SIDE,UPLO
* ..
* .. Array Arguments ..
* DOUBLE PRECISION A(LDA,*),B(LDB,*),C(LDC,*)
* ..
*
*
*> \par Purpose:
* =============
*>
*> \verbatim
*>
*> DSKEWSYMM performs one of the matrix-matrix operations
*>
*> C := alpha*A*B + beta*C,
*>
*> or
*>
*> C := alpha*B*A + beta*C,
*>
*> where alpha and beta are scalars, A is a skew-symmetric matrix and B and
*> C are m by n matrices.
*> \endverbatim
*
* Arguments:
* ==========
*
*> \param[in] SIDE
*> \verbatim
*> SIDE is CHARACTER*1
*> On entry, SIDE specifies whether the skew-symmetric matrix A
*> appears on the left or right in the operation as follows:
*>
*> SIDE = 'L' or 'l' C := alpha*A*B + beta*C,
*>
*> SIDE = 'R' or 'r' C := alpha*B*A + beta*C,
*> \endverbatim
*>
*> \param[in] UPLO
*> \verbatim
*> UPLO is CHARACTER*1
*> On entry, UPLO specifies whether the upper or lower
*> triangular part of the skew-symmetric matrix A is to be
*> referenced as follows:
*>
*> UPLO = 'U' or 'u' Only the upper triangular part of the
*> skew-symmetric matrix is to be referenced.
*>
*> UPLO = 'L' or 'l' Only the lower triangular part of the
*> skew-symmetric matrix is to be referenced.
*> \endverbatim
*>
*> \param[in] M
*> \verbatim
*> M is INTEGER
*> On entry, M specifies the number of rows of the matrix C.
*> M must be at least zero.
*> \endverbatim
*>
*> \param[in] N
*> \verbatim
*> N is INTEGER
*> On entry, N specifies the number of columns of the matrix C.
*> N must be at least zero.
*> \endverbatim
*>
*> \param[in] ALPHA
*> \verbatim
*> ALPHA is DOUBLE PRECISION
*> On entry, ALPHA specifies the scalar alpha.
*> \endverbatim
*>
*> \param[in] A
*> \verbatim
*> A is DOUBLE PRECISION array, dimension ( LDA, ka ), where ka is
*> m when SIDE = 'L' or 'l' and is n otherwise.
*> Before entry with SIDE = 'L' or 'l', the m by m part of
*> the array A must contain the skew-symmetric matrix, such that
*> when UPLO = 'U' or 'u', the strictly m by m upper triangular
*> part of the array A must contain the upper triangular part
*> of the skew-symmetric matrix and the leading lower triangular
*> part of A is not referenced, and when UPLO = 'L' or 'l',
*> the strictly m by m lower triangular part of the array A
*> must contain the lower triangular part of the skew-symmetric
*> matrix and the leading upper triangular part of A is not
*> referenced.
*> Before entry with SIDE = 'R' or 'r', the n by n part of
*> the array A must contain the skew-symmetric matrix, such that
*> when UPLO = 'U' or 'u', the strictly n by n upper triangular
*> part of the array A must contain the upper triangular part
*> of the skew-symmetric matrix and the leading lower triangular
*> part of A is not referenced, and when UPLO = 'L' or 'l',
*> the strictly n by n lower triangular part of the array A
*> must contain the lower triangular part of the skew-symmetric
*> matrix and the leading upper triangular part of A is not
*> referenced.
*> \endverbatim
*>
*> \param[in] LDA
*> \verbatim
*> LDA is INTEGER
*> On entry, LDA specifies the first dimension of A as declared
*> in the calling (sub) program. When SIDE = 'L' or 'l' then
*> LDA must be at least max( 1, m ), otherwise LDA must be at
*> least max( 1, n ).
*> \endverbatim
*>
*> \param[in] B
*> \verbatim
*> B is DOUBLE PRECISION array, dimension ( LDB, N )
*> Before entry, the leading m by n part of the array B must
*> contain the matrix B.
*> \endverbatim
*>
*> \param[in] LDB
*> \verbatim
*> LDB is INTEGER
*> On entry, LDB specifies the first dimension of B as declared
*> in the calling (sub) program. LDB must be at least
*> max( 1, m ).
*> \endverbatim
*>
*> \param[in] BETA
*> \verbatim
*> BETA is DOUBLE PRECISION.
*> On entry, BETA specifies the scalar beta. When BETA is
*> supplied as zero then C need not be set on input.
*> \endverbatim
*>
*> \param[in,out] C
*> \verbatim
*> C is DOUBLE PRECISION array, dimension ( LDC, N )
*> Before entry, the leading m by n part of the array C must
*> contain the matrix C, except when beta is zero, in which
*> case C need not be set on entry.
*> On exit, the array C is overwritten by the m by n updated
*> matrix.
*> \endverbatim
*>
*> \param[in] LDC
*> \verbatim
*> LDC is INTEGER
*> On entry, LDC specifies the first dimension of C as declared
*> in the calling (sub) program. LDC must be at least
*> max( 1, m ).
*> \endverbatim
*
* Authors:
* ========
*
*> \author Univ. of Tennessee
*> \author Univ. of California Berkeley
*> \author Univ. of Colorado Denver
*> \author NAG Ltd.
*
*> \ingroup skewhemm
*
*> \par Further Details:
* =====================
*>
*> \verbatim
*>
*> Level 3 Blas routine.
*> Derived from subroutine dsymm.
*>
*> -- Written on 6-Jul-2025.
*> Shuo Zheng, China.
*> \endverbatim
*>
* =====================================================================
SUBROUTINE DSKEWSYMM(SIDE,UPLO,M,N,ALPHA,A,LDA,B,
+ LDB,BETA,C,LDC)
IMPLICIT NONE
*
* -- Reference BLAS level3 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
*
* .. Scalar Arguments ..
DOUBLE PRECISION ALPHA,BETA
INTEGER LDA,LDB,LDC,M,N
CHARACTER SIDE,UPLO
* ..
* .. Array Arguments ..
DOUBLE PRECISION A(LDA,*),B(LDB,*),C(LDC,*)
* ..
*
* =====================================================================
*
* .. External Functions ..
LOGICAL LSAME
EXTERNAL LSAME
* ..
* .. External Subroutines ..
EXTERNAL XERBLA
* ..
* .. Intrinsic Functions ..
INTRINSIC MAX
* ..
* .. Local Scalars ..
DOUBLE PRECISION TEMP1,TEMP2
INTEGER I,INFO,J,K,NROWA
LOGICAL UPPER
* ..
* .. Parameters ..
DOUBLE PRECISION ONE,ZERO
PARAMETER (ONE=1.0D+0,ZERO=0.0D+0)
* ..
*
* Set NROWA as the number of rows of A.
*
IF (LSAME(SIDE,'L')) THEN
NROWA = M
ELSE
NROWA = N
END IF
UPPER = LSAME(UPLO,'U')
*
* Test the input parameters.
*
INFO = 0
IF ((.NOT.LSAME(SIDE,'L')) .AND.
+ (.NOT.LSAME(SIDE,'R'))) THEN
INFO = 1
ELSE IF ((.NOT.UPPER) .AND.
+ (.NOT.LSAME(UPLO,'L'))) THEN
INFO = 2
ELSE IF (M.LT.0) THEN
INFO = 3
ELSE IF (N.LT.0) THEN
INFO = 4
ELSE IF (LDA.LT.MAX(1,NROWA)) THEN
INFO = 7
ELSE IF (LDB.LT.MAX(1,M)) THEN
INFO = 9
ELSE IF (LDC.LT.MAX(1,M)) THEN
INFO = 12
END IF
IF (INFO.NE.0) THEN
CALL XERBLA('DSKEWSYMM ',INFO)
RETURN
END IF
*
* Quick return if possible.
*
IF ((M.EQ.0) .OR. (N.EQ.0) .OR.
+ ((ALPHA.EQ.ZERO).AND. (BETA.EQ.ONE))) RETURN
*
* And when alpha.eq.zero.
*
IF (ALPHA.EQ.ZERO) THEN
IF (BETA.EQ.ZERO) THEN
DO 20 J = 1,N
DO 10 I = 1,M
C(I,J) = ZERO
10 CONTINUE
20 CONTINUE
ELSE
DO 40 J = 1,N
DO 30 I = 1,M
C(I,J) = BETA*C(I,J)
30 CONTINUE
40 CONTINUE
END IF
RETURN
END IF
*
* Start the operations.
*
IF (LSAME(SIDE,'L')) THEN
*
* Form C := alpha*A*B + beta*C.
*
IF (UPPER) THEN
DO 70 J = 1,N
DO 60 I = 1,M
TEMP1 = ALPHA*B(I,J)
TEMP2 = ZERO
DO 50 K = 1,I - 1
C(K,J) = C(K,J) + TEMP1*A(K,I)
TEMP2 = TEMP2 - B(K,J)*A(K,I)
50 CONTINUE
IF (BETA.EQ.ZERO) THEN
C(I,J) = ALPHA*TEMP2
ELSE
C(I,J) = BETA*C(I,J) +
+ ALPHA*TEMP2
END IF
60 CONTINUE
70 CONTINUE
ELSE
DO 100 J = 1,N
DO 90 I = M,1,-1
TEMP1 = ALPHA*B(I,J)
TEMP2 = ZERO
DO 80 K = I + 1,M
C(K,J) = C(K,J) + TEMP1*A(K,I)
TEMP2 = TEMP2 - B(K,J)*A(K,I)
80 CONTINUE
IF (BETA.EQ.ZERO) THEN
C(I,J) = ALPHA*TEMP2
ELSE
C(I,J) = BETA*C(I,J) +
+ ALPHA*TEMP2
END IF
90 CONTINUE
100 CONTINUE
END IF
ELSE
*
* Form C := alpha*B*A + beta*C.
*
DO 170 J = 1,N
IF (BETA.EQ.ZERO) THEN
DO 110 I = 1,M
C(I,J) = ZERO
110 CONTINUE
ELSE
DO 120 I = 1,M
C(I,J) = BETA*C(I,J)
120 CONTINUE
END IF
DO 140 K = 1,J - 1
IF (UPPER) THEN
TEMP1 = ALPHA*A(K,J)
ELSE
TEMP1 = -ALPHA*A(J,K)
END IF
DO 130 I = 1,M
C(I,J) = C(I,J) + TEMP1*B(I,K)
130 CONTINUE
140 CONTINUE
DO 160 K = J + 1,N
IF (UPPER) THEN
TEMP1 = -ALPHA*A(J,K)
ELSE
TEMP1 = ALPHA*A(K,J)
END IF
DO 150 I = 1,M
C(I,J) = C(I,J) + TEMP1*B(I,K)
150 CONTINUE
160 CONTINUE
170 CONTINUE
END IF
*
RETURN
*
* End of DSKEWSYMM
*
END
+327
View File
@@ -0,0 +1,327 @@
*> \brief \b DSKEWSYMV
*
* =========== DOCUMENTATION ===========
*
* Online html documentation available at
* http://www.netlib.org/lapack/explore-html/
*
* Definition:
* ===========
*
* SUBROUTINE DSKEWSYMV(UPLO,N,ALPHA,A,LDA,X,INCX,BETA,Y,INCY)
*
* .. Scalar Arguments ..
* DOUBLE PRECISION ALPHA,BETA
* INTEGER INCX,INCY,LDA,N
* CHARACTER UPLO
* ..
* .. Array Arguments ..
* DOUBLE PRECISION A(LDA,*),X(*),Y(*)
* ..
*
*
*> \par Purpose:
* =============
*>
*> \verbatim
*>
*> DSKEWSYMV performs the matrix-vector operation
*>
*> y := alpha*A*x + beta*y,
*>
*> where alpha and beta are scalars, x and y are n element vectors and
*> A is an n by n skew-symmetric matrix.
*> \endverbatim
*
* Arguments:
* ==========
*
*> \param[in] UPLO
*> \verbatim
*> UPLO is CHARACTER*1
*> On entry, UPLO specifies whether the upper or lower
*> triangular part of the array A is to be referenced as
*> follows:
*>
*> UPLO = 'U' or 'u' Only the upper triangular part of A
*> is to be referenced.
*>
*> UPLO = 'L' or 'l' Only the lower triangular part of A
*> is to be referenced.
*> \endverbatim
*>
*> \param[in] N
*> \verbatim
*> N is INTEGER
*> On entry, N specifies the order of the matrix A.
*> N must be at least zero.
*> \endverbatim
*>
*> \param[in] ALPHA
*> \verbatim
*> ALPHA is DOUBLE PRECISION
*> On entry, ALPHA specifies the scalar alpha.
*> \endverbatim
*>
*> \param[in] A
*> \verbatim
*> A is DOUBLE PRECISION array, dimension ( LDA, N )
*> Before entry with UPLO = 'U' or 'u', the strictly n by n
*> upper triangular part of the array A must contain the upper
*> triangular part of the skew-symmetric matrix and the leading
*> lower triangular part of A is not referenced.
*> Before entry with UPLO = 'L' or 'l', the strictly n by n
*> lower triangular part of the array A must contain the lower
*> triangular part of the skew-symmetric matrix and the leading
*> upper triangular part of A is not referenced.
*> \endverbatim
*>
*> \param[in] LDA
*> \verbatim
*> LDA is INTEGER
*> On entry, LDA specifies the first dimension of A as declared
*> in the calling (sub) program. LDA must be at least
*> max( 1, n ).
*> \endverbatim
*>
*> \param[in] X
*> \verbatim
*> X is DOUBLE PRECISION array, dimension at least
*> ( 1 + ( n - 1 )*abs( INCX ) ).
*> Before entry, the incremented array X must contain the n
*> element vector x.
*> \endverbatim
*>
*> \param[in] INCX
*> \verbatim
*> INCX is INTEGER
*> On entry, INCX specifies the increment for the elements of
*> X. INCX must not be zero.
*> \endverbatim
*>
*> \param[in] BETA
*> \verbatim
*> BETA is DOUBLE PRECISION.
*> On entry, BETA specifies the scalar beta. When BETA is
*> supplied as zero then Y need not be set on input.
*> \endverbatim
*>
*> \param[in,out] Y
*> \verbatim
*> Y is DOUBLE PRECISION array, dimension at least
*> ( 1 + ( n - 1 )*abs( INCY ) ).
*> Before entry, the incremented array Y must contain the n
*> element vector y. On exit, Y is overwritten by the updated
*> vector y.
*> \endverbatim
*>
*> \param[in] INCY
*> \verbatim
*> INCY is INTEGER
*> On entry, INCY specifies the increment for the elements of
*> Y. INCY must not be zero.
*> \endverbatim
*
* Authors:
* ========
*
*> \author Univ. of Tennessee
*> \author Univ. of California Berkeley
*> \author Univ. of Colorado Denver
*> \author NAG Ltd.
*
*> \ingroup skewhemv
*
*> \par Further Details:
* =====================
*>
*> \verbatim
*>
*> Level 2 Blas routine.
*> The vector and matrix arguments are not referenced when N = 0, or M = 0
*> Derived from subroutine dsymv.
*>
*> -- Written on 6-Jul-2025.
*> Shuo Zheng, China.
*> \endverbatim
*>
* =====================================================================
SUBROUTINE DSKEWSYMV(UPLO,N,ALPHA,A,LDA,X,INCX,BETA,Y,INCY)
IMPLICIT NONE
*
* -- Reference BLAS level2 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
*
* .. Scalar Arguments ..
DOUBLE PRECISION ALPHA,BETA
INTEGER INCX,INCY,LDA,N
CHARACTER UPLO
* ..
* .. Array Arguments ..
DOUBLE PRECISION A(LDA,*),X(*),Y(*)
* ..
*
* =====================================================================
*
* .. Parameters ..
DOUBLE PRECISION ONE,ZERO
PARAMETER (ONE=1.0D+0,ZERO=0.0D+0)
* ..
* .. Local Scalars ..
DOUBLE PRECISION TEMP1,TEMP2
INTEGER I,INFO,IX,IY,J,JX,JY,KX,KY
* ..
* .. External Functions ..
LOGICAL LSAME
EXTERNAL LSAME
* ..
* .. External Subroutines ..
EXTERNAL XERBLA
* ..
* .. Intrinsic Functions ..
INTRINSIC MAX
* ..
*
* Test the input parameters.
*
INFO = 0
IF (.NOT.LSAME(UPLO,'U') .AND. .NOT.LSAME(UPLO,'L')) THEN
INFO = 1
ELSE IF (N.LT.0) THEN
INFO = 2
ELSE IF (LDA.LT.MAX(1,N)) THEN
INFO = 5
ELSE IF (INCX.EQ.0) THEN
INFO = 7
ELSE IF (INCY.EQ.0) THEN
INFO = 10
END IF
IF (INFO.NE.0) THEN
CALL XERBLA('DSKEWSYMV ',INFO)
RETURN
END IF
*
* Quick return if possible.
*
IF ((N.EQ.0) .OR. ((ALPHA.EQ.ZERO).AND. (BETA.EQ.ONE))) RETURN
*
* Set up the start points in X and Y.
*
IF (INCX.GT.0) THEN
KX = 1
ELSE
KX = 1 - (N-1)*INCX
END IF
IF (INCY.GT.0) THEN
KY = 1
ELSE
KY = 1 - (N-1)*INCY
END IF
*
* Start the operations. In this version the elements of A are
* accessed sequentially with one pass through the triangular part
* of A.
*
* First form y := beta*y.
*
IF (BETA.NE.ONE) THEN
IF (INCY.EQ.1) THEN
IF (BETA.EQ.ZERO) THEN
DO 10 I = 1,N
Y(I) = ZERO
10 CONTINUE
ELSE
DO 20 I = 1,N
Y(I) = BETA*Y(I)
20 CONTINUE
END IF
ELSE
IY = KY
IF (BETA.EQ.ZERO) THEN
DO 30 I = 1,N
Y(IY) = ZERO
IY = IY + INCY
30 CONTINUE
ELSE
DO 40 I = 1,N
Y(IY) = BETA*Y(IY)
IY = IY + INCY
40 CONTINUE
END IF
END IF
END IF
IF (ALPHA.EQ.ZERO) RETURN
IF (LSAME(UPLO,'U')) THEN
*
* Form y when A is stored in upper triangle.
*
IF ((INCX.EQ.1) .AND. (INCY.EQ.1)) THEN
DO 60 J = 1,N
TEMP1 = ALPHA*X(J)
TEMP2 = ZERO
DO 50 I = 1,J - 1
Y(I) = Y(I) + TEMP1*A(I,J)
TEMP2 = TEMP2 - A(I,J)*X(I)
50 CONTINUE
Y(J) = Y(J) + ALPHA*TEMP2
60 CONTINUE
ELSE
JX = KX
JY = KY
DO 80 J = 1,N
TEMP1 = ALPHA*X(JX)
TEMP2 = ZERO
IX = KX
IY = KY
DO 70 I = 1,J - 1
Y(IY) = Y(IY) + TEMP1*A(I,J)
TEMP2 = TEMP2 - A(I,J)*X(IX)
IX = IX + INCX
IY = IY + INCY
70 CONTINUE
Y(JY) = Y(JY) + ALPHA*TEMP2
JX = JX + INCX
JY = JY + INCY
80 CONTINUE
END IF
ELSE
*
* Form y when A is stored in lower triangle.
*
IF ((INCX.EQ.1) .AND. (INCY.EQ.1)) THEN
DO 100 J = 1,N
TEMP1 = ALPHA*X(J)
TEMP2 = ZERO
DO 90 I = J + 1,N
Y(I) = Y(I) + TEMP1*A(I,J)
TEMP2 = TEMP2 - A(I,J)*X(I)
90 CONTINUE
Y(J) = Y(J) + ALPHA*TEMP2
100 CONTINUE
ELSE
JX = KX
JY = KY
DO 120 J = 1,N
TEMP1 = ALPHA*X(JX)
TEMP2 = ZERO
IX = JX
IY = JY
DO 110 I = J + 1,N
IX = IX + INCX
IY = IY + INCY
Y(IY) = Y(IY) + TEMP1*A(I,J)
TEMP2 = TEMP2 - A(I,J)*X(IX)
110 CONTINUE
Y(JY) = Y(JY) + ALPHA*TEMP2
JX = JX + INCX
JY = JY + INCY
120 CONTINUE
END IF
END IF
*
RETURN
*
* End of DSKEWSYMV
*
END
+294
View File
@@ -0,0 +1,294 @@
*> \brief \b DSKEWSYR2
*
* =========== DOCUMENTATION ===========
*
* Online html documentation available at
* http://www.netlib.org/lapack/explore-html/
*
* Definition:
* ===========
*
* SUBROUTINE DSKEWSYR2(UPLO,N,ALPHA,X,INCX,Y,INCY,A,LDA)
*
* .. Scalar Arguments ..
* DOUBLE PRECISION ALPHA
* INTEGER INCX,INCY,LDA,N
* CHARACTER UPLO
* ..
* .. Array Arguments ..
* DOUBLE PRECISION A(LDA,*),X(*),Y(*)
* ..
*
*
*> \par Purpose:
* =============
*>
*> \verbatim
*>
*> DSKEWSYR2 performs the skew-symmetric rank 2 operation
*>
*> A := -alpha*x*y**T + alpha*y*x**T + A,
*>
*> where alpha is a scalar, x and y are n element vectors and A is an n
*> by n skew-symmetric matrix.
*> \endverbatim
*
* Arguments:
* ==========
*
*> \param[in] UPLO
*> \verbatim
*> UPLO is CHARACTER*1
*> On entry, UPLO specifies whether the upper or lower
*> triangular part of the array A is to be referenced as
*> follows:
*>
*> UPLO = 'U' or 'u' Only the upper triangular part of A
*> is to be referenced.
*>
*> UPLO = 'L' or 'l' Only the lower triangular part of A
*> is to be referenced.
*> \endverbatim
*>
*> \param[in] N
*> \verbatim
*> N is INTEGER
*> On entry, N specifies the order of the matrix A.
*> N must be at least zero.
*> \endverbatim
*>
*> \param[in] ALPHA
*> \verbatim
*> ALPHA is DOUBLE PRECISION
*> On entry, ALPHA specifies the scalar alpha.
*> \endverbatim
*>
*> \param[in] X
*> \verbatim
*> X is DOUBLE PRECISION array, dimension at least
*> ( 1 + ( n - 1 )*abs( INCX ) ).
*> Before entry, the incremented array X must contain the n
*> element vector x.
*> \endverbatim
*>
*> \param[in] INCX
*> \verbatim
*> INCX is INTEGER
*> On entry, INCX specifies the increment for the elements of
*> X. INCX must not be zero.
*> \endverbatim
*>
*> \param[in] Y
*> \verbatim
*> Y is DOUBLE PRECISION array, dimension at least
*> ( 1 + ( n - 1 )*abs( INCY ) ).
*> Before entry, the incremented array Y must contain the n
*> element vector y.
*> \endverbatim
*>
*> \param[in] INCY
*> \verbatim
*> INCY is INTEGER
*> On entry, INCY specifies the increment for the elements of
*> Y. INCY must not be zero.
*> \endverbatim
*>
*> \param[in,out] A
*> \verbatim
*> A is DOUBLE PRECISION array, dimension ( LDA, N )
*> Before entry with UPLO = 'U' or 'u', the strictly n by n
*> upper triangular part of the array A must contain the upper
*> triangular part of the skew-symmetric matrix and the leading
*> lower triangular part of A is not referenced. On exit, the
*> upper triangular part of the array A is overwritten by the
*> upper triangular part of the updated matrix.
*> Before entry with UPLO = 'L' or 'l', the strictly n by n
*> lower triangular part of the array A must contain the lower
*> triangular part of the skew-symmetric matrix and the leading
*> upper triangular part of A is not referenced. On exit, the
*> lower triangular part of the array A is overwritten by the
*> lower triangular part of the updated matrix.
*> \endverbatim
*>
*> \param[in] LDA
*> \verbatim
*> LDA is INTEGER
*> On entry, LDA specifies the first dimension of A as declared
*> in the calling (sub) program. LDA must be at least
*> max( 1, n ).
*> \endverbatim
*
* Authors:
* ========
*
*> \author Univ. of Tennessee
*> \author Univ. of California Berkeley
*> \author Univ. of Colorado Denver
*> \author NAG Ltd.
*
*> \ingroup skewher2
*
*> \par Further Details:
* =====================
*>
*> \verbatim
*>
*> Level 2 Blas routine.
*> Derived from subroutine dsyr2.
*>
*> -- Written on 6-Jul-2025.
*> Shuo Zheng, China.
*> \endverbatim
*>
* =====================================================================
SUBROUTINE DSKEWSYR2(UPLO,N,ALPHA,X,INCX,Y,INCY,A,LDA)
IMPLICIT NONE
*
* -- Reference BLAS level2 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
*
* .. Scalar Arguments ..
DOUBLE PRECISION ALPHA
INTEGER INCX,INCY,LDA,N
CHARACTER UPLO
* ..
* .. Array Arguments ..
DOUBLE PRECISION A(LDA,*),X(*),Y(*)
* ..
*
* =====================================================================
*
* .. Parameters ..
DOUBLE PRECISION ZERO
PARAMETER (ZERO=0.0D+0)
* ..
* .. Local Scalars ..
DOUBLE PRECISION TEMP1,TEMP2
INTEGER I,INFO,IX,IY,J,JX,JY,KX,KY
* ..
* .. External Functions ..
LOGICAL LSAME
EXTERNAL LSAME
* ..
* .. External Subroutines ..
EXTERNAL XERBLA
* ..
* .. Intrinsic Functions ..
INTRINSIC MAX
* ..
*
* Test the input parameters.
*
INFO = 0
IF (.NOT.LSAME(UPLO,'U') .AND. .NOT.LSAME(UPLO,'L')) THEN
INFO = 1
ELSE IF (N.LT.0) THEN
INFO = 2
ELSE IF (INCX.EQ.0) THEN
INFO = 5
ELSE IF (INCY.EQ.0) THEN
INFO = 7
ELSE IF (LDA.LT.MAX(1,N)) THEN
INFO = 9
END IF
IF (INFO.NE.0) THEN
CALL XERBLA('DSKEWSYR2 ',INFO)
RETURN
END IF
*
* Quick return if possible.
*
IF ((N.EQ.0) .OR. (ALPHA.EQ.ZERO)) RETURN
*
* Set up the start points in X and Y if the increments are not both
* unity.
*
IF ((INCX.NE.1) .OR. (INCY.NE.1)) THEN
IF (INCX.GT.0) THEN
KX = 1
ELSE
KX = 1 - (N-1)*INCX
END IF
IF (INCY.GT.0) THEN
KY = 1
ELSE
KY = 1 - (N-1)*INCY
END IF
JX = KX
JY = KY
END IF
*
* Start the operations. In this version the elements of A are
* accessed sequentially with one pass through the triangular part
* of A.
*
IF (LSAME(UPLO,'U')) THEN
*
* Form A when A is stored in the upper triangle.
*
IF ((INCX.EQ.1) .AND. (INCY.EQ.1)) THEN
DO 20 J = 1,N
IF ((X(J).NE.ZERO) .OR. (Y(J).NE.ZERO)) THEN
TEMP1 = ALPHA*Y(J)
TEMP2 = ALPHA*X(J)
DO 10 I = 1,J-1
A(I,J) = A(I,J) - X(I)*TEMP1 + Y(I)*TEMP2
10 CONTINUE
END IF
20 CONTINUE
ELSE
DO 40 J = 1,N
IF ((X(JX).NE.ZERO) .OR. (Y(JY).NE.ZERO)) THEN
TEMP1 = ALPHA*Y(JY)
TEMP2 = ALPHA*X(JX)
IX = KX
IY = KY
DO 30 I = 1,J-1
A(I,J) = A(I,J) - X(IX)*TEMP1 + Y(IY)*TEMP2
IX = IX + INCX
IY = IY + INCY
30 CONTINUE
END IF
JX = JX + INCX
JY = JY + INCY
40 CONTINUE
END IF
ELSE
*
* Form A when A is stored in the lower triangle.
*
IF ((INCX.EQ.1) .AND. (INCY.EQ.1)) THEN
DO 60 J = 1,N
IF ((X(J).NE.ZERO) .OR. (Y(J).NE.ZERO)) THEN
TEMP1 = ALPHA*Y(J)
TEMP2 = ALPHA*X(J)
DO 50 I = J+1,N
A(I,J) = A(I,J) - X(I)*TEMP1 + Y(I)*TEMP2
50 CONTINUE
END IF
60 CONTINUE
ELSE
DO 80 J = 1,N
IF ((X(JX).NE.ZERO) .OR. (Y(JY).NE.ZERO)) THEN
TEMP1 = ALPHA*Y(JY)
TEMP2 = ALPHA*X(JX)
IX = JX + INCX
IY = JY + INCY
DO 70 I = J+1,N
A(I,J) = A(I,J) - X(IX)*TEMP1 + Y(IY)*TEMP2
IX = IX + INCX
IY = IY + INCY
70 CONTINUE
END IF
JX = JX + INCX
JY = JY + INCY
80 CONTINUE
END IF
END IF
*
RETURN
*
* End of DSKEWSYR2
*
END
+395
View File
@@ -0,0 +1,395 @@
*> \brief \b DSKEWSYR2K
*
* =========== DOCUMENTATION ===========
*
* Online html documentation available at
* http://www.netlib.org/lapack/explore-html/
*
* Definition:
* ===========
*
* SUBROUTINE DSKEWSYR2K(UPLO,TRANS,N,K,ALPHA,A,LDA,B,LDB,BETA,C,LDC)
*
* .. Scalar Arguments ..
* DOUBLE PRECISION ALPHA,BETA
* INTEGER K,LDA,LDB,LDC,N
* CHARACTER TRANS,UPLO
* ..
* .. Array Arguments ..
* DOUBLE PRECISION A(LDA,*),B(LDB,*),C(LDC,*)
* ..
*
*
*> \par Purpose:
* =============
*>
*> \verbatim
*>
*> DSKEWSYR2K performs one of the skew-symmetric rank 2k operations
*>
*> C := -alpha*A*B**T + alpha*B*A**T + beta*C,
*>
*> or
*>
*> C := -alpha*A**T*B + alpha*B**T*A + beta*C,
*>
*> where alpha and beta are scalars, C is an n by n skew-symmetric matrix
*> and A and B are n by k matrices in the first case and k by n
*> matrices in the second case.
*> \endverbatim
*
* Arguments:
* ==========
*
*> \param[in] UPLO
*> \verbatim
*> UPLO is CHARACTER*1
*> On entry, UPLO specifies whether the upper or lower
*> triangular part of the array C is to be referenced as
*> follows:
*>
*> UPLO = 'U' or 'u' Only the upper triangular part of C
*> is to be referenced.
*>
*> UPLO = 'L' or 'l' Only the lower triangular part of C
*> is to be referenced.
*> \endverbatim
*>
*> \param[in] TRANS
*> \verbatim
*> TRANS is CHARACTER*1
*> On entry, TRANS specifies the operation to be performed as
*> follows:
*>
*> TRANS = 'N' or 'n' C := -alpha*A*B**T + alpha*B*A**T +
*> beta*C.
*>
*> TRANS = 'T' or 't' C := -alpha*A**T*B + alpha*B**T*A +
*> beta*C.
*>
*> TRANS = 'C' or 'c' C := -alpha*A**T*B + alpha*B**T*A +
*> beta*C.
*> \endverbatim
*>
*> \param[in] N
*> \verbatim
*> N is INTEGER
*> On entry, N specifies the order of the matrix C. N must be
*> at least zero.
*> \endverbatim
*>
*> \param[in] K
*> \verbatim
*> K is INTEGER
*> On entry with TRANS = 'N' or 'n', K specifies the number
*> of columns of the matrices A and B, and on entry with
*> TRANS = 'T' or 't' or 'C' or 'c', K specifies the number
*> of rows of the matrices A and B. K must be at least zero.
*> \endverbatim
*>
*> \param[in] ALPHA
*> \verbatim
*> ALPHA is DOUBLE PRECISION.
*> On entry, ALPHA specifies the scalar alpha.
*> \endverbatim
*>
*> \param[in] A
*> \verbatim
*> A is DOUBLE PRECISION array, dimension ( LDA, ka ), where ka is
*> k when TRANS = 'N' or 'n', and is n otherwise.
*> Before entry with TRANS = 'N' or 'n', the leading n by k
*> part of the array A must contain the matrix A, otherwise
*> the leading k by n part of the array A must contain the
*> matrix A.
*> \endverbatim
*>
*> \param[in] LDA
*> \verbatim
*> LDA is INTEGER
*> On entry, LDA specifies the first dimension of A as declared
*> in the calling (sub) program. When TRANS = 'N' or 'n'
*> then LDA must be at least max( 1, n ), otherwise LDA must
*> be at least max( 1, k ).
*> \endverbatim
*>
*> \param[in] B
*> \verbatim
*> B is DOUBLE PRECISION array, dimension ( LDB, kb ), where kb is
*> k when TRANS = 'N' or 'n', and is n otherwise.
*> Before entry with TRANS = 'N' or 'n', the leading n by k
*> part of the array B must contain the matrix B, otherwise
*> the leading k by n part of the array B must contain the
*> matrix B.
*> \endverbatim
*>
*> \param[in] LDB
*> \verbatim
*> LDB is INTEGER
*> On entry, LDB specifies the first dimension of B as declared
*> in the calling (sub) program. When TRANS = 'N' or 'n'
*> then LDB must be at least max( 1, n ), otherwise LDB must
*> be at least max( 1, k ).
*> \endverbatim
*>
*> \param[in] BETA
*> \verbatim
*> BETA is DOUBLE PRECISION.
*> On entry, BETA specifies the scalar beta.
*> \endverbatim
*>
*> \param[in,out] C
*> \verbatim
*> C is DOUBLE PRECISION array, dimension ( LDC, N )
*> Before entry with UPLO = 'U' or 'u', the strictly n by n
*> upper triangular part of the array C must contain the upper
*> triangular part of the skew-symmetric matrix and the leading
*> lower triangular part of C is not referenced. On exit, the
*> upper triangular part of the array C is overwritten by the
*> upper triangular part of the updated matrix.
*> Before entry with UPLO = 'L' or 'l', the strictly n by n
*> lower triangular part of the array C must contain the lower
*> triangular part of the skew-symmetric matrix and the leading
*> upper triangular part of C is not referenced. On exit, the
*> lower triangular part of the array C is overwritten by the
*> lower triangular part of the updated matrix.
*> \endverbatim
*>
*> \param[in] LDC
*> \verbatim
*> LDC is INTEGER
*> On entry, LDC specifies the first dimension of C as declared
*> in the calling (sub) program. LDC must be at least
*> max( 1, n ).
*> \endverbatim
*
* Authors:
* ========
*
*> \author Univ. of Tennessee
*> \author Univ. of California Berkeley
*> \author Univ. of Colorado Denver
*> \author NAG Ltd.
*
*> \ingroup skewher2k
*
*> \par Further Details:
* =====================
*>
*> \verbatim
*>
*> Level 3 Blas routine.
*> Derived from subroutine dsyr2k.
*>
*> -- Written on 6-Jul-2025.
*> Shuo Zheng, China.
*> \endverbatim
*>
* =====================================================================
SUBROUTINE DSKEWSYR2K(UPLO,TRANS,N,K,ALPHA,A,LDA,B,
+ LDB,BETA,C,LDC)
IMPLICIT NONE
*
* -- Reference BLAS level3 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
*
* .. Scalar Arguments ..
DOUBLE PRECISION ALPHA,BETA
INTEGER K,LDA,LDB,LDC,N
CHARACTER TRANS,UPLO
* ..
* .. Array Arguments ..
DOUBLE PRECISION A(LDA,*),B(LDB,*),C(LDC,*)
* ..
*
* =====================================================================
*
* .. External Functions ..
LOGICAL LSAME
EXTERNAL LSAME
* ..
* .. External Subroutines ..
EXTERNAL XERBLA
* ..
* .. Intrinsic Functions ..
INTRINSIC MAX
* ..
* .. Local Scalars ..
DOUBLE PRECISION TEMP1,TEMP2
INTEGER I,INFO,J,L,NROWA
LOGICAL UPPER
* ..
* .. Parameters ..
DOUBLE PRECISION ONE,ZERO
PARAMETER (ONE=1.0D+0,ZERO=0.0D+0)
* ..
*
* Test the input parameters.
*
IF (LSAME(TRANS,'N')) THEN
NROWA = N
ELSE
NROWA = K
END IF
UPPER = LSAME(UPLO,'U')
*
INFO = 0
IF ((.NOT.UPPER) .AND. (.NOT.LSAME(UPLO,'L'))) THEN
INFO = 1
ELSE IF ((.NOT.LSAME(TRANS,'N')) .AND.
+ (.NOT.LSAME(TRANS,'T')) .AND.
+ (.NOT.LSAME(TRANS,'C'))) THEN
INFO = 2
ELSE IF (N.LT.0) THEN
INFO = 3
ELSE IF (K.LT.0) THEN
INFO = 4
ELSE IF (LDA.LT.MAX(1,NROWA)) THEN
INFO = 7
ELSE IF (LDB.LT.MAX(1,NROWA)) THEN
INFO = 9
ELSE IF (LDC.LT.MAX(1,N)) THEN
INFO = 12
END IF
IF (INFO.NE.0) THEN
CALL XERBLA('DSKEWSYR2K',INFO)
RETURN
END IF
*
* Quick return if possible.
*
IF ((N.EQ.0) .OR. (((ALPHA.EQ.ZERO).OR.
+ (K.EQ.0)).AND. (BETA.EQ.ONE))) RETURN
*
* And when alpha.eq.zero.
*
IF (ALPHA.EQ.ZERO) THEN
IF (UPPER) THEN
IF (BETA.EQ.ZERO) THEN
DO 20 J = 1,N
DO 10 I = 1,J-1
C(I,J) = ZERO
10 CONTINUE
20 CONTINUE
ELSE
DO 40 J = 1,N
DO 30 I = 1,J-1
C(I,J) = BETA*C(I,J)
30 CONTINUE
40 CONTINUE
END IF
ELSE
IF (BETA.EQ.ZERO) THEN
DO 60 J = 1,N
DO 50 I = J+1,N
C(I,J) = ZERO
50 CONTINUE
60 CONTINUE
ELSE
DO 80 J = 1,N
DO 70 I = J+1,N
C(I,J) = BETA*C(I,J)
70 CONTINUE
80 CONTINUE
END IF
END IF
RETURN
END IF
*
* Start the operations.
*
IF (LSAME(TRANS,'N')) THEN
*
* Form C := alpha*A*B**T + alpha*B*A**T + C.
*
IF (UPPER) THEN
DO 130 J = 1,N
IF (BETA.EQ.ZERO) THEN
DO 90 I = 1,J-1
C(I,J) = ZERO
90 CONTINUE
ELSE IF (BETA.NE.ONE) THEN
DO 100 I = 1,J-1
C(I,J) = BETA*C(I,J)
100 CONTINUE
END IF
DO 120 L = 1,K
IF ((A(J,L).NE.ZERO) .OR. (B(J,L).NE.ZERO)) THEN
TEMP1 = ALPHA*B(J,L)
TEMP2 = ALPHA*A(J,L)
DO 110 I = 1,J-1
C(I,J) = C(I,J) - A(I,L)*TEMP1 +
+ B(I,L)*TEMP2
110 CONTINUE
END IF
120 CONTINUE
130 CONTINUE
ELSE
DO 180 J = 1,N
IF (BETA.EQ.ZERO) THEN
DO 140 I = J+1,N
C(I,J) = ZERO
140 CONTINUE
ELSE IF (BETA.NE.ONE) THEN
DO 150 I = J+1,N
C(I,J) = BETA*C(I,J)
150 CONTINUE
END IF
DO 170 L = 1,K
IF ((A(J,L).NE.ZERO) .OR. (B(J,L).NE.ZERO)) THEN
TEMP1 = ALPHA*B(J,L)
TEMP2 = ALPHA*A(J,L)
DO 160 I = J+1,N
C(I,J) = C(I,J) - A(I,L)*TEMP1 +
+ B(I,L)*TEMP2
160 CONTINUE
END IF
170 CONTINUE
180 CONTINUE
END IF
ELSE
*
* Form C := alpha*A**T*B + alpha*B**T*A + C.
*
IF (UPPER) THEN
DO 210 J = 1,N
DO 200 I = 1,J-1
TEMP1 = ZERO
TEMP2 = ZERO
DO 190 L = 1,K
TEMP1 = TEMP1 + A(L,I)*B(L,J)
TEMP2 = TEMP2 + B(L,I)*A(L,J)
190 CONTINUE
IF (BETA.EQ.ZERO) THEN
C(I,J) = -ALPHA*TEMP1 + ALPHA*TEMP2
ELSE
C(I,J) = BETA*C(I,J) - ALPHA*TEMP1 +
+ ALPHA*TEMP2
END IF
200 CONTINUE
210 CONTINUE
ELSE
DO 240 J = 1,N
DO 230 I = J+1,N
TEMP1 = ZERO
TEMP2 = ZERO
DO 220 L = 1,K
TEMP1 = TEMP1 + A(L,I)*B(L,J)
TEMP2 = TEMP2 + B(L,I)*A(L,J)
220 CONTINUE
IF (BETA.EQ.ZERO) THEN
C(I,J) = -ALPHA*TEMP1 + ALPHA*TEMP2
ELSE
C(I,J) = BETA*C(I,J) - ALPHA*TEMP1 +
+ ALPHA*TEMP2
END IF
230 CONTINUE
240 CONTINUE
END IF
END IF
*
RETURN
*
* End of DSKEWSYR2K
*
END
+1
View File
@@ -144,6 +144,7 @@
*>
* =====================================================================
SUBROUTINE DSPMV(UPLO,N,ALPHA,AP,X,INCX,BETA,Y,INCY)
IMPLICIT NONE
*
* -- Reference BLAS level2 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -124,6 +124,7 @@
*>
* =====================================================================
SUBROUTINE DSPR(UPLO,N,ALPHA,X,INCX,AP)
IMPLICIT NONE
*
* -- Reference BLAS level2 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -139,6 +139,7 @@
*>
* =====================================================================
SUBROUTINE DSPR2(UPLO,N,ALPHA,X,INCX,Y,INCY,AP)
IMPLICIT NONE
*
* -- Reference BLAS level2 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -79,6 +79,7 @@
*>
* =====================================================================
SUBROUTINE DSWAP(N,DX,INCX,DY,INCY)
IMPLICIT NONE
*
* -- Reference BLAS level1 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -186,6 +186,7 @@
*>
* =====================================================================
SUBROUTINE DSYMM(SIDE,UPLO,M,N,ALPHA,A,LDA,B,LDB,BETA,C,LDC)
IMPLICIT NONE
*
* -- Reference BLAS level3 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -149,6 +149,7 @@
*>
* =====================================================================
SUBROUTINE DSYMV(UPLO,N,ALPHA,A,LDA,X,INCX,BETA,Y,INCY)
IMPLICIT NONE
*
* -- Reference BLAS level2 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -129,6 +129,7 @@
*>
* =====================================================================
SUBROUTINE DSYR(UPLO,N,ALPHA,X,INCX,A,LDA)
IMPLICIT NONE
*
* -- Reference BLAS level2 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -144,6 +144,7 @@
*>
* =====================================================================
SUBROUTINE DSYR2(UPLO,N,ALPHA,X,INCX,Y,INCY,A,LDA)
IMPLICIT NONE
*
* -- Reference BLAS level2 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -189,6 +189,7 @@
*>
* =====================================================================
SUBROUTINE DSYR2K(UPLO,TRANS,N,K,ALPHA,A,LDA,B,LDB,BETA,C,LDC)
IMPLICIT NONE
*
* -- Reference BLAS level3 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -166,6 +166,7 @@
*>
* =====================================================================
SUBROUTINE DSYRK(UPLO,TRANS,N,K,ALPHA,A,LDA,BETA,C,LDC)
IMPLICIT NONE
*
* -- Reference BLAS level3 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -183,6 +183,7 @@
*>
* =====================================================================
SUBROUTINE DTBMV(UPLO,TRANS,DIAG,N,K,A,LDA,X,INCX)
IMPLICIT NONE
*
* -- Reference BLAS level2 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -186,6 +186,7 @@
*>
* =====================================================================
SUBROUTINE DTBSV(UPLO,TRANS,DIAG,N,K,A,LDA,X,INCX)
IMPLICIT NONE
*
* -- Reference BLAS level2 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -139,6 +139,7 @@
*>
* =====================================================================
SUBROUTINE DTPMV(UPLO,TRANS,DIAG,N,AP,X,INCX)
IMPLICIT NONE
*
* -- Reference BLAS level2 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -141,6 +141,7 @@
*>
* =====================================================================
SUBROUTINE DTPSV(UPLO,TRANS,DIAG,N,AP,X,INCX)
IMPLICIT NONE
*
* -- Reference BLAS level2 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -174,6 +174,7 @@
*>
* =====================================================================
SUBROUTINE DTRMM(SIDE,UPLO,TRANSA,DIAG,M,N,ALPHA,A,LDA,B,LDB)
IMPLICIT NONE
*
* -- Reference BLAS level3 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -144,6 +144,7 @@
*>
* =====================================================================
SUBROUTINE DTRMV(UPLO,TRANS,DIAG,N,A,LDA,X,INCX)
IMPLICIT NONE
*
* -- Reference BLAS level2 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -178,6 +178,7 @@
*>
* =====================================================================
SUBROUTINE DTRSM(SIDE,UPLO,TRANSA,DIAG,M,N,ALPHA,A,LDA,B,LDB)
IMPLICIT NONE
*
* -- Reference BLAS level3 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -140,6 +140,7 @@
*
* =====================================================================
SUBROUTINE DTRSV(UPLO,TRANS,DIAG,N,A,LDA,X,INCX)
IMPLICIT NONE
*
* -- Reference BLAS level2 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -69,6 +69,7 @@
*>
* =====================================================================
DOUBLE PRECISION FUNCTION DZASUM(N,ZX,INCX)
IMPLICIT NONE
*
* -- Reference BLAS level1 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+3 -2
View File
@@ -86,11 +86,12 @@
!> \endverbatim
!>
! =====================================================================
function DZNRM2( n, x, incx )
function DZNRM2( n, x, incx )
implicit none
integer, parameter :: wp = kind(1.d0)
real(wp) :: DZNRM2
!
! -- Reference BLAS level1 routine (version 3.9.1) --
! -- Reference BLAS level1 routine --
! -- Reference BLAS is a software package provided by Univ. of Tennessee, --
! -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
! March 2021
+193
View File
@@ -0,0 +1,193 @@
!> \brief \b ICAMAX
!
! =========== DOCUMENTATION ===========
!
! Online html documentation available at
! http://www.netlib.org/lapack/explore-html/
!
! Definition:
! ===========
!
! INTEGER FUNCTION ICAMAX(N,X,INCX)
!
! .. Scalar Arguments ..
! INTEGER INCX,N
! ..
! .. Array Arguments ..
! COMPLEX X(*)
! ..
!
!
!> \par Purpose:
! =============
!>
!> \verbatim
!>
!> ICAMAX finds the index of the first element having maximum |Re(.)| + |Im(.)|
!> \endverbatim
!
! Arguments:
! ==========
!
!> \param[in] N
!> \verbatim
!> N is INTEGER
!> number of elements in input vector(s)
!> \endverbatim
!>
!> \param[in] X
!> \verbatim
!> X is COMPLEX array, dimension ( 1 + ( N - 1 )*abs( INCX ) )
!> \endverbatim
!>
!> \param[in] INCX
!> \verbatim
!> INCX is INTEGER
!> storage spacing between elements of X
!> \endverbatim
!
! Authors:
! ========
!
!> James Demmel, University of California Berkeley, USA
!> Weslley Pereira, National Renewable Energy Laboratory, USA
!
!> \ingroup iamax
!
!> \par Further Details:
! =====================
!>
!> \verbatim
!>
!> James Demmel et al. Proposed Consistent Exception Handling for the BLAS and
!> LAPACK, 2022 (https://arxiv.org/abs/2207.09281).
!>
!> \endverbatim
!>
! =====================================================================
integer function icamax(n, x, incx)
implicit none
integer, parameter :: wp = kind(1.e0)
!
! -- Reference BLAS level1 routine --
! -- Reference BLAS is a software package provided by Univ. of Tennessee, --
! -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
!
! .. Constants ..
real(wp), parameter :: hugeval = huge(0.0_wp)
!
! .. Scalar Arguments ..
integer :: n, incx
!
! .. Array Arguments ..
complex(wp) :: x(*)
! ..
! .. Local Scalars ..
integer :: i, j, ix, jx
real(wp) :: val, smax
logical :: scaledsmax
! ..
! .. Intrinsic Functions ..
intrinsic :: abs, aimag, huge, real
!
! Quick return if possible
!
icamax = 0
if (n < 1 .or. incx < 1) return
!
icamax = 1
if (n == 1) return
!
icamax = 0
scaledsmax = .false.
smax = -1
!
! scaledsmax = .true. indicates that x(icamax) is finite but
! abs(real(x(icamax))) + abs(aimag(x(icamax))) overflows
!
if (incx == 1) then
! code for increment equal to 1
do i = 1, n
if (x(i) /= x(i)) then
! return when first NaN found
icamax = i
return
elseif (abs(real(x(i))) > hugeval .or. abs(aimag(x(i))) > hugeval) then
! keep looking for first NaN
do j = i+1, n
if (x(j) /= x(j)) then
! return when first NaN found
icamax = j
return
endif
enddo
! record location of first Inf
icamax = i
return
else ! still no Inf found yet
if (.not. scaledsmax) then
! no abs(real(x(i))) + abs(aimag(x(i))) = Inf yet
val = abs(real(x(i))) + abs(aimag(x(i)))
if (val > hugeval) then
scaledsmax = .true.
smax = 0.25*abs(real(x(i))) + 0.25*abs(aimag(x(i)))
icamax = i
elseif (val > smax) then ! everything finite so far
smax = val
icamax = i
endif
else ! scaledsmax
val = 0.25*abs(real(x(i))) + 0.25*abs(aimag(x(i)))
if (val > smax) then
smax = val
icamax = i
endif
endif
endif
end do
else
! code for increment not equal to 1
ix = 1
do i = 1, n
if (x(ix) /= x(ix)) then
! return when first NaN found
icamax = i
return
elseif (abs(real(x(ix))) > hugeval .or. abs(aimag(x(ix))) > hugeval) then
! keep looking for first NaN
jx = ix + incx
do j = i+1, n
if (x(jx) /= x(jx)) then
! return when first NaN found
icamax = j
return
endif
jx = jx + incx
enddo
! record location of first Inf
icamax = i
return
else ! still no Inf found yet
if (.not. scaledsmax) then
! no abs(real(x(ix))) + abs(aimag(x(ix))) = Inf yet
val = abs(real(x(ix))) + abs(aimag(x(ix)))
if (val > hugeval) then
scaledsmax = .true.
smax = 0.25*abs(real(x(ix))) + 0.25*abs(aimag(x(ix)))
icamax = i
elseif (val > smax) then ! everything finite so far
smax = val
icamax = i
endif
else ! scaledsmax
val = 0.25*abs(real(x(ix))) + 0.25*abs(aimag(x(ix)))
if (val > smax) then
smax = val
icamax = i
endif
endif
endif
ix = ix + incx
end do
endif
end
+1
View File
@@ -68,6 +68,7 @@
*>
* =====================================================================
INTEGER FUNCTION IDAMAX(N,DX,INCX)
IMPLICIT NONE
*
* -- Reference BLAS level1 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -68,6 +68,7 @@
*>
* =====================================================================
INTEGER FUNCTION ISAMAX(N,SX,INCX)
IMPLICIT NONE
*
* -- Reference BLAS level1 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+193
View File
@@ -0,0 +1,193 @@
!> \brief \b IZAMAX
!
! =========== DOCUMENTATION ===========
!
! Online html documentation available at
! http://www.netlib.org/lapack/explore-html/
!
! Definition:
! ===========
!
! INTEGER FUNCTION IZAMAX(N,X,INCX)
!
! .. Scalar Arguments ..
! INTEGER INCX,N
! ..
! .. Array Arguments ..
! DOUBLE COMPLEX X(*)
! ..
!
!
!> \par Purpose:
! =============
!>
!> \verbatim
!>
!> IZAMAX finds the index of the first element having maximum |Re(.)| + |Im(.)|
!> \endverbatim
!
! Arguments:
! ==========
!
!> \param[in] N
!> \verbatim
!> N is INTEGER
!> number of elements in input vector(s)
!> \endverbatim
!>
!> \param[in] X
!> \verbatim
!> X is DOUBLE COMPLEX array, dimension ( 1 + ( N - 1 )*abs( INCX ) )
!> \endverbatim
!>
!> \param[in] INCX
!> \verbatim
!> INCX is INTEGER
!> storage spacing between elements of X
!> \endverbatim
!
! Authors:
! ========
!
!> James Demmel, University of California Berkeley, USA
!> Weslley Pereira, National Renewable Energy Laboratory, USA
!
!> \ingroup iamax
!
!> \par Further Details:
! =====================
!>
!> \verbatim
!>
!> James Demmel et al. Proposed Consistent Exception Handling for the BLAS and
!> LAPACK, 2022 (https://arxiv.org/abs/2207.09281).
!>
!> \endverbatim
!>
! =====================================================================
integer function izamax(n, x, incx)
implicit none
integer, parameter :: wp = kind(1.d0)
!
! -- Reference BLAS level1 routine --
! -- Reference BLAS is a software package provided by Univ. of Tennessee, --
! -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
!
! .. Constants ..
real(wp), parameter :: hugeval = huge(0.0_wp)
!
! .. Scalar Arguments ..
integer :: n, incx
!
! .. Array Arguments ..
complex(wp) :: x(*)
! ..
! .. Local Scalars ..
integer :: i, j, ix, jx
real(wp) :: val, smax
logical :: scaledsmax
! ..
! .. Intrinsic Functions ..
intrinsic :: abs, dimag, huge, real
!
! Quick return if possible
!
izamax = 0
if (n < 1 .or. incx < 1) return
!
izamax = 1
if (n == 1) return
!
izamax = 0
scaledsmax = .false.
smax = -1
!
! scaledsmax = .true. indicates that x(izamax) is finite but
! abs(real(x(izamax))) + abs(dimag(x(izamax))) overflows
!
if (incx == 1) then
! code for increment equal to 1
do i = 1, n
if (x(i) /= x(i)) then
! return when first NaN found
izamax = i
return
elseif (abs(real(x(i))) > hugeval .or. abs(dimag(x(i))) > hugeval) then
! keep looking for first NaN
do j = i+1, n
if (x(j) /= x(j)) then
! return when first NaN found
izamax = j
return
endif
enddo
! record location of first Inf
izamax = i
return
else ! still no Inf found yet
if (.not. scaledsmax) then
! no abs(real(x(i))) + abs(dimag(x(i))) = Inf yet
val = abs(real(x(i))) + abs(dimag(x(i)))
if (val > hugeval) then
scaledsmax = .true.
smax = 0.25*abs(real(x(i))) + 0.25*abs(dimag(x(i)))
izamax = i
elseif (val > smax) then ! everything finite so far
smax = val
izamax = i
endif
else ! scaledsmax
val = 0.25*abs(real(x(i))) + 0.25*abs(dimag(x(i)))
if (val > smax) then
smax = val
izamax = i
endif
endif
endif
end do
else
! code for increment not equal to 1
ix = 1
do i = 1, n
if (x(ix) /= x(ix)) then
! return when first NaN found
izamax = i
return
elseif (abs(real(x(ix))) > hugeval .or. abs(dimag(x(ix))) > hugeval) then
! keep looking for first NaN
jx = ix + incx
do j = i+1, n
if (x(jx) /= x(jx)) then
! return when first NaN found
izamax = j
return
endif
jx = jx + incx
enddo
! record location of first Inf
izamax = i
return
else ! still no Inf found yet
if (.not. scaledsmax) then
! no abs(real(x(ix))) + abs(dimag(x(ix))) = Inf yet
val = abs(real(x(ix))) + abs(dimag(x(ix)))
if (val > hugeval) then
scaledsmax = .true.
smax = 0.25*abs(real(x(ix))) + 0.25*abs(dimag(x(ix)))
izamax = i
elseif (val > smax) then ! everything finite so far
smax = val
izamax = i
endif
else ! scaledsmax
val = 0.25*abs(real(x(ix))) + 0.25*abs(dimag(x(ix)))
if (val > smax) then
smax = val
izamax = i
endif
endif
endif
ix = ix + incx
end do
endif
end
+1
View File
@@ -50,6 +50,7 @@
*
* =====================================================================
LOGICAL FUNCTION LSAME(CA,CB)
IMPLICIT NONE
*
* -- Reference BLAS level1 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+1
View File
@@ -69,6 +69,7 @@
*>
* =====================================================================
REAL FUNCTION SASUM(N,SX,INCX)
IMPLICIT NONE
*
* -- Reference BLAS level1 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
+148
View File
@@ -0,0 +1,148 @@
*> \brief \b SAXPBY
*
* =========== DOCUMENTATION ===========
*
* Online html documentation available at
* http://www.netlib.org/lapack/explore-html/
*
* Definition:
* ===========
*
* SUBROUTINE SAXPBY(N,SA,SX,INCX,SB,SY,INCY)
*
* .. Scalar Arguments ..
* REAL SA,SB
* INTEGER INCX,INCY,N
* ..
* .. Array Arguments ..
* REAL SX(*),SY(*)
* ..
*
*
*> \par Purpose:
* =============
*>
*> \verbatim
*>
*> SAXPBY constant times a vector plus constant times a vector.
*>
*> Y = ALPHA * X + BETA * Y
*>
*> \endverbatim
*
* Arguments:
* ==========
*
*> \param[in] N
*> \verbatim
*> N is INTEGER
*> number of elements in input vector(s)
*> \endverbatim
*>
*> \param[in] SA
*> \verbatim
*> SA is REAL
*> On entry, SA specifies the scalar alpha.
*> \endverbatim
*>
*> \param[in] SX
*> \verbatim
*> SX is REAL array, dimension ( 1 + ( N - 1 )*abs( INCX ) )
*> \endverbatim
*>
*> \param[in] INCX
*> \verbatim
*> INCX is INTEGER
*> storage spacing between elements of SX
*> \endverbatim
*>
*> \param[in] SB
*> \verbatim
*> SB is REAL
*> On entry, SB specifies the scalar beta.
*> \endverbatim
*>
*> \param[in,out] SY
*> \verbatim
*> SY is REAL array, dimension ( 1 + ( N - 1 )*abs( INCY ) )
*> \endverbatim
*>
*> \param[in] INCY
*> \verbatim
*> INCY is INTEGER
*> storage spacing between elements of SY
*> \endverbatim
*
* Authors:
* ========
*
*> \author Univ. of Tennessee
*> \author Univ. of California Berkeley
*> \author Univ. of Colorado Denver
*> \author NAG Ltd.
*> \author Martin Koehler, MPI Magdeburg
*
*> \ingroup axpby
*
* =====================================================================
SUBROUTINE SAXPBY(N,SA,SX,INCX,SB,SY,INCY)
IMPLICIT NONE
*
* -- Reference BLAS level1 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
*
* .. Scalar Arguments ..
REAL SA,SB
INTEGER INCX,INCY,N
* ..
* .. Array Arguments ..
REAL SX(*),SY(*)
* ..
* .. External Subroutines ..
EXTERNAL SSCAL
*
* =====================================================================
*
* .. Local Scalars ..
INTEGER I,IX,IY,M,MP1
* ..
* .. Intrinsic Functions ..
INTRINSIC MOD
* ..
IF (N.LE.0) RETURN
* Scale if SA.EQ.0
IF (SA.EQ.0.0E0 .AND. SB.NE.0.0E0) THEN
CALL SSCAL(N, SB, SY, INCY)
RETURN
END IF
IF (INCX.EQ.1 .AND. INCY.EQ.1) THEN
*
* code for both increments equal to 1
*
DO I = 1,N
SY(I) = SB*SY(I) + SA*SX(I)
END DO
ELSE
*
* code for unequal increments or equal increments
* not equal to 1
*
IX = 1
IY = 1
IF (INCX.LT.0) IX = (-N+1)*INCX + 1
IF (INCY.LT.0) IY = (-N+1)*INCY + 1
DO I = 1,N
SY(IY) = SB*SY(IY) + SA*SX(IX)
IX = IX + INCX
IY = IY + INCY
END DO
END IF
RETURN
*
* End of SAXPBY
*
END
+1
View File
@@ -86,6 +86,7 @@
*>
* =====================================================================
SUBROUTINE SAXPY(N,SA,SX,INCX,SY,INCY)
IMPLICIT NONE
*
* -- Reference BLAS level1 routine --
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --

Some files were not shown because too many files have changed in this diff Show More