TABLE OF CONTENTS
ABINIT/mlwfovlp_qp [ Functions ]
NAME
mlwfovlp_qp
FUNCTION
Routine which computes replaces LDA wave functions and eigenvalues with GW quasiparticle ones using previously computed qp wave functions in LDA bloch function representation for Wannier code (www.wannier.org f90 version).
COPYRIGHT
Copyright (C) 2008-2018 ABINIT group (DRH) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
dtset <type(dataset_type)>=all input variables for this dataset dtfil <type(datafiles_type)>=variables related to files mband=maximum number of bands mcg=size of wave-functions array (cg) =mpw*nspinor*mband*mkmem*nsppol mcprj=size of projected wave-functions array (cprj) =nspinor*mband*mkmem*nsppol mkmem =number of k points treated by this node. mpw=maximum dimensioned size of npw. natom=number of atoms in cell. nkpt=number of k points. npwarr(nkpt)=number of planewaves in basis at this k point nspden=number of spin-density components nsppol=1 for unpolarized, 2 for spin-polarized rprimd(3,3)=dimensional primitive translations for real space (bohr) Hdr<Hdr_type>=The abinit header. MPI_enreg=information about MPI parallelization Cprj_BZ(natom,nspinor*mband*mkmem*nsppol)= <p_lmn|Cnk> coefficients for each WF |Cnk> and each |p_lmn> non-local projector
OUTPUT
SIDE EFFECTS
cg(2,mcg)=planewave coefficients of wavefunctions replaced by quasiparticle wavefunctions eigen(mband*nkpt*nsppol)=array for holding eigenvalues replaced by qp eigenvalues(hartree)
NOTES
Number of bands for wannier calculation must be identical to number used for gw calculation. Bands not wanted for wannier calculation must be excluded in exclude_band statement in wannier90.win file. Full plane-wave basis for LDA wavefunctions must be used in GW calculation, or inaccuracies may result. This is at best a beta version of this code, with little consistency checking, so the user must be very careful or the results may be invalid.
PARENTS
outscfcv
CHILDREN
crystal_free,crystal_from_hdr,destroy_mpi_enreg,ebands_free,ebands_init initmpi_seq,kmesh_free,kmesh_init,listkk,metric,pawcprj_getdim,rdgw rdqps,sort_dp,update_cprj,wrtout,zgemm
SOURCE
61 #if defined HAVE_CONFIG_H 62 #include "config.h" 63 #endif 64 65 #include "abi_common.h" 66 67 subroutine mlwfovlp_qp(cg,Cprj_BZ,dtset,dtfil,eigen,mband,mcg,mcprj,mkmem,mpw,natom,& 68 & nkpt,npwarr,nspden,nsppol,ntypat,Hdr,Pawtab,rprimd,MPI_enreg) 69 70 use defs_basis 71 use defs_datatypes 72 use defs_abitypes 73 use defs_wannier90 74 use m_errors 75 use m_profiling_abi 76 77 use m_mpinfo, only : destroy_mpi_enreg 78 use m_pawtab, only : pawtab_type 79 use m_pawcprj, only : pawcprj_type, paw_overlap, pawcprj_getdim 80 use m_numeric_tools, only : isordered 81 use m_geometry, only : metric 82 use m_crystal, only : crystal_t, crystal_free 83 use m_crystal_io, only : crystal_from_hdr 84 use m_bz_mesh, only : kmesh_t, kmesh_init, kmesh_free 85 use m_ebands, only : ebands_init, ebands_free 86 use m_qparticles, only : rdqps, rdgw 87 use m_sort, only : sort_dp 88 89 !This section has been created automatically by the script Abilint (TD). 90 !Do not modify the following lines by hand. 91 #undef ABI_FUNC 92 #define ABI_FUNC 'mlwfovlp_qp' 93 use interfaces_14_hidewrite 94 use interfaces_51_manage_mpi 95 use interfaces_56_recipspace 96 use interfaces_70_gw, except_this_one => mlwfovlp_qp 97 !End of the abilint section 98 99 implicit none 100 !Arguments ------------------------------------ 101 !scalars 102 integer,intent(in) :: mband,mcg,mcprj,mkmem,mpw,nkpt,nspden,natom,ntypat 103 integer,intent(in) :: nsppol 104 type(dataset_type),intent(in) :: dtset 105 type(datafiles_type),intent(in) :: dtfil 106 type(Hdr_type),intent(in) :: Hdr 107 type(MPI_type),intent(in) :: MPI_enreg 108 type(pawcprj_type),target,intent(inout) :: Cprj_BZ(natom,mcprj) 109 type(Pawtab_type),intent(in) :: Pawtab(ntypat*Dtset%usepaw) 110 !arrays 111 integer,intent(in) :: npwarr(nkpt) 112 real(dp),intent(inout) :: cg(2,mcg) 113 real(dp),intent(inout) :: eigen(mband*nkpt*nsppol) 114 real(dp),intent(in) :: rprimd(3,3) 115 116 !Local variables------------------------------- 117 !scalars 118 integer,parameter :: from_QPS_FILE=1,from_GW_FILE=2 119 integer :: sppoldbl,timrev,bantot_ibz,ikibz,ikbz,dimrho 120 integer :: iband,icg,icg_shift,ii,ipw,isppol,my_nspinor,nband_k,ord_iband 121 integer :: nfftot,ikpt,irzkpt,npw_k,ikg 122 integer :: nscf,nbsc,itimrev,band_index,nkibz,nkbz 123 integer :: gw_timrev,input !,jb_idx,ib_idx,ijpack, jband, 124 integer :: nprocs,ios 125 real(dp) :: TOL_SORT=tol12 126 real(dp) :: dksqmax,ucvol !ortho_err, 127 logical :: ltest,qpenek_is_ordered,g0w0_exists 128 character(len=500) :: msg 129 character(len=fnlen) :: gw_fname 130 type(ebands_t) :: QP_bst 131 type(crystal_t) :: Cryst 132 type(kmesh_t) :: Kibz_mesh 133 type(MPI_type) :: MPI_enreg_seq 134 !arrays 135 integer :: indkk(nkpt,6),my_ngfft(18) 136 integer,allocatable :: npwarr_ibz(:),nband_ibz(:),ibz2bz(:,:),istwfk_ibz(:) 137 integer,allocatable :: dimlmn(:),iord(:),nattyp_dum(:) 138 real(dp) :: gmet(3,3),gprimd(3,3),rmet(3,3) !,paw_ovlp(2) 139 real(dp),allocatable :: qp_rhor(:,:),sorted_qpene(:) 140 real(dp),allocatable :: kibz(:,:),wtk_ibz(:) 141 real(dp),allocatable :: doccde_ibz(:),occfact_ibz(:),eigen_ibz(:) 142 real(dp),allocatable :: igwene(:,:,:) 143 complex(dpc),allocatable :: m_lda_to_qp(:,:,:,:),m_lda_to_qp_BZ(:,:,:,:) !,ortho(:) 144 complex(dpc),allocatable :: m_tmp(:,:),cg_k(:,:),cg_qpk(:,:) 145 type(Pawrhoij_type),allocatable :: prev_Pawrhoij(:) 146 !type(pawcprj_type),pointer :: Cp1(:,:),Cp2(:,:) 147 148 !************************************************************************ 149 150 ABI_UNUSED(mkmem) 151 152 DBG_ENTER("COLL") 153 154 write(msg,'(17a)')ch10,& 155 & ' mlwfovlp_qp: WARNING',ch10,& 156 & ' The input *_WFK file of LDA wavefunctions to be converted',ch10,& 157 & ' to GW quasiparticle wavefunctions MUST have been written in',ch10,& 158 & ' the run that produced the GW *_KSS file using kssform 3,',ch10,& 159 & ' the ONLY value of kssform permitted for GW Wannier functions.',ch10,& 160 & ' Otherwise, the *_QPS file needed here will be inconsistent,',ch10,& 161 & ' and the output quasiparticle wavefunctions will be garbage.',ch10,& 162 & ' No internal check that can verify this is presently possible.',ch10 163 call wrtout(std_out,msg,'COLL') 164 165 ! === Some features are not implemented yet === 166 ABI_CHECK(Dtset%nspinor==1,'nspinor==2 not implemented') 167 ABI_CHECK(Dtset%nsppol==1,'nsppol==2 not implemented, check wannier90') 168 ltest=ALL(Dtset%nband(1:Dtset%nkpt*Dtset%nsppol)==Dtset%nband(1)) 169 ABI_CHECK(ltest,'nband(:) should be constant') 170 ! 171 ! MPI initialization 172 nprocs=MPI_enreg%nproc_cell 173 174 if (nprocs/=1) then 175 MSG_ERROR("mlwfovlp_qp not programmed for parallel execution") 176 end if 177 178 ! Compute reciprocal space metric gmet for unit cell of disk wf 179 call metric(gmet,gprimd,-1,rmet,rprimd,ucvol) 180 181 ! Compute k points from gw irreducible set equivalent to full-zone wannier set 182 sppoldbl=1 ; timrev=1 ; my_nspinor=max(1,Dtset%nspinor/MPI_enreg%nproc_spinor) 183 call listkk(dksqmax,gmet,indkk,dtset%kptgw,dtset%kpt,dtset%nkptgw,nkpt,& 184 & dtset%nsym,sppoldbl,dtset%symafm,dtset%symrel,timrev) 185 186 if (dksqmax>tol8) then 187 write(msg,'(5a)')& 188 & 'Set of GW irreducible-zone kptgw in input file is inconsistent',ch10,& 189 & 'with full-zone set being used for wannier90 setup.',ch10,& 190 & 'Action : correct input data' 191 MSG_ERROR(msg) 192 end if 193 ! 194 ! === Initialize object defining the Band strucuture === 195 ! * Initialize with KS results using IBZ indexing. 196 ! * After rdqps, QP_bst will contain the QP amplitudes. 197 nkibz=Dtset%nkptgw 198 ABI_MALLOC(kibz,(3,nkibz)) 199 ABI_MALLOC(wtk_ibz,(nkibz)) 200 kibz=Dtset%kptgw(:,1:Dtset%nkptgw) 201 202 ! MG: This part is needed to get the IBZ weight that will be reported 203 ! on ab_out thus we should be consistent. Ideally Cryst should be 204 ! one of the basic abinit objects and it should be passed to this routine. 205 206 gw_timrev=1; if (timrev==1) gw_timrev=2 !different conventions are used in GW and abinit!! 207 208 call crystal_from_hdr(Cryst,Hdr,gw_timrev) 209 call kmesh_init(Kibz_mesh,Cryst,nkibz,kibz,Dtset%kptopt) 210 wtk_ibz=Kibz_mesh%wt 211 call crystal_free(Cryst) 212 call kmesh_free(Kibz_mesh) 213 214 ABI_MALLOC(ibz2bz,(nkibz,6)) 215 call listkk(dksqmax,gmet,ibz2bz,dtset%kpt,dtset%kptgw,nkpt,dtset%nkptgw,& 216 & dtset%nsym,sppoldbl,dtset%symafm,dtset%symrel,timrev) 217 218 ltest=ALL(ibz2bz(:,2)==1) 219 ABI_CHECK(ltest,'Not able to found irreducible points in the BZ set!') 220 221 if (dksqmax>tol8) then 222 write(msg,'(5a)')& 223 & 'Set of GW irreducible-zone kptgw in input file is inconsistent',ch10,& 224 & 'with full-zone set being used for wannier90 setup.',ch10,& 225 & 'Action: correct input data' 226 MSG_ERROR(msg) 227 end if 228 229 ABI_MALLOC(npwarr_ibz,(nkibz)) 230 ABI_MALLOC(istwfk_ibz,(nkibz)) 231 ABI_MALLOC(nband_ibz,(nkibz*nsppol)) 232 233 do isppol=1,nsppol 234 do ikibz=1,nkibz 235 ikbz=ibz2bz(ikibz+(sppoldbl-1)*(isppol-1)*nkibz,1) 236 npwarr_ibz(ikibz)= npwarr(ikbz) 237 istwfk_ibz(ikibz)=Dtset%istwfk(ikbz) 238 nband_ibz(ikibz+(isppol-1)*nkibz)=Dtset%nband(ikbz+(isppol-1)*nkpt) 239 end do 240 end do 241 242 bantot_ibz=SUM(nband_ibz) 243 ABI_MALLOC(doccde_ibz,(bantot_ibz)) 244 ABI_MALLOC(eigen_ibz,(bantot_ibz)) 245 ABI_MALLOC(occfact_ibz,(bantot_ibz)) 246 doccde_ibz(:)=zero ; eigen_ibz(:)=zero ; occfact_ibz(:)=zero 247 248 band_index=0 249 do isppol=1,nsppol 250 do ikibz=1,nkibz 251 ikbz=ibz2bz(ikibz+(sppoldbl-1)*(isppol-1)*nkibz,1) 252 nband_k=nband_ibz(ikibz+(isppol-1)*nkibz) 253 ii=SUM(Dtset%nband(1:ikbz+(isppol-1)*nkpt))-nband_k 254 eigen_ibz(band_index+1:band_index+nband_k)=eigen(ii+1:ii+nband_k) 255 band_index=band_index+nband_k 256 end do 257 end do 258 259 call ebands_init(bantot_ibz,QP_bst,Dtset%nelect,doccde_ibz,eigen_ibz,istwfk_ibz,kibz,nband_ibz,& 260 & nkibz,npwarr_ibz,nsppol,Dtset%nspinor,Dtset%tphysel,Dtset%tsmear,Dtset%occopt,occfact_ibz,wtk_ibz,& 261 & dtset%charge,dtset%kptopt,dtset%kptrlatt_orig,dtset%nshiftk_orig,dtset%shiftk_orig,& 262 & dtset%kptrlatt,dtset%nshiftk,dtset%shiftk) 263 264 ABI_FREE(kibz) 265 ABI_FREE(wtk_ibz) 266 ABI_FREE(ibz2bz) 267 ABI_FREE(npwarr_ibz) 268 ABI_FREE(istwfk_ibz) 269 ABI_FREE(nband_ibz) 270 ABI_FREE(doccde_ibz) 271 ABI_FREE(eigen_ibz) 272 ABI_FREE(occfact_ibz) 273 274 ! === Read in quasiparticle information === 275 ! * Initialize QP amplitudes with KS, QP_bst% presently contains KS energies. 276 ! * If file not found return, everything has been already initialized with KS values 277 ! Here qp_rhor is not needed thus dimrho=0 278 ! TODO paral_kgb not implemented but this is the last problem to be solved 279 280 ABI_MALLOC(m_lda_to_qp,(mband,mband,dtset%nkptgw,nsppol)) 281 m_lda_to_qp=czero 282 do iband=1,mband 283 m_lda_to_qp(iband,iband,:,:)=cone 284 end do 285 286 ! * Fake MPI_type for rdqps 287 call initmpi_seq(MPI_enreg_seq) 288 289 my_ngfft=Dtset%ngfft; if (Dtset%usepaw==1.and.ALL(Dtset%ngfftdg(1:3)/=0)) my_ngfft=Dtset%ngfftdg 290 nfftot=PRODUCT(my_ngfft(1:3)); dimrho=0 291 292 ! Change gw_fname to read a GW file instead of the QPS file. 293 ! TODO not so sure that wannier90 can handle G0W0 eigenvalues that are not ordered, though! 294 gw_fname = "g0w0" 295 g0w0_exists = .FALSE. 296 inquire(file=gw_fname,iostat=ios,exist=g0w0_exists) 297 if (ios/=0) then 298 MSG_ERROR('File g0w0 exists but iostat returns nonzero value!') 299 end if 300 301 if (.not.g0w0_exists) then ! read QPS file (default behavior). 302 input = from_QPS_FILE 303 ABI_DT_MALLOC(prev_Pawrhoij,(Cryst%natom*Dtset%usepaw)) 304 ABI_MALLOC(qp_rhor,(nfftot,nspden*dimrho)) 305 306 call rdqps(QP_bst,Dtfil%fnameabi_qps,Dtset%usepaw,Dtset%nspden,dimrho,nscf,& 307 & nfftot,my_ngfft,ucvol,Dtset%paral_kgb,Cryst,Pawtab,MPI_enreg_seq,nbsc,m_lda_to_qp,qp_rhor,prev_Pawrhoij) 308 309 ABI_FREE(qp_rhor) 310 ABI_DT_FREE(prev_Pawrhoij) 311 312 else ! Read GW file (m_lda_to_qp has been already set to 1, no extrapolation is performed) 313 write(msg,*) ' READING GW CORRECTIONS FROM FILE g0w0 !' 314 MSG_WARNING(msg) 315 input = from_GW_FILE 316 ABI_MALLOC(igwene,(QP_bst%mband,QP_bst%nkpt,QP_bst%nsppol)) 317 call rdgw(QP_bst,gw_fname,igwene,extrapolate=.FALSE.) 318 ABI_FREE(igwene) 319 end if 320 321 ! === Begin big loop over full-zone k points and spin (not implemented) === 322 ! * Wannier90 treats only a single spin, changes in wannier90 are needed 323 ABI_MALLOC(cg_k,(mpw,mband)) 324 ABI_MALLOC(cg_qpk,(mpw,mband)) 325 ABI_MALLOC(m_tmp,(mband,mband)) 326 327 band_index=0 ; icg=0 ; ikg=0 328 do isppol=1,nsppol 329 do ikpt=1,nkpt 330 331 irzkpt =indkk(ikpt+(sppoldbl-1)*(isppol-1)*nkpt,1) 332 itimrev=indkk(ikpt+(sppoldbl-1)*(isppol-1)*nkpt,6) 333 npw_k=npwarr(ikpt) 334 nband_k=dtset%nband(ikpt+(isppol-1)*nkpt) 335 336 if (nband_k/=mband) then 337 write(msg,'(a,i0,7a)')& 338 & 'Number of bands for k point ',ikpt,' is inconsistent with number',ch10,& 339 & 'specified for wannier90 calculation',ch10,& 340 & 'Action : correct input so all band numbers are equal for GW',ch10,& 341 & 'and wannier90 datasets.' 342 MSG_ERROR(msg) 343 end if 344 345 ! === Load KS states for this kbz and spin === 346 do iband=1,nband_k 347 icg_shift=npw_k*my_nspinor*(iband-1)+icg 348 do ipw=1,npw_k 349 cg_k(ipw,iband)=DCMPLX(cg(1,ipw+icg_shift),cg(2,ipw+icg_shift)) 350 end do 351 end do 352 353 ! If time reversal is used for relating ikpt to irzkpt, then multiply by 354 ! the complex conjugate of the lda-to-qp transformation matrix 355 if (itimrev==0) then 356 m_tmp(:,:)=m_lda_to_qp(:,:,irzkpt,isppol) 357 else if (itimrev==1) then 358 m_tmp(:,:)=conjg(m_lda_to_qp(:,:,irzkpt,isppol)) 359 else 360 write(msg,'(2(a,i0))')'Invalid indkk(ikpt,6) ',itimrev,'from routine listkk for k-point ',ikpt 361 MSG_BUG(msg) 362 end if 363 364 call ZGEMM('N','N',npw_k,mband,mband,cone,cg_k,mpw,m_tmp,mband,czero,cg_qpk,mpw) 365 366 ! === Orthonormality test === 367 ! * nband >= maxval(bndgw) for this to pass, but may be less than nband used in GW. 368 ! * Unfortunately, does not test WFK and QPS consistency. 369 !allocate(ortho(nband_k*(nband_k+1)/2)) 370 !ortho=czero; ijpack=0 371 !do jband=1,nband_k 372 ! jb_idx=band_index+jband 373 ! if (dtset%usepaw==1) Cp2 => Cprj_BZ(:,jband:jband+(my_nspinor-1)) 374 ! do iband=1,jband 375 ! ib_idx=band_index+iband 376 ! ijpack=ijpack+1 377 ! ortho(ijpack)=sum(conjg(cg_qpk(1:npw_k,iband))*cg_qpk(1:npw_k,jband)) 378 ! if (dtset%usepaw==1) then 379 ! Cp1 => Cprj_BZ(:,iband:iband+(my_nspinor-1)) 380 ! paw_ovlp = paw_overlap(Cp2,Cp1,Cryst%typat,Pawtab) 381 ! ortho(ijpack) = ortho(ijpack) + CMPLX(paw_ovlp(1),paw_ovlp(2)) 382 ! end if 383 ! if (jband==iband) ortho(ijpack)=ortho(ijpack)-cone 384 ! end do 385 !end do 386 !ortho_err=maxval(abs(ortho)) 387 388 !write(std_out,*)' drh - mlwfovlp_qp: ikpt,ortho_err',ikpt,ortho_err 389 !if (ortho_err>tol6) then 390 ! write(msg, '(3a,i4,a,i6,a,1p,e8.1,3a)' )& 391 !& ' orthonormality error for quasiparticle wave functions.',ch10,& 392 !& ' spin=',isppol,' k point=',ikpt,' ortho_err=',ortho_err,' >1E-6',ch10,& 393 !& ' Action : Be sure input nband>=maxval(bndgw)' 394 ! MSG_ERROR(msg) 395 !end if 396 !deallocate(ortho) 397 398 ! Replace lda wave functions and eigenvalues with quasiparticle ones. 399 qpenek_is_ordered = isordered(nband_k,QP_bst%eig(:,irzkpt,isppol),">",TOL_SORT) 400 401 if (input==from_QPS_FILE .and. .not.qpenek_is_ordered) then 402 write(msg,'(3a)')& 403 & " QP energies read from QPS file are not ordered, likely nband_k>nbdgw. ",ch10,& 404 & " Change nband in the input file so that it equals the number of GW states calculated" 405 MSG_WARNING(msg) 406 end if 407 408 if ( .TRUE. ) then 409 do iband=1,nband_k 410 icg_shift=npw_k*my_nspinor*(iband-1)+icg 411 eigen(iband+band_index)=QP_bst%eig(iband,irzkpt,isppol) 412 do ipw=1,npw_k 413 cg(1,ipw+icg_shift)= real(cg_qpk(ipw,iband)) 414 cg(2,ipw+icg_shift)=aimag(cg_qpk(ipw,iband)) 415 end do 416 end do 417 else 418 ! FIXME There's a problem in twannier90 since nband_k > nbdgw and therefore we also read KS states from the QPS file! 419 ! Automatic test has to be changed! 420 write(msg,'(2a,3f8.4,3a)')ch10,& 421 & "QP energies at k-point ",QP_bst%kptns(:,irzkpt)," are not sorted in ascending numerical order!",ch10,& 422 & "Performing reordering of energies and wavefunctions to be written on the final WKF file." 423 MSG_ERROR(msg) 424 !write(std_out,*)"eig",(QP_bst%eig(ii,irzkpt,isppol),ii=1,nband_k) 425 ABI_MALLOC(sorted_qpene,(nband_k)) 426 ABI_MALLOC(iord,(nband_k)) 427 sorted_qpene = QP_bst%eig(1:nband_k,irzkpt,isppol) 428 iord = (/(ii, ii=1,nband_k)/) 429 430 call sort_dp(nband_k,sorted_qpene,iord,TOL_SORT) 431 do ii=1,nband_k 432 write(std_out,*)"%eig, sorted_qpene, iord",QP_bst%eig(ii,irzkpt,isppol)*Ha_eV,sorted_qpene(ii)*Ha_eV,iord(ii) 433 end do 434 435 do iband=1,nband_k 436 ord_iband = iord(iband) 437 icg_shift=npw_k*my_nspinor*(iband-1)+icg 438 !eigen(iband+band_index)=QP_bst%eig(iband,irzkpt,isppol) 439 eigen(iband+band_index)=QP_bst%eig(ord_iband,irzkpt,isppol) 440 do ipw=1,npw_k 441 !cg(1,ipw+icg_shift)= real(cg_qpk(ipw,iband)) 442 !cg(2,ipw+icg_shift)=aimag(cg_qpk(ipw,iband)) 443 cg(1,ipw+icg_shift)= real(cg_qpk(ipw,ord_iband)) 444 cg(2,ipw+icg_shift)=aimag(cg_qpk(ipw,ord_iband)) 445 end do 446 end do 447 ABI_FREE(sorted_qpene) 448 ABI_FREE(iord) 449 end if 450 451 band_index=band_index+nband_k 452 icg=icg+npw_k*my_nspinor*nband_k 453 ikg=ikg+npw_k 454 end do !ikpt 455 end do !isppol 456 457 ABI_FREE(cg_k) 458 ABI_FREE(cg_qpk) 459 ABI_FREE(m_tmp) 460 461 ! === If PAW, update projections in BZ === 462 ! * Since I am lazy and here I do not care about memory, I just reconstruct m_lda_to_qp in the BZ. 463 ! * update_cprj will take care of updating the PAW projections to get <p_lmn|QP_{nks]> 464 ! This allows some CPU saving, no need to call ctocprj. 465 ! FIXME this part should be tested, automatic test to be provided 466 if (Dtset%usepaw==1) then 467 ABI_MALLOC(dimlmn,(natom)) 468 call pawcprj_getdim(dimlmn,dtset%natom,nattyp_dum,ntypat,Dtset%typat,pawtab,'R') 469 470 nkbz=nkpt 471 ABI_MALLOC(m_lda_to_qp_BZ,(mband,mband,nkbz,nsppol)) 472 do isppol=1,nsppol 473 do ikbz=1,nkbz 474 ikibz =indkk(ikibz+(sppoldbl-1)*(isppol-1)*nkbz,1) 475 itimrev=indkk(ikibz+(sppoldbl-1)*(isppol-1)*nkbz,6) 476 select case (itimrev) 477 case (0) 478 m_lda_to_qp_BZ(:,:,ikbz,isppol)=m_lda_to_qp(:,:,ikibz,isppol) 479 case (1) 480 m_lda_to_qp_BZ(:,:,ikbz,isppol)=CONJG(m_lda_to_qp(:,:,ikibz,isppol)) 481 case default 482 write(msg,'(a,i3)')"Wrong itimrev= ",itimrev 483 MSG_BUG(msg) 484 end select 485 end do 486 end do 487 488 call update_cprj(natom,nkbz,mband,nsppol,my_nspinor,m_lda_to_qp_BZ,dimlmn,Cprj_BZ) 489 ABI_FREE(dimlmn) 490 ABI_FREE(m_lda_to_qp_BZ) 491 end if !PAW 492 493 write(msg,'(6a)')ch10,& 494 & ' mlwfovlp_qp: Input KS wavefuctions have been converted',ch10,& 495 & ' to GW quasiparticle wavefunctions for maximally localized wannier',ch10,& 496 & ' function construction by wannier90.' 497 call wrtout(ab_out,msg,'COLL') 498 call wrtout(std_out,msg,'COLL') 499 500 ABI_FREE(m_lda_to_qp) 501 call ebands_free(QP_bst) 502 call destroy_mpi_enreg(MPI_enreg_seq) 503 504 DBG_EXIT("COLL") 505 506 end subroutine mlwfovlp_qp