TABLE OF CONTENTS
ABINIT/mlwfovlp_proj [ Functions ]
NAME
mlwfovlp_proj
FUNCTION
Routine which computes projection A_{mn}(k) for Wannier code (www.wannier.org f90 version).
COPYRIGHT
Copyright (C) 2005-2018 ABINIT group (BAmadon,FJollet,TRangel,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
cg(2,mpw*nspinor*mband*mkmem*nsppol)=planewave coefficients of wavefunctions cprj(natom,nspinor*mband*mkmem*nsppol)= <p_lmn|Cnk> coefficients for each WF |Cnk> and each |p_lmn> non-local projector dtset <type(dataset_type)>=all input variables for this dataset filew90_win = secondary input file for wannier90 (WAS NOT USED IN v6.7.1 - so has been temporarily removed) kg(3,mpw*mkmem)=reduced planewave coordinates. lproj= flag 0: no projections, 1: random projections, 2: projections on atomic orbitals 3: projections on projectors mband=maximum number of bands mkmem =number of k points treated by this node. npwarr(nkpt)=number of planewaves in basis at this k point mpi_enreg=information about MPI parallelization mpw=maximum dimensioned size of npw. natom=number of atoms in cell. nattyp(ntypat)= # atoms of each type. nkpt=number of k points. nspinor=number of spinorial components of the wavefunctions nsppol=1 for unpolarized, 2 for spin-polarized ntypat=number of types of atoms in unit cell. num_bands=number of bands actually used to construct the wannier function nwan= number of wannier fonctions (read in wannier90.win). proj_l(mband)= angular part of the projection function (quantum number l) proj_m(mband)= angular part of the projection function (quantum number m) proj_radial(mband)= radial part of the projection. proj_site(3,mband)= site of the projection. proj_x(3,mband)= x axis for the projection. proj_z(3,mband)= z axis for the projection. proj_zona(mband)= extension of the radial part. psps <type(pseudopotential_type)>=variables related to pseudopotentials spin = just used for nsppol>1 ; 0 both, 1 just spin up, 2 just spin down
OUTPUT
A_matrix(num_bands,nwan,nkpt,nsppol)= Matrix of projections needed by wannier_run ( also wannier90random.amn is written)
SIDE EFFECTS
(only writing, printing)
NOTES
PARENTS
mlwfovlp
CHILDREN
mlwfovlp_radial,mlwfovlp_ylmfac,wrtout,ylm_cmplx
SOURCE
66 #if defined HAVE_CONFIG_H 67 #include "config.h" 68 #endif 69 70 #include "abi_common.h" 71 72 subroutine mlwfovlp_proj(A_matrix,band_in,cg,cprj,dtset,gprimd,just_augmentation,kg,& 73 &lproj,max_num_bands,mband,mkmem,mpi_enreg,mpw,mwan,natom,nattyp,& 74 &nkpt,npwarr,nspinor,& 75 &nsppol,ntypat,num_bands,nwan,pawtab,proj_l,proj_m,proj_radial,& 76 &proj_site,proj_x,proj_z,proj_zona,psps,spin,ucvol) 77 78 use defs_basis 79 use defs_datatypes 80 use defs_abitypes 81 use defs_wannier90 82 use m_errors 83 use m_profiling_abi 84 use m_xmpi 85 use m_sort 86 87 use m_numeric_tools, only : uniformrandom 88 use m_pawtab, only : pawtab_type 89 use m_pawcprj, only : pawcprj_type 90 use m_paw_sphharm, only : ylm_cmplx 91 92 !This section has been created automatically by the script Abilint (TD). 93 !Do not modify the following lines by hand. 94 #undef ABI_FUNC 95 #define ABI_FUNC 'mlwfovlp_proj' 96 use interfaces_14_hidewrite 97 use interfaces_67_common, except_this_one => mlwfovlp_proj 98 !End of the abilint section 99 100 implicit none 101 102 !Arguments ------------------------------------ 103 !scalars 104 complex(dpc),parameter :: c1=(1._dp,0._dp) 105 integer,intent(in) :: lproj,max_num_bands,mband,mkmem,mpw,mwan,natom,nkpt,nspinor,nsppol 106 integer,intent(in) :: ntypat,spin 107 type(MPI_type),intent(in) :: mpi_enreg 108 type(dataset_type),intent(in) :: dtset 109 type(pseudopotential_type),intent(in) :: psps 110 !arrays 111 integer ::nattyp(ntypat) 112 integer,intent(in) :: kg(3,mpw*mkmem),npwarr(nkpt),num_bands(nsppol),nwan(nsppol),proj_l(mband,nsppol) 113 integer,intent(in) :: proj_m(mband,nsppol) 114 integer,intent(inout)::proj_radial(mband,nsppol) 115 real(dp),intent(in) :: cg(2,mpw*nspinor*mband*mkmem*nsppol) 116 real(dp),intent(in) :: gprimd(3,3),proj_site(3,mband,nsppol) 117 real(dp),intent(in) :: proj_x(3,mband,nsppol),proj_z(3,mband,nsppol),proj_zona(mband,nsppol) 118 complex(dpc),intent(out) :: A_matrix(max_num_bands,mwan,nkpt,nsppol) 119 !character(len=fnlen),intent(in) :: filew90_win(nsppol) 120 logical,intent(in) :: band_in(mband,nsppol) 121 logical,intent(in)::just_augmentation(mwan,nsppol) 122 type(pawcprj_type) :: cprj(natom,nspinor*mband*mkmem*nsppol) 123 type(pawtab_type),intent(in) :: pawtab(psps%ntypat*psps%usepaw) 124 125 !Local variables------------------------------- 126 !scalars 127 integer :: iatom,iatprjn,iband,iband1,iband2,ibg,icat,icg,icg_shift 128 integer :: idum,idx,ikg,ikpt,ilmn,ipw,iproj 129 integer :: ispinor,isppol,itypat,iwan,jband,jj1,libprjn 130 integer :: lmn_size,natprjn,nband_k,nbprjn,npw_k 131 integer :: sumtmp 132 integer :: max_lmax,max_lmax2,mproj,nprocs,spaceComm,rank 133 real(dp),parameter :: qtol=2.0d-8 134 real(dp) :: arg,norm_error,norm_error_bar 135 real(dp) :: ucvol,x1,x2,xnorm,xnormb,xx,yy,zz 136 complex(dpc) :: amn_tmp(nspinor) 137 complex(dpc) :: cstr_fact 138 character(len=500) :: message 139 !arrays 140 integer :: kg_k(3,mpw),lmax(nsppol),lmax2(nsppol),nproj(nsppol) 141 integer,allocatable :: lprjn(:),npprjn(:) 142 real(dp) :: kpg(3),kpt(3) 143 real(dp),allocatable :: amn(:,:,:,:,:),amn2(:,:,:,:,:,:,:) 144 real(dp),allocatable :: gsum2(:),kpg2(:),radial(:) 145 complex(dpc),allocatable :: gf(:,:),gft_lm(:) 146 complex(dpc),allocatable :: ylmc_fac(:,:,:),ylmcp(:) 147 148 !no_abirules 149 !Tables 3.1 & 3.2, User guide 150 integer,save :: orb_l_defs(-5:3)=(/2,2,1,1,1,0,1,2,3/) 151 ! integer,parameter :: mtransfo(0:3,7)=& 152 !& reshape((/1,0,0,0,0,0,0,1,1,1,0,0,0,0,0,-2,-1,2,1,0,0,0,-1,1,2,-2,-3,3/),(/4,7/)) 153 154 !************************************************************************ 155 156 157 !mpi initialization 158 spaceComm=MPI_enreg%comm_cell 159 nprocs=xmpi_comm_size(spaceComm) 160 rank=MPI_enreg%me_kpt 161 162 !Check input variables 163 if ((lproj/=1).and.(lproj/=2).and.(lproj/=5)) then 164 write(message, '(3a)' )& 165 & ' Value of lproj no allowed ',ch10,& 166 & ' Action : change lproj.' 167 MSG_ERROR(message) 168 end if 169 170 write(message, '(a,a)' )ch10,& 171 & '** mlwfovlp_proj: compute A_matrix of initial guess for wannier functions' 172 call wrtout(std_out,message,'COLL') 173 174 !Initialize to 0.d0 175 A_matrix(:,:,:,:)=cmplx(0.d0,0.d0) 176 177 178 ! 179 !End of preliminarities 180 ! 181 182 !********************* Write Random projectors 183 if(lproj==1) then 184 idum=123456 185 ! Compute random projections 186 ABI_ALLOCATE(amn,(2,mband,mwan,nkpt,nsppol)) 187 amn=zero 188 do isppol=1,nsppol 189 if(spin.ne.0 .and. spin.ne.isppol) cycle 190 do ikpt=1,nkpt 191 ! 192 ! MPI: cycle over kpts not treated by this node 193 ! 194 if (ABS(MPI_enreg%proc_distrb(ikpt,1,isppol)-rank)/=0) CYCLE 195 ! write(std_out,'("kpt loop2: ikpt",i3," rank ",i3)') ikpt,rank 196 197 ! 198 do iband1=1,mband 199 xnormb=0.d0 200 do iband2=1,nwan(isppol) 201 x1=uniformrandom(idum) 202 x2=uniformrandom(idum) 203 xnorm=sqrt(x1**2+x2**2) 204 xnormb=xnormb+xnorm 205 amn(1,iband1,iband2,ikpt,isppol)=x1 206 amn(2,iband1,iband2,ikpt,isppol)=x2 207 end do 208 do iband2=1,nwan(isppol) 209 amn(1,iband1,iband2,ikpt,isppol)=amn(1,iband1,iband2,ikpt,isppol)/xnormb 210 amn(2,iband1,iband2,ikpt,isppol)=amn(2,iband1,iband2,ikpt,isppol)/xnormb 211 end do !iband2 212 end do !iband1 213 end do !ikpt 214 end do !isppol 215 do isppol=1,nsppol 216 if(spin.ne.0 .and. spin.ne.isppol) cycle 217 do ikpt=1,nkpt 218 ! 219 ! MPI: cycle over kpts not treated by this node 220 ! 221 if (ABS(MPI_enreg%proc_distrb(ikpt,1,isppol)-rank)/=0) CYCLE 222 ! 223 do iband2=1,nwan(isppol) 224 jband=0 225 do iband1=1,mband 226 if(band_in(iband1,isppol)) then 227 jband=jband+1 228 if(jband.gt.num_bands(isppol)) then 229 write(message, '(3a)' )& 230 & ' Value of jband is above num_bands ',ch10,& 231 & ' Action : contact Abinit group' 232 MSG_ERROR(message) 233 end if 234 A_matrix(jband,iband2,ikpt,isppol)=cmplx(amn(1,iband1,iband2,ikpt,isppol),amn(2,iband1,iband2,ikpt,isppol)) 235 end if 236 end do !iband1 237 end do !iband2 238 end do !ikpt 239 end do !isppol 240 ABI_DEALLOCATE(amn) 241 end if 242 243 !********************* Projection on atomic orbitals based on .win file 244 if( lproj==2) then !based on .win file 245 nproj(:)=nwan(:)/nspinor !if spinors, then the number of projections are 246 mproj=maxval(nproj(:)) 247 ! half the total of wannier functions 248 ! 249 ! obtain lmax and lmax2 250 lmax(:)=0 251 lmax2(:)=0 252 ! 253 do isppol=1,nsppol 254 if(spin.ne.0 .and. spin.ne.isppol) cycle 255 do iproj=1,nproj(isppol) 256 lmax(isppol)=max(lmax(isppol),orb_l_defs(proj_l(iproj,isppol))) 257 end do !iproj 258 lmax2(isppol)=(lmax(isppol)+1)**2 259 end do !isppol 260 max_lmax=maxval(lmax(:)) 261 max_lmax2=maxval(lmax2(:)) 262 ! Allocate arrays 263 ABI_ALLOCATE(ylmc_fac,(max_lmax2,mproj,nsppol)) 264 ! 265 ! get ylmfac, factor used for rotations and hybrid orbitals 266 do isppol=1,nsppol 267 if(spin.ne.0 .and. spin.ne.isppol) cycle 268 call mlwfovlp_ylmfac(ylmc_fac(1:lmax2(isppol),1:nproj(isppol),isppol),lmax(isppol),lmax2(isppol),& 269 & mband,nproj(isppol),proj_l(:,isppol),proj_m(:,isppol),proj_x(:,:,isppol),proj_z(:,:,isppol)) 270 end do 271 ! 272 norm_error=zero 273 norm_error_bar=zero 274 icg=0 275 ! 276 do isppol=1,nsppol 277 ! Allocate arrays 278 if(spin.eq.0 .or. spin.eq.isppol) then 279 ! this has to be done this way because the variable icg changes at the end of the 280 ! cycle. We cannot just skip the hole cycle. 281 ABI_ALLOCATE(gf,(mpw,nproj(isppol))) 282 ABI_ALLOCATE(gft_lm,(lmax2(isppol))) 283 ABI_ALLOCATE(gsum2,(nproj(isppol))) 284 ABI_ALLOCATE(kpg2,(mpw)) 285 ABI_ALLOCATE(radial,(lmax2(isppol))) 286 ABI_ALLOCATE(ylmcp,(lmax2(isppol))) 287 end if 288 ! 289 ikg=0 290 do ikpt=1, nkpt 291 ! 292 ! MPI: cycle over kpts not treated by this node 293 ! 294 if (ABS(MPI_enreg%proc_distrb(ikpt,1,isppol)-rank)/=0) CYCLE 295 ! 296 if(spin.eq.0 .or. spin.eq.isppol) then 297 write(message, '(a,i6,a,2i6)' ) & 298 & ' processor',rank,' will compute k-point,spin=',ikpt,isppol 299 call wrtout(std_out, message,'COLL') 300 end if 301 ! 302 ! Initialize variables 303 npw_k=npwarr(ikpt) 304 gsum2(:)=0.d0 305 gf(:,:) = (0.d0,0.d0) 306 kpt(:)=dtset%kpt(:,ikpt) 307 kg_k(:,1:npw_k)=kg(:,1+ikg:npw_k+ikg) 308 309 do ipw=1, npw_k 310 kpg(1)= (kpt(1) + real(kg_k(1,ipw),dp)) !k+G 311 kpg(2)= (kpt(2) + real(kg_k(2,ipw),dp)) 312 kpg(3)= (kpt(3) + real(kg_k(3,ipw),dp)) 313 ! 314 ! Calculate modulus of k+G 315 xx=gprimd(1,1)*kpg(1)+gprimd(1,2)*kpg(2)+gprimd(1,3)*kpg(3) 316 yy=gprimd(2,1)*kpg(1)+gprimd(2,2)*kpg(2)+gprimd(2,3)*kpg(3) 317 zz=gprimd(3,1)*kpg(1)+gprimd(3,2)*kpg(2)+gprimd(3,3)*kpg(3) 318 kpg2(ipw)= two_pi*sqrt(xx**2+yy**2+zz**2) 319 ! 320 ! Complex Y_lm for k+G 321 if(lmax(isppol)==0) then 322 ylmcp(1)=c1/sqrt(four_pi) 323 else 324 call ylm_cmplx(lmax(isppol),ylmcp,xx,yy,zz) 325 end if 326 ! 327 if(spin.eq.0 .or. spin.eq.isppol) then 328 ! ! 329 do iproj=1,nproj(isppol) 330 ! 331 ! In PAW, we can use proj_radial > 4 to indicate that we just 332 ! want the in-sphere contribution 333 ! 334 if( psps%usepaw==1) then 335 if( just_augmentation(iproj,isppol)) cycle 336 end if 337 ! 338 ! obtain radial part 339 call mlwfovlp_radial(proj_zona(iproj,isppol),lmax(isppol),lmax2(isppol)& 340 & ,radial,proj_radial(iproj,isppol),kpg2(ipw)) 341 ! 342 ! scale complex representation of projector orbital with radial functions 343 ! of appropriate l 344 gft_lm(:)=radial(:)*ylmc_fac(1:lmax2(isppol),iproj,isppol) 345 ! 346 ! complex structure factor for projector orbital position 347 arg = ( kpg(1)*proj_site(1,iproj,isppol) + & 348 & kpg(2)*proj_site(2,iproj,isppol) + & 349 & kpg(3)*proj_site(3,iproj,isppol) ) * 2*pi 350 cstr_fact = cmplx(cos(arg), -sin(arg) ) 351 ! 352 ! obtain guiding functions 353 gf(ipw,iproj)=cstr_fact*dot_product(ylmcp,gft_lm) 354 ! 355 gsum2(iproj)=gsum2(iproj)+real(gf(ipw,iproj))**2+aimag(gf(ipw,iproj))**2 356 end do !iproj 357 end if !spin 358 end do !ipw 359 ! 360 if(spin.eq.0 .or. spin.eq.isppol) then 361 do iproj=1,nproj(isppol) 362 ! 363 ! In PAW, we can use proj_radial > 4 to indicate that we just 364 ! want the in-sphere contribution 365 ! 366 if(psps%usepaw==1 ) then 367 if (just_augmentation(iproj,isppol)) cycle 368 end if 369 ! 370 gsum2(iproj)=16._dp*pi**2*gsum2(iproj)/ucvol 371 gf(:,iproj)=gf(:,iproj)/sqrt(gsum2(iproj)) 372 norm_error=max(abs(gsum2(iproj)-one),norm_error) 373 norm_error_bar=norm_error_bar+(gsum2(iproj)-one)**2 374 end do !iproj 375 ! 376 ! Guiding functions are computed. 377 ! compute overlaps of gaussian projectors and wave functions 378 do iproj=1,nproj(isppol) 379 ! 380 ! In PAW, we can use proj_radial > 4 to indicate that we just 381 ! want the in-sphere contribution 382 ! 383 if(psps%usepaw==1 ) then 384 if ( just_augmentation(iproj,isppol)) cycle 385 end if 386 ! 387 jband=0 388 do iband=1,mband 389 if(band_in(iband,isppol)) then 390 icg_shift=npw_k*nspinor*(iband-1)+icg 391 jband=jband+1 392 amn_tmp(:)=cmplx(0.d0,0.d0) 393 do ispinor=1,nspinor 394 do ipw=1,npw_k 395 ! 396 ! The case of spinors is tricky, we have nproj = nwan/2 397 ! so we project to spin up and spin down separately, to have at 398 ! the end an amn matrix with nwan projections. 399 ! 400 ! idx=ipw*nspinor - (nspinor-ispinor) 401 idx=ipw+(ispinor-1)*npw_k 402 amn_tmp(ispinor)=amn_tmp(ispinor)+gf(ipw,iproj)*cmplx(cg(1,idx+icg_shift),-cg(2,idx+icg_shift)) 403 end do !ipw 404 end do !ispinor 405 do ispinor=1,nspinor 406 iwan=(iproj*nspinor)- (nspinor-ispinor) 407 A_matrix(jband,iwan,ikpt,isppol)=amn_tmp(ispinor) 408 end do 409 end if !band_in 410 end do !iband 411 end do !iproj 412 end if !spin==isppol 413 icg=icg+npw_k*nspinor*mband 414 ikg=ikg+npw_k 415 end do !ikpt 416 ! Deallocations 417 if(spin.eq.0 .or. spin.eq.isppol) then 418 ABI_DEALLOCATE(gf) 419 ABI_DEALLOCATE(gft_lm) 420 ABI_DEALLOCATE(gsum2) 421 ABI_DEALLOCATE(kpg2) 422 ABI_DEALLOCATE(radial) 423 ABI_DEALLOCATE(ylmcp) 424 end if 425 end do !isppol 426 ! 427 ! if(isppol==1) then 428 ! norm_error_bar=sqrt(norm_error_bar/real(nkpt*(nwan(1)),dp)) 429 ! else 430 ! if(spin==0) norm_error_bar=sqrt(norm_error_bar/real(nkpt*(nwan(1)+nwan(2)),dp)) 431 ! if(spin==1) norm_error_bar=sqrt(norm_error_bar/real(nkpt*nwan(1),dp)) 432 ! if(spin==2) norm_error_bar=sqrt(norm_error_bar/real(nkpt*nwan(2),dp)) 433 ! end if 434 ! if(norm_error>0.05_dp) then 435 ! write(message, '(6a,f6.3,a,f6.3,12a)' )ch10,& 436 ! & ' mlwfovlp_proj : WARNING',ch10,& 437 ! & ' normalization error for wannier projectors',ch10,& 438 ! & ' is',norm_error_bar,' (average) and',norm_error,' (max).',ch10,& 439 ! & ' this may indicate more cell-to-cell overlap of the radial functions',ch10,& 440 ! & ' than you want.',ch10,& 441 ! & ' Action : modify zona (inverse range of radial functions)',ch10,& 442 ! ' under "begin projectors" in ',trim(filew90_win),' file',ch10 443 ! call wrtout(std_out,message,'COLL') 444 ! end if 445 ! 446 ! !Deallocate 447 ! deallocate(ylmc_fac) 448 ! 449 ABI_DEALLOCATE(ylmc_fac) 450 end if !lproj==2 451 452 453 !*************** computes projection from PROJECTORS ******************** 454 if(lproj==3) then !! if LPROJPRJ 455 ! ----- set values for projections --------------------- ! INPUT 456 ! nbprjn:number of different l-values for projectors 457 ! lprjn: value of l for each projectors par ordre croissant 458 ! npprjn: number of projectors for each lprjn 459 natprjn=1 ! atoms with wannier functions are first 460 if(natprjn/=1) then ! in this case lprjn should depend on iatprjn 461 MSG_ERROR("natprjn/=1") 462 end if 463 nbprjn=2 464 ABI_ALLOCATE(lprjn,(nbprjn)) 465 lprjn(1)=0 466 lprjn(2)=1 467 ABI_ALLOCATE(npprjn,(0:lprjn(nbprjn))) 468 npprjn(0)=1 469 npprjn(1)=1 470 ! --- test coherence of nbprjn and nwan 471 sumtmp=0 472 do iatprjn=1,natprjn 473 do libprjn=0,lprjn(nbprjn) 474 sumtmp=sumtmp+(2*libprjn+1)*npprjn(libprjn) 475 end do 476 end do 477 if(sumtmp/=nwan(1)) then 478 write(std_out,*) "Number of Wannier orbitals is not equal to number of projections" 479 write(std_out,*) "Action: check values of lprjn,npprjn % nwan" 480 write(std_out,*) "nwan, sumtmp=",nwan,sumtmp 481 MSG_ERROR("Aborting now") 482 end if 483 ! --- end test of coherence 484 ABI_ALLOCATE(amn2,(2,natom,nsppol,nkpt,mband,nspinor,nwan(1))) 485 if(psps%usepaw==1) then 486 amn2=zero 487 ibg=0 488 do isppol=1,nsppol 489 do ikpt=1,nkpt 490 nband_k=dtset%nband(ikpt+(isppol-1)*nkpt) 491 do iband=1,nband_k 492 ! write(std_out,*)"amn2",iband,ibg,ikpt 493 do ispinor=1,nspinor 494 icat=1 495 do itypat=1,dtset%ntypat 496 lmn_size=pawtab(itypat)%lmn_size 497 do iatom=icat,icat+nattyp(itypat)-1 498 jj1=0 499 do ilmn=1,lmn_size 500 if(iatom.le.natprjn) then 501 ! do iwan=1,nwan 502 do libprjn=0,lprjn(nbprjn) 503 ! if (psps%indlmn(1,ilmn,itypat)==proj_l(iwan)) then 504 ! if (psps%indlmn(2,ilmn,itypat)==mtransfo(proj_l(iwan),proj_m(iwan))) then 505 if (psps%indlmn(1,ilmn,itypat)==libprjn) then 506 if (psps%indlmn(3,ilmn,itypat)<=npprjn(libprjn)) then 507 if(band_in(iband,isppol)) then 508 jj1=jj1+1 509 if(jj1>nwan(isppol)) then 510 write(std_out,*) "number of wannier orbitals is lower than lmn_size" 511 write(std_out,*) jj1,nwan(isppol) 512 MSG_ERROR("Aborting now") 513 end if 514 amn2(1,iatom,isppol,ikpt,iband,ispinor,jj1)=cprj(iatom,iband+ibg)%cp(1,ilmn) 515 amn2(2,iatom,isppol,ikpt,iband,ispinor,jj1)=cprj(iatom,iband+ibg)%cp(2,ilmn) 516 end if 517 end if 518 end if 519 end do ! libprjn 520 ! endif 521 ! endif 522 ! enddo ! iwan 523 end if ! natprjn 524 end do !ilmn 525 end do ! iatom 526 icat=icat+nattyp(itypat) 527 end do ! itypat 528 end do ! ispinor 529 end do !iband 530 ibg=ibg+nband_k*nspinor 531 ! write(std_out,*)'amn2b',iband,ibg,ikpt 532 end do !ikpt 533 end do ! isppol 534 535 ! ----------------------- Save Amn -------------------- 536 do isppol=1,nsppol 537 do ikpt=1,nkpt 538 do iband2=1,nwan(isppol) 539 jband=0 540 do iband1=1,mband 541 if(band_in(iband1,isppol)) then 542 jband=jband+1 543 A_matrix(jband,iband2,ikpt,isppol)=& 544 & cmplx(amn2(1,1,1,ikpt,iband1,1,iband2),amn2(2,1,1,ikpt,iband1,1,iband2)) 545 end if 546 end do 547 end do 548 end do 549 end do 550 end if !usepaw 551 ABI_DEALLOCATE(amn2) 552 ABI_DEALLOCATE(npprjn) 553 ABI_DEALLOCATE(lprjn) 554 555 end if ! lproj==3 556 557 558 end subroutine mlwfovlp_proj