TABLE OF CONTENTS
- ABINIT/m_blas
- m_blas/blas_cholesky_ortho_dpc
- m_blas/blas_cholesky_ortho_spc
- m_blas/sqmat_iconjgtrans_dpc
- m_blas/sqmat_iconjgtrans_spc
- m_blas/sqmat_itranspose_dp
- m_blas/sqmat_itranspose_dpc
- m_blas/sqmat_itranspose_sp
- m_blas/sqmat_itranspose_spc
- m_blas/sqmat_oconjgtrans_dpc
- m_blas/sqmat_oconjgtrans_spc
- m_blas/sqmat_otranspose_dp
- m_blas/sqmat_otranspose_dpc
- m_blas/sqmat_otranspose_sp
- m_blas/sqmat_otranspose_spc
ABINIT/m_blas [ Modules ]
NAME
m_blas
FUNCTION
This module defines interfaces for overloading BLAS routines. whose goal is twofold. On one hand, using generic interfaces renders the code more readable, especially when the routine can be compiled with different precision type (single-precision or double precision as done for example in the GW code) On the other hand, the generic interfaces defined here introduce a programming layer that can be exploited for interfacing non-standard libraries such as for example CUBLAS routines for GPU computations.
COPYRIGHT
Copyright (C) 1992-2018 ABINIT group (MG) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt . For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
NOTES
The convention about names of interfaced routine is: x<name>, where <name> is usually equal to the name of the standard routine without the first character specifying the type and kind. The full list of names is reported below. BLAS procedures interfaced in this module are marked with an asterisk. A complete list of possible overloaded interfaces is provided as guide for future additions. ================ ==== BLAS 1 ==== ================ FUNCTION idamax isamax icamax izamax ---> XIAMAX(n,dx,incx) * FUNCTION snrm2 dnrm2 scnrm2 dznmr2 ---> XNRM2(n,x,incx) FUNCTION sasum dasum scasum dzasum ---> XASUM(n,x,incx) * FUNCTION cdotu zdotu ---> XDOTU(n,x,incx,y,incy) * FUNCTION cdotc zdotc ---> XDOTC(n,x,incx,y,incy) FUNCTION sdot ddot ---> XDOT(n,x,incx,y,incy) FUNCTION sdsdot sdot ---> XDSDOT(n,x,incx,y,incy) SUBROUTINE saxpy daxpy caxpy zaxpy ---> XAXPY(n,ca,cx,incx,cy,incy) * SUBROUTINE scopy dcopy ccopy zcopy ---> XCOPY(n,cx,incx,cy,incy) SUBROUTINE srotg drotg crotg zrotg ---> XROTG(a,b,c,s) SUBROUTINE srot drot csrot zdrot ---> XROT(n,x,incx,y,incy,c,s) * SUBROUTINE sscal dscal cscal zscal csscal zdscal ---> XSCAL(n,a,x,incx) SUBROUTINE sswap dswap cswap zswap ---> XSWAP(n,x,incx,y,incy) ================ ==== BLAS 2 ==== ================ SUBROUTINE sgbmv dgbmv cgbmv zgbmv ---> XGBMV(trans,m,kl,ku,n,alpha,A,lda,X,incx,beta,Y,incy) * SUBROUTINE sgemv dgemv cgemv zgemv ---> XGEMV(trans,m,n,alpha,A,lda,X,incx,beta,Y,incy) * SUBROUTINE cgerc zgerc ---> XGERC(m,n,alpha,x,incx,y,incy,A,lda) SUBROUTINE cgeru zgeru ---> XGERU(m,n,alpha,x,incx,y,incy,A,lda) SUBROUTINE chbmv zhbmv ---> XHBMV(uplo,n,k,alpha,A,lda,X,incx,beta,Y,incy) SUBROUTINE chemv zhemv ---> XHEMV(uplo,n,alpha,A,lda,X,incx,beta,Y,incy) SUBROUTINE cher zher ---> XHER(uplo,n,alpha,x,incx,A,lda) SUBROUTINE cher2 zher2 ---> XHER2(uplo,n,alpha,x,incx,y,incy,A,lda) SUBROUTINE chpr zhpr ---> XHPR(uplo,n,alpha,x,incx,AP) SUBROUTINE chpr2 zhpr2 ---> XHPR2(uplo,n,alpha,x,incx,y,incy,AP) SUBROUTINE chpmv zhpmv ---> XHPMV(uplo,n,alpha,AP,X,incx,beta,Y,incy) SUBROUTINE stbmv dtbmv ctbmv ztbmv ---> XTBMV(uplo,trans,diag,n,k,A,lda,X,incx) SUBROUTINE stpmv dtpmv ctpmv ztpmv ---> XTPMV(uplo,trans,diag,n,AP,X,incx) SUBROUTINE strmv dtrmv ctrmv ztrmv ---> XTRMV(uplo,trans,diag,n,A,lda,X,incx) SUBROUTINE ssymv dsymv ---> XSYMV(uplo,n,alpha,A,lda,X,incx,beta,Y,incy) SUBROUTINE ssbmv dsbmv ---> XSBMV(uplo,n,k,alpha,A,lda,X,incx,beta,Y,incy) SUBROUTINE sspmv dspmv ---> XSPMV(uplo,n,alpha,AP,X,incx,beta,Y,incy) SUBROUTINE stbsv dtbsv ctbsv ztbsv ---> XTBSV(uplo,trans,diag,n,k,A,lda,X,incx) SUBROUTINE stpsv dtpsv ctpsv ztpsv ---> XTPSV(uplo,trans,diag,n,AP,X,incx) SUBROUTINE strsv dtrsv ctrsv ztrsv ---> XTRSV(uplo,trans,diag,n,A,lda,X,incx) SUBROUTINE sger dger ---> XGER(m,n,alpha,x,incx,y,incy,A,lda) SUBROUTINE sspr dspr ---> XSPR(uplo,n,alpha,x,incx,AP) SUBROUTINE sspr2 dspr2 ---> XSPR2(uplo,n,alpha,x,incx,y,incy,AP) SUBROUTINE ssyr dsyr ---> XSYR(uplo,n,alpha,x,incx,A,lda) SUBROUTINE ssyr2 dsyr2 ---> XSYR2(uplo,n,alpha,x,incx,y,incy,A,lda) ================ ==== BLAS 3 ==== ================ * SUBROUTINE sgemm dgemm cgemm zgemm ---> XGEMM(transA,transB,m,n,k,alpha,A,lda,B,ldb,beta,C,ldc) SUBROUTINE chemm zhemm ---> XHEMM(side,uplo,m,n,alpha,A,lda,B,ldb,beta,C,ldc) SUBROUTINE cher2k zher2k ---> XHER2K(uplo,trans,n,k,alpha,A,lda,B,ldb,beta,C,ldc) * SUBROUTINE cherk zherk ---> XHERK(uplo,trans,n,k,alpha,A,lda,beta,C,ldc) SUBROUTINE ssymm dsymm csymm zsymm ---> XSYMM(side,uplo,m,n,alpha,A,lda,B,ldb,beta,C,ldc) SUBROUTINE ssyr2k dsyr2k csyr2k zsyr2k ---> XSYR2K(uplo,trans,n,k,alpha,A,lda,B,ldb,beta,C,ldc) SUBROUTINE ssyrk dsyrk csyrk zsyrk ---> XSYRK(uplo,trans,n,k,alpha,A,lda,beta,C,ldc) SUBROUTINE strmm dtrmm ctrmm ztrmm ---> XTRMM(side,uplo,transa,diag,m,n,alpha,A,lda,B,ldb) SUBROUTINE strsm dtrsm ctrsm ztrsm ---> XTRSM(side,uplo,transa,diag,m,n,alpha,A,lda,B,ldb)
SOURCE
94 #if defined HAVE_CONFIG_H 95 #include "config.h" 96 #endif 97 98 #include "abi_common.h" 99 100 MODULE m_blas 101 102 use defs_basis 103 use m_profiling_abi 104 use m_errors 105 106 implicit none 107 108 private 109 110 !BLAS1 111 public :: xnrm2 112 public :: xscal 113 public :: xdotu 114 public :: xdotc 115 public :: xcopy 116 117 !BLAS2 118 public :: xgemv 119 public :: xgerc 120 121 !BLAS3 122 public :: xgemm 123 public :: xherk 124 125 ! Helper functions 126 public :: blas_cholesky_ortho ! Cholesky orthogonalization. 127 128 public :: sqmat_itranspose ! In-place transposition of a square matrix. 129 public :: sqmat_otranspose ! out-of-place transposition of a square matrix. 130 131 public :: sqmat_iconjgtrans ! in-place conjugate transpose of a square matrix. 132 public :: sqmat_oconjgtrans ! out-of-place conjugate transpose of a square matrix. 133 134 !---------------------------------------------------------------------- 135 136 interface xnrm2 137 ! 138 function snrm2 ( n, x, incx ) 139 use defs_basis 140 real(sp) :: snrm2 141 integer,intent(in) :: incx, n 142 real(sp),intent(in) :: x( * ) 143 end function snrm2 144 ! 145 function dnrm2 ( n, x, incx ) 146 use defs_basis 147 real(dp) :: dnrm2 148 integer,intent(in) :: incx, n 149 real(dp),intent(in) :: x( * ) 150 end function dnrm2 151 ! 152 function scnrm2( n, x, incx ) 153 use defs_basis 154 real(sp) :: scnrm2 155 integer,intent(in) :: incx, n 156 complex(spc),intent(in) :: x( * ) 157 end function scnrm2 158 ! 159 function dznrm2( n, x, incx ) 160 use defs_basis 161 real(dp) :: dznrm2 162 integer,intent(in) :: incx, n 163 complex(dpc),intent(in) :: x( * ) 164 end function dznrm2 165 ! 166 end interface xnrm2 167 168 !------------------------------------------------------------------------------- 169 170 interface xscal 171 ! 172 subroutine sscal(n,sa,sx,incx) 173 use defs_basis 174 implicit none 175 integer :: incx 176 integer :: n 177 real(sp) :: sa 178 real(sp) :: sx(*) 179 end subroutine sscal 180 ! 181 subroutine dscal(n,da,dx,incx) 182 use defs_basis 183 implicit none 184 integer :: incx 185 integer :: n 186 real(dp):: da 187 real(dp):: dx(*) 188 end subroutine dscal 189 ! 190 subroutine cscal(n,ca,cx,incx) 191 use defs_basis 192 implicit none 193 integer :: incx 194 integer :: n 195 complex(spc) :: ca 196 complex(spc) :: cx(*) 197 end subroutine cscal 198 ! 199 subroutine zscal(n,za,zx,incx) 200 use defs_basis 201 implicit none 202 integer :: incx 203 integer :: n 204 complex(dpc) :: za 205 complex(dpc) :: zx(*) 206 end subroutine zscal 207 ! 208 subroutine csscal(n,sa,cx,incx) 209 use defs_basis 210 implicit none 211 integer :: incx 212 integer :: n 213 real(sp) :: sa 214 complex(spc) :: cx(*) 215 end subroutine csscal 216 ! 217 subroutine zdscal(n,da,zx,incx) 218 use defs_basis 219 implicit none 220 integer :: incx 221 integer :: n 222 real(dp) :: da 223 complex(dpc) :: zx(*) 224 end subroutine zdscal 225 ! 226 end interface xscal 227 228 !------------------------------------------------------------------------------- 229 230 interface xdotu 231 ! 232 #ifdef HAVE_LINALG_ZDOTU_BUG 233 module procedure cdotu 234 module procedure zdotu 235 #else 236 function cdotu(n,cx,incx,cy,incy) 237 use defs_basis 238 complex(spc) :: cdotu 239 complex(spc),intent(in) :: cx(*),cy(*) 240 integer,intent(in) :: incx,incy,n 241 end function cdotu 242 ! 243 function zdotu(n,zx,incx,zy,incy) 244 use defs_basis 245 complex(dpc) :: zdotu 246 complex(dpc),intent(in) :: zx(*),zy(*) 247 integer,intent(in) :: incx,incy,n 248 end function zdotu 249 #endif 250 ! 251 end interface xdotu 252 253 !------------------------------------------------------------------------------- 254 255 256 ! CDOTC, CDOTU, ZDOTC, and ZDOTU are problematic if Mac OS X's Vec lib is used. 257 ! See http://developer.apple.com/hardwaredrivers/ve/errata.html. 258 ! If needed, we replace them with plain Fortran code. 259 260 interface xdotc 261 ! 262 #ifdef HAVE_LINALG_ZDOTC_BUG 263 module procedure cdotc 264 module procedure zdotc 265 #else 266 function cdotc(n,cx,incx,cy,incy) 267 use defs_basis 268 complex(spc) :: cdotc 269 complex(spc),intent(in) :: cx(*),cy(*) 270 integer,intent(in) :: incx,incy,n 271 end function cdotc 272 ! 273 function zdotc(n,zx,incx,zy,incy) 274 use defs_basis 275 complex(dpc) :: zdotc 276 complex(dpc),intent(in) :: zx(*),zy(*) 277 integer,intent(in) :: incx,incy,n 278 end function zdotc 279 #endif 280 ! 281 end interface xdotc 282 283 !------------------------------------------------------------------------------- 284 285 interface xcopy 286 !module procedure ABI_xcopy 287 ! 288 subroutine scopy(n,sx,incx,sy,incy) 289 use defs_basis 290 implicit none 291 integer,intent(in) :: incx 292 integer,intent(in) :: incy 293 integer,intent(in) :: n 294 real(sp),intent(in) :: sx(*) 295 real(sp),intent(inout) :: sy(*) 296 end subroutine scopy 297 ! 298 subroutine dcopy(n,dx,incx,dy,incy) 299 use defs_basis 300 implicit none 301 integer,intent(in) :: incx 302 integer,intent(in) :: incy 303 integer,intent(in) :: n 304 real(dp),intent(in) :: dx(*) 305 real(dp),intent(inout) :: dy(*) 306 end subroutine dcopy 307 ! 308 subroutine ccopy(n,cx,incx,cy,incy) 309 use defs_basis 310 implicit none 311 integer,intent(in) :: incx 312 integer,intent(in) :: incy 313 integer,intent(in) :: n 314 complex(spc),intent(in) :: cx(*) 315 complex(spc),intent(inout) :: cy(*) 316 end subroutine ccopy 317 ! 318 subroutine zcopy(n,cx,incx,cy,incy) 319 use defs_basis 320 implicit none 321 integer,intent(in) :: incx 322 integer,intent(in) :: incy 323 integer,intent(in) :: n 324 complex(dpc),intent(in) :: cx(*) 325 complex(dpc),intent(inout) :: cy(*) 326 end subroutine zcopy 327 ! 328 end interface xcopy 329 330 !------------------------------------------------------------------------------- 331 332 interface xgemv 333 ! 334 subroutine sgemv ( trans, m, n, alpha, a, lda, x, incx, beta, y, incy ) 335 use defs_basis 336 real(sp),intent(in) :: alpha, beta 337 integer,intent(in) :: incx, incy, lda, m, n 338 character(len=1),intent(in) :: trans 339 real(sp),intent(in) :: a( lda, * ), x( * ) 340 real(sp),intent(inout) :: y( * ) 341 end subroutine sgemv 342 ! 343 subroutine dgemv ( trans, m, n, alpha, a, lda, x, incx, beta, y, incy ) 344 use defs_basis 345 real(dp),intent(in) :: alpha, beta 346 integer,intent(in) :: incx, incy, lda, m, n 347 character(len=1),intent(in) :: trans 348 real(dp),intent(in) :: a( lda, * ), x( * ) 349 real(dp),intent(inout) :: y( * ) 350 end subroutine dgemv 351 ! 352 subroutine cgemv ( trans, m, n, alpha, a, lda, x, incx, beta, y, incy ) 353 use defs_basis 354 complex(spc),intent(in) :: alpha, beta 355 integer,intent(in) :: incx, incy, lda, m, n 356 character(len=1),intent(in) :: trans 357 complex(spc),intent(in) :: a( lda, * ), x( * ) 358 complex(spc),intent(inout) :: y( * ) 359 end subroutine cgemv 360 ! 361 subroutine zgemv ( trans, m, n, alpha, a, lda, x, incx, beta, y, incy ) 362 use defs_basis 363 complex(dpc),intent(in) :: alpha, beta 364 integer,intent(in) :: incx, incy, lda, m, n 365 character(len=1),intent(in) :: trans 366 complex(dpc),intent(in) :: a( lda, * ), x( * ) 367 complex(dpc),intent(inout) :: y( * ) 368 end subroutine zgemv 369 ! 370 end interface xgemv 371 372 !------------------------------------------------------------------------------- 373 374 interface xgerc 375 ! 376 subroutine cgerc ( m, n, alpha, x, incx, y, incy, a, lda ) 377 use defs_basis 378 complex(spc),intent(in) :: alpha 379 integer,intent(in) :: incx, incy, lda, m, n 380 complex(spc),intent(inout) :: a( lda, * ) 381 complex(spc),intent(in) :: x( * ), y( * ) 382 end subroutine cgerc 383 ! 384 subroutine zgerc ( m, n, alpha, x, incx, y, incy, a, lda ) 385 use defs_basis 386 complex(dpc),intent(in) :: alpha 387 integer,intent(in) :: incx, incy, lda, m, n 388 complex(dpc),intent(inout) :: a( lda, * ) 389 complex(dpc),intent(in) :: x( * ), y( * ) 390 end subroutine zgerc 391 ! 392 end interface xgerc 393 394 !------------------------------------------------------------------------------- 395 396 interface xgemm 397 ! 398 subroutine sgemm ( transa, transb, m, n, k, alpha, a, lda, b, ldb, beta, c, ldc ) 399 use defs_basis 400 character(len=1),intent(in) :: transa, transb 401 integer,intent(in) :: m, n, k, lda, ldb, ldc 402 real(sp),intent(in) :: alpha, beta 403 real(sp),intent(in) :: a( lda, * ), b( ldb, * ) 404 real(sp),intent(inout) :: c( ldc, * ) 405 end subroutine sgemm 406 ! 407 subroutine dgemm ( transa, transb, m, n, k, alpha, a, lda, b, ldb, beta, c, ldc ) 408 use defs_basis 409 character(len=1),intent(in) :: transa, transb 410 integer,intent(in) :: m, n, k, lda, ldb, ldc 411 real(dp),intent(in) :: alpha, beta 412 real(dp),intent(in) :: a( lda, * ), b( ldb, * ) 413 real(dp),intent(inout) :: c( ldc, * ) 414 end subroutine dgemm 415 ! 416 subroutine cgemm ( transa, transb, m, n, k, alpha, a, lda, b, ldb, beta, c, ldc ) 417 use defs_basis 418 character(len=1),intent(in) :: transa, transb 419 integer,intent(in) :: m, n, k, lda, ldb, ldc 420 complex(spc),intent(in) :: alpha, beta 421 complex(spc),intent(in) :: a( lda, * ), b( ldb, * ) 422 complex(spc),intent(inout) :: c( ldc, * ) 423 end subroutine cgemm 424 ! 425 subroutine zgemm ( transa, transb, m, n, k, alpha, a, lda, b, ldb, beta, c, ldc ) 426 use defs_basis 427 character(len=1),intent(in) :: transa, transb 428 integer,intent(in) :: m, n, k, lda, ldb, ldc 429 complex(dpc),intent(in) :: alpha, beta 430 complex(dpc),intent(in) :: a( lda, * ), b( ldb, * ) 431 complex(dpc),intent(inout) :: c( ldc, * ) 432 end subroutine zgemm 433 ! 434 end interface xgemm 435 436 interface xherk 437 ! 438 subroutine cherk( uplo, trans, n, k, alpha, a, lda, beta, c, ldc ) 439 use defs_basis 440 character(len=1),intent(in) :: uplo 441 character(len=1),intent(in) :: trans 442 integer,intent(in) :: n,k,lda,ldc 443 real(sp),intent(in) :: alpha 444 complex(spc),intent(in) :: a( lda, * ) 445 real(sp),intent(in) :: beta 446 complex(spc),intent(inout) :: c( ldc, * ) 447 end subroutine cherk 448 ! 449 subroutine zherk( uplo, trans, n, k, alpha, a, lda, beta, c, ldc ) 450 use defs_basis 451 character(len=1),intent(in) :: uplo 452 character(len=1),intent(in) :: trans 453 integer,intent(in) :: n,k,lda,ldc 454 real(dp),intent(in) :: alpha 455 complex(dpc),intent(in) :: a( lda, * ) 456 real(dp),intent(in) :: beta 457 complex(dpc),intent(inout) :: c( ldc, * ) 458 end subroutine zherk 459 ! 460 end interface xherk 461 462 !------------------------------------------------------------------------------- 463 464 interface blas_cholesky_ortho 465 module procedure blas_cholesky_ortho_spc 466 module procedure blas_cholesky_ortho_dpc 467 end interface blas_cholesky_ortho 468 469 interface sqmat_itranspose 470 module procedure sqmat_itranspose_sp 471 module procedure sqmat_itranspose_dp 472 module procedure sqmat_itranspose_spc 473 module procedure sqmat_itranspose_dpc 474 end interface sqmat_itranspose 475 476 interface sqmat_otranspose 477 module procedure sqmat_otranspose_sp 478 module procedure sqmat_otranspose_dp 479 module procedure sqmat_otranspose_spc 480 module procedure sqmat_otranspose_dpc 481 end interface sqmat_otranspose 482 483 interface sqmat_iconjgtrans 484 module procedure sqmat_iconjgtrans_spc 485 module procedure sqmat_iconjgtrans_dpc 486 end interface sqmat_iconjgtrans 487 488 interface sqmat_oconjgtrans 489 module procedure sqmat_oconjgtrans_spc 490 module procedure sqmat_oconjgtrans_dpc 491 end interface sqmat_oconjgtrans 492 493 real(sp),private,parameter :: zero_sp = 0._sp 494 real(sp),private,parameter :: one_sp = 1._sp 495 496 real(dp),private,parameter :: zero_dp = 0._dp 497 real(dp),private,parameter :: one_dp = 1._dp 498 499 complex(spc),private,parameter :: czero_spc = (0._sp,0._sp) 500 complex(spc),private,parameter :: cone_spc = (1._sp,0._sp) 501 502 complex(dpc),private,parameter :: czero_dpc = (0._dp,0._dp) 503 complex(dpc),private,parameter :: cone_dpc = (1._dp,0._dp) 504 505 CONTAINS !======================================================================================== 506 507 ! CDOTC, CDOTU, ZDOTC, and ZDOTU are problematic if Mac OS X's Vec lib is used. 508 ! See http://developer.apple.com/hardwaredrivers/ve/errata.html. 509 ! Here we replace them with plain Fortran code. 510 511 #ifdef HAVE_LINALG_ZDOTC_BUG 512 !#warning "Using internal replacement for zdotc. External library cannot be used" 513 #include "replacements/cdotc.f" 514 #include "replacements/zdotc.f" 515 #endif 516 517 #ifdef HAVE_LINALG_ZDOTU_BUG 518 !#warning "Using internal replacement for zdotu. External library cannot be used" 519 #include "replacements/cdotu.f" 520 #include "replacements/zdotu.f" 521 #endif 522 523 !----------------------------------------------------------------------
m_blas/blas_cholesky_ortho_dpc [ Functions ]
[ Top ] [ m_blas ] [ Functions ]
NAME
blas_cholesky_ortho_dpc
FUNCTION
Performs the Cholesky orthonormalization of the vectors stored in iomat.
INPUTS
vec_size=Size of each vector. nvec=Number of vectors in iomat
OUTPUT
cf_ovlp=Cholesky factorization of the overlap matrix. ovlp = U^H U with U upper triangle matrix returned in cf_ovlp
SIDE EFFECTS
iomat(vec_size,nvec) input: Input set of vectors. output: Orthonormalized set.
PARENTS
CHILDREN
mkl_zomatcopy
SOURCE
629 subroutine blas_cholesky_ortho_dpc(vec_size,nvec,iomat,cf_ovlp,use_gemm) 630 631 632 !This section has been created automatically by the script Abilint (TD). 633 !Do not modify the following lines by hand. 634 #undef ABI_FUNC 635 #define ABI_FUNC 'blas_cholesky_ortho_dpc' 636 !End of the abilint section 637 638 implicit none 639 640 !Arguments ------------------------------------ 641 integer,intent(in) :: vec_size,nvec 642 logical,optional,intent(in) :: use_gemm 643 complex(dpc),intent(inout) :: iomat(vec_size,nvec) 644 complex(dpc),intent(out) :: cf_ovlp(nvec,nvec) 645 646 !Local variables ------------------------------ 647 !scalars 648 integer :: ierr 649 logical :: my_usegemm 650 character(len=500) :: msg 651 652 ! ************************************************************************* 653 654 ! 1) Calculate overlap_ij = <phi_i|phi_j> 655 my_usegemm = .FALSE.; if (PRESENT(use_gemm)) my_usegemm = use_gemm 656 657 if (my_usegemm) then 658 call xgemm("Conjugate","Normal",nvec,nvec,vec_size,cone_dpc,iomat,vec_size,iomat,vec_size,czero_dpc,cf_ovlp,nvec) 659 else 660 call xherk("U","C", nvec, vec_size, one_dp, iomat, vec_size, zero_dp, cf_ovlp, nvec) 661 end if 662 ! 663 ! 2) Cholesky factorization: ovlp = U^H U with U upper triangle matrix. 664 call ZPOTRF('U',nvec,cf_ovlp,nvec,ierr) 665 if (ierr/=0) then 666 write(msg,'(a,i0)')' ZPOTRF returned info= ',ierr 667 MSG_ERROR(msg) 668 end if 669 ! 670 ! 3) Solve X U = io_mat. On exit io_mat is orthonormalized. 671 call ZTRSM('Right','Upper','Normal','Normal',vec_size,nvec,cone_dpc,cf_ovlp,nvec,iomat,vec_size) 672 673 end subroutine blas_cholesky_ortho_dpc
m_blas/blas_cholesky_ortho_spc [ Functions ]
[ Top ] [ m_blas ] [ Functions ]
NAME
blas_cholesky_ortho_spc
FUNCTION
Performs the Cholesky orthonormalization of the vectors stored in iomat.
INPUTS
vec_size=Size of each vector. nvec=Number of vectors in iomat
OUTPUT
cf_ovlp=Cholesky factorization of the overlap matrix. ovlp = U^H U with U upper triangle matrix returned in cf_ovlp
SIDE EFFECTS
iomat(vec_size,nvec) input: Input set of vectors. output: Orthonormalized set.
PARENTS
CHILDREN
mkl_zomatcopy
SOURCE
553 subroutine blas_cholesky_ortho_spc(vec_size,nvec,iomat,cf_ovlp,use_gemm) 554 555 556 !This section has been created automatically by the script Abilint (TD). 557 !Do not modify the following lines by hand. 558 #undef ABI_FUNC 559 #define ABI_FUNC 'blas_cholesky_ortho_spc' 560 !End of the abilint section 561 562 implicit none 563 564 !Arguments ------------------------------------ 565 integer,intent(in) :: vec_size,nvec 566 logical,optional,intent(in) :: use_gemm 567 complex(spc),intent(inout) :: iomat(vec_size,nvec) 568 complex(spc),intent(out) :: cf_ovlp(nvec,nvec) 569 570 !Local variables ------------------------------ 571 !scalars 572 integer :: ierr 573 logical :: my_usegemm 574 character(len=500) :: msg 575 576 ! ************************************************************************* 577 578 ! 1) Calculate overlap_ij = <phi_i|phi_j> 579 ! TODO: use dsyrk 580 my_usegemm = .FALSE.; if (PRESENT(use_gemm)) my_usegemm = use_gemm 581 582 if (my_usegemm) then 583 call xgemm("Conjugate","Normal",nvec,nvec,vec_size,cone_spc,iomat,vec_size,iomat,vec_size,czero_spc,cf_ovlp,nvec) 584 else 585 call xherk("U","C", nvec, vec_size, one_sp, iomat, vec_size, zero_sp, cf_ovlp, nvec) 586 end if 587 ! 588 ! 2) Cholesky factorization: ovlp = U^H U with U upper triangle matrix. 589 call CPOTRF('U',nvec,cf_ovlp,nvec,ierr) 590 if (ierr/=0) then 591 write(msg,'(a,i0)')' ZPOTRF returned info= ',ierr 592 MSG_ERROR(msg) 593 end if 594 ! 595 ! 3) Solve X U = io_mat. On exit iomat is orthonormalized. 596 call CTRSM('Right','Upper','Normal','Normal',vec_size,nvec,cone_spc,cf_ovlp,nvec,iomat,vec_size) 597 598 end subroutine blas_cholesky_ortho_spc
m_blas/sqmat_iconjgtrans_dpc [ Functions ]
[ Top ] [ m_blas ] [ Functions ]
NAME
sqmat_iconjgtrans_dpc
FUNCTION
Compute alpha * CONJG(mat^T) in place. target: double precision complex matrix.
INPUTS
n=size of the matrix [alpha]=scalar, set to 1.0 if not present
SIDE EFFECTS
mat(n,n)=in output, it contains alpha * CONJG(mat^T).
PARENTS
CHILDREN
mkl_zomatcopy
SOURCE
1292 subroutine sqmat_iconjgtrans_dpc(n,mat,alpha) 1293 1294 1295 !This section has been created automatically by the script Abilint (TD). 1296 !Do not modify the following lines by hand. 1297 #undef ABI_FUNC 1298 #define ABI_FUNC 'sqmat_iconjgtrans_dpc' 1299 !End of the abilint section 1300 1301 implicit none 1302 1303 !Arguments ------------------------------------ 1304 !scalars 1305 integer,intent(in) :: n 1306 complex(dpc),optional,intent(in) :: alpha 1307 !arrays 1308 complex(dpc),intent(inout) :: mat(n,n) 1309 1310 ! ************************************************************************* 1311 1312 #if defined HAVE_LINALG_MKL_IMATCOPY 1313 if (PRESENT(alpha)) then 1314 call mkl_zimatcopy("Column", "C", n, n, alpha, mat, n, n) 1315 else 1316 call mkl_zimatcopy("Column", "C", n, n, cone_dpc, mat, n, n) 1317 end if 1318 #else 1319 ! Fallback to Fortran. 1320 if (PRESENT(alpha)) then 1321 !$OMP PARALLEL WORKSHARE 1322 mat = alpha * TRANSPOSE(CONJG(mat)) 1323 !$OMP END PARALLEL WORKSHARE 1324 else 1325 !$OMP PARALLEL WORKSHARE 1326 mat = TRANSPOSE(CONJG(mat)) 1327 !$OMP END PARALLEL WORKSHARE 1328 end if 1329 #endif 1330 1331 end subroutine sqmat_iconjgtrans_dpc
m_blas/sqmat_iconjgtrans_spc [ Functions ]
[ Top ] [ m_blas ] [ Functions ]
NAME
sqmat_iconjgtrans_spc
FUNCTION
Compute alpha * CONJG(mat^T) in place. target: single precision complex matrix.
INPUTS
n=size of the matrix [alpha]=scalar, set to 1.0 if not present
SIDE EFFECTS
mat(n,n)=in output, it contains alpha * CONJG(mat^T).
PARENTS
CHILDREN
mkl_zomatcopy
SOURCE
1227 subroutine sqmat_iconjgtrans_spc(n,mat,alpha) 1228 1229 1230 !This section has been created automatically by the script Abilint (TD). 1231 !Do not modify the following lines by hand. 1232 #undef ABI_FUNC 1233 #define ABI_FUNC 'sqmat_iconjgtrans_spc' 1234 !End of the abilint section 1235 1236 implicit none 1237 1238 !Arguments ------------------------------------ 1239 !scalars 1240 integer,intent(in) :: n 1241 complex(spc),optional,intent(in) :: alpha 1242 !arrays 1243 complex(spc),intent(inout) :: mat(n,n) 1244 1245 ! ************************************************************************* 1246 1247 #if defined HAVE_LINALG_MKL_IMATCOPY 1248 if (PRESENT(alpha)) then 1249 call mkl_cimatcopy("Column", "C", n, n, alpha, mat, n, n) 1250 else 1251 call mkl_cimatcopy("Column", "C", n, n, cone_spc, mat, n, n) 1252 end if 1253 #else 1254 ! Fallback to Fortran. 1255 if (PRESENT(alpha)) then 1256 !$OMP PARALLEL WORKSHARE 1257 mat = alpha * TRANSPOSE(CONJG(mat)) 1258 !$OMP END PARALLEL WORKSHARE 1259 else 1260 !$OMP PARALLEL WORKSHARE 1261 mat = TRANSPOSE(CONJG(mat)) 1262 !$OMP END PARALLEL WORKSHARE 1263 end if 1264 #endif 1265 1266 end subroutine sqmat_iconjgtrans_spc
m_blas/sqmat_itranspose_dp [ Functions ]
[ Top ] [ m_blas ] [ Functions ]
NAME
sqmat_itranspose_dp
FUNCTION
Compute alpha * mat^T in place. target: double precision real matrix.
INPUTS
n=size of the matrix [alpha]=scalar, set to 1.0 if not present
SIDE EFFECTS
mat(n,n)=in output, it contains alpha * mat^T.
PARENTS
CHILDREN
mkl_zomatcopy
SOURCE
764 subroutine sqmat_itranspose_dp(n,mat,alpha) 765 766 767 !This section has been created automatically by the script Abilint (TD). 768 !Do not modify the following lines by hand. 769 #undef ABI_FUNC 770 #define ABI_FUNC 'sqmat_itranspose_dp' 771 !End of the abilint section 772 773 implicit none 774 775 !Arguments ------------------------------------ 776 !scalars 777 integer,intent(in) :: n 778 real(dp),optional,intent(in) :: alpha 779 !arrays 780 real(dp),intent(inout) :: mat(n,n) 781 782 ! ************************************************************************* 783 784 #if defined HAVE_LINALG_MKL_IMATCOPY 785 if (PRESENT(alpha)) then 786 call mkl_dimatcopy("Column", "Trans", n, n, alpha, mat, n, n) 787 else 788 call mkl_dimatcopy("Column", "Trans", n, n, one_dp, mat, n, n) 789 end if 790 #else 791 ! Fallback to Fortran. 792 if (PRESENT(alpha)) then 793 !$OMP PARALLEL WORKSHARE 794 mat = alpha * TRANSPOSE(mat) 795 !$OMP END PARALLEL WORKSHARE 796 else 797 !$OMP PARALLEL WORKSHARE 798 mat = TRANSPOSE(mat) 799 !$OMP END PARALLEL WORKSHARE 800 end if 801 #endif 802 803 end subroutine sqmat_itranspose_dp
m_blas/sqmat_itranspose_dpc [ Functions ]
[ Top ] [ m_blas ] [ Functions ]
NAME
sqmat_itranspose_dpc
FUNCTION
Compute alpha * mat^T in place. target: double precision complex matrix.
INPUTS
n=size of the matrix [alpha]=scalar, set to 1.0 if not present
SIDE EFFECTS
mat(n,n)=in output, it contains alpha * mat^T.
PARENTS
CHILDREN
mkl_zomatcopy
SOURCE
894 subroutine sqmat_itranspose_dpc(n,mat,alpha) 895 896 897 !This section has been created automatically by the script Abilint (TD). 898 !Do not modify the following lines by hand. 899 #undef ABI_FUNC 900 #define ABI_FUNC 'sqmat_itranspose_dpc' 901 !End of the abilint section 902 903 implicit none 904 905 !Arguments ------------------------------------ 906 !scalars 907 integer,intent(in) :: n 908 complex(dpc),optional,intent(in) :: alpha 909 !arrays 910 complex(dpc),intent(inout) :: mat(n,n) 911 912 ! ************************************************************************* 913 914 #if defined HAVE_LINALG_MKL_IMATCOPY 915 if (PRESENT(alpha)) then 916 call mkl_zimatcopy("Column", "Trans", n, n, alpha, mat, n, n) 917 else 918 call mkl_zimatcopy("Column", "Trans", n, n, cone_dpc, mat, n, n) 919 end if 920 #else 921 ! Fallback to Fortran. 922 if (PRESENT(alpha)) then 923 !$OMP PARALLEL WORKSHARE 924 mat = alpha * TRANSPOSE(mat) 925 !$OMP END PARALLEL WORKSHARE 926 else 927 !$OMP PARALLEL WORKSHARE 928 mat = TRANSPOSE(mat) 929 !$OMP END PARALLEL WORKSHARE 930 end if 931 #endif 932 933 end subroutine sqmat_itranspose_dpc
m_blas/sqmat_itranspose_sp [ Functions ]
[ Top ] [ m_blas ] [ Functions ]
NAME
sqmat_itranspose_sp
FUNCTION
Compute alpha * mat^T in place. target: single precision real matrix.
INPUTS
n=size of the matrix [alpha]=scalar, set to 1.0 if not present
SIDE EFFECTS
mat(n,n)=in output, it contains alpha * mat^T.
PARENTS
CHILDREN
mkl_zomatcopy
SOURCE
699 subroutine sqmat_itranspose_sp(n,mat,alpha) 700 701 702 !This section has been created automatically by the script Abilint (TD). 703 !Do not modify the following lines by hand. 704 #undef ABI_FUNC 705 #define ABI_FUNC 'sqmat_itranspose_sp' 706 !End of the abilint section 707 708 implicit none 709 710 !Arguments ------------------------------------ 711 !scalars 712 integer,intent(in) :: n 713 real(sp),optional,intent(in) :: alpha 714 !arrays 715 real(sp),intent(inout) :: mat(n,n) 716 717 ! ************************************************************************* 718 719 #if defined HAVE_LINALG_MKL_IMATCOPY 720 if (PRESENT(alpha)) then 721 call mkl_simatcopy("Column", "Trans", n, n, alpha, mat, n, n) 722 else 723 call mkl_simatcopy("Column", "Trans", n, n, one_sp, mat, n, n) 724 end if 725 #else 726 ! Fallback to Fortran. 727 if (PRESENT(alpha)) then 728 !$OMP PARALLEL WORKSHARE 729 mat = alpha * TRANSPOSE(mat) 730 !$OMP END PARALLEL WORKSHARE 731 else 732 !$OMP PARALLEL WORKSHARE 733 mat = TRANSPOSE(mat) 734 !$OMP END PARALLEL WORKSHARE 735 end if 736 #endif 737 738 end subroutine sqmat_itranspose_sp
m_blas/sqmat_itranspose_spc [ Functions ]
[ Top ] [ m_blas ] [ Functions ]
NAME
sqmat_itranspose_spc
FUNCTION
Compute alpha * mat^T in place. target: single precision complex matrix.
INPUTS
n=size of the matrix [alpha]=scalar, set to 1.0 if not present
SIDE EFFECTS
mat(n,n)=in output, it contains alpha * mat^T.
PARENTS
CHILDREN
mkl_zomatcopy
SOURCE
829 subroutine sqmat_itranspose_spc(n,mat,alpha) 830 831 832 !This section has been created automatically by the script Abilint (TD). 833 !Do not modify the following lines by hand. 834 #undef ABI_FUNC 835 #define ABI_FUNC 'sqmat_itranspose_spc' 836 !End of the abilint section 837 838 implicit none 839 840 !Arguments ------------------------------------ 841 !scalars 842 integer,intent(in) :: n 843 complex(spc),optional,intent(in) :: alpha 844 !arrays 845 complex(spc),intent(inout) :: mat(n,n) 846 847 ! ************************************************************************* 848 849 #if defined HAVE_LINALG_MKL_IMATCOPY 850 if (PRESENT(alpha)) then 851 call mkl_cimatcopy("Column", "Trans", n, n, alpha, mat, n, n) 852 else 853 call mkl_cimatcopy("Column", "Trans", n, n, cone_spc, mat, n, n) 854 end if 855 #else 856 ! Fallback to Fortran. 857 if (PRESENT(alpha)) then 858 !$OMP PARALLEL WORKSHARE 859 mat = alpha * TRANSPOSE(mat) 860 !$OMP END PARALLEL WORKSHARE 861 else 862 !$OMP PARALLEL WORKSHARE 863 mat = TRANSPOSE(mat) 864 !$OMP END PARALLEL WORKSHARE 865 end if 866 #endif 867 868 end subroutine sqmat_itranspose_spc
m_blas/sqmat_oconjgtrans_dpc [ Functions ]
[ Top ] [ m_blas ] [ Functions ]
NAME
sqmat_oconjgtrans_dpc
FUNCTION
Compute alpha * CONJG(mat^T) out-of-place. target: double precision complex matrix.
INPUTS
n=size of the matrix [alpha]=scalar, set to 1.0 if not present imat(n,n)=Input matrix.
OUTPUT
omat(n,n)=contains alpha * CONJG(imat^T).
PARENTS
CHILDREN
mkl_zomatcopy
SOURCE
1425 subroutine sqmat_oconjgtrans_dpc(n,imat,omat,alpha) 1426 1427 1428 !This section has been created automatically by the script Abilint (TD). 1429 !Do not modify the following lines by hand. 1430 #undef ABI_FUNC 1431 #define ABI_FUNC 'sqmat_oconjgtrans_dpc' 1432 !End of the abilint section 1433 1434 implicit none 1435 1436 !Arguments ------------------------------------ 1437 !scalars 1438 integer,intent(in) :: n 1439 complex(dpc),optional,intent(in) :: alpha 1440 !arrays 1441 complex(dpc),intent(in) :: imat(n,n) 1442 complex(dpc),intent(out) :: omat(n,n) 1443 1444 ! ************************************************************************* 1445 1446 #if defined HAVE_LINALG_MKL_OMATCOPY 1447 if (PRESENT(alpha)) then 1448 call mkl_zomatcopy("Column", "C", n, n, alpha, imat, n, omat, n) 1449 else 1450 call mkl_zomatcopy("Column", "C", n, n, cone_dpc, imat, n, omat, n) 1451 end if 1452 #else 1453 ! Fallback to Fortran. 1454 if (PRESENT(alpha)) then 1455 !$OMP PARALLEL WORKSHARE 1456 omat = alpha * TRANSPOSE(CONJG(imat)) 1457 !$OMP END PARALLEL WORKSHARE 1458 else 1459 !$OMP PARALLEL WORKSHARE 1460 omat = TRANSPOSE(CONJG(imat)) 1461 !$OMP END PARALLEL WORKSHARE 1462 end if 1463 #endif 1464 1465 end subroutine sqmat_oconjgtrans_dpc
m_blas/sqmat_oconjgtrans_spc [ Functions ]
[ Top ] [ m_blas ] [ Functions ]
NAME
sqmat_oconjgtrans_spc
FUNCTION
Compute alpha * CONJG(mat^T) out-of-place. target: single precision complex matrix.
INPUTS
n=size of the matrix [alpha]=scalar, set to 1.0 if not present imat(n,n)=Input matrix.
OUTPUT
omat(n,n)=contains alpha * CONJG(imat^T).
PARENTS
CHILDREN
mkl_zomatcopy
SOURCE
1358 subroutine sqmat_oconjgtrans_spc(n,imat,omat,alpha) 1359 1360 1361 !This section has been created automatically by the script Abilint (TD). 1362 !Do not modify the following lines by hand. 1363 #undef ABI_FUNC 1364 #define ABI_FUNC 'sqmat_oconjgtrans_spc' 1365 !End of the abilint section 1366 1367 implicit none 1368 1369 !Arguments ------------------------------------ 1370 !scalars 1371 integer,intent(in) :: n 1372 complex(spc),optional,intent(in) :: alpha 1373 !arrays 1374 complex(spc),intent(in) :: imat(n,n) 1375 complex(spc),intent(out) :: omat(n,n) 1376 1377 ! ************************************************************************* 1378 1379 #if defined HAVE_LINALG_MKL_OMATCOPY 1380 if (PRESENT(alpha)) then 1381 call mkl_comatcopy("Column", "C", n, n, alpha, imat, n, omat, n) 1382 else 1383 call mkl_comatcopy("Column", "C", n, n, cone_spc, imat, n, omat, n) 1384 end if 1385 #else 1386 ! Fallback to Fortran. 1387 if (PRESENT(alpha)) then 1388 !$OMP PARALLEL WORKSHARE 1389 omat = alpha * TRANSPOSE(CONJG(imat)) 1390 !$OMP END PARALLEL WORKSHARE 1391 else 1392 !$OMP PARALLEL WORKSHARE 1393 omat = TRANSPOSE(CONJG(imat)) 1394 !$OMP END PARALLEL WORKSHARE 1395 end if 1396 #endif 1397 1398 end subroutine sqmat_oconjgtrans_spc
m_blas/sqmat_otranspose_dp [ Functions ]
[ Top ] [ m_blas ] [ Functions ]
NAME
sqmat_otranspose_dp
FUNCTION
Compute alpha * mat^T out-of-place. target: double precision real matrix.
INPUTS
n=size of the matrix [alpha]=scalar, set to 1.0 if not present imat(n,n)=Input matrix.
OUTPUT
omat(n,n)=contains alpha * imat^T.
PARENTS
CHILDREN
mkl_zomatcopy
SOURCE
1027 subroutine sqmat_otranspose_dp(n,imat,omat,alpha) 1028 1029 1030 !This section has been created automatically by the script Abilint (TD). 1031 !Do not modify the following lines by hand. 1032 #undef ABI_FUNC 1033 #define ABI_FUNC 'sqmat_otranspose_dp' 1034 !End of the abilint section 1035 1036 implicit none 1037 1038 !Arguments ------------------------------------ 1039 !scalars 1040 integer,intent(in) :: n 1041 real(dp),optional,intent(in) :: alpha 1042 !arrays 1043 real(dp),intent(in) :: imat(n,n) 1044 real(dp),intent(out) :: omat(n,n) 1045 1046 ! ************************************************************************* 1047 1048 #if defined HAVE_LINALG_MKL_OMATCOPY 1049 if (PRESENT(alpha)) then 1050 call mkl_domatcopy("Column", "Transpose", n, n, alpha, imat, n, omat, n) 1051 else 1052 call mkl_domatcopy("Column", "Transpose", n, n, one_dp, imat, n, omat, n) 1053 end if 1054 #else 1055 ! Fallback to Fortran. 1056 if (PRESENT(alpha)) then 1057 !$OMP PARALLEL WORKSHARE 1058 omat = alpha * TRANSPOSE(imat) 1059 !$OMP END PARALLEL WORKSHARE 1060 else 1061 !$OMP PARALLEL WORKSHARE 1062 omat = TRANSPOSE(imat) 1063 !$OMP END PARALLEL WORKSHARE 1064 end if 1065 #endif 1066 1067 end subroutine sqmat_otranspose_dp
m_blas/sqmat_otranspose_dpc [ Functions ]
[ Top ] [ m_blas ] [ Functions ]
NAME
sqmat_otranspose_dpc
FUNCTION
Compute alpha * mat^T out-of-place. target: double precision complex matrix.
INPUTS
n=size of the matrix [alpha]=scalar, set to 1.0 if not present imat(n,n)=Input matrix.
OUTPUT
omat(n,n)=contains alpha * imat^T.
PARENTS
CHILDREN
mkl_zomatcopy
SOURCE
1161 subroutine sqmat_otranspose_dpc(n,imat,omat,alpha) 1162 1163 1164 !This section has been created automatically by the script Abilint (TD). 1165 !Do not modify the following lines by hand. 1166 #undef ABI_FUNC 1167 #define ABI_FUNC 'sqmat_otranspose_dpc' 1168 !End of the abilint section 1169 1170 implicit none 1171 1172 !Arguments ------------------------------------ 1173 !scalars 1174 integer,intent(in) :: n 1175 complex(dpc),optional,intent(in) :: alpha 1176 !arrays 1177 complex(dpc),intent(in) :: imat(n,n) 1178 complex(dpc),intent(out) :: omat(n,n) 1179 1180 ! ************************************************************************* 1181 1182 #if defined HAVE_LINALG_MKL_OMATCOPY 1183 if (PRESENT(alpha)) then 1184 call mkl_zomatcopy("Column", "Transpose", n, n, alpha, imat, n, omat, n) 1185 else 1186 call mkl_zomatcopy("Column", "Transpose", n, n, cone_dpc, imat, n, omat, n) 1187 end if 1188 #else 1189 ! Fallback to Fortran. 1190 if (PRESENT(alpha)) then 1191 !$OMP PARALLEL WORKSHARE 1192 omat = alpha * TRANSPOSE(imat) 1193 !$OMP END PARALLEL WORKSHARE 1194 else 1195 !$OMP PARALLEL WORKSHARE 1196 omat = TRANSPOSE(imat) 1197 !$OMP END PARALLEL WORKSHARE 1198 end if 1199 #endif 1200 1201 end subroutine sqmat_otranspose_dpc
m_blas/sqmat_otranspose_sp [ Functions ]
[ Top ] [ m_blas ] [ Functions ]
NAME
sqmat_otranspose_sp
FUNCTION
Compute alpha * mat^T out-of-place. target: single precision real matrix.
INPUTS
n=size of the matrix [alpha]=scalar, set to 1.0 if not present imat(n,n)=Input matrix.
OUTPUT
omat(n,n)=contains alpha * imat^T.
PARENTS
CHILDREN
mkl_zomatcopy
SOURCE
960 subroutine sqmat_otranspose_sp(n,imat,omat,alpha) 961 962 963 !This section has been created automatically by the script Abilint (TD). 964 !Do not modify the following lines by hand. 965 #undef ABI_FUNC 966 #define ABI_FUNC 'sqmat_otranspose_sp' 967 !End of the abilint section 968 969 implicit none 970 971 !Arguments ------------------------------------ 972 !scalars 973 integer,intent(in) :: n 974 real(sp),optional,intent(in) :: alpha 975 !arrays 976 real(sp),intent(in) :: imat(n,n) 977 real(sp),intent(out) :: omat(n,n) 978 979 ! ************************************************************************* 980 981 #if defined HAVE_LINALG_MKL_OMATCOPY 982 if (PRESENT(alpha)) then 983 call mkl_somatcopy("Column", "Transpose", n, n, alpha, imat, n, omat, n) 984 else 985 call mkl_somatcopy("Column", "Transpose", n, n, one_sp, imat, n, omat, n) 986 end if 987 #else 988 ! Fallback to Fortran. 989 if (PRESENT(alpha)) then 990 !$OMP PARALLEL WORKSHARE 991 omat = alpha * TRANSPOSE(imat) 992 !$OMP END PARALLEL WORKSHARE 993 else 994 !$OMP PARALLEL WORKSHARE 995 omat = TRANSPOSE(imat) 996 !$OMP END PARALLEL WORKSHARE 997 end if 998 #endif 999 1000 end subroutine sqmat_otranspose_sp
m_blas/sqmat_otranspose_spc [ Functions ]
[ Top ] [ m_blas ] [ Functions ]
NAME
sqmat_otranspose_spc
FUNCTION
Compute alpha * mat^T out-of-place. target: single precision complex matrix.
INPUTS
n=size of the matrix [alpha]=scalar, set to 1.0 if not present imat(n,n)=Input matrix.
OUTPUT
omat(n,n)=contains alpha * imat^T.
PARENTS
CHILDREN
mkl_zomatcopy
SOURCE
1094 subroutine sqmat_otranspose_spc(n,imat,omat,alpha) 1095 1096 1097 !This section has been created automatically by the script Abilint (TD). 1098 !Do not modify the following lines by hand. 1099 #undef ABI_FUNC 1100 #define ABI_FUNC 'sqmat_otranspose_spc' 1101 !End of the abilint section 1102 1103 implicit none 1104 1105 !Arguments ------------------------------------ 1106 !scalars 1107 integer,intent(in) :: n 1108 complex(spc),optional,intent(in) :: alpha 1109 !arrays 1110 complex(spc),intent(in) :: imat(n,n) 1111 complex(spc),intent(out) :: omat(n,n) 1112 1113 ! ************************************************************************* 1114 1115 #if defined HAVE_LINALG_MKL_OMATCOPY 1116 if (PRESENT(alpha)) then 1117 call mkl_comatcopy("Column", "Transpose", n, n, alpha, imat, n, omat, n) 1118 else 1119 call mkl_comatcopy("Column", "Transpose", n, n, cone_spc, imat, n, omat, n) 1120 end if 1121 #else 1122 ! Fallback to Fortran. 1123 if (PRESENT(alpha)) then 1124 !$OMP PARALLEL WORKSHARE 1125 omat = alpha * TRANSPOSE(imat) 1126 !$OMP END PARALLEL WORKSHARE 1127 else 1128 !$OMP PARALLEL WORKSHARE 1129 omat = TRANSPOSE(imat) 1130 !$OMP END PARALLEL WORKSHARE 1131 end if 1132 #endif 1133 1134 end subroutine sqmat_otranspose_spc