TABLE OF CONTENTS
ABINIT/getchc [ Functions ]
NAME
getchc
FUNCTION
Compute <C_left|H|C> for input vectors |C> and |C_left>. Note that |C_left> can be an array of "ndat" wavefunctions if only the non-local part is computed
INPUTS
cpopt=flag defining the status of cwaveprj%cp(:)=<Proj_i|Cnk> scalars (PAW only) (same meaning as in nonlop.F90 routine) if cpopt=-1, <p_lmn|in> (and derivatives) are computed here (and not saved) if cpopt= 0, <p_lmn|in> are computed here and saved if cpopt= 1, <p_lmn|in> and first derivatives are computed here and saved if cpopt= 2 <p_lmn|in> are already in memory; if cpopt= 3 <p_lmn|in> are already in memory; first derivatives are computed here and saved if cpopt= 4 <p_lmn|in> and first derivatives are already in memory; cwavef(2,npw*my_nspinor)=planewave coefficients of wavefunction. cwavef_left(2,npw*my_nspinor)=planewave coefficients of wavefunction left. cwaveprj(natom,my_nspinor*(1+cpopt))= wave function projected on nl projectors cwaveprj_left(natom,my_nspinor*(1+cpopt))= wave function projected on nl projectors (for left WF) cwavef_r(2,n4,n5,n6,nspinor) = wave function in real space cwavef_left_r(2,n4,n5,n6,nspinor) = wave function in real space (for left WF) gs_ham <type(gs_hamiltonian_type)>=all data for the Hamiltonian to be applied lambda=factor to be used when computing <G|H-lambda.S|C> - only for sij_opt=-1 Typically lambda is the eigenvalue (or its guess) mpi_enreg=information about MPI parallelization ndat=number of left wavefunctions sij_opt= -PAW ONLY- if 0, only matrix elements <G|H|C> have to be computed (S=overlap) if 1, matrix elements <G|S|C> have to be computed in gsc in addition to ghc if -1, matrix elements <G|H-lambda.S|C> have to be computed in ghc (gsc not used) type_calc= option governing which part of Hamitonian is to be applied: 1: local part only 2: non-local+Fock+kinetic only (added to the existing Hamiltonian) 3: local + kinetic only (added to the existing Hamiltonian) ===== Optional inputs ===== [kg_fft_k(3,:)]=optional, (k+G) vector coordinates to be used for the FFT tranformation instead of the one contained in gs_ham datastructure. Typically used for real WF (in parallel) which are FFT-transformed 2 by 2. [kg_fft_kp(3,:)]=optional, (k^prime+G) vector coordinates to be used for the FFT tranformation [select_k]=optional, option governing the choice of k points to be used. gs_ham datastructure contains quantities needed to apply Hamiltonian in reciprocal space between 2 kpoints, k and k^prime (equal in most cases); if select_k=1, <k^prime|H|k> is applied [default] if select_k=2, <k|H|k^prime> is applied if select_k=3, <k|H|k> is applied if select_k=4, <k^prime|H|k^prime> is applied
OUTPUT
chc(2*ndat)=matrix elements <C_left|H|C> (if sij_opt>=0) or <C_left|H-lambda.S|C> (if sij_opt=-1)
SIDE EFFECTS
SOURCE
108 subroutine getchc(chc,cpopt,cwavef,cwavef_left,cwaveprj,cwaveprj_left,cwavef_r,cwavef_left_r,& 109 & gs_ham,lambda,mpi_enreg,ndat,& 110 & sij_opt,type_calc,& 111 & kg_fft_k,kg_fft_kp,select_k) ! optional arguments 112 113 !Arguments ------------------------------------ 114 !scalars 115 integer,intent(in) :: cpopt,ndat 116 integer,intent(in) :: sij_opt,type_calc 117 integer,intent(in),optional :: select_k 118 real(dp),intent(in) :: lambda 119 real(dp),intent(inout) :: chc(2*ndat) 120 type(MPI_type),intent(in) :: mpi_enreg 121 type(gs_hamiltonian_type),intent(inout),target :: gs_ham 122 !arrays 123 integer,intent(in),optional,target :: kg_fft_k(:,:),kg_fft_kp(:,:) 124 real(dp),intent(inout) :: cwavef(:,:),cwavef_left(:,:),cwavef_r(:,:,:,:,:),cwavef_left_r(:,:,:,:,:) 125 type(pawcprj_type),intent(inout),target :: cwaveprj(:,:),cwaveprj_left(:,:) 126 127 !Local variables------------------------------- 128 !scalars 129 integer,parameter :: re=1,im=2 130 integer :: choice,cpopt_here,i1,i2,i3,idat,idir 131 integer :: ig,igspinor,ispinor,ispinor_left,my_nspinor 132 integer :: nnlout,nffttot,npw,npw_k1,npw_k2,nspinortot,n1,n2,n3 133 integer :: paw_opt,select_k_,shift1,shift2,signs,tim_nonlop 134 logical :: k1_eq_k2,has_fock 135 logical :: nspinor1TreatedByThisProc,nspinor2TreatedByThisProc 136 character(len=500) :: msg 137 !arrays 138 integer, pointer :: gbound_k1(:,:),gbound_k2(:,:),kg_k1(:,:),kg_k2(:,:) 139 real(dp) :: enlout(ndat),enlout_im(ndat),lambda_ndat(ndat),tsec(2),z_tmp(2) 140 real(dp),allocatable :: gsc(:,:),gvnlxc(:,:) 141 real(dp), pointer :: kinpw_k1(:),kinpw_k2(:),kpt_k1(:),kpt_k2(:) 142 143 ! ********************************************************************* 144 145 DBG_ENTER("COLL") 146 147 !Keep track of total time spent in getchc: 148 call timab(1370,1,tsec) 149 150 !Select k-dependent objects according to select_k input parameter 151 select_k_=1;if (present(select_k)) select_k_=select_k 152 if (select_k_==KPRIME_H_K) then 153 ! <k^prime|H|k> 154 npw_k1 = gs_ham%npw_fft_k ; npw_k2 = gs_ham%npw_fft_kp 155 kpt_k1 => gs_ham%kpt_k ; kpt_k2 => gs_ham%kpt_kp 156 kg_k1 => gs_ham%kg_k ; kg_k2 => gs_ham%kg_kp 157 gbound_k1 => gs_ham%gbound_k ; gbound_k2 => gs_ham%gbound_kp 158 kinpw_k1 => gs_ham%kinpw_k ; kinpw_k2 => gs_ham%kinpw_kp 159 else if (select_k_==K_H_KPRIME) then 160 ! <k|H|k^prime> 161 npw_k1 = gs_ham%npw_fft_kp; npw_k2 = gs_ham%npw_fft_k 162 kpt_k1 => gs_ham%kpt_kp ; kpt_k2 => gs_ham%kpt_k 163 kg_k1 => gs_ham%kg_kp ; kg_k2 => gs_ham%kg_k 164 gbound_k1 => gs_ham%gbound_kp ; gbound_k2 => gs_ham%gbound_k 165 kinpw_k1 => gs_ham%kinpw_kp ; kinpw_k2 => gs_ham%kinpw_k 166 else if (select_k_==K_H_K) then 167 ! <k|H|k> 168 npw_k1 = gs_ham%npw_fft_k ; npw_k2 = gs_ham%npw_fft_k 169 kpt_k1 => gs_ham%kpt_k ; kpt_k2 => gs_ham%kpt_k 170 kg_k1 => gs_ham%kg_k ; kg_k2 => gs_ham%kg_k 171 gbound_k1 => gs_ham%gbound_k ; gbound_k2 => gs_ham%gbound_k 172 kinpw_k1 => gs_ham%kinpw_k ; kinpw_k2 => gs_ham%kinpw_k 173 else if (select_k_==KPRIME_H_KPRIME) then 174 ! <k^prime|H|k^prime> 175 npw_k1 = gs_ham%npw_fft_kp; npw_k2 = gs_ham%npw_fft_kp 176 kpt_k1 => gs_ham%kpt_kp ; kpt_k2 => gs_ham%kpt_kp 177 kg_k1 => gs_ham%kg_kp ; kg_k2 => gs_ham%kg_kp 178 gbound_k1 => gs_ham%gbound_kp ; gbound_k2 => gs_ham%gbound_kp 179 kinpw_k1 => gs_ham%kinpw_kp ; kinpw_k2 => gs_ham%kinpw_kp 180 end if 181 k1_eq_k2=(all(abs(kpt_k1(:)-kpt_k2(:))<tol8)) 182 if (.not.k1_eq_k2) then 183 ABI_ERROR('getchc is not implemented yet for k1/=k2') 184 end if 185 186 !Check sizes 187 my_nspinor=max(1,gs_ham%nspinor/mpi_enreg%nproc_spinor) 188 if (size(cwavef)<2*npw_k1*my_nspinor) then 189 msg='wrong size for cwavef!' 190 ABI_BUG(msg) 191 end if 192 if (gs_ham%usepaw==1.and.cpopt>=0) then 193 if (size(cwaveprj)<gs_ham%natom*my_nspinor) then 194 msg='wrong size for cwaveprj!' 195 ABI_BUG(msg) 196 end if 197 end if 198 if (gs_ham%usepaw==1) then 199 if (size(cwaveprj_left)<gs_ham%natom*my_nspinor*ndat) then 200 msg='wrong size for cwaveprj_left!' 201 ABI_BUG(msg) 202 end if 203 end if 204 205 !Eventually overwrite plane waves data for FFT 206 if (present(kg_fft_k)) then 207 kg_k1 => kg_fft_k ; kg_k2 => kg_fft_k 208 npw_k1=size(kg_k1,2) ; npw_k2=size(kg_k2,2) 209 end if 210 if (present(kg_fft_kp)) then 211 kg_k2 => kg_fft_kp ; npw_k2=size(kg_k2,2) 212 end if 213 214 !paral_kgb constraint 215 if (mpi_enreg%paral_kgb==1.and.(.not.k1_eq_k2)) then 216 msg='paral_kgb=1 not allowed for k/=k_^prime!' 217 ABI_BUG(msg) 218 end if 219 220 !Do we add Fock exchange term ? 221 has_fock=(associated(gs_ham%fockcommon)) 222 if (has_fock) then 223 ABI_BUG('Fock not implemented yet') 224 end if 225 ! if (has_fock) fock => gs_ham%fockcommon 226 227 !Parallelization over spinors management 228 if (mpi_enreg%paral_spinor==0) then 229 shift1=npw_k1;shift2=npw_k2 230 nspinor1TreatedByThisProc=.true. 231 nspinor2TreatedByThisProc=(nspinortot==2) 232 else 233 shift1=0;shift2=0 234 nspinor1TreatedByThisProc=(mpi_enreg%me_spinor==0) 235 nspinor2TreatedByThisProc=(mpi_enreg%me_spinor==1) 236 end if 237 238 npw=gs_ham%npw_k 239 nspinortot=gs_ham%nspinor 240 241 !============================================================ 242 ! Application of the local potential 243 !============================================================ 244 245 if ((type_calc==0).or.(type_calc==1).or.(type_calc==3)) then 246 247 call timab(1371,1,tsec) 248 ! Need a Vlocal 249 if (.not.associated(gs_ham%vlocal)) then 250 ABI_BUG("We need vlocal in gs_ham!") 251 end if 252 if (ndat>1) then 253 ABI_ERROR("ndat should be 1 for the local part") 254 end if 255 256 n1=gs_ham%ngfft(1) 257 n2=gs_ham%ngfft(2) 258 n3=gs_ham%ngfft(3) 259 nffttot=n1*n2*n3 260 chc = zero 261 ! Treat scalar local potentials 262 if (gs_ham%nvloc==1) then 263 if (gs_ham%istwf_k==2) then 264 do i3=1,n3 265 do i2=1,n2 266 do i1=1,n1 267 chc(1) = chc(1) + gs_ham%vlocal(i1,i2,i3,1)*cwavef_r(1,i1,i2,i3,1)*cwavef_left_r(1,i1,i2,i3,1) 268 end do 269 end do 270 end do 271 chc(2)=zero 272 else 273 do ispinor=1,my_nspinor 274 do i3=1,n3 275 do i2=1,n2 276 do i1=1,n1 277 z_tmp(1) = cwavef_r(1,i1,i2,i3,ispinor)*cwavef_left_r(1,i1,i2,i3,ispinor) & 278 & +cwavef_r(2,i1,i2,i3,ispinor)*cwavef_left_r(2,i1,i2,i3,ispinor) 279 z_tmp(2) = cwavef_r(2,i1,i2,i3,ispinor)*cwavef_left_r(1,i1,i2,i3,ispinor) & 280 & -cwavef_r(1,i1,i2,i3,ispinor)*cwavef_left_r(2,i1,i2,i3,ispinor) 281 chc(1) = chc(1) + gs_ham%vlocal(i1,i2,i3,1)*z_tmp(1) 282 chc(2) = chc(2) + gs_ham%vlocal(i1,i2,i3,1)*z_tmp(2) 283 end do 284 end do 285 end do 286 end do 287 end if 288 else ! nvloc = 4 289 do ispinor=1,my_nspinor 290 do ispinor_left=1,my_nspinor 291 do i3=1,n3 292 do i2=1,n2 293 do i1=1,n1 294 z_tmp(1) = cwavef_r(1,i1,i2,i3,ispinor)*cwavef_left_r(1,i1,i2,i3,ispinor_left) & 295 & +cwavef_r(2,i1,i2,i3,ispinor)*cwavef_left_r(2,i1,i2,i3,ispinor_left) 296 z_tmp(2) = cwavef_r(2,i1,i2,i3,ispinor)*cwavef_left_r(1,i1,i2,i3,ispinor_left) & 297 & -cwavef_r(1,i1,i2,i3,ispinor)*cwavef_left_r(2,i1,i2,i3,ispinor_left) 298 if (ispinor==ispinor_left) then 299 ! Then vloc is real : vloc_uu = vloc(1) and vloc_dd = vloc(2) 300 chc(1) = chc(1) + gs_ham%vlocal(i1,i2,i3,ispinor)*z_tmp(1) 301 chc(2) = chc(2) + gs_ham%vlocal(i1,i2,i3,ispinor)*z_tmp(2) 302 else if (ispinor==1.and.ispinor_left==2) then ! Psi(left)_d^* Psi_u vloc_ud 303 ! Otherwise vloc is complex Re(vloc_ud) = vloc(3) 304 ! Im(vloc_ud) = vloc(4) 305 ! vloc_du = (vloc_ud)^* 306 chc(1) = chc(1) + gs_ham%vlocal(i1,i2,i3,3)*z_tmp(1) + gs_ham%vlocal(i1,i2,i3,4)*z_tmp(2) 307 chc(2) = chc(2) + gs_ham%vlocal(i1,i2,i3,3)*z_tmp(1) - gs_ham%vlocal(i1,i2,i3,4)*z_tmp(2) 308 else if (ispinor==2.and.ispinor_left==1) then ! Psi(left)_u^* Psi_d vloc_ud 309 chc(1) = chc(1) + gs_ham%vlocal(i1,i2,i3,3)*z_tmp(1) - gs_ham%vlocal(i1,i2,i3,4)*z_tmp(2) 310 chc(2) = chc(2) + gs_ham%vlocal(i1,i2,i3,3)*z_tmp(1) + gs_ham%vlocal(i1,i2,i3,4)*z_tmp(2) 311 end if 312 end do 313 end do 314 end do 315 end do 316 end do 317 end if 318 chc = chc / dble(nffttot) 319 320 call timab(1371,2,tsec) 321 322 end if ! type_calc 323 324 if ((type_calc==0).or.(type_calc==2).or.(type_calc==3).or.(type_calc==4)) then 325 326 !============================================================ 327 ! Application of the non-local potential and the Fock potential 328 !============================================================ 329 330 if ((type_calc==0).or.(type_calc==2).or.(type_calc==4)) then 331 332 signs=1 ; choice=1 ; nnlout=1 ; idir=0 ; tim_nonlop=15 333 cpopt_here=-1;if (gs_ham%usepaw==1) cpopt_here=cpopt 334 ABI_MALLOC(gvnlxc,(0,0)) 335 ABI_MALLOC(gsc,(0,0)) 336 ! if (has_fock) then 337 ! if (gs_ham%usepaw==1) then 338 ! cpopt_here=max(cpopt,0) 339 ! if (cpopt<2) then 340 ! ABI_DATATYPE_ALLOCATE(cwaveprj_fock,(gs_ham%natom,my_nspinor*ndat)) 341 ! ABI_ALLOCATE(dimcprj,(gs_ham%natom)) 342 ! call pawcprj_getdim(dimcprj,gs_ham%natom,gs_ham%nattyp,gs_ham%ntypat,& 343 !& gs_ham%typat,fock%pawtab,'O') 344 ! call pawcprj_alloc(cwaveprj_fock,0,dimcprj) 345 ! ABI_DEALLOCATE(dimcprj) 346 ! else 347 ! cwaveprj_fock=>cwaveprj 348 ! end if 349 ! cwaveprj_nonlop=>cwaveprj_fock 350 ! else 351 ! cwaveprj_nonlop=>cwaveprj 352 ! cwaveprj_fock=>cwaveprj 353 ! end if 354 ! else 355 ! end if 356 paw_opt=gs_ham%usepaw ; if (sij_opt/=0) paw_opt=sij_opt+3 357 lambda_ndat = lambda 358 359 enlout=zero 360 enlout_im=zero 361 call nonlop(choice,cpopt_here,cwaveprj,enlout,gs_ham,idir,lambda_ndat,mpi_enreg,1,& 362 & nnlout,paw_opt,signs,gsc,tim_nonlop,cwavef,gvnlxc,select_k=select_k_,& 363 & cprjin_left=cwaveprj_left,enlout_im=enlout_im,ndat_left=ndat) 364 365 do idat=1,ndat 366 chc(2*idat-1) = chc(2*idat-1) + enlout(idat) 367 chc(2*idat ) = chc(2*idat ) + enlout_im(idat) 368 end do 369 370 ABI_FREE(gvnlxc) 371 ABI_FREE(gsc) 372 373 end if ! if(type_calc... 374 375 !============================================================ 376 ! Assemble kinetic, local, nonlocal and Fock contributions 377 !============================================================ 378 379 if (type_calc==0.or.type_calc==2.or.type_calc==3) then 380 381 if (ndat>1) then 382 ABI_ERROR("ndat should be 1 for the kinetic part") 383 end if 384 385 call timab(1372,1,tsec) 386 ! Add modified kinetic contributions 387 ! to <CP|H|C(n,k)>. 388 if (gs_ham%istwf_k==1) then 389 ! !!$OMP PARALLEL DO PRIVATE(igspinor) COLLAPSE(2) 390 do ispinor=1,my_nspinor 391 do ig=1,npw_k2 392 igspinor=ig+npw_k2*(ispinor-1) 393 if (kinpw_k2(ig)<huge(zero)*1.d-11) then 394 chc(1) = chc(1) + kinpw_k2(ig)*cwavef(re,igspinor)*cwavef_left(re,igspinor) 395 chc(1) = chc(1) + kinpw_k2(ig)*cwavef(im,igspinor)*cwavef_left(im,igspinor) 396 chc(2) = chc(2) + kinpw_k2(ig)*cwavef(im,igspinor)*cwavef_left(re,igspinor) 397 chc(2) = chc(2) - kinpw_k2(ig)*cwavef(re,igspinor)*cwavef_left(im,igspinor) 398 end if 399 end do ! ig 400 end do ! ispinor 401 else if (gs_ham%istwf_k==2.and.mpi_enreg%me_g0==1) then 402 if (kinpw_k2(1)<huge(zero)*1.d-11) then 403 chc(1) = chc(1) + kinpw_k2(1)*cwavef(re,1)*cwavef_left(re,1) 404 end if 405 do ig=2,npw_k2 406 if (kinpw_k2(ig)<huge(zero)*1.d-11) then 407 chc(1) = chc(1) + 2*kinpw_k2(ig)*cwavef(re,ig)*cwavef_left(re,ig) 408 chc(1) = chc(1) + 2*kinpw_k2(ig)*cwavef(im,ig)*cwavef_left(im,ig) 409 end if 410 end do 411 else 412 do ig=1,npw_k2 413 if (kinpw_k2(ig)<huge(zero)*1.d-11) then 414 chc(1) = chc(1) + 2*kinpw_k2(ig)*cwavef(re,ig)*cwavef_left(re,ig) 415 chc(1) = chc(1) + 2*kinpw_k2(ig)*cwavef(im,ig)*cwavef_left(im,ig) 416 end if 417 end do ! ig 418 end if 419 ! Special case of PAW + Fock : only return Fock operator contribution in gvnlxc 420 ! if (gs_ham%usepaw==1 .and. has_fock)then 421 ! gvnlxc=gvnlxc-gvnlc 422 ! ABI_DEALLOCATE(gvnlc) 423 ! endif 424 ! 425 ! if ((type_calc==0).or.(type_calc==2)) then 426 ! if (has_fock.and.gs_ham%usepaw==1.and.cpopt<2) then 427 ! call pawcprj_free(cwaveprj_fock) 428 ! ABI_DATATYPE_DEALLOCATE(cwaveprj_fock) 429 ! end if 430 ! end if 431 call timab(1372,2,tsec) 432 end if 433 434 end if ! type_calc 435 436 call timab(1370,2,tsec) 437 438 DBG_EXIT("COLL") 439 440 end subroutine getchc
ABINIT/getcsc [ Functions ]
NAME
getcsc
FUNCTION
Compute <C_left|S|C> for input vectors |C> and |C_left>. Note that |C_left> can be an array of "ndat" wavefunctions
INPUTS
cpopt=flag defining the status of cwaveprj%cp(:)=<Proj_i|Cnk> scalars (PAW only) (same meaning as in nonlop.F90 routine) if cpopt=-1, <p_lmn|in> (and derivatives) are computed here (and not saved) if cpopt= 0, <p_lmn|in> are computed here and saved if cpopt= 1, <p_lmn|in> and first derivatives are computed here and saved if cpopt= 2 <p_lmn|in> are already in memory; if cpopt= 3 <p_lmn|in> are already in memory; first derivatives are computed here and saved if cpopt= 4 <p_lmn|in> and first derivatives are already in memory; cwavef(2,npw*my_nspinor)=planewave coefficients of wavefunction. cwavef_left(2,npw*my_nspinor)=planewave coefficients of wavefunction left. cprj(natom,my_nspinor*(1+cpopt))= wave function projected on nl projectors cprj_left(natom,my_nspinor*(1+cpopt))= wave function projected on nl projectors (for left WF) gs_ham <type(gs_hamiltonian_type)>=all data for the Hamiltonian to be applied lambda=factor to be used when computing <G|H-lambda.S|C> - only for sij_opt=-1 Typically lambda is the eigenvalue (or its guess) mpi_enreg=information about MPI parallelization ndat=number of left wavefunctions mpi_enreg=information about MPI parallelization
OUTPUT
csc(2*ndat)=matrix elements <C_left|S|C>
SOURCE
478 subroutine getcsc(csc,cpopt,cwavef,cwavef_left,cprj,cprj_left,gs_ham,mpi_enreg,ndat,& 479 & tim_getcsc,& 480 & select_k) ! optional arguments 481 482 !Arguments ------------------------------------ 483 !scalars 484 integer,intent(in) :: cpopt,ndat,tim_getcsc 485 integer,intent(in),optional :: select_k 486 real(dp),intent(out) :: csc(2*ndat) 487 type(MPI_type),intent(in) :: mpi_enreg 488 type(gs_hamiltonian_type),intent(inout) :: gs_ham 489 !arrays 490 real(dp),intent(inout) :: cwavef(:,:) 491 real(dp),intent(inout),target :: cwavef_left(:,:) 492 type(pawcprj_type),intent(inout) :: cprj(:,:) 493 type(pawcprj_type),intent(inout),target :: cprj_left(:,:) 494 495 !Local variables------------------------------- 496 !scalars 497 integer :: choice,idat,idir,istwf_k 498 integer :: npw,nspinor,paw_opt,select_k_,signs,tim_nonlop,nnlout 499 !character(len=500) :: msg 500 !arrays 501 real(dp) :: tsec(2),real_csc(ndat) 502 real(dp),allocatable :: gsc(:,:),gvnlxc(:,:) 503 real(dp),allocatable :: enlout(:),enlout_im(:) 504 ! real(dp) :: dum 505 ! real(dp),pointer :: cwavef_left_idat(:,:) 506 ! ********************************************************************* 507 508 DBG_ENTER("COLL") 509 510 call timab(1360+tim_getcsc,1,tsec) 511 512 istwf_k = gs_ham%istwf_k 513 npw = gs_ham%npw_k 514 nspinor = gs_ham%nspinor 515 516 if (size(cwavef,2)/=npw*nspinor) then 517 ABI_BUG('Wrong size for cwavef') 518 end if 519 if (size(cwavef_left,2)/=npw*nspinor*ndat) then 520 ABI_BUG('Wrong size for cwavef_left') 521 end if 522 523 call timab(1361,1,tsec) 524 if (istwf_k==1) then 525 call zgemv('C',npw*nspinor,ndat,cone,cwavef_left,npw*nspinor,cwavef,1,czero,csc,1) 526 else ! nspinor==1 in that case 527 call dgemv('C',2*npw,ndat,one,cwavef_left,2*npw,cwavef,1,zero,real_csc,1) 528 do idat=1,ndat 529 csc(2*idat-1) = two*real_csc(idat) 530 csc(2*idat ) = zero 531 end do 532 if (istwf_k==2 .and. mpi_enreg%me_g0==1) then ! Gamma k-point and I have G=0 533 do idat=1,ndat 534 csc(2*idat-1) = csc(2*idat-1) - cwavef_left(1,1+npw*(idat-1))*cwavef(1,1) 535 end do 536 end if 537 end if 538 call timab(1361,2,tsec) 539 540 541 if (gs_ham%usepaw==1) then 542 543 if (size(cprj,2)/=nspinor) then 544 ABI_BUG('Wrong size for cprj') 545 end if 546 if (size(cprj_left,2)/=nspinor*ndat) then 547 ABI_BUG('Wrong size for cprj_left') 548 end if 549 550 select_k_=1;if (present(select_k)) select_k_=select_k 551 choice=1 ; nnlout=1 ; idir=0 ; tim_nonlop=16 ; paw_opt=3 552 ABI_MALLOC(gsc,(0,0)) 553 ABI_MALLOC(gvnlxc,(0,0)) 554 ABI_MALLOC(enlout ,(ndat)) 555 ABI_MALLOC(enlout_im,(ndat)) 556 enlout=zero 557 enlout_im=zero 558 signs=1 559 call nonlop(choice,cpopt,cprj,enlout,gs_ham,idir,(/zero/),mpi_enreg,1,& 560 & nnlout,paw_opt,signs,gsc,tim_nonlop,cwavef,gvnlxc,select_k=select_k_,& 561 & cprjin_left=cprj_left,enlout_im=enlout_im,ndat_left=ndat) 562 do idat=1,ndat 563 csc(2*idat-1) = csc(2*idat-1) + enlout(idat) 564 csc(2*idat ) = csc(2*idat ) + enlout_im(idat) 565 end do 566 ABI_FREE(gsc) 567 ABI_FREE(gvnlxc) 568 ABI_FREE(enlout ) 569 ABI_FREE(enlout_im) 570 571 end if 572 573 call timab(1360+tim_getcsc,2,tsec) 574 575 DBG_EXIT("COLL") 576 577 end subroutine getcsc
ABINIT/m_getchc [ Modules ]
NAME
m_getchc
FUNCTION
Compute <G|H|C> for input vector |C> expressed in reciprocal space;
COPYRIGHT
Copyright (C) 1998-2024 ABINIT group (DCA, XG, GMR, LSI, MT) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
SOURCE
16 #if defined HAVE_CONFIG_H 17 #include "config.h" 18 #endif 19 20 #include "abi_common.h" 21 22 module m_getchc 23 24 use defs_basis 25 use m_errors 26 use m_abicore 27 use m_xmpi 28 29 use defs_abitypes, only : mpi_type 30 use m_time, only : timab 31 use m_pawcprj, only : pawcprj_type!, pawcprj_alloc, pawcprj_free, pawcprj_getdim 32 ! use m_bandfft_kpt, only : bandfft_kpt, bandfft_kpt_get_ikpt 33 use m_hamiltonian, only : gs_hamiltonian_type, KPRIME_H_K, K_H_KPRIME, K_H_K, KPRIME_H_KPRIME 34 use m_nonlop, only : nonlop 35 ! use m_fock, only : fock_common_type, fock_get_getghc_call 36 ! use m_fock_getghc, only : fock_getghc, fock_ACE_getghc 37 use m_cgtools, only : dotprod_g 38 39 implicit none 40 41 private