TABLE OF CONTENTS


ABINIT/ortho_reim [ Functions ]

[ Top ] [ Functions ]

NAME

 ortho_reim

FUNCTION

 This routine computes the overlap of two wavefunctions (for a given number of bands)
 and orthonormalizes it:
      - Computes the products of two rectangular matrices
         containing the wavefunctions psi and S.psi (where S is the
         overlap (with the PAW terms if necessary)).
      - Does a Cholesky decomposition of this overlap
      - rotates the initial matrix blockvectorx by the triangular matrix to
         have an orthonormal set of wavefunctions

 This version operates on arrays in which the real and the imaginary part
 are packed together (real parts first, them imaginary parts), used when istwfk=2

INPUTS

  blockvectorbx = matrix of dimension (blocksize,vectsize)
                  (e.g. block of overlap*wavefunction)
  blocksize     = dimension of matrices (e.g number of bands)
  spaceComm     = communicator used for MPI parallelization
  vectsize      = dimension of matrices (e.g number of G vector)

OUTPUT

  sqgram        = Choleski decomposition of transpose(blockvector)*blockvectorx

SIDE EFFECTS

  blockvectorx  = on input, matrix of dimension (vectsize,blocksize)
                  (e.g block of wavefunction)
  blockvectorx  = on output, orthonormalized wavefunction.

PARENTS

      lobpcgIIwf,m_lobpcg,m_lobpcgIIIwf,pw_orthon

CHILDREN

      abi_xgemm,abi_xpotrf,abi_xtrsm,wrtout,xmpi_sum

SOURCE

161 subroutine ortho_reim(blockvectorx,blockvectorbx,blocksize,spaceComm,sqgram,vectsize)
162 
163 !This section has been created automatically by the script Abilint (TD).
164 !Do not modify the following lines by hand.
165 #undef ABI_FUNC
166 #define ABI_FUNC 'ortho_reim'
167 !End of the abilint section
168 
169  implicit none
170 
171 !Arguments ------------------------------------
172 !scalars
173  integer,intent(in) :: blocksize,vectsize,spaceComm
174 !arrays
175  real(dp),intent(in) :: blockvectorbx(vectsize,blocksize)
176  real(dp),intent(inout) :: blockvectorx(vectsize,blocksize)
177  real(dp),intent(out) :: sqgram(blocksize,blocksize)
178 
179 !Local variables-------------------------------
180 !scalars
181  integer :: ierr,info
182  character(len=500) :: message
183 
184 ! *********************************************************************
185 
186  call abi_xgemm('t','n',blocksize,blocksize,vectsize,cone,blockvectorx,&
187 &   vectsize,blockvectorbx,vectsize,czero,sqgram,blocksize)
188 
189  call xmpi_sum(sqgram,spaceComm,ierr)
190 
191  !Cholesky factorization of sqgram (ouside upper Triangular of sqgram)
192  call abi_d2zpotrf('u',blocksize,sqgram,blocksize,info) !vz_d
193 
194  if (info /= 0 )  then
195    write(message,'(a,i0)')'dpotrf, info=',info
196    MSG_ERROR(message)
197  end if
198 
199 !Find X  X*sqgram=blockvectorx
200  call abi_xtrsm('r','u','n','n',vectsize,blocksize,one,sqgram,blocksize,blockvectorx,vectsize)
201 
202 end subroutine ortho_reim

ABINIT/zorthonormalize [ Functions ]

[ Top ] [ Functions ]

NAME

 zorthonormalize

FUNCTION

 This routine computes the overlap of two complex wavefunctions (for a given number of bands)
 and orthonormalizes it:
      - Computes the products of two rectangular matrices
         containing the wavefunctions psi and S.psi (where S is the
         overlap (with the PAW terms if necessary)).
      - Does a Cholesky decomposition of this overlap
      - rotates the initial matrix blockvectorx by the triangular matrix to
         have an orthonormal set of wavefunctions

INPUTS

  blockvectorbx = matrix of dimension (blocksize,vectsize)
                  (e.g. block of overlap*wavefunction)
  blocksize     = dimension of matrices (e.g number of bands)
  spaceComm     = communicator used for MPI parallelization
  vectsize      = dimension of matrices (e.g number of G vector)

OUTPUT

  sqgram        = Choleski decomposition of transpose(blockvector)*blockvectorx

SIDE EFFECTS

  blockvectorx  = on input, matrix of dimension (vectsize,blocksize)
                  (e.g block of wavefunction)
  blockvectorx  = on output, orthonormalized wavefunction.

PARENTS

      lobpcgccIIIwf,lobpcgccIIwf,m_lobpcg,pw_orthon

CHILDREN

      wrtout,xmpi_sum,zgemm,zpotrf,ztrsm

SOURCE

