TABLE OF CONTENTS
ABINIT/pawprt [ Functions ]
NAME
pawprt
FUNCTION
Print out data concerning PAW formalism (pseudopotential strength, augmentation occupancies...) To be called at the end of the SCF cycle
COPYRIGHT
Copyright (C) 1998-2018 ABINIT group (FJ,MT,BA) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt . For the initials of contributors, see ~abinit/doc/developers/contributors.txt.
INPUTS
dtset <type(dataset_type)>=all input variables in this dataset | enunit=parameter determining units of output energies | kptopt=option for the generation of k points | natom=number of atoms in cell | ntypat = number of atom types | pawprtvol= printing volume | pawspnorb=flag: 1 if spin-orbit coupling is activated | typat(natom)=type of each atom electronpositron <type(electronpositron_type)>=quantities for the electron-positron annihilation (optional argument) mpi_atmtab(:)=--optional-- indexes of the atoms treated by current proc comm_atom=--optional-- MPI communicator over atoms my_natom=number of atoms treated by current processor paw_ij(my_natom) <type(paw_ij_type)>=paw arrays given on (i,j) channels pawrhoij(my_natom) <type(pawrhoij_type)>= paw rhoij occupancies and related data pawtab(ntypat*usepaw) <type(pawtab_type)>=paw tabulated starting data
OUTPUT
(only printing)
PARENTS
bethe_salpeter,outscfcv,screening,sigma
CHILDREN
free_my_atmtab,get_my_atmtab,mat_mlms2jmj,mat_slm2ylm,paw_ij_free paw_ij_gather,paw_ij_nullify,pawio_print_ij,pawrhoij_free pawrhoij_gather,pawrhoij_nullify,setnoccmmp,wrtout,xmpi_comm_group xmpi_group_free,xmpi_group_translate_ranks
SOURCE
49 #if defined HAVE_CONFIG_H 50 #include "config.h" 51 #endif 52 53 #include "abi_common.h" 54 55 subroutine pawprt(dtset,my_natom,paw_ij,pawrhoij,pawtab,& 56 & electronpositron,& ! optional argument 57 & mpi_atmtab,comm_atom) ! optional arguments (parallelism) 58 59 60 use defs_basis 61 use defs_abitypes 62 use defs_parameters 63 use m_profiling_abi 64 use m_errors 65 use m_xmpi 66 67 use m_paral_atom, only : get_my_atmtab, free_my_atmtab 68 use m_electronpositron, only : electronpositron_type,electronpositron_calctype,EP_POSITRON 69 use m_pawang, only : pawang_type, mat_slm2ylm, mat_mlms2jmj 70 use m_pawtab, only : pawtab_type 71 use m_paw_ij, only : paw_ij_type, paw_ij_free, paw_ij_nullify, paw_ij_gather 72 use m_pawrhoij, only : pawrhoij_type, pawrhoij_free, pawrhoij_gather, pawrhoij_nullify 73 use m_paw_io, only : pawio_print_ij 74 use m_paw_sphharm, only : mat_mlms2jmj, mat_slm2ylm 75 76 !This section has been created automatically by the script Abilint (TD). 77 !Do not modify the following lines by hand. 78 #undef ABI_FUNC 79 #define ABI_FUNC 'pawprt' 80 use interfaces_14_hidewrite 81 use interfaces_65_paw, except_this_one => pawprt 82 !End of the abilint section 83 84 implicit none 85 86 !Arguments ------------------------------------ 87 !scalars 88 integer,intent(in) :: my_natom 89 integer,optional,intent(in) :: comm_atom 90 type(dataset_type),intent(in) :: dtset 91 type(electronpositron_type),pointer,optional :: electronpositron 92 !arrays 93 integer,optional,target,intent(in) :: mpi_atmtab(:) 94 type(paw_ij_type),target,intent(inout) :: paw_ij(my_natom) 95 type(pawrhoij_type),target,intent(inout) :: pawrhoij(my_natom) 96 type(pawtab_type),target,intent(in) :: pawtab(dtset%ntypat) 97 98 !Local variables------------------------------- 99 !scalars 100 integer,parameter :: natmax=2 101 integer :: cplex_dij,group1,group2,iat,iatom,ierr,ii,im1,im2,ipositron,irhoij,ispden 102 integer :: i_unitfi,itypat,jrhoij,ll,llp,me_atom,my_comm_atom,natprt,ndij,nspden,nsppol 103 integer :: optsym,sz1,sz2,unitfi,unt 104 real(dp) :: mnorm,mx,my,mz,ntot,valmx,localm 105 logical,parameter :: debug=.false. 106 logical :: my_atmtab_allocated,paral_atom,useexexch,usepawu 107 type(pawang_type):: pawang_dum 108 character(len=7),parameter :: dspin1(6)=(/"up ","down ","up-up ","dwn-dwn","up-dwn ","dwn-up "/) 109 character(len=8),parameter :: dspin2(6)=(/"up ","down ","dens (n)","magn (x)","magn (y)","magn (z)"/) 110 character(len=9),parameter :: dspin3(6)=(/"up ","down ","up-up ","down-down","Re[up-dn]","Im[up-dn]"/) 111 character(len=500) :: msg0,msg 112 !arrays 113 integer :: idum(1) 114 integer :: idum1(0),idum3(0,0,0) 115 integer,allocatable :: jatom(:),opt_l_index(:) 116 integer,pointer :: my_atmtab(:) 117 real(dp) :: rdum2(0,0),rdum4(0,0,0,0) 118 real(dp),allocatable :: rhoijs(:,:) 119 complex(dpc),allocatable :: noccmmp_ylm(:,:,:),noccmmp_jmj(:,:),noccmmp_slm(:,:,:) 120 type(paw_ij_type), ABI_CONTIGUOUS pointer :: paw_ij_all(:) 121 type(pawrhoij_type),ABI_CONTIGUOUS pointer :: pawrhoij_all(:) 122 123 ! ********************************************************************* 124 125 DBG_ENTER("COLL") 126 127 !Set up parallelism over atoms 128 paral_atom=(present(comm_atom).and.(my_natom/=dtset%natom)) 129 nullify(my_atmtab);if (present(mpi_atmtab)) my_atmtab => mpi_atmtab 130 my_comm_atom=xmpi_comm_self;if (present(comm_atom)) my_comm_atom=comm_atom 131 call get_my_atmtab(my_comm_atom,my_atmtab,my_atmtab_allocated,paral_atom,dtset%natom,my_natom_ref=my_natom) 132 me_atom=xmpi_comm_rank(my_comm_atom) 133 134 !Continue only if comm_atom contains the master of the output comm 135 if (paral_atom) then 136 call xmpi_comm_group(abinit_comm_output,group1,ierr) 137 call xmpi_comm_group(my_comm_atom,group2,ierr) 138 call xmpi_group_translate_ranks(group1,1,(/0/),group2,idum,ierr) 139 call xmpi_group_free(group1) 140 call xmpi_group_free(group2) 141 if (idum(1)==xmpi_undefined) then 142 call free_my_atmtab(my_atmtab,my_atmtab_allocated) 143 return 144 end if 145 end if 146 147 !Initializations 148 natprt=natmax;if (dtset%natom==1) natprt=1 149 if (dtset%pawprtvol<0) natprt=dtset%natom 150 ABI_ALLOCATE(jatom,(natprt)) 151 if (natprt==1) then 152 jatom(1)=1 153 else if (natprt==2) then 154 jatom(1)=1;jatom(2)=dtset%natom 155 else if (natprt==dtset%natom) then 156 do iat=1,dtset%natom 157 jatom(iat)=iat 158 end do 159 else 160 MSG_BUG("invalid value of natprt!") 161 end if 162 usepawu=(count(pawtab(:)%usepawu>0)>0) 163 useexexch=(count(pawtab(:)%useexexch>0)>0) 164 ipositron=0 165 if (present(electronpositron)) then 166 if (associated(electronpositron)) ipositron=electronpositron%calctype 167 end if 168 169 !Main title 170 write(msg, '(2a)' ) ch10,& 171 & ' ==== Results concerning PAW augmentation regions ====' 172 call wrtout(ab_out,msg,'COLL') 173 call wrtout(std_out,msg,'COLL') 174 msg=' ' 175 call wrtout(ab_out,msg,'COLL') 176 call wrtout(std_out,msg,'COLL') 177 178 !If atomic data are distributed, retrieve all Dij on master proc 179 if (paral_atom) then 180 if (me_atom==0) then 181 ABI_DATATYPE_ALLOCATE(paw_ij_all,(dtset%natom)) 182 call paw_ij_nullify(paw_ij_all) 183 else 184 ABI_DATATYPE_ALLOCATE(paw_ij_all,(0)) 185 end if 186 call paw_ij_gather(paw_ij,paw_ij_all,0,my_comm_atom) 187 else 188 paw_ij_all => paw_ij 189 end if 190 191 !Print out pseudopotential strength 192 !---------------------------------- 193 if (me_atom==0) then 194 do i_unitfi=1,2 195 unitfi=ab_out;if (i_unitfi==2) unitfi=std_out 196 do unt=1,2 197 if ((unt==1).and.(dtset%enunit==0.or.dtset%enunit==2)) then 198 write(msg,'(a)') ' Total pseudopotential strength Dij (hartree):' 199 call wrtout(unitfi,msg,'COLL') 200 else if ((unt==2).and.(dtset%enunit==1.or.dtset%enunit==2)) then 201 write(msg,'(a)') ' Total pseudopotential strength Dij (eV):' 202 call wrtout(unitfi,msg,'COLL') 203 end if 204 if (ipositron>0) then 205 if (((unt==1).and.(dtset%enunit==0.or.dtset%enunit==2)).or.& 206 & ((unt==2).and.(dtset%enunit==1.or.dtset%enunit==2))) then 207 if (electronpositron%has_pos_ham==0) then 208 write(msg,'(a)') ' -Note: these are the electronic Dij' 209 else 210 write(msg,'(a)') ' -Note: these are the positronic Dij' 211 end if 212 call wrtout(unitfi,msg,'COLL') 213 end if 214 end if 215 if (((unt==1).and.(dtset%enunit==0.or.dtset%enunit==2)).or.& 216 & ((unt==2).and.(dtset%enunit==1.or.dtset%enunit==2))) then 217 do iat=1,natprt 218 iatom=jatom(iat);nspden=paw_ij_all(iatom)%nspden;ndij=paw_ij_all(iatom)%ndij 219 cplex_dij=paw_ij_all(iatom)%cplex_dij 220 optsym=2;if (cplex_dij==2.and.dtset%nspinor==1) optsym=1 221 do ispden=1,ndij 222 valmx=100._dp;if (ispden==1) valmx=-1._dp 223 msg='' ; msg0='' 224 if (dtset%natom>1.or.nspden>1.or.ndij==4) write(msg0, '(a,i3)' ) ' Atom #',iatom 225 if (nspden==1.and.ndij/=4) write(msg,'(a)') trim(msg0) 226 if (nspden==2) write(msg,'(2a,i1)') trim(msg0),' - Spin component ',ispden 227 if (ndij==4) write(msg,'(3a)') trim(msg0),' - Component ',trim(dspin1(ispden+2*(ndij/4))) 228 if (dtset%natom>1.or.nspden>1.or.ndij==4) then 229 call wrtout(unitfi,msg,'COLL') 230 end if 231 if (ndij/=4.or.ispden<=2) then 232 call pawio_print_ij(unitfi,paw_ij_all(iatom)%dij(:,ispden),paw_ij_all(iatom)%lmn2_size,& 233 & cplex_dij,paw_ij_all(iatom)%lmn_size,-1,idum,0,dtset%pawprtvol,& 234 & idum,valmx,unt,opt_sym=optsym) 235 else 236 call pawio_print_ij(unitfi,paw_ij_all(iatom)%dij(:,ispden),paw_ij_all(iatom)%lmn2_size,& 237 & cplex_dij,paw_ij_all(iatom)%lmn_size,-1,idum,0,dtset%pawprtvol,& 238 & idum,valmx,unt,asym_ij=paw_ij_all(iatom)%dij(:,7-ispden),opt_sym=optsym) 239 end if 240 end do 241 end do 242 end if 243 msg=' ' 244 call wrtout(unitfi,msg,'COLL') 245 end do 246 end do 247 end if 248 if (paral_atom.and.(.not.usepawu).and.(.not.useexexch)) then 249 call paw_ij_free(paw_ij_all) 250 ABI_DATATYPE_DEALLOCATE(paw_ij_all) 251 end if 252 253 !If atomic data are distributed, retrieve all Rhoij on master proc 254 if (paral_atom) then 255 if (me_atom==0) then 256 ABI_DATATYPE_ALLOCATE(pawrhoij_all,(dtset%natom)) 257 else 258 ABI_DATATYPE_ALLOCATE(pawrhoij_all,(0)) 259 end if 260 call pawrhoij_nullify(pawrhoij_all) 261 call pawrhoij_gather(pawrhoij,pawrhoij_all,0,my_comm_atom,& 262 & with_grhoij=.false.,with_lmnmix=.false.,& 263 & with_rhoij_=.false.,with_rhoijres=.false.) 264 else 265 pawrhoij_all => pawrhoij 266 end if 267 268 !Print out SYMMETRIZED occupancies of the partial waves 269 !------------------------------------------------------ 270 if (me_atom==0) then 271 do i_unitfi=1,2 272 unitfi=ab_out;if (i_unitfi==2) unitfi=std_out 273 write(msg,'(a)') ' Augmentation waves occupancies Rhoij:' 274 call wrtout(unitfi,msg,'COLL') 275 if (ipositron>0) then 276 if (electronpositron%particle==EP_POSITRON) then 277 write(msg,'(a)') ' -Note: these are the electronic Rhoij' 278 else 279 write(msg,'(a)') ' -Note: these are the positronic Rhoij' 280 end if 281 call wrtout(unitfi,msg,'COLL') 282 end if 283 if (dtset%pawspnorb>0.and.pawrhoij_all(1)%cplex==1.and.dtset%kptopt/=1.and.dtset%kptopt/=2) then 284 write(msg,'(6a)') ' pawprt: - WARNING:',ch10,& 285 & ' Spin-orbit coupling is activated but only real part of Rhoij occupancies',ch10,& 286 & ' has been computed; they could have an imaginary part (not printed here).' 287 call wrtout(unitfi,msg,'COLL') 288 end if 289 do iat=1,natprt 290 iatom=jatom(iat);nspden=pawrhoij_all(iatom)%nspden 291 optsym=2;if (pawrhoij_all(iatom)%cplex==2.and.dtset%nspinor==1) optsym=1 292 do ispden=1,nspden 293 valmx=25._dp;if (ispden==1) valmx=-1._dp 294 msg='' ; msg0='' 295 if (dtset%natom>1.or.nspden>1) write(msg0, '(a,i3)' ) ' Atom #',iatom 296 if (nspden==1) write(msg,'(a)') trim(msg0) 297 if (nspden==2) write(msg,'(2a,i1)') trim(msg0),' - Spin component ',ispden 298 if (nspden==4) write(msg,'(3a)') trim(msg0),' - Component ',dspin2(ispden+2*(nspden/4)) 299 if (dtset%natom>1.or.nspden>1) then 300 call wrtout(unitfi,msg,'COLL') 301 end if 302 call pawio_print_ij(unitfi,pawrhoij_all(iatom)%rhoijp(:,ispden),pawrhoij_all(iatom)%nrhoijsel,& 303 & pawrhoij_all(iatom)%cplex,pawrhoij_all(iatom)%lmn_size,-1,idum,1,dtset%pawprtvol,& 304 & pawrhoij_all(iatom)%rhoijselect(:),valmx,1,opt_sym=optsym) 305 end do 306 end do 307 msg=' ' 308 call wrtout(unitfi,msg,'COLL') 309 end do 310 end if 311 312 !PAW+U or local exact-exchange: print out +U components of occupancies 313 !------------------------------------------------------------------------------- 314 if ((usepawu.or.useexexch).and.ipositron/=1.and.me_atom==0) then 315 do i_unitfi=1,2 316 unitfi=ab_out;if (i_unitfi==2) unitfi=std_out 317 if(useexexch) write(msg,'(a)') & 318 & ' "Local exact-exchange" part of augmentation waves occupancies Rhoij:' 319 if(usepawu) write(msg,'(a)') & 320 & ' "PAW+U" part of augmentation waves occupancies Rhoij:' 321 call wrtout(unitfi,msg,'COLL') 322 valmx=-1._dp 323 do iatom=1,dtset%natom 324 nspden=pawrhoij_all(iatom)%nspden;itypat=dtset%typat(iatom) 325 cplex_dij=paw_ij_all(iatom)%cplex_dij 326 ll=-1;llp=-1 327 if (pawtab(itypat)%usepawu>0) ll=pawtab(itypat)%lpawu 328 if (pawtab(itypat)%useexexch>0) llp=pawtab(itypat)%lexexch 329 if (ll/=llp.and.ll/=-1.and.llp/=-1) then 330 MSG_BUG(" lpawu/=lexexch forbidden!") 331 end if 332 ABI_ALLOCATE(opt_l_index, (pawtab(itypat)%lmn_size)) 333 ll=max(ll,llp) 334 if (ll>=0) then 335 optsym=2;if (pawrhoij_all(iatom)%cplex==2.and.dtset%nspinor==1) optsym=1 336 do ispden=1,nspden 337 msg='' ; msg0='' 338 write(msg0,'(a,i3,a,i1,a)') ' Atom #',iatom,' - L=',ll,' ONLY' 339 if (nspden==1) write(msg,'(a)') trim(msg0) 340 if (nspden==2) write(msg,'(2a,i1)') trim(msg0),' - Spin component ',ispden 341 if (nspden==4) write(msg,'(3a)') trim(msg0),' - Component ',dspin2(ispden+2*(nspden/4)) 342 call wrtout(unitfi,msg,'COLL') 343 344 opt_l_index = pawtab(itypat)%indlmn(1,1:pawtab(itypat)%lmn_size) 345 call pawio_print_ij(unitfi,pawrhoij_all(iatom)%rhoijp(:,ispden),pawrhoij_all(iatom)%nrhoijsel,& 346 & pawrhoij_all(iatom)%cplex,pawrhoij_all(iatom)%lmn_size,ll,& 347 & opt_l_index,1,dtset%pawprtvol,& 348 & pawrhoij_all(iatom)%rhoijselect(:),valmx,1,opt_sym=optsym) 349 350 end do 351 if (debug.and.nspden==4) then 352 sz1=paw_ij_all(iatom)%lmn2_size*cplex_dij 353 sz2=paw_ij_all(iatom)%ndij 354 ABI_ALLOCATE(rhoijs,(sz1,sz2)) 355 do irhoij=1,pawrhoij_all(iatom)%nrhoijsel 356 jrhoij=cplex_dij*(irhoij-1)+1 357 rhoijs(jrhoij,1)=pawrhoij_all(iatom)%rhoijp(jrhoij,1)+pawrhoij_all(iatom)%rhoijp(jrhoij,4) 358 rhoijs(jrhoij+1,1)=pawrhoij_all(iatom)%rhoijp(jrhoij+1,1)+pawrhoij_all(iatom)%rhoijp(jrhoij+1,4) 359 rhoijs(jrhoij,2)=pawrhoij_all(iatom)%rhoijp(jrhoij,1)-pawrhoij_all(iatom)%rhoijp(jrhoij,4) 360 rhoijs(jrhoij+1,2)=pawrhoij_all(iatom)%rhoijp(jrhoij+1,1)-pawrhoij_all(iatom)%rhoijp(jrhoij+1,4) 361 rhoijs(jrhoij,3)=pawrhoij_all(iatom)%rhoijp(jrhoij,2)+pawrhoij_all(iatom)%rhoijp(jrhoij+1,3) 362 rhoijs(jrhoij+1,3)=pawrhoij_all(iatom)%rhoijp(jrhoij+1,2)-pawrhoij_all(iatom)%rhoijp(jrhoij,3) 363 rhoijs(jrhoij,4)=pawrhoij_all(iatom)%rhoijp(jrhoij ,2)-pawrhoij_all(iatom)%rhoijp(jrhoij+1,3) 364 rhoijs(jrhoij+1,4)=pawrhoij_all(iatom)%rhoijp(jrhoij+1,2)+pawrhoij_all(iatom)%rhoijp(jrhoij,3) 365 end do 366 do ispden=1,nspden 367 write(msg,'(3a)') trim(msg0),' - Component ',dspin1(ispden+2*(nspden/4)) 368 opt_l_index = pawtab(itypat)%indlmn(1,1:pawtab(itypat)%lmn_size) 369 370 call pawio_print_ij(unitfi,rhoijs(:,ispden),pawrhoij_all(iatom)%nrhoijsel,& 371 & pawrhoij_all(iatom)%cplex,pawrhoij_all(iatom)%lmn_size,ll,& 372 & opt_l_index,1,dtset%pawprtvol,& 373 & pawrhoij_all(iatom)%rhoijselect(:),valmx,1,opt_sym=optsym) 374 375 call wrtout(unitfi,"WARNING: half of the array is not correct",'COLL') 376 end do 377 ABI_DEALLOCATE(rhoijs) 378 end if 379 end if 380 ABI_DEALLOCATE(opt_l_index) 381 end do ! iatom 382 call wrtout(unitfi,' ','COLL') 383 end do 384 end if 385 386 !PAW+U: print out occupations for correlated orbitals 387 !---------------------------------------------------- 388 if (usepawu.and.ipositron/=1.and.me_atom==0) then 389 do i_unitfi=1,2 390 unitfi=ab_out;if (i_unitfi==2) unitfi=std_out 391 write(msg,'(3a)') & 392 & ' ---------- LDA+U DATA --------------------------------------------------- ',ch10 393 call wrtout(unitfi,msg,'COLL') 394 do iatom=1,dtset%natom 395 itypat=dtset%typat(iatom);ll=pawtab(itypat)%lpawu 396 nspden=paw_ij_all(iatom)%nspden;ndij=paw_ij_all(iatom)%ndij 397 cplex_dij=paw_ij_all(iatom)%cplex_dij 398 if ((ll>=0).and.(pawtab(itypat)%usepawu>0)) then 399 write(msg,fmt='(a,i5,a,i4,a)') " ====== For Atom ", iatom,& 400 & ", occupations for correlated orbitals. lpawu =",ll,ch10 401 call wrtout(unitfi,msg,'COLL') 402 if(pawtab(itypat)%usepawu>=10) then 403 write(msg,fmt='(a)') " (This is PAW atomic orbital occupations)" 404 call wrtout(unitfi,msg,'COLL') 405 write(msg,fmt='(a)') " (For Wannier orbital occupations, refer to DFT+DMFT occupations above)" 406 call wrtout(unitfi,msg,'COLL') 407 end if 408 if(nspden==2) then 409 do ispden=1,nspden 410 write(msg,fmt='(a,i4,a,i3,a,f10.5)') " Atom", iatom,& 411 & ". Occ. for lpawu and for spin",ispden," =",paw_ij_all(iatom)%nocctot(ispden) 412 call wrtout(unitfi,msg,'COLL') 413 end do 414 localm=paw_ij_all(iatom)%nocctot(2)-paw_ij_all(iatom)%nocctot(1) 415 write(msg,fmt='(a,i4,a,2x,f12.6)') " => On atom",iatom,& 416 & ", local Mag. for lpawu is ",localm 417 call wrtout(unitfi,msg,'COLL') 418 end if 419 if(ndij==4) then 420 ntot=paw_ij_all(iatom)%nocctot(1) 421 mx=paw_ij_all(iatom)%nocctot(2) 422 my=paw_ij_all(iatom)%nocctot(3) 423 mz=paw_ij_all(iatom)%nocctot(4) 424 mnorm=sqrt(mx*mx+my*my+mz*mz) 425 write(msg,'(a,i4,a,2x,e15.8)') " => On atom",iatom,", for lpawu, local Mag. x is ",mx 426 call wrtout(unitfi,msg,'COLL') 427 write(msg,'(14x,a,2x,e15.8)') " local Mag. y is ",my 428 call wrtout(unitfi,msg,'COLL') 429 write(msg,'(14x,a,2x,e15.8)') " local Mag. z is ",mz 430 call wrtout(unitfi,msg,'COLL') 431 write(msg,'(14x,a,2x,e15.8)') " norm of Mag. is ",mnorm 432 call wrtout(unitfi,msg,'COLL') 433 write(msg,fmt='(8x,a,2x,f10.5)') " (along mag axis) occ. for majority spin is = ",& 434 & half*(ntot+mnorm) 435 call wrtout(unitfi,msg,'COLL') 436 write(msg,fmt='(8x,a,2x,f10.5)') " (along mag axis) occ. for minority spin is = ",& 437 & half*(ntot-mnorm) 438 call wrtout(unitfi,msg,'COLL') 439 end if 440 write(msg,'(3a)') ch10," == Occupation matrix for correlated orbitals:",ch10 441 call wrtout(unitfi,msg,'COLL') 442 do ispden=1,ndij 443 if (nspden==1.and.ndij/=4.and.(cplex_dij==1)) write(msg,fmt='(a)') " Up component only..." 444 if (nspden==2) write(msg,fmt='(a,i3)')" Occupation matrix for spin",ispden 445 if (ndij==4.or.(cplex_dij==2)) & 446 & write(msg,fmt='(2a)') " Occupation matrix for component ",trim(dspin1(ispden+2*(ndij/4))) 447 call wrtout(unitfi,msg,'COLL') 448 do im1=1,ll*2+1 449 if(cplex_dij==1)& 450 & write(msg,'(12(1x,9(1x,f10.5)))') (paw_ij_all(iatom)%noccmmp(1,im1,im2,ispden),im2=1,ll*2+1) 451 if(cplex_dij==2)& 452 & write(msg,'(12(1x,9(1x,"(",f9.5,",",f9.5,")")))')& 453 & (paw_ij_all(iatom)%noccmmp(:,im1,im2,ispden),im2=1,ll*2+1) 454 call wrtout(unitfi,msg,'COLL') 455 end do 456 write(msg,'(2a)') ch10,' ' 457 call wrtout(unitfi,msg,'COLL') 458 end do 459 ! Transformation matrices: real->complex spherical harmonics 460 if(paw_ij_all(iatom)%ndij==4) then 461 ABI_ALLOCATE(noccmmp_ylm,(2*ll+1,2*ll+1,paw_ij_all(iatom)%ndij)) 462 noccmmp_ylm=czero 463 ABI_ALLOCATE(noccmmp_slm,(2*ll+1,2*ll+1,paw_ij_all(iatom)%ndij)) 464 noccmmp_slm=czero 465 ! Go from real notation for complex noccmmp to complex notation in noccmmp_slm 466 noccmmp_slm(:,:,:)=cmplx(paw_ij_all(iatom)%noccmmp(1,:,:,:),paw_ij_all(iatom)%noccmmp(2,:,:,:)) 467 ii=std_out;if (unitfi==ab_out) ii=-1 468 call mat_slm2ylm(ll,noccmmp_slm,noccmmp_ylm,paw_ij_all(iatom)%ndij,& 469 & 1,1,dtset%pawprtvol,ii,'COLL') ! optspin=1 because up spin are first 470 do ispden=1,paw_ij_all(iatom)%ndij 471 write(msg,'(3a)') ch10,& 472 & "== Occupation matrix in the complex harmonics basis for component ",& 473 & trim(dspin1(ispden+2*(paw_ij_all(iatom)%ndij/4))) 474 call wrtout(unitfi,msg,'COLL') 475 do im1=1,ll*2+1 476 write(msg,'(12(1x,9(1x,"(",f9.5,",",f9.5,")")))') & 477 & (noccmmp_ylm(im1,im2,ispden),im2=1,ll*2+1) 478 call wrtout(unitfi,msg,'COLL') 479 end do 480 end do 481 write(msg,'(a)') ch10 482 call wrtout(unitfi,msg,'COLL') 483 if (dtset%pawspnorb>0) then 484 ABI_ALLOCATE(noccmmp_jmj,(2*(2*ll+1),2*(2*ll+1))) 485 noccmmp_jmj=czero 486 ii=std_out;if (unitfi==ab_out) ii=-1 487 call mat_mlms2jmj(ll,noccmmp_ylm,noccmmp_jmj,paw_ij_all(iatom)%ndij,& 488 & 1,1,dtset%pawprtvol,-1,'COLL') ! optspin=1: up spin are first 489 write(msg,'(3a)') ch10,"== Occupation matrix in the J (= L-1/2, L+1/2) and M_J basis" 490 call wrtout(unitfi,msg,'COLL') 491 do im1=1,2*(ll*2+1) 492 write(msg,'(12(1x,18(1x,"(",f7.3,",",f7.3,")")))') & 493 & (noccmmp_jmj(im1,im2),im2=1,2*(ll*2+1)) 494 call wrtout(unitfi,msg,'COLL') 495 end do 496 write(msg,'(a)') ch10 497 call wrtout(unitfi,msg,'COLL') 498 ABI_DEALLOCATE(noccmmp_jmj) 499 end if ! pawspnorb 500 ABI_DEALLOCATE(noccmmp_ylm) 501 ABI_DEALLOCATE(noccmmp_slm) 502 end if ! ndij==4 503 end if ! ((ll>=0).and.(pawtab(itypat)%usepawu>0)) 504 end do 505 end do 506 end if 507 508 !Exact exchange: print out occupations for correlated orbitals 509 !------------------------------------------------------------- 510 if (useexexch.and.ipositron/=1.and.me_atom==0) then 511 nspden=paw_ij_all(1)%nspden;nsppol=paw_ij_all(1)%nsppol;ndij=paw_ij_all(1)%ndij 512 do iatom=1,dtset%natom 513 itypat=dtset%typat(iatom);ll=pawtab(itypat)%lexexch 514 cplex_dij=paw_ij_all(iatom)%cplex_dij 515 if (ll>=0.and.pawtab(itypat)%useexexch>0) then 516 ABI_ALLOCATE(paw_ij_all(iatom)%noccmmp,(cplex_dij,2*ll+1,2*ll+1,ndij)) 517 ABI_ALLOCATE(paw_ij_all(iatom)%nocctot,(nspden)) 518 end if 519 end do 520 call setnoccmmp(1,0,rdum4,0,0,idum3,dtset%natom,dtset%natom,0,1,nsppol,0,dtset%ntypat,& 521 & paw_ij_all,pawang_dum,dtset%pawprtvol,pawrhoij_all,pawtab,rdum2,idum1,dtset%typat,1,0) 522 do i_unitfi=1,2 523 unitfi=ab_out;if (i_unitfi==2) unitfi=std_out 524 write(msg, '(3a)' ) & 525 & ' ---------- Exact Exchange --------------------------------------------------- ',ch10 526 call wrtout(unitfi,msg,'COLL') 527 do iatom=1,dtset%natom 528 itypat=dtset%typat(iatom);ll=pawtab(itypat)%lexexch 529 cplex_dij=paw_ij_all(iatom)%cplex_dij 530 if ((ll>=0).and.(pawtab(itypat)%useexexch>0)) then 531 write(msg,fmt='(a,i5,a,i4,a)') " ====== For Atom",iatom,& 532 & ", occupations for correlated orbitals. l =",ll,ch10 533 call wrtout(unitfi,msg,'COLL') 534 do ispden=1,ndij 535 if (nspden==1.and.ndij/=4) write(msg,fmt='(a)') " Up component only..." 536 if (nspden==2) write(msg,fmt='(a,i3)')" Occupation matrix for spin",ispden 537 if (ndij==4) write(msg,fmt='(2a)') " Occupation matrix for component ",& 538 & trim(dspin2(ispden+2*(ndij/4))) 539 call wrtout(unitfi,msg,'COLL') 540 do im1=1,ll*2+1 541 if(cplex_dij==1)& 542 & write(msg,'(12(1x,9(1x,f10.5)))')& 543 & (paw_ij_all(iatom)%noccmmp(1,im1,im2,ispden),im2=1,ll*2+1) 544 if(cplex_dij==2)& 545 & write(msg,'(12(1x,9(1x,"(",f7.3,",",f7.3,")")))') & 546 & (paw_ij_all(iatom)%noccmmp(:,im1,im2,ispden),im2=1,ll*2+1) 547 call wrtout(unitfi,msg,'COLL') 548 end do 549 call wrtout(unitfi,' ','COLL') 550 end do 551 end if 552 end do 553 end do 554 do iatom=1,dtset%natom 555 if (allocated(paw_ij_all(iatom)%noccmmp)) then 556 ABI_DEALLOCATE(paw_ij_all(iatom)%noccmmp) 557 end if 558 if (allocated(paw_ij_all(iatom)%nocctot)) then 559 ABI_DEALLOCATE(paw_ij_all(iatom)%nocctot) 560 end if 561 end do 562 end if 563 564 msg=' ' 565 call wrtout(ab_out,msg,'COLL') 566 call wrtout(std_out,msg,'COLL') 567 568 !Destroy temporary stored atomic data 569 ABI_DEALLOCATE(jatom) 570 call free_my_atmtab(my_atmtab,my_atmtab_allocated) 571 if (paral_atom) then 572 call pawrhoij_free(pawrhoij_all) 573 ABI_DATATYPE_DEALLOCATE(pawrhoij_all) 574 if (usepawu.or.useexexch) then 575 call paw_ij_free(paw_ij_all) 576 ABI_DATATYPE_DEALLOCATE(paw_ij_all) 577 end if 578 end if 579 580 DBG_EXIT("COLL") 581 582 end subroutine pawprt