TABLE OF CONTENTS
ABINIT/pawuj_red [ Functions ]
NAME
pawuj_red
FUNCTION
Store atomic occupancies, potential shift, positions in dtpawuj datastructure.
COPYRIGHT
Copyright (C) 1998-2018 ABINIT group (DJA) 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
fatvshift=factor that multiplies atvshift 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 natom=number of atoms in cell ntypat = number of atom types paw_ij(my_natom) <type(paw_ij_type)>=paw arrays given on (i,j) channels pawprtvol= printing volume pawtab(ntypat*usepaw) <type(pawtab_type)>=paw tabulated starting data
OUTPUT
dtpawuj(0:ndtpawuj) (initialization of fields vsh, occ, iuj,nnat)
PARENTS
scfcv
CHILDREN
free_my_atmtab,get_my_atmtab,linvmat,prmat,wrtout,xmpi_lor,xmpi_sum
SOURCE
39 #if defined HAVE_CONFIG_H 40 #include "config.h" 41 #endif 42 43 #include "abi_common.h" 44 45 subroutine pawuj_red(dtset,dtpawuj,fatvshift,my_natom,natom,ntypat,paw_ij,pawrad,pawtab,ndtpawuj,& 46 & mpi_atmtab,comm_atom) ! optional arguments (parallelism) 47 48 use defs_basis 49 use defs_abitypes 50 use defs_parameters 51 use m_profiling_abi 52 use m_xmpi 53 54 use m_pptools, only : prmat 55 use m_pawrad, only : pawrad_type 56 use m_pawtab, only : pawtab_type 57 use m_paw_ij, only : paw_ij_type 58 use m_paral_atom, only : get_my_atmtab, free_my_atmtab 59 60 !This section has been created automatically by the script Abilint (TD). 61 !Do not modify the following lines by hand. 62 #undef ABI_FUNC 63 #define ABI_FUNC 'pawuj_red' 64 use interfaces_14_hidewrite 65 use interfaces_65_paw, except_this_one => pawuj_red 66 !End of the abilint section 67 68 implicit none 69 70 !Arguments ------------------------------------ 71 !scalars 72 integer,intent(in) :: my_natom,natom,ntypat,ndtpawuj 73 integer,optional,intent(in) :: comm_atom 74 real(dp),intent(in) :: fatvshift 75 !arrays 76 integer,optional,target,intent(in) :: mpi_atmtab(:) 77 type(paw_ij_type),intent(in) :: paw_ij(my_natom) 78 type(pawtab_type),intent(in) :: pawtab(ntypat) 79 type(pawrad_type),intent(in) :: pawrad(ntypat) 80 type(dataset_type),intent(in) :: dtset 81 type(macro_uj_type),intent(inout) :: dtpawuj(0:ndtpawuj) 82 83 !Local variables------------------------------- 84 !scalars 85 integer,parameter :: natmax=2,ncoeff=3 86 integer :: iatom,iatom_tot,ierr,im1,im2,ispden,itypat,ll,nspden,nsppol,iuj 87 integer :: my_comm_atom,nnat,natpawu,natvshift,pawujat,ndtset,typawujat 88 logical :: usepawu !antiferro, 89 logical :: my_atmtab_allocated,paral_atom 90 character(len=1000) :: message,hstr 91 character(len=500) :: messg 92 !arrays 93 logical :: musk(3,natom) 94 integer,pointer :: my_atmtab(:) 95 real(dp) :: rrtab(ncoeff),wftab(ncoeff),a(ncoeff,ncoeff),b(ncoeff,ncoeff)! ,a(ncoeff,ncoeff) 96 real(dp),allocatable :: nnocctot(:,:) !,magv(:) 97 real(dp),allocatable :: atvshift(:,:,:) ! atvshift(natvshift,2,natom) 98 logical,allocatable :: dmusk(:,:),atvshmusk(:,:,:) !atvshmusk(natvshift,2,natom) 99 100 ! ********************************************************************* 101 102 !Initializations 103 nspden=1;nsppol=1 104 if (my_natom>0) then 105 nspden=paw_ij(1)%nspden ; nsppol=paw_ij(1)%nsppol 106 end if 107 108 natvshift=dtset%natvshift 109 pawujat=dtset%pawujat 110 natpawu=dtset%natpawu ; ndtset=dtset%ndtset 111 ABI_ALLOCATE(atvshift,(natvshift,nspden,natom)) 112 ABI_ALLOCATE(atvshmusk,(natvshift,nspden,natom)) 113 ABI_ALLOCATE(dmusk,(nspden,natom)) 114 ABI_ALLOCATE(nnocctot,(nspden,natom)) 115 musk=.false.; dmusk=.false. 116 atvshift=fatvshift*dtset%atvshift 117 atvshmusk=.false. 118 nnocctot=zero 119 typawujat=dtset%typat(pawujat) 120 usepawu=(count(pawtab(:)%usepawu>0)>0) 121 122 !Set up parallelism over atoms 123 paral_atom=(present(comm_atom).and.(my_natom/=natom)) 124 nullify(my_atmtab);if (present(mpi_atmtab)) my_atmtab => mpi_atmtab 125 my_comm_atom=xmpi_comm_self;if (present(comm_atom)) my_comm_atom=comm_atom 126 call get_my_atmtab(my_comm_atom,my_atmtab,my_atmtab_allocated,paral_atom,natom,my_natom_ref=my_natom) 127 128 nnat=0 129 if (usepawu) then 130 write(message,'(3a)') ch10, '---------- pawuj_red ------ ',ch10 131 call wrtout(std_out, message,'COLL'); 132 do iatom=1,my_natom 133 iatom_tot=iatom;if (paral_atom) iatom_tot=my_atmtab(iatom) 134 itypat=paw_ij(iatom)%itypat;ll=pawtab(itypat)%lpawu 135 if ((ll>=0).and.(pawtab(itypat)%usepawu>0).and.itypat==typawujat) then 136 musk(:,iatom_tot)=(/.true., .true., .true. /) 137 atvshmusk(:,:,iatom_tot)=reshape((/ (( (im1==1), im1=1,natvshift) ,im2=1,nspden ) /),(/natvshift,nspden/)) 138 do ispden=1,nspden 139 nnocctot(ispden,iatom_tot)=paw_ij(iatom)%nocctot(ispden) 140 dmusk(ispden,iatom_tot)=.true. 141 end do 142 nnat=nnat+1 143 end if 144 end do 145 146 ! Reduction in case of parallelism 147 if (paral_atom) then 148 call xmpi_sum(nnocctot ,my_comm_atom,ierr) 149 call xmpi_lor(atvshmusk,my_comm_atom) ! dim=natpawu ??? 150 call xmpi_lor(dmusk ,my_comm_atom) 151 call xmpi_lor(musk ,my_comm_atom) 152 call xmpi_sum(nnat ,my_comm_atom,ierr) 153 end if 154 155 iuj=maxval(dtpawuj(:)%iuj) 156 !write(std_out,*)' pawuj_red: iuj',iuj 157 !write(std_out,*)' pawuj_red: dtpawuj(:)%iuj ',dtpawuj(:)%iuj 158 159 if (iuj==1.or.iuj==3) then ! 1 and 3: non-scf steps 160 dtpawuj(iuj+1)%iuj=iuj+1 161 end if 162 163 ! TODO: check that this is correct: this point is passed several times 164 ! for a given value of iuj - should the stuff be accumulated instead of replaced? 165 if(allocated(dtpawuj(iuj)%vsh)) then 166 ABI_DEALLOCATE(dtpawuj(iuj)%vsh) 167 end if 168 ABI_ALLOCATE(dtpawuj(iuj)%vsh,(nspden,nnat)) 169 if(allocated(dtpawuj(iuj)%occ)) then 170 ABI_DEALLOCATE(dtpawuj(iuj)%occ) 171 end if 172 ABI_ALLOCATE(dtpawuj(iuj)%occ,(nspden,nnat)) 173 if(allocated(dtpawuj(iuj)%xred)) then 174 ABI_DEALLOCATE(dtpawuj(iuj)%xred) 175 end if 176 ABI_ALLOCATE(dtpawuj(iuj)%xred,(3,nnat)) 177 178 if (iuj==1) then 179 ABI_ALLOCATE(dtpawuj(0)%vsh,(nspden,nnat)) 180 ABI_ALLOCATE(dtpawuj(0)%occ,(nspden,nnat)) 181 ABI_ALLOCATE(dtpawuj(0)%xred,(3,nnat)) 182 dtpawuj(0)%vsh=0 183 dtpawuj(0)%occ=0 184 dtpawuj(0)%xred=0 185 end if 186 187 rrtab=(/0.75_dp,0.815_dp,1.0_dp/)*pawtab(typawujat)%rpaw 188 wftab=pawtab(typawujat)%phi(pawtab(typawujat)%mesh_size,pawtab(typawujat)%lnproju(1)) 189 190 ! DEBUG 191 ! write(std_out,*)' pawuj_red: rrtab ',rrtab 192 ! write(std_out,*)' pawuj_red: wftab ',wftab 193 ! END DEBUG 194 195 do im1=1,ncoeff 196 ! DEBUG 197 ! write(std_out,*)' pawuj_red: ncoeff ',ncoeff,' im1 ',im1 198 ! END DEBUG 199 if (pawrad(typawujat)%mesh_type==1) then 200 im2=nint(rrtab(im1)/pawrad(typawujat)%rstep+1) 201 else if (pawrad(typawujat)%mesh_type==2) then 202 im2=nint(log(rrtab(im1)/pawrad(typawujat)%rstep+1)/pawrad(typawujat)%lstep+1) 203 else if (pawrad(typawujat)%mesh_type==3) then 204 im2=nint(log(rrtab(im1)/pawrad(typawujat)%rstep)/pawrad(typawujat)%lstep+1) 205 else if (pawrad(typawujat)%mesh_type==4) then 206 im2=nint(pawtab(typawujat)%mesh_size*(1-exp((-one)*rrtab(im1)/pawrad(typawujat)%rstep))+1) 207 end if 208 209 ! DEBUG 210 ! write(std_out,*)' pawuj_red: im2 ',im2 211 ! END DEBUG 212 213 rrtab(im1)=pawrad(typawujat)%rad(im2) 214 wftab(im1)=pawtab(typawujat)%phi(im2,pawtab(typawujat)%lnproju(1)) 215 end do 216 write(message,fmt='(a,i3,a,10f10.5)')' pawuj_red: mesh_type',pawrad(typawujat)%mesh_type,' rrtab:', rrtab 217 call wrtout(std_out,message,'COLL') 218 write(message,fmt='(a,10f10.5)')' pawuj_red: wftab', wftab 219 call wrtout(std_out,message,'COLL') 220 a=reshape((/ (( rrtab(im2)**(im1-1), im1=1,3) ,im2=1,3 )/),(/ncoeff,ncoeff/)) 221 write(messg,fmt='(a)')'A' 222 call linvmat(a,b,ncoeff,messg,2,0.0_dp,3) ! linvmat(inmat,oumat,nat,nam,option,gam,prtvol) 223 write(std_out,*) 'pawuj_red: a,b ', a,b 224 wftab=matmul(wftab,b) 225 write(std_out,*) 'pawuj_red: wftab ', wftab 226 dtpawuj(iuj)%wfchr(4:6)=wftab 227 228 dtpawuj(iuj)%nat=nnat 229 write(std_out,*) 'pawuj_red: m1' 230 dtpawuj(iuj)%vsh=reshape(pack(atvshift,atvshmusk),(/ nspden,nnat /)) 231 ! factor in next line to compensate nocctot contains just occ of 1 spin channel for nspden=1 232 write(std_out,*) 'pawuj_red: m2' 233 dtpawuj(iuj)%occ=reshape(pack(nnocctot,dmusk),(/nspden,nnat/))*(3-nspden) 234 write(std_out,*) 'pawuj_red: m3' 235 ! dtpawuj(iuj)%occ=dtpawuj(iuj)%occ/pawtab(typawujat)%ph0phiint(1) 236 237 write(std_out,*) 'pawuj_red: occ ', dtpawuj(iuj)%occ 238 239 dtpawuj(iuj)%xred=reshape(pack(dtset%xred_orig(:,:,1),musk),(/3,nnat/)) 240 dtpawuj(iuj)%ph0phiint=pawtab(typawujat)%ph0phiint(1) 241 dtpawuj(iuj)%wfchr(1:3)=(/ pawtab(typawujat)%zioneff(1)*(dtset%lpawu(typawujat)+2),& 242 & one*(dtset%lpawu(typawujat)+1),one*(dtset%lpawu(typawujat))/) 243 dtpawuj(iuj)%pawrad=pawtab(typawujat)%rpaw 244 245 write(std_out,*) 'pawuj_red: wfchr ',dtpawuj(iuj)%wfchr 246 247 248 write (hstr,'(I0)') iuj 249 write(message,'(a,a,I3,I3,a)') ch10, '---------- MARK ------ ',iuj,maxval(dtpawuj(:)%iuj) ,ch10 250 call wrtout(std_out,message,'COLL') 251 write(message,fmt='(a)') 'vsh'//trim(hstr) 252 call wrtout(std_out,message,'COLL') 253 call prmat(dtpawuj(iuj)%vsh(:,:),1,nnat*nspden,1) 254 write(message,fmt='(a)') 'occ'//trim(hstr) 255 call wrtout(std_out,message,'COLL') 256 call prmat(dtpawuj(iuj)%occ(:,:),1,nnat*nspden,1) 257 write(message, '(3a)' )'---------- MARK ---------- ',ch10 258 call wrtout(std_out,message,'COLL') 259 end if !usepawu 260 261 ABI_DEALLOCATE(nnocctot) 262 ABI_DEALLOCATE(dmusk) 263 ABI_DEALLOCATE(atvshift) 264 ABI_DEALLOCATE(atvshmusk) 265 266 !Destroy atom table used for parallelism 267 call free_my_atmtab(my_atmtab,my_atmtab_allocated) 268 269 end subroutine pawuj_red