TABLE OF CONTENTS


m_abi_linalg/abi_d2zgemm [ Functions ]

[ Top ] [ m_abi_linalg ] [ Functions ]

NAME

 abi_d2zgemm

FUNCTION

INPUTS

PARENTS

CHILDREN

SOURCE

115 subroutine abi_d2zgemm(TRANSA,TRANSB,M,N,K,ALPHA,A,LDA,B,LDB,BETA,C,LDC,&
116 &       x_cplx)
117 
118 !This section has been created automatically by the script Abilint (TD).
119 !Do not modify the following lines by hand.
120 #undef ABI_FUNC
121 #define ABI_FUNC 'abi_d2zgemm'
122 !End of the abilint section
123 
124  implicit none
125 
126 !Arguments ------------------------------------
127  character(len=1), intent(in) :: transa
128  character(len=1), intent(in) :: transb
129  integer, intent(in) :: lda
130  integer, intent(in) :: ldb
131  integer, intent(in) :: ldc
132  integer, intent(in) :: m
133  integer, intent(in) :: n
134  integer, intent(in) :: k
135  complex(dpc),intent(in) :: alpha
136  complex(dpc),intent(in) :: beta
137  real(dp),target, intent(in)    :: a(lda,*)    ! FIXME should be lda * x_cplx
138  real(dp),target, intent(in)    :: b(ldb,*)
139  real(dp),target, intent(inout) :: c(ldc,*)
140 !Only for lobpcgwf
141  integer, intent(in), optional :: x_cplx
142 
143 !Local variables-------------------------------
144  integer  :: cplx_,info
145 #ifdef DEV_LINALG_TIMING
146  real(dp) :: tsec(2)
147  call timab(TIMAB_XGEMM,1,tsec)
148 #endif
149 
150  cplx_=1 ; if(PRESENT(x_cplx)) cplx_ = x_cplx
151 
152  if (XPLASMA_ISON) then
153    info = -1
154 #ifdef HAVE_LINALG_PLASMA
155    !write(std_out,*) "Will call plasma_[zd]gemm"
156    if (cplx_ == 2) then
157       info = PLASMA_zgemm_c(trans_plasma(TRANSA),trans_plasma(TRANSB),M,N,K,&
158 &      ALPHA,c_loc(A),LDA,c_loc(B),LDB,BETA,c_loc(C),LDC)
159    else
160       info = PLASMA_dgemm_c(trans_plasma(TRANSA),trans_plasma(TRANSB),M,N,K,&
161 &      real(ALPHA,dp),c_loc(A),LDA,c_loc(B),LDB,real(BETA,dp),c_loc(C),LDC)
162    end if
163 #endif
164    ABI_CHECK(info==0,"PLASMA_[z,d]gemm_c returned info !=0")
165  else
166    if (cplx_ == 2) then
167      if (use_zgemm3m(m,n,k)) then
168        call _ZGEMM3M(TRANSA,TRANSB,M,N,K,ALPHA,A,LDA,B,LDB,BETA,C,LDC)
169      else
170        call zgemm(TRANSA,TRANSB,M,N,K,ALPHA,A,LDA,B,LDB,BETA,C,LDC)
171      end if
172    else
173      call dgemm(TRANSA,TRANSB,M,N,K,real(ALPHA,dp),A,LDA,B,LDB,real(BETA,dp),C,LDC)
174    end if
175  end if
176 
177 #ifdef DEV_LINALG_TIMING
178  call timab(TIMAB_XGEMM,2,tsec)
179 #endif
180 
181 end subroutine abi_d2zgemm

m_abi_linalg/abi_xgemm [ Functions ]

[ Top ] [ m_abi_linalg ] [ Functions ]

NAME

  abi_xgemm

FUNCTION

  abi_xgemm is the generic function that solve :
 *     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 m by k matrix,  op( B )  a  k by n matrix and  C an m by n matrix.

COPYRIGHT

  Copyright (C) 2001-2018 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_zgemm_2d [ Functions ]

[ Top ] [ m_abi_linalg ] [ Functions ]

NAME

 abi_zgemm_2d

FUNCTION

INPUTS

PARENTS

SOURCE

42  subroutine abi_zgemm_2d(TRANSA,TRANSB,M,N,K,ALPHA,A,LDA,B,LDB,BETA,C,LDC)
43 
44 !This section has been created automatically by the script Abilint (TD).
45 !Do not modify the following lines by hand.
46 #undef ABI_FUNC
47 #define ABI_FUNC 'abi_zgemm_2d'
48 
49 !End of the abilint section
50 
51  implicit none
52 
53  !Arguments------------------------------------
54  character(len=1),intent(in) :: TRANSA
55  character(len=1),intent(in) :: TRANSB
56  integer,intent(in) :: K
57  integer,intent(in) :: LDA
58  integer,intent(in) :: LDB
59  integer,intent(in) :: LDC
60  integer,intent(in) :: M
61  integer,intent(in) :: N
62  complex(dpc),intent(in) :: ALPHA
63  complex(dpc),intent(in) :: BETA
64  complex(dpc),target,intent(in) :: A(lda,*)
65  complex(dpc),target,intent(in) :: B(ldb,*)
66  complex(dpc),target,intent(inout) :: C(ldc,*)
67 
68  integer :: info
69 #ifdef DEV_LINALG_TIMING
70  real(dp) :: tsec(2)
71  call timab(TIMAB_XGEMM,1,tsec)
72 #endif
73 
74  if (XPLASMA_ISON) then
75    info = -1
76 #ifdef HAVE_LINALG_PLASMA
77    !write(std_out,*)"Will call PLASMA_zgemm_c"
78    info = PLASMA_zgemm_c(trans_plasma(TRANSA),trans_plasma(TRANSB),M,N,K,ALPHA,&
79 &    c_loc(A),LDA,c_loc(B),LDB,BETA,c_loc(C),LDC)
80 #endif
81    ABI_CHECK(info==0,"PLASMA_zgemm_c returned info !=0")
82  else
83    if (use_zgemm3m(m,n,k)) then
84      call _ZGEMM3M(TRANSA,TRANSB,M,N,K,ALPHA,A,LDA,B,LDB,BETA,C,LDC)
85    else
86      call zgemm(TRANSA,TRANSB,M,N,K,ALPHA,A,LDA,B,LDB,BETA,C,LDC)
87    end if
88  end if
89 
90 #ifdef DEV_LINALG_TIMING
91  call timab(TIMAB_XGEMM,2,tsec)
92 #endif
93 
94 end subroutine abi_zgemm_2d