TABLE OF CONTENTS
ABINIT/psp2nl [ Functions ]
NAME
psp2nl
FUNCTION
Goedecker-Teter-Hutter nonlocal pseudopotential (from preprint of 1996). Uses Gaussians for fully nonlocal form, analytic expressions.
COPYRIGHT
Copyright (C) 1998-2017 ABINIT group (DCA, XG, GMR, 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 . For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
INPUTS
h1p=factor defining strength of 1st projector for l=1 channel h1s=factor defining strength of 1st projector for l=0 channel h2s=factor defining strength of 2nd projector for l=0 channel lnmax=max. number of (l,n) components over all type of psps mqgrid=number of grid points for qgrid qgrid(mqgrid)=array of |G| values rrp=core radius for p channel (bohr) rrs=core radius for s channel (bohr)
OUTPUT
ekb(lnmax)=Kleinman-Bylander energy ffspl(mqgrid,2,lnmax)=Kleinman-Bylander form factor f_l(q) and second derivative from spline fit for each angular momentum and each projector
PARENTS
psp2in
CHILDREN
spline
SOURCE
41 #if defined HAVE_CONFIG_H 42 #include "config.h" 43 #endif 44 45 #include "abi_common.h" 46 47 subroutine psp2nl(ekb,ffspl,h1p,h1s,h2s,lnmax,mqgrid,qgrid,rrp,rrs) 48 49 use defs_basis 50 use m_splines 51 use m_profiling_abi 52 53 !This section has been created automatically by the script Abilint (TD). 54 !Do not modify the following lines by hand. 55 #undef ABI_FUNC 56 #define ABI_FUNC 'psp2nl' 57 !End of the abilint section 58 59 implicit none 60 61 !Arguments ------------------------------------ 62 !scalars 63 integer,intent(in) :: lnmax,mqgrid 64 real(dp),intent(in) :: h1p,h1s,h2s,rrp,rrs 65 !arrays 66 real(dp),intent(in) :: qgrid(mqgrid) 67 real(dp),intent(inout) :: ekb(lnmax),ffspl(mqgrid,2,lnmax) !vz_i 68 69 !Local variables------------------------------- 70 !scalars 71 integer :: iln,iqgrid 72 real(dp) :: qmax,yp1,ypn 73 !arrays 74 real(dp),allocatable :: work(:) 75 76 ! ************************************************************************* 77 78 ABI_ALLOCATE(work,(mqgrid)) 79 80 !Kleinman-Bylander energies ekb were set to zero in calling program 81 82 !Compute KB energies 83 iln=0 84 if (abs(h1s)>1.d-12) then 85 iln=iln+1 86 ekb(iln)=h1s*32.d0*rrs**3*(pi**(2.5d0)/(4.d0*pi)**2) 87 end if 88 if (abs(h2s)>1.d-12) then 89 iln=iln+1 90 ekb(iln) =h2s*(128.d0/15.d0)*rrs**3*(pi**(2.5d0)/(4.d0*pi)**2) 91 end if 92 if (abs(h1p)>1.d-12) then 93 iln=iln+1 94 ekb(iln)=h1p*(64.d0/3.d0)*rrp**5*(pi**(2.5d0)/(4.d0*pi)**2) 95 end if 96 97 !Compute KB form factor 98 iln=0 99 100 !l=0 first projector 101 if (abs(h1s)>1.d-12) then 102 iln=iln+1 103 do iqgrid=1,mqgrid 104 ffspl(iqgrid,1,iln)=exp(-0.5d0*(two_pi*qgrid(iqgrid)*rrs)**2) 105 end do 106 ! Compute yp1,ypn=derivatives of f(q) at q=0, q=qgrid(mqgrid) 107 yp1=0.d0 108 qmax=qgrid(mqgrid) 109 ypn=-4.d0*pi**2*qmax*rrs**2*exp(-0.5d0*(two_pi*qmax*rrs)**2) 110 ! Fit spline to get second derivatives by spline fit 111 call spline(qgrid,ffspl(:,1,iln),mqgrid,yp1,ypn,ffspl(:,2,iln)) 112 ! else 113 ! or else put first projector nonlocal correction at l=0 to 0 114 ! ffspl(:,:,iln)=0.0d0 115 end if 116 117 !l=0 second projector 118 if (abs(h2s)>1.d-12) then 119 iln=iln+1 120 do iqgrid=1,mqgrid 121 ffspl(iqgrid,1,iln)=exp(-0.5d0*(two_pi*qgrid(iqgrid)*rrs)**2) * & 122 & (3.d0-(two_pi*qgrid(iqgrid)*rrs)**2) 123 end do 124 ! Compute yp1,ypn=derivatives of f(q) at q=0, q=qgrid(mqgrid) 125 yp1=0.d0 126 qmax=qgrid(mqgrid) 127 ypn=4.d0*pi**2*qmax*rrs**2*exp(-0.5d0*(two_pi*qmax*rrs)**2) * & 128 & (-5.d0+(two_pi*qmax*rrs)**2) 129 ! Fit spline to get second derivatives by spline fit 130 call spline(qgrid,ffspl(:,1,iln),mqgrid,yp1,ypn,ffspl(:,2,iln)) 131 ! else if(mproj>=2)then 132 ! or else put second projector nonlocal correction at l=0 to 0 133 ! ffspl(:,:,iln)=0.0d0 134 end if 135 136 !l=1 first projector 137 if (abs(h1p)>1.d-12) then 138 iln=iln+1 139 do iqgrid=1,mqgrid 140 ffspl(iqgrid,1,iln)=exp(-0.5d0*(two_pi*qgrid(iqgrid)*rrp)**2) * & 141 & (two_pi*qgrid(iqgrid)) 142 end do 143 ! Compute yp1,ypn=derivatives of f(q) at q=0, q=qgrid(mqgrid) 144 yp1=two_pi 145 qmax=qgrid(mqgrid) 146 ypn=-two_pi*((two_pi*qmax*rrp)**2-1.d0) * exp(-0.5d0*(two_pi*qmax*rrp)**2) 147 ! Fit spline to get second derivatives by spline fit 148 call spline(qgrid,ffspl(:,1,iln),mqgrid,yp1,ypn,ffspl(:,2,iln)) 149 ! else if(mpsang>=2)then 150 ! or else put first projector l=1 nonlocal correction to 0 151 ! ffspl(:,:,iln)=0.0d0 152 end if 153 154 ABI_DEALLOCATE(work) 155 156 end subroutine psp2nl