TABLE OF CONTENTS
- ABINIT/m_paw_pwaves_lmn
- m_paw_pwaves_lmn/paw_pwaves_lmn_free
- m_paw_pwaves_lmn/paw_pwaves_lmn_init
- m_paw_pwaves_lmn/paw_pwaves_lmn_t
ABINIT/m_paw_pwaves_lmn [ Modules ]
NAME
m_paw_pwaves_lmn
FUNCTION
This module contains the definition of the paw_pwaves_lmn structured datatype, as well as related functions and methods. paw_pwaves_lmn is used to store the 3D values of the all-electron and of the pseudized part of the PAW partial waves on the set of FFT points falling inside the spheres around each atom.
COPYRIGHT
Copyright (C) 2008-2024 ABINIT group (MG,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 .
NOTES
* Routines tagged with "@type_name" are strongly connected to the definition of the data type. Strongly connected means that the proper functioning of the implementation relies on the assumption that the tagged procedure is consistent with the type declaration. Every time a developer changes the structure "type_name" adding new entries, he/she has to make sure that all the strongly connected routines are changed accordingly to accommodate the modification of the data type. Typical examples of strongly connected routines are creation, destruction or reset methods.
SOURCE
28 #if defined HAVE_CONFIG_H 29 #include "config.h" 30 #endif 31 32 #include "abi_common.h" 33 34 MODULE m_paw_pwaves_lmn 35 36 use defs_basis 37 use m_errors 38 use m_abicore 39 use m_xmpi 40 use m_sort 41 42 use m_numeric_tools, only : wrap2_zero_one 43 use m_geometry, only : xcart2xred 44 use m_pawrad, only : pawrad_type, nderiv_gen, pawrad_deducer0 45 use m_pawtab, only : pawtab_type 46 use m_pawfgrtab, only : pawfgrtab_type 47 use m_paw_numeric, only : paw_spline, paw_splint 48 use m_paw_sphharm, only : initylmr 49 use m_paral_atom, only : get_my_atmtab, free_my_atmtab 50 51 implicit none 52 53 private
m_paw_pwaves_lmn/paw_pwaves_lmn_free [ Functions ]
[ Top ] [ m_paw_pwaves_lmn ] [ Functions ]
NAME
paw_pwaves_lmn_free
FUNCTION
INPUTS
OUTPUT
SOURCE
429 subroutine paw_pwaves_lmn_free(Paw_onsite) 430 431 implicit none 432 433 !Arguments ------------------------------------ 434 !arrays 435 type(paw_pwaves_lmn_t),intent(inout) :: Paw_onsite(:) 436 437 !Local variables------------------------------- 438 !scalars 439 integer :: iatom 440 441 ! ************************************************************************* 442 443 !@paw_pwaves_lmn_t 444 445 do iatom=LBOUND(Paw_onsite,DIM=1),UBOUND(Paw_onsite,DIM=1) 446 if (allocated(Paw_onsite(iatom)%phi)) then 447 ABI_FREE(Paw_onsite(iatom)%phi) 448 end if 449 if (allocated(Paw_onsite(iatom)%tphi)) then 450 ABI_FREE(Paw_onsite(iatom)%tphi) 451 end if 452 if (allocated(Paw_onsite(iatom)%r0shift)) then 453 ABI_FREE(Paw_onsite(iatom)%r0shift) 454 end if 455 if (allocated(Paw_onsite(iatom)%phi_gr )) then 456 ABI_FREE(Paw_onsite(iatom)%phi_gr) 457 end if 458 if (allocated(Paw_onsite(iatom)%tphi_gr)) then 459 ABI_FREE(Paw_onsite(iatom)%tphi_gr) 460 end if 461 end do 462 463 end subroutine paw_pwaves_lmn_free
m_paw_pwaves_lmn/paw_pwaves_lmn_init [ Functions ]
[ Top ] [ m_paw_pwaves_lmn ] [ Functions ]
NAME
paw_pwaves_lmn_init
FUNCTION
INPUTS
mpi_atmtab(:)=--optional-- indexes of the atoms treated by current proc comm_atom=--optional-- MPI communicator over atoms
OUTPUT
SOURCE
126 subroutine paw_pwaves_lmn_init(Paw_onsite,my_natom,natom,ntypat,rprimd,xcart,Pawtab, & 127 & Pawrad,local_pawfgrtab,optgrad,& 128 & mpi_atmtab,comm_atom) ! optional arguments (parallelism) 129 130 implicit none 131 132 !Arguments ------------------------------------ 133 !scalars 134 integer,intent(in) :: my_natom,natom,ntypat 135 integer,optional,intent(in) :: optgrad,comm_atom 136 !arrays 137 integer,optional,target,intent(in) :: mpi_atmtab(:) 138 real(dp),intent(in) :: xcart(3,natom),rprimd(3,3) 139 type(pawtab_type),target,intent(in) :: Pawtab(ntypat) 140 type(pawrad_type),intent(in) :: Pawrad(ntypat) 141 type(pawfgrtab_type),intent(in) :: local_pawfgrtab(my_natom) 142 type(paw_pwaves_lmn_t),intent(out) :: Paw_onsite(my_natom) 143 144 !Local variables------------------------------- 145 !scalars 146 integer :: itypat,ln_size,lmn_size,mesh_size,inl,iatom,iatom1,my_comm_atom,my_optgrad 147 integer :: nfgd,ifgd,ipsang,option_ylmr,normchoice,ii,jlmn,jl,jm,jlm,jln 148 logical :: my_atmtab_allocated,paral_atom 149 real(dp) :: phj,rR,tphj,ybcbeg,ybcend 150 !arrays 151 integer, allocatable :: iperm(:) 152 integer,ABI_CONTIGUOUS pointer :: indlmn(:,:) 153 integer,pointer :: my_atmtab(:) 154 real(dp) :: yvals(4),red(3),shift(3) 155 real(dp),allocatable :: ff(:),nrm(:),nrm_sort(:),phigrd(:,:),tphigrd(:,:),ylm_tmp(:,:),ylm(:,:),ylm_gr(:,:,:) 156 real(dp),allocatable :: rsph_red(:,:),rsph_cart(:,:),phigrd_gr(:,:),tphigrd_gr(:,:),gg(:) 157 type(paw_pwaves_lmn_t),allocatable :: Paw_lmn_spline(:) 158 159 ! ************************************************************************* 160 161 !@paw_pwaves_lmn_t 162 163 if (my_natom==0) return 164 my_optgrad = -1; if (present(optgrad)) my_optgrad = optgrad 165 166 ABI_CHECK(all(local_pawfgrtab(:)%rfgd_allocated==1),"R vectors not allocated in pawfgrtab") 167 168 ! Set up parallelism over atoms 169 paral_atom=(present(comm_atom).and.(my_natom/=natom)) 170 nullify(my_atmtab);if (present(mpi_atmtab)) my_atmtab => mpi_atmtab 171 my_comm_atom=xmpi_comm_self;if (present(comm_atom)) my_comm_atom=comm_atom 172 call get_my_atmtab(my_comm_atom,my_atmtab,my_atmtab_allocated,paral_atom,natom,my_natom_ref=my_natom) 173 174 ! Prepare the spline. Calculate 2nd derivatives of partial waves for each atom type. 175 ABI_MALLOC(Paw_lmn_spline,(ntypat)) 176 177 do itypat=1,ntypat 178 ln_size =Pawtab(itypat)%basis_size 179 mesh_size=Pawtab(itypat)%mesh_size 180 181 ABI_MALLOC(Paw_lmn_spline(itypat)%phi ,(mesh_size,ln_size)) 182 ABI_MALLOC(Paw_lmn_spline(itypat)%tphi,(mesh_size,ln_size)) 183 184 do inl=1,ln_size ! Calculate 2nd derivatives of %phi and %tphi for each ln component. 185 ybcbeg=zero; ybcend=zero 186 call paw_spline(Pawrad(itypat)%rad,Pawtab(itypat)%phi(:,inl), mesh_size,& 187 & ybcbeg,ybcend,Paw_lmn_spline(itypat)%phi(:,inl)) 188 189 ybcbeg=zero; ybcend=zero 190 call paw_spline(Pawrad(itypat)%rad,Pawtab(itypat)%tphi(:,inl),mesh_size,& 191 & ybcbeg,ybcend,Paw_lmn_spline(itypat)%tphi(:,inl)) 192 end do 193 end do 194 ! 195 ! === spline for each atom === 196 ! * FFT points within PAW sphere depend on the atom site. 197 do iatom1=1,my_natom 198 iatom=iatom1;if (paral_atom) iatom=my_atmtab(iatom1) 199 200 itypat = local_pawfgrtab(iatom1)%itypat 201 ln_size = Pawtab(itypat)%basis_size 202 lmn_size = Pawtab(itypat)%lmn_size 203 mesh_size= Pawrad(itypat)%mesh_size 204 nfgd = local_pawfgrtab(iatom1)%nfgd ! no. of points in the fine grid for this PAW sphere 205 indlmn => Pawtab(itypat)%indlmn 206 207 ! The points in the PAW sphere might belong to a different unit cell, in this case one has to 208 ! reconstruct the contribution to the AE psi_k given by the lattice-symmetric atom in another sphere. 209 ! Here I wrap rr back into the first unit cell keeping trace of the lattice vector L 210 ! needed so that rr = rr_first_ucell + L 211 ! The contribution to the AE u(r) in the first unit cell has to be multiplied by e^{-ikL}. 212 213 Paw_onsite(iatom1)%nfgd = nfgd 214 Paw_onsite(iatom1)%lmn_size = lmn_size 215 ! 216 ABI_MALLOC(Paw_onsite(iatom1)%r0shift,(3,nfgd)) 217 218 ABI_MALLOC(rsph_red,(3,nfgd)) 219 ABI_MALLOC(rsph_cart,(3,nfgd)) 220 do ifgd=1,nfgd 221 rsph_cart(:,ifgd) = local_pawfgrtab(iatom1)%rfgd(:,ifgd) + xcart(:,iatom) 222 end do 223 call xcart2xred(nfgd,rprimd,rsph_cart,rsph_red) ! go to reduced coordinates. 224 ABI_FREE(rsph_cart) 225 226 do ifgd=1,nfgd 227 call wrap2_zero_one(rsph_red(1,ifgd),red(1),shift(1)) ! rr = r_cell + shift 228 call wrap2_zero_one(rsph_red(2,ifgd),red(2),shift(2)) 229 call wrap2_zero_one(rsph_red(3,ifgd),red(3),shift(3)) 230 Paw_onsite(iatom1)%r0shift(:,ifgd) = NINT(shift) 231 !if (ANY( ABS(shift) > tol12)) then 232 !ABI_WARNING("rmR_red is outside the first unit cell.") 233 !write(ab_out,*)rsph_red(:,ifgd),shift 234 !end if 235 end do 236 ABI_FREE(rsph_red) 237 ! 238 ! * Obtain |r-R| on fine grid, note that rfgd is given in Cartesian coordinates. 239 ABI_MALLOC(nrm,(nfgd)) 240 do ifgd=1,nfgd 241 nrm(ifgd) = sqrt(dot_product(local_pawfgrtab(iatom1)%rfgd(:,ifgd),local_pawfgrtab(iatom1)%rfgd(:,ifgd))) 242 end do 243 ! 244 ! * Compute Ylm for each r-R vector. 245 ipsang = 1 + (Pawtab(itypat)%l_size-1)/2 ! recall l_size=2*l_max-1 where l_max is shifted by 1. 246 ABI_MALLOC(ylm_tmp,(ipsang**2,nfgd)) 247 normchoice = 1 ! Use computed norms of input vectors. 248 if (my_optgrad==1) then 249 option_ylmr=2 ! Compute Ylm(r-R) and its gradient 250 ABI_MALLOC(ylm_gr,(3,ipsang**2,nfgd)) 251 else 252 option_ylmr= 1 ! To compute Ylm(r-R). 253 ABI_MALLOC(ylm_gr,(3,3,0)) 254 end if 255 call initylmr(ipsang,normchoice,nfgd,nrm,option_ylmr,local_pawfgrtab(iatom1)%rfgd,ylm_tmp,ylm_gr) 256 ! 257 ! Exchange dimensions for better memory access. 258 ABI_MALLOC(ylm,(nfgd,ipsang**2)) 259 do ii=1,ipsang**2 260 ylm(:,ii) = ylm_tmp(ii,:) 261 end do 262 ABI_FREE(ylm_tmp) 263 ! 264 ! In order to do spline fits, the |r-R| data must be sorted 265 ! Here we sort the nrm points, keeping track of which goes where 266 ABI_MALLOC(nrm_sort,(nfgd)) 267 nrm_sort = nrm 268 269 ABI_MALLOC(iperm,(nfgd)) 270 do ifgd=1,nfgd 271 iperm(ifgd)=ifgd 272 end do 273 274 call sort_dp(nfgd,nrm_sort,iperm,tol8) 275 276 ! Now make spline fits of phi and tphi onto the fine grid around the atom 277 ABI_MALLOC(phigrd,(nfgd,ln_size)) 278 ABI_MALLOC(tphigrd,(nfgd,ln_size)) 279 ABI_MALLOC(ff,(nfgd)) 280 if (my_optgrad==1) then 281 ABI_MALLOC(phigrd_gr,(nfgd,ln_size)) 282 ABI_MALLOC(tphigrd_gr,(nfgd,ln_size)) 283 ABI_MALLOC(gg,(mesh_size)) 284 end if 285 286 do inl=1,ln_size 287 ! 288 ! * splint phi onto points and reorder indices. 289 call paw_splint(mesh_size,Pawrad(itypat)%rad,Pawtab(itypat)%phi(:,inl),& 290 & Paw_lmn_spline(itypat)%phi(:,inl),nfgd,nrm_sort,ff) 291 do ifgd=1,nfgd 292 ii=iperm(ifgd) 293 phigrd(ii,inl) = ff(ifgd) 294 end do 295 ! 296 ! * compute d phi/dr, interpolate onto points and reorder indices. 297 if (my_optgrad==1) then 298 ybcbeg=zero; ybcend=zero 299 call nderiv_gen(gg,Pawtab(itypat)%phi(:,inl),Pawrad(itypat)) 300 call paw_spline(Pawrad(itypat)%rad,gg,mesh_size,ybcbeg,ybcend,& 301 & Paw_lmn_spline(itypat)%phi(:,inl)) 302 call paw_splint(mesh_size,Pawrad(itypat)%rad,Paw_lmn_spline(itypat)%phi(:,inl),& 303 & Paw_lmn_spline(itypat)%phi(:,inl),nfgd,nrm_sort,ff) 304 do ifgd=1,nfgd 305 ii=iperm(ifgd) 306 phigrd_gr(ii,inl) = ff(ifgd) 307 end do 308 end if 309 ! 310 ! * compute d tphi/dr, interpolate onto points and reorder indices. 311 call paw_splint(mesh_size,Pawrad(itypat)%rad,Pawtab(itypat)%tphi(:,inl),& 312 & Paw_lmn_spline(itypat)%tphi(:,inl),nfgd,nrm_sort,ff) 313 do ifgd=1,nfgd 314 ii=iperm(ifgd) 315 tphigrd(ii,inl) = ff(ifgd) 316 end do 317 if (my_optgrad==1) then 318 ybcbeg=zero; ybcend=zero 319 call nderiv_gen(gg,Pawtab(itypat)%tphi(:,inl),Pawrad(itypat)) 320 call paw_spline(Pawrad(itypat)%rad,gg,mesh_size,ybcbeg,ybcend,& 321 & Paw_lmn_spline(itypat)%tphi(:,inl)) 322 call paw_splint(mesh_size,Pawrad(itypat)%rad,Paw_lmn_spline(itypat)%tphi(:,inl),& 323 & Paw_lmn_spline(itypat)%phi(:,inl),nfgd,nrm_sort,ff) 324 do ifgd=1,nfgd 325 ii=iperm(ifgd) 326 tphigrd_gr(ii,inl) = ff(ifgd) 327 end do 328 end if 329 end do !inl 330 331 ABI_FREE(ff) 332 if (my_optgrad==1) then 333 ABI_FREE(gg) 334 end if 335 ! 336 ! === Calculate AE and PS partial waves inside the sphere === 337 ! * recall that <r|phi>=u(r)*Slm(r^)/r, hence avoid division by zero except for s-waves. 338 ABI_MALLOC(Paw_onsite(iatom1)%phi ,(nfgd,lmn_size)) 339 ABI_MALLOC(Paw_onsite(iatom1)%tphi,(nfgd,lmn_size)) 340 341 if (my_optgrad==1) then 342 ABI_MALLOC(Paw_onsite(iatom1)%phi_gr ,(3,nfgd,lmn_size)) 343 ABI_MALLOC(Paw_onsite(iatom1)%tphi_gr,(3,nfgd,lmn_size)) 344 end if 345 346 do jlmn=1,lmn_size 347 jl = indlmn(1,jlmn) 348 jm = indlmn(2,jlmn) 349 jlm = indlmn(4,jlmn) 350 jln = indlmn(5,jlmn) 351 352 do ifgd=1,nfgd ! loop over fine grid points in current PAW sphere 353 !if (nrm(ifgd)>tol16) then 354 if (nrm(ifgd)>tol10) then ! tol10 to be consistent with initylmr. 355 rR = nrm(ifgd) ! value of |r-R| 356 !write(ab_out,*) 'rR:',rR,' phigrd:',phigrd(ifgd,jln),' tphigrd:',tphigrd(ifgd,jln),' ylm:',ylm(ifgd,jlm) 357 phj = phigrd (ifgd,jln)*ylm(ifgd,jlm)/rR 358 tphj= tphigrd(ifgd,jln)*ylm(ifgd,jlm)/rR 359 Paw_onsite(iatom1)%phi (ifgd,jlmn) = phj 360 Paw_onsite(iatom1)%tphi(ifgd,jlmn) = tphj 361 362 if (my_optgrad==1) then 363 Paw_onsite(iatom1)%phi_gr (1:3,ifgd,jlmn) = phigrd (ifgd,jln)*ylm_gr(1:3,jlm,ifgd) & 364 & + phigrd_gr(ifgd,jln)*local_pawfgrtab(iatom1)%rfgd(1:3,ifgd)*ylm(ifgd,jlm) 365 Paw_onsite(iatom1)%tphi_gr (1:3,ifgd,jlmn) = tphigrd (ifgd,jln)*ylm_gr(1:3,jlm,ifgd) & 366 & + tphigrd_gr(ifgd,jln)*local_pawfgrtab(iatom1)%rfgd(1:3,ifgd)*ylm(ifgd,jlm) 367 end if 368 369 else ! Extrapolate if the point is at the origin 370 yvals(1) = zero 371 if (jl==0) then 372 yvals(2:4) = Pawtab(itypat)%phi(2:4,jln)/Pawrad(itypat)%rad(2:4) 373 call pawrad_deducer0(yvals,4,Pawrad(itypat)) 374 end if 375 Paw_onsite(iatom1)%phi(ifgd,jlmn) = yvals(1) * ylm(ifgd,jlm) 376 yvals(1) = zero 377 if (jl==0) then 378 yvals(2:4) = Pawtab(itypat)%tphi(2:4,jln)/Pawrad(itypat)%rad(2:4) 379 call pawrad_deducer0(yvals,4,pawrad(itypat)) 380 end if 381 Paw_onsite(iatom1)%tphi(ifgd,jlmn) = yvals(1) * ylm(ifgd,jlm) 382 ! The gradient is expected to go to zero at the origin 383 if (my_optgrad==1) then 384 Paw_onsite(iatom1)%phi_gr (1:3,ifgd,jlmn) = zero 385 Paw_onsite(iatom1)%tphi_gr (1:3,ifgd,jlmn) = zero 386 end if 387 end if 388 389 end do !nfgd 390 end do !jlmn 391 392 ABI_FREE(nrm) 393 ABI_FREE(nrm_sort) 394 ABI_FREE(iperm) 395 ABI_FREE(phigrd) 396 ABI_FREE(tphigrd) 397 ABI_FREE(ylm) 398 ABI_FREE(ylm_gr) 399 if (my_optgrad==1) then 400 ABI_FREE(phigrd_gr) 401 ABI_FREE(tphigrd_gr) 402 end if 403 end do !iatom 404 ! 405 !* Free 2nd derivates used for spline. 406 call paw_pwaves_lmn_free(Paw_lmn_spline) 407 ABI_FREE(Paw_lmn_spline) 408 409 ! Destroy atom table used for parallelism 410 call free_my_atmtab(my_atmtab,my_atmtab_allocated) 411 412 end subroutine paw_pwaves_lmn_init
m_paw_pwaves_lmn/paw_pwaves_lmn_t [ Types ]
[ Top ] [ m_paw_pwaves_lmn ] [ Types ]
NAME
paw_pwaves_lmn_t
FUNCTION
Datatype used to store the 3D values of the all-electron and of the pseudized part of the PAW partial waves on the set of FFT points falling inside the spheres around each atom. The data is mainly used for plotting the true PAW wavefunctions in real space.
SOURCE
69 type,public :: paw_pwaves_lmn_t 70 71 integer :: nfgd 72 73 integer :: lmn_size 74 75 !$integer :: ngfft(18) 76 77 integer,allocatable :: r0shift(:,:) 78 ! r0shift(3,nfgd) 79 80 !real(dp),allocatable:: phk_atm(:,:) 81 ! phk_atmt(2,nfgd) 82 83 real(dp),allocatable :: phi(:,:) 84 ! phi (nfgd,lmn_size) 85 ! \phi_{nlm}(ifgd) for each point of the FFT mesh located inside the PAW sphere (see pawfgrtab_type). 86 87 real(dp),allocatable :: tphi(:,:) 88 ! tphi (nfgd,lmn_size) 89 ! \tphi_{nlm}(ifgd) for each point of the FFT mesh located inside the PAW sphere (see pawfgrtab_type). 90 91 real(dp),allocatable :: phi_gr(:,:,:) 92 ! phi_gr (3,nfgd,lmn_size) 93 ! gradient, in cartesian coordinates, of \phi_{nlm}(ifgd) for each point of the FFT mesh 94 ! located inside the PAW sphere (see pawfgrtab_type). 95 96 real(dp),allocatable :: tphi_gr(:,:,:) 97 ! tphi_gr (3,nfgd,lmn_size) 98 ! gradient, in cartesian coordinates, of \tphi_{nlm}(ifgd) for each point of the FFT mesh 99 ! located inside the PAW sphere (see pawfgrtab_type). 100 101 end type paw_pwaves_lmn_t 102 103 public :: paw_pwaves_lmn_init 104 public :: paw_pwaves_lmn_free