243 subroutine zorthonormalize(blockvectorx,blockvectorbx,blocksize,spaceComm,sqgram,vectsize)
244 
245 !This section has been created automatically by the script Abilint (TD).
246 !Do not modify the following lines by hand.
247 #undef ABI_FUNC
248 #define ABI_FUNC 'zorthonormalize'
249 !End of the abilint section
250 
251  implicit none
252 
253 !Arguments ------------------------------------
254 !scalars
255  integer,intent(in) :: blocksize,spaceComm,vectsize
256 !arrays
257  complex(dpc),intent(in) :: blockvectorbx(vectsize,blocksize)
258  complex(dpc),intent(inout) :: blockvectorx(vectsize,blocksize)
259  complex(dpc),intent(out) :: sqgram(blocksize,blocksize)
260 
261 !Local variables-------------------------------
262 !scalars
263  integer :: ierr,info
264  character(len=500) :: message
265 
266 ! *********************************************************************
267 
268  call abi_xgemm('c','n',blocksize,blocksize,vectsize,cone,blockvectorx,&
269 & vectsize,blockvectorbx,vectsize,czero,sqgram,blocksize)
270 
271  call xmpi_sum(sqgram,spaceComm,ierr)
272 
273  call abi_xpotrf('u',blocksize,sqgram,blocksize,info)
274 
275  if (info /= 0 )  then
276    write(message,'(a,i0)')'zpotrf, info=',info
277    MSG_ERROR(message)
278  end if
279 
280  call abi_xtrsm('r','u','n','n',vectsize,blocksize,cone,sqgram,blocksize,blockvectorx,vectsize)
281 
282 end subroutine zorthonormalize

m_abi_linalg/abi_xorthonormalize [ Functions ]

[ Top ] [ m_abi_linalg ] [ Functions ]

NAME

  abi_xorthonormalize

FUNCTION

  abi_xorthonormalize is the generic function for computing the
  overlap of two complex wavefunctions (for a given number of bands)
  and orthonormalizes it:

COPYRIGHT

  Copyright (C) 2001-2018 ABINIT group (LNguyen,FDahm (CS), FBottin, GZ, AR, MT)
  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/xorthonormalize [ Functions ]

[ Top ] [ m_abi_linalg ] [ Functions ]

NAME

  xorthonormalize

FUNCTION

  This routine computes the overlap of two complex wavefunctions (for a given number of bands)
  and orthonormalizes it:
      - Computes the products of two rectangular matrices
         containing the wavefunctions psi and S.psi (where S is the
         overlap (with the PAW terms if necessary)).
      - Does a Cholesky decomposition of this overlap
      - rotates the initial matrix blockvectorx by the triangular matrix to
         have an orthonormal set of wavefunctions

INPUTS

  blockvectorbx = matrix of dimension (blocksize,vectsize)
                  (e.g. block of overlap*wavefunction)
  blocksize     = dimension of matrices (e.g number of bands)
  spaceComm     = communicator used for MPI parallelization
  vectsize      = dimension of matrices (e.g number of G vector)

OUTPUT

  sqgram        = Choleski decomposition of transpose(blockvector)*blockvectorx

SIDE EFFECTS

  blockvectorx  = on input, matrix of dimension (vectsize,blocksize)
                  (e.g block of wavefunction)
  blockvectorx  = on output, orthonormalized wavefunction.

PARENTS

      lobpcgwf

CHILDREN

SOURCE

 58 subroutine xorthonormalize(blockvectorx,blockvectorbx,blocksize,spaceComm,sqgram,vectsize,&
 59 &                          x_cplx,timopt,tim_xortho) ! optional arguments
 60 
 61 
 62 !This section has been created automatically by the script Abilint (TD).
 63 !Do not modify the following lines by hand.
 64 #undef ABI_FUNC
 65 #define ABI_FUNC 'xorthonormalize'
 66 !End of the abilint section
 67 
 68  implicit none
 69 !Arguments ------------------------------------
 70 !scalars
 71  integer,intent(in) :: blocksize,vectsize,spaceComm,x_cplx
 72  integer, intent(in), optional :: timopt,tim_xortho
 73  !arrays
 74  real(dp),intent(in) :: blockvectorbx(vectsize,blocksize)
 75  real(dp),intent(inout) :: blockvectorx(vectsize,blocksize)
 76  real(dp),intent(out) :: sqgram(x_cplx*blocksize,blocksize)
 77 
 78 !Local variables-------------------------------
 79  real(dp) :: tsec(2)
 80  integer  :: ierr,info
 81  character(len=500) :: message
 82  character, dimension(2) :: cparam
 83 
 84  ! *********************************************************************
 85 
 86  if (present(tim_xortho).and.present(timopt)) then
 87    if(abs(timopt)==3) then
 88      call timab(tim_xortho,1,tsec)
 89    end if
 90  end if
 91 
 92  cparam(1)='t'
 93  cparam(2)='c'
 94 
 95  call abi_xgemm(cparam(x_cplx),'n',blocksize,blocksize,vectsize,cone,blockvectorx,&
 96 &   vectsize,blockvectorbx,vectsize,czero,sqgram,blocksize,x_cplx=x_cplx)
 97 
 98  call xmpi_sum(sqgram,spaceComm,ierr)
 99 
100  !Cholesky factorization of sqgram (ouside upper Triangular of sqgram)
101  call abi_xpotrf('u',blocksize,sqgram,blocksize,info,x_cplx=x_cplx)
102 
103  if (info /= 0 )  then
104    write(message,'(a,i0)')'abi_xpotrf, info=',info
105    MSG_ERROR(message)
106  end if
107 
108  !Find X  X*sqgram=blockvectorx
109  call abi_xtrsm('r','u','n','n',vectsize,blocksize,cone,sqgram,blocksize,&
110 &   blockvectorx,vectsize,x_cplx=x_cplx)
111 
112  if (present(tim_xortho).and.present(timopt)) then
113    if(abs(timopt)==3) then
114      call timab(tim_xortho,2,tsec)
115    end if
116  end if
117 
118 end subroutine xorthonormalize