TABLE OF CONTENTS
- m_abi_linalg/abi_d2ztrsm
- m_abi_linalg/abi_d2ztrsm_3d
- m_abi_linalg/abi_dtrsm
- m_abi_linalg/abi_xtrsm
- m_abi_linalg/abi_ztrsm
m_abi_linalg/abi_d2ztrsm [ Functions ]
[ Top ] [ m_abi_linalg ] [ Functions ]
NAME
abi_d2ztrsm
FUNCTION
INPUTS
SOURCE
149 subroutine abi_d2ztrsm(side,uplo,transa,diag,m,n,alpha,a,lda,b,ldb,& 150 & x_cplx) 151 152 !Arguments------------------------------------- 153 character(len=1), intent(in) :: side,uplo,transa,diag 154 integer, intent(in) :: m,n,lda,ldb 155 complex(dpc), intent(in) :: alpha 156 real(dp),target, intent(in) :: a(lda,*) ! FIXME should be lda * x_cplx 157 real(dp),target, intent(inout) :: b(ldb,*) 158 !Only for lobpcgwf 159 integer, intent(in), optional :: x_cplx 160 161 !Local variables------------------------------- 162 integer :: cplx_ 163 #ifdef HAVE_LINALG_PLASMA 164 integer :: info 165 #endif 166 167 #ifdef DEV_LINALG_TIMING 168 real(dp) :: tsec(2) 169 call timab(TIMAB_XTRSM,1,tsec) 170 #endif 171 172 cplx_=1 ; if(PRESENT(x_cplx)) cplx_ = x_cplx 173 174 if (ABI_LINALG_PLASMA_ISON) then 175 #ifdef HAVE_LINALG_PLASMA 176 if(cplx_ == 2) then 177 info = PLASMA_ztrsm_c(side_plasma(side),uplo_plasma(uplo),trans_plasma(TRANSA),diag_plasma(diag),& 178 & m,n,alpha,c_loc(a),lda,c_loc(b),ldb) 179 else 180 info = PLASMA_dtrsm_c(side_plasma(side),uplo_plasma(uplo),trans_plasma(TRANSA),diag_plasma(diag),& 181 & m,n,real(alpha,dp),c_loc(a),lda,c_loc(b),ldb) 182 end if 183 #endif 184 else 185 if(cplx_ == 2) then 186 call ztrsm(side,uplo,transa,diag,m,n,alpha,a,lda,b,ldb) 187 else 188 call dtrsm(side,uplo,transa,diag,m,n,real(alpha,dp),a,lda,b,ldb) 189 end if 190 end if 191 192 #ifdef DEV_LINALG_TIMING 193 call timab(TIMAB_XTRSM,2,tsec) 194 #endif 195 196 end subroutine abi_d2ztrsm
m_abi_linalg/abi_d2ztrsm_3d [ Functions ]
[ Top ] [ m_abi_linalg ] [ Functions ]
NAME
abi_d2ztrsm_3d
FUNCTION
INPUTS
SOURCE
212 subroutine abi_d2ztrsm_3d(side,uplo,transa,diag,m,n,alpha,a,lda,b,ldb) 213 214 !Arguments------------------------------------- 215 character(len=1), intent(in) :: side,uplo,transa,diag 216 integer, intent(in) :: m,n,lda,ldb 217 complex(dpc), intent(in) :: alpha 218 real(dp), target,intent(in) :: a(2,lda,*) 219 real(dp), target,intent(inout) :: b(2,ldb,*) 220 221 !Local variables------------------------------- 222 #ifdef HAVE_LINALG_PLASMA 223 integer :: info 224 #endif 225 226 #ifdef DEV_LINALG_TIMING 227 real(dp) :: tsec(2) 228 call timab(TIMAB_XTRSM,1,tsec) 229 #endif 230 231 if (ABI_LINALG_PLASMA_ISON) then 232 #ifdef HAVE_LINALG_PLASMA 233 info = PLASMA_ztrsm_c(side_plasma(side),uplo_plasma(uplo),trans_plasma(TRANSA),diag_plasma(diag),& 234 & m,n,alpha,c_loc(a),lda,c_loc(b),ldb) 235 #endif 236 else 237 call ztrsm(side,uplo,transa,diag,m,n,alpha,a,lda,b,ldb) 238 end if 239 240 #ifdef DEV_LINALG_TIMING 241 call timab(TIMAB_XTRSM,2,tsec) 242 #endif 243 244 end subroutine abi_d2ztrsm_3d
m_abi_linalg/abi_dtrsm [ Functions ]
[ Top ] [ m_abi_linalg ] [ Functions ]
NAME
abi_dtrsm
FUNCTION
INPUTS
SOURCE
87 subroutine abi_dtrsm(side,uplo,transa,diag,m,n,alpha,a,lda,b,ldb,& 88 & x_cplx) 89 90 !Arguments------------------------------------- 91 character(len=1), intent(in) :: side,uplo,transa,diag 92 integer, intent(in) :: m,n,lda,ldb 93 real(dp), intent(in) :: alpha 94 real(dp),target, intent(in) :: a(lda,*) ! FIXME should be lda * x_cplx 95 real(dp),target, intent(inout) :: b(ldb,*) 96 !Only for lobpcgwf 97 integer, intent(in), optional :: x_cplx 98 99 !Local variables------------------------------- 100 integer :: cplx_ 101 #ifdef HAVE_LINALG_PLASMA 102 integer :: info 103 #endif 104 105 #ifdef DEV_LINALG_TIMING 106 real(dp) :: tsec(2) 107 call timab(TIMAB_XTRSM,1,tsec) 108 #endif 109 110 cplx_=1 ; if(PRESENT(x_cplx)) cplx_ = x_cplx 111 112 if (ABI_LINALG_PLASMA_ISON) then 113 #ifdef HAVE_LINALG_PLASMA 114 if(cplx_ == 2) then 115 info = PLASMA_ztrsm_c(side_plasma(side),uplo_plasma(uplo),trans_plasma(TRANSA),diag_plasma(diag),& 116 & m,n,cmplx(alpha,0.d0,dpc),c_loc(a),lda,c_loc(b),ldb) 117 else 118 info = PLASMA_dtrsm_c(side_plasma(side),uplo_plasma(uplo),trans_plasma(TRANSA),diag_plasma(diag),& 119 & m,n,alpha,c_loc(a),lda,c_loc(b),ldb) 120 end if 121 #endif 122 else 123 if(cplx_ == 2) then 124 call ztrsm(side,uplo,transa,diag,m,n,cmplx(alpha,0.d0,dpc),a,lda,b,ldb) 125 else 126 call dtrsm(side,uplo,transa,diag,m,n,alpha,a,lda,b,ldb) 127 end if 128 end if 129 130 #ifdef DEV_LINALG_TIMING 131 call timab(TIMAB_XTRSM,2,tsec) 132 #endif 133 134 end subroutine abi_dtrsm
m_abi_linalg/abi_xtrsm [ Functions ]
[ Top ] [ m_abi_linalg ] [ Functions ]
NAME
abi_xtrsm
FUNCTION
abi_xtrsm is the generic function that solve : * op( A )*X = alpha*B, or X*op( A ) = alpha*B, * * where alpha is a scalar, X and B are m by n matrices, A is a unit, or * non-unit, upper or lower triangular matrix and op( A ) is one of * * op( A ) = A or op( A ) = A**T. * * The matrix X is overwritten on B.
COPYRIGHT
Copyright (C) 2001-2024 ABINIT group (LNguyen,FDahm (CS)) This file is distributed under the terms of the GNU General Public License, see ~ABINIT/Infos/copyright or http://www.gnu.org/copyleft/gpl.txt .
SOURCE
m_abi_linalg/abi_ztrsm [ Functions ]
[ Top ] [ m_abi_linalg ] [ Functions ]
NAME
abi_ztrsm
FUNCTION
INPUTS
SOURCE
37 subroutine abi_ztrsm(side,uplo,transa,diag,m,n,alpha,a,lda,b,ldb) 38 39 !Arguments------------------------------------- 40 character(len=1), intent(in) :: side 41 character(len=1), intent(in) :: uplo 42 character(len=1), intent(in) :: transa 43 character(len=1), intent(in) :: diag 44 integer, intent(in) :: m,n,ldb,lda 45 complex(dpc), intent(in) :: alpha 46 complex(dpc),target,intent(in) :: a(lda,*) 47 complex(dpc),target,intent(inout) :: b(ldb,*) 48 49 !Local variables------------------------------- 50 #ifdef HAVE_LINALG_PLASMA 51 integer :: info 52 #endif 53 54 #ifdef DEV_LINALG_TIMING 55 real(dp) :: tsec(2) 56 call timab(TIMAB_XTRSM,1,tsec) 57 #endif 58 59 if (ABI_LINALG_PLASMA_ISON) then 60 #ifdef HAVE_LINALG_PLASMA 61 info = PLASMA_ztrsm_c(side_plasma(side),uplo_plasma(uplo),trans_plasma(transa),diag_plasma(diag),& 62 & m,n,alpha,c_loc(a),lda,c_loc(b),ldb) 63 #endif 64 else 65 call ztrsm(side,uplo,transa,diag,m,n,alpha,a,lda,b,ldb) 66 end if 67 68 #ifdef DEV_LINALG_TIMING 69 call timab(TIMAB_XTRSM,2,tsec) 70 #endif 71 72 end subroutine abi_ztrsm