TABLE OF CONTENTS
ABINIT/m_paw_onsite [ Modules ]
NAME
m_paw_onsite
FUNCTION
This module contains a set of routines to compute various PAW on-site quantities i.e. quantities expressed with <Phi_i|...|Phi_j> and/or <tild_Phi_i|...|tild_Phi_j>.
COPYRIGHT
Copyright (C) 2013-2018 ABINIT group (MT,FJ) 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
FOR DEVELOPPERS: in order to preserve the portability of libPAW library, please consult ~abinit/src/??_libpaw/libpaw-coding-rules.txt
SOURCE
22 #include "libpaw.h" 23 24 MODULE m_paw_onsite 25 26 USE_DEFS 27 USE_MSG_HANDLING 28 USE_MEMORY_PROFILING 29 30 use m_pawrad, only : pawrad_type, pawrad_deducer0, simp_gen, nderiv_gen 31 use m_pawtab, only : pawtab_type 32 use m_paw_sphharm, only : setnabla_ylm 33 34 implicit none 35 36 private 37 38 !public procedures. 39 public :: pawnabla_init ! Evaluate valence-valence on-site contribs of the nabla operator in cart. coord. 40 public :: pawnabla_core_init ! Evaluate core-valence on-site contribs of the nabla operator in cart. coord.
m_paw_onsite/pawnabla_core_init [ Functions ]
[ Top ] [ m_paw_onsite ] [ Functions ]
NAME
pawnabla_core_init
FUNCTION
Evaluate core-valence onsite contributions of the nabla operator in cartesian coordinates. Core wave-functions are only given for one atom type. Store values in the pawtab% data structure.
INPUTS
mpsang=1+maximum angular momentum
ntypat=Number of types of atoms in cell
Pawrad(ntypat)<Pawrad_type>=PAW radial mesh and related data:
%rad(mesh_size)=The coordinates of all the points of the radial mesh
Pawtab(ntypat) <type(pawtab_type>=PAW tabulated starting data:
%mesh_size=Dimension of radial mesh
%lmn_size=Number of (l,m,n) elements for the PAW basis
phi_cor(mesh_size,nphicor)=core wave-functions for the current type of atoms.
indlmn_cor(6,nlmn_core)= array giving l,m,n,lm,ln,s for i=lmn, for the core wave functions.
OUTPUT
See side effects
SIDE EFFECTS
Pawtab(ntypat) <type(pawtab_type>=PAW tabulated starting data: %has_nabla=set to 1 in matrix elements are calculated and stored %nabla_ij(3,lmn_size,lmn_size)= <phi_i|nabla|phi_j>-<tphi_i|nabla|tphi_j>
NOTES
MG extracted this piece of code from optics_paw.F90 in order to have something more reusable! Note however the storage mode of nabla_ij differs from optics_paw (here Cartesian coordinates run faster). Besides nabla_ij contains the matrix elements of \nabla instead of the elements of the momentum operator p.
PARENTS
CHILDREN
setnabla_ylm,nderiv_gen,pawrad_deducer0,simp_gen
SOURCE
264 subroutine pawnabla_core_init(mpsang,ntypat,pawrad,pawtab,phi_cor,indlmn_cor) 265 266 267 !This section has been created automatically by the script Abilint (TD). 268 !Do not modify the following lines by hand. 269 #undef ABI_FUNC 270 #define ABI_FUNC 'pawnabla_core_init' 271 !End of the abilint section 272 273 implicit none 274 275 !Arguments ------------------------------------ 276 !scalars 277 integer,intent(in) :: mpsang,ntypat 278 !arrays 279 integer,intent(in) :: indlmn_cor(:,:) 280 real(dp) :: phi_cor(:,:) 281 type(pawtab_type),target,intent(inout) :: pawtab(ntypat) 282 type(pawrad_type),intent(in) :: pawrad(ntypat) 283 284 !Local variables------------------------------- 285 !scalars 286 integer :: nln,nln_cor,il,ilm,ilmn,iln,itypat 287 integer :: jl,jlm,jlmn,jln,lmn_size,lmncmax,mesh_size,mesh_size_cor 288 real(dp) :: intg 289 character(len=500) :: msg 290 !arrays 291 integer, LIBPAW_CONTIGUOUS pointer :: indlmn(:,:) 292 real(dp) :: ang_phipphj(mpsang**2,mpsang**2,8) 293 real(dp),allocatable :: dphi(:),ff(:),int1(:,:),int2(:,:),rad(:) 294 295 ! ************************************************************************* 296 297 mesh_size_cor=size(phi_cor,1) 298 nln_cor=size(phi_cor,2) 299 lmncmax=size(indlmn_cor,2) 300 301 if (mpsang>4)then 302 write(msg,'(3a)')& 303 & 'Not designed for angular momentum greater than 3 ',ch10,& 304 & 'Modification in the table defined in routine setnabla_ylm is required.' 305 MSG_BUG(msg) 306 end if 307 308 if (mesh_size_cor/=pawrad(1)%mesh_size) then 309 write(msg,'(a)') 'Wrong mesh_size_cor value (1)!' 310 MSG_BUG(msg) 311 end if 312 313 if (any(mesh_size_cor/=pawtab(:)%mesh_size)) then 314 write(msg,'(3a)') 'Wrong mesh_size_cor value (2)!',ch10,& 315 & 'Should have only one type of atom.' 316 MSG_ERROR(msg) 317 end if 318 319 !Integration of the angular part: all angular integrals have been computed 320 !outside Abinit and tabulated for each (l,m) value up to l=3 321 call setnabla_ylm(ang_phipphj,mpsang) 322 323 do itypat=1,ntypat 324 325 ! COMPUTE nabla_ij := <phi_i|nabla|phi_cor> for this type 326 mesh_size=pawtab(itypat)%mesh_size 327 lmn_size=pawtab(itypat)%lmn_size 328 nln=pawtab(itypat)%basis_size 329 330 if (allocated(pawtab(itypat)%nabla_ij)) then 331 LIBPAW_DEALLOCATE(pawtab(itypat)%nabla_ij) 332 end if 333 LIBPAW_ALLOCATE(pawtab(itypat)%nabla_ij,(3,lmn_size,lmncmax)) 334 pawtab(itypat)%has_nabla=1 335 336 LIBPAW_ALLOCATE(ff,(mesh_size)) 337 LIBPAW_ALLOCATE(rad,(mesh_size)) 338 LIBPAW_ALLOCATE(int2,(lmn_size,lmncmax)) 339 LIBPAW_ALLOCATE(int1,(lmn_size,lmncmax)) 340 LIBPAW_ALLOCATE(dphi,(mesh_size)) 341 342 indlmn => pawtab(itypat)%indlmn 343 rad(1:mesh_size)=pawrad(itypat)%rad(1:mesh_size) 344 345 ! int1=\int phi phi_core/r dr 346 do jln=1,nln_cor 347 do iln=1,nln 348 ff(2:mesh_size)=(pawtab(itypat)%phi(2:mesh_size,iln)*phi_cor(2:mesh_size,jln))/rad(2:mesh_size) 349 call pawrad_deducer0(ff,mesh_size,pawrad(itypat)) 350 call simp_gen(intg,ff,pawrad(itypat)) 351 int1(iln,jln)=intg 352 end do 353 end do 354 355 ! int2=\int phi/r d/dr(phi_core/r) r^2dr - phi phi_core/r dr 356 do jln=1,nln_cor 357 ff(1:mesh_size)=phi_cor(1:mesh_size,jln) 358 call nderiv_gen(dphi,ff,pawrad(itypat)) 359 do iln=1,nln 360 ff(2:mesh_size)=pawtab(itypat)%phi(2:mesh_size,iln)*dphi(2:mesh_size) & 361 & -pawtab(itypat)%phi(2:mesh_size,iln)*phi_cor(2:mesh_size,jln)/rad(2:mesh_size) 362 call pawrad_deducer0(ff,mesh_size,pawrad(itypat)) 363 call simp_gen(intg,ff,pawrad(itypat)) 364 int2(iln,jln)=intg 365 end do 366 end do 367 368 ! Integration of the radial part, Note unpacked loop 369 do jlmn=1,lmncmax 370 jlm=indlmn_cor(4,jlmn) 371 jl =indlmn_cor(5,jlmn) 372 do ilmn=1,lmn_size 373 ilm=indlmn(4,ilmn) 374 il =indlmn(5,ilmn) 375 376 pawtab(itypat)%nabla_ij(1,ilmn,jlmn)= & 377 & int2(il,jl)* ang_phipphj(ilm,jlm,1) & 378 & +int1(il,jl)*(ang_phipphj(ilm,jlm,2)+ang_phipphj(ilm,jlm,3)) 379 380 pawtab(itypat)%nabla_ij(2,ilmn,jlmn)= & 381 & int2(il,jl)* ang_phipphj(ilm,jlm,4) & 382 & +int1(il,jl)*(ang_phipphj(ilm,jlm,5)+ang_phipphj(ilm,jlm,6)) 383 384 pawtab(itypat)%nabla_ij(3,ilmn,jlmn)= & 385 & int2(il,jl)* ang_phipphj(ilm,jlm,7) & 386 & +int1(il,jl)* ang_phipphj(ilm,jlm,8) 387 388 end do !ilmn 389 end do !jlmn 390 391 pawtab(itypat)%has_nabla=3 392 LIBPAW_DEALLOCATE(ff) 393 LIBPAW_DEALLOCATE(rad) 394 LIBPAW_DEALLOCATE(int2) 395 LIBPAW_DEALLOCATE(int1) 396 LIBPAW_DEALLOCATE(dphi) 397 398 end do !itypat 399 400 end subroutine pawnabla_core_init
m_paw_onsite/pawnabla_init [ Functions ]
[ Top ] [ m_paw_onsite ] [ Functions ]
NAME
pawnabla_init
FUNCTION
Evaluate all valence-valence onsite contributions of the nabla operator in cartesian coordinates. Store values in the pawtab% data structure.
INPUTS
mpsang=1+maximum angular momentum
ntypat=Number of types of atoms in cell
Pawrad(ntypat)<Pawrad_type>=PAW radial mesh and related data:
%rad(mesh_size)=The coordinates of all the points of the radial mesh
Pawtab(ntypat) <type(pawtab_type>=PAW tabulated starting data:
%mesh_size=Dimension of radial mesh
%lmn_size=Number of (l,m,n) elements for the PAW basis
OUTPUT
See side effects
SIDE EFFECTS
Pawtab(ntypat) <type(pawtab_type>=PAW tabulated starting data: %has_nabla=set to 1 in matrix elements are calculated and stored %nabla_ij(3,lmn_size,lmn_size)= <phi_i|nabla|phi_j>-<tphi_i|nabla|tphi_j>
NOTES
MG extracted this piece of code from optics_paw.F90 in order to have something more reusable! Note however the storage mode of nabla_ij differs from optics_paw (here Cartesian coordinates run faster). Besides nabla_ij contains the matrix elements of \nabla instead of the elements of the momentum operator p.
PARENTS
CHILDREN
setnabla_ylm,nderiv_gen,pawrad_deducer0,simp_gen
SOURCE
89 subroutine pawnabla_init(mpsang,ntypat,pawrad,pawtab) 90 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 'pawnabla_init' 96 !End of the abilint section 97 98 implicit none 99 100 !Arguments ------------------------------------ 101 !scalars 102 integer,intent(in) :: mpsang,ntypat 103 !arrays 104 type(pawtab_type),target,intent(inout) :: pawtab(ntypat) 105 type(pawrad_type),intent(in) :: pawrad(ntypat) 106 107 !Local variables------------------------------- 108 !scalars 109 integer :: nln,il,ilm,ilmn,iln,itypat 110 integer :: jl,jlm,jlmn,jln,lmn_size,mesh_size 111 real(dp) :: intg 112 character(len=500) :: msg 113 !arrays 114 integer, LIBPAW_CONTIGUOUS pointer :: indlmn(:,:) 115 real(dp) :: ang_phipphj(mpsang**2,mpsang**2,8) 116 real(dp),allocatable :: dphi(:),dtphi(:),ff(:),int1(:,:),int2(:,:),rad(:) 117 118 ! ************************************************************************* 119 120 if (mpsang>4)then 121 write(msg,'(3a)')& 122 & 'Not designed for angular momentum greater than 3 ',ch10,& 123 & 'Modification in the table defined in routine setnabla_ylm is required.' 124 MSG_BUG(msg) 125 end if 126 127 !Integration of the angular part: all angular integrals have been computed 128 !outside Abinit and tabulated for each (l,m) value up to l=3 129 call setnabla_ylm(ang_phipphj,mpsang) 130 131 do itypat=1,ntypat 132 133 ! COMPUTE nabla_ij := <phi_i|nabla|phi_j>-<tphi_i|nabla|tphi_j> for this type 134 mesh_size=pawtab(itypat)%mesh_size 135 lmn_size=pawtab(itypat)%lmn_size 136 nln=pawtab(itypat)%basis_size 137 138 if (allocated(pawtab(itypat)%nabla_ij)) then 139 LIBPAW_DEALLOCATE(pawtab(itypat)%nabla_ij) 140 end if 141 LIBPAW_ALLOCATE(pawtab(itypat)%nabla_ij,(3,lmn_size,lmn_size)) 142 pawtab(itypat)%has_nabla=1 143 144 LIBPAW_ALLOCATE(ff,(mesh_size)) 145 LIBPAW_ALLOCATE(rad,(mesh_size)) 146 LIBPAW_ALLOCATE(int2,(lmn_size,lmn_size)) 147 LIBPAW_ALLOCATE(int1,(lmn_size,lmn_size)) 148 LIBPAW_ALLOCATE(dphi,(mesh_size)) 149 LIBPAW_ALLOCATE(dtphi,(mesh_size)) 150 151 indlmn => pawtab(itypat)%indlmn 152 rad(1:mesh_size)=pawrad(itypat)%rad(1:mesh_size) 153 154 ! int1=\int phi phj/r dr - \int tphi tphj /r dr 155 do jln=1,nln 156 do iln=1,nln 157 ff(2:mesh_size)= ( & 158 & pawtab(itypat)%phi (2:mesh_size,iln)*pawtab(itypat)%phi (2:mesh_size,jln) & 159 & -pawtab(itypat)%tphi(2:mesh_size,iln)*pawtab(itypat)%tphi(2:mesh_size,jln) ) /rad(2:mesh_size) 160 call pawrad_deducer0(ff,mesh_size,pawrad(itypat)) 161 call simp_gen(intg,ff,pawrad(itypat)) 162 int1(iln,jln)=intg 163 end do 164 end do 165 166 ! int2=\int phi/r d/dr(phj/r) r^2dr - \int tphi/r d/dr(tphj/r)r^2 dr 167 do jln=1,nln 168 ff(1:mesh_size)=pawtab(itypat)%phi(1:mesh_size,jln) 169 call nderiv_gen(dphi,ff,pawrad(itypat)) 170 ff(1:mesh_size)=pawtab(itypat)%tphi(1:mesh_size,jln) 171 call nderiv_gen(dtphi,ff,pawrad(itypat)) 172 173 do iln=1,nln 174 ff(2:mesh_size)= & 175 & pawtab(itypat)%phi (2:mesh_size,iln)*dphi (2:mesh_size) & 176 & -pawtab(itypat)%phi (2:mesh_size,iln)*pawtab(itypat)%phi (2:mesh_size,jln)/rad(2:mesh_size) & 177 & -( pawtab(itypat)%tphi(2:mesh_size,iln)*dtphi(2:mesh_size) & 178 & -pawtab(itypat)%tphi(2:mesh_size,iln)*pawtab(itypat)%tphi(2:mesh_size,jln)/rad(2:mesh_size) ) 179 call pawrad_deducer0(ff,mesh_size,pawrad(itypat)) 180 call simp_gen(intg,ff,pawrad(itypat)) 181 int2(iln,jln)=intg 182 end do 183 end do 184 185 ! Integration of the radial part, Note unpacked loop 186 do jlmn=1,lmn_size 187 jlm=indlmn(4,jlmn) 188 jl =indlmn(5,jlmn) 189 do ilmn=1,lmn_size 190 ilm=indlmn(4,ilmn) 191 il =indlmn(5,ilmn) 192 193 pawtab(itypat)%nabla_ij(1,ilmn,jlmn)= & 194 & int2(il,jl)* ang_phipphj(ilm,jlm,1) & 195 & +int1(il,jl)*(ang_phipphj(ilm,jlm,2)+ang_phipphj(ilm,jlm,3)) 196 197 pawtab(itypat)%nabla_ij(2,ilmn,jlmn)= & 198 & int2(il,jl)* ang_phipphj(ilm,jlm,4) & 199 & +int1(il,jl)*(ang_phipphj(ilm,jlm,5)+ang_phipphj(ilm,jlm,6)) 200 201 pawtab(itypat)%nabla_ij(3,ilmn,jlmn)= & 202 & int2(il,jl)* ang_phipphj(ilm,jlm,7) & 203 & +int1(il,jl)* ang_phipphj(ilm,jlm,8) 204 205 end do !ilmn 206 end do !jlmn 207 208 pawtab(itypat)%has_nabla=2 209 LIBPAW_DEALLOCATE(ff) 210 LIBPAW_DEALLOCATE(rad) 211 LIBPAW_DEALLOCATE(int2) 212 LIBPAW_DEALLOCATE(int1) 213 LIBPAW_DEALLOCATE(dphi) 214 LIBPAW_DEALLOCATE(dtphi) 215 216 end do !itypat 217 218 end subroutine pawnabla_init