scnc
/
quantum-espresso
mirror of https://gitlab.com/QEF/q-e.git
1
0
Fork 0
Code Issues Packages Projects Releases Wiki Activity

The following pseudopotential-related variables in module uspp_param: 

zp, psd, dion, betar, jjj, qqq, qfunc, qfcoef, vloc_at, rinner,
   nbeta, kkbeta, nqf, nqlc, lll, tvanp
have been replaced by the corresponding variables in structure 'upf'.
There shouldn't be any side effects, but who knows. There is still a
copy of the above variables that will be removed sooner or later.
Basically : variable([i,j,k,..,]n) => upf(n)%variable [(i,j,k,..)]
Note that upf%qfunc has for the time being three indices instead of two,
and that upf%kkbeta is the analogous of kkbeta and not what it used to be.
The logic of this operation will be clearer when it will be completed


git-svn-id: http://qeforge.qe-forge.org/svn/q-e/trunk/espresso@4308 c92efa57-630b-4861-b058-cf58834340f0
This commit is contained in:
giannozz 2007-10-05 09:26:23 +00:00
parent d9af5862d6
commit 216c32ccd3
44 changed files with 520 additions and 528 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

1
CPV/make.depend
View File

@ -929,7 +929,6 @@ pseudopot_sub.o : ../Modules/kind.o
pseudopot_sub.o : ../Modules/mp.o
pseudopot_sub.o : ../Modules/mp_global.o
pseudopot_sub.o : ../Modules/parameters.o
pseudopot_sub.o : ../Modules/pseudo_types.o
pseudopot_sub.o : ../Modules/recvec.o
pseudopot_sub.o : ../Modules/smallbox.o
pseudopot_sub.o : ../Modules/uspp.o

180
CPV/pseudopot_sub.f90
View File

@ -19,7 +19,7 @@
USE kinds, ONLY: DP
use uspp, only: dvan, nhtolm, indv
use uspp_param, only: nhm, nh, dion
use uspp_param, only: upf, nhm, nh
use ions_base, only: nsp
use atom, only: numeric
!
@ -42,7 +42,7 @@
do iv=1,nh(is)
do jv=1,nh(is)
if ( nhtolm(iv,is) == nhtolm(jv,is) ) then
dvan( iv, jv, is ) = fac * dion( indv(iv,is), indv(jv,is), is )
dvan( iv, jv, is ) = fac * upf(is)%dion( indv(iv,is), indv(jv,is) )
endif
end do
end do
@ -64,14 +64,11 @@
use cvan, only: ish !
use uspp, only: nkb, & !
nkbus !
use uspp_param, only: nbeta, &!
use uspp_param, only: upf, &!
lmaxkb, &!
lll, &!
nhm, &!
nbetam, &!
nh, &!
tvanp, &!
nqlc, &!
lmaxq !
use uspp, only: nhtol, &!
nhtolm, &!
@ -97,18 +94,18 @@
!
do is = 1, nsp
ind = 0
do iv = 1, nbeta(is)
lmaxkb = max( lmaxkb, lll( iv, is ) )
ind = ind + 2 * lll( iv, is ) + 1
do iv = 1, upf(is)%nbeta
lmaxkb = max( lmaxkb, upf(is)%lll( iv ) )
ind = ind + 2 * upf(is)%lll( iv ) + 1
end do
nh(is) = ind
ish(is)=nkb
nkb = nkb + na(is) * nh(is)
if( tvanp(is) ) nkbus = nkbus + na(is) * nh(is)
if( upf(is)%tvanp ) nkbus = nkbus + na(is) * nh(is)
nlcc_any = nlcc_any .OR. nlcc(is)
end do
nhm = MAXVAL( nh(1:nsp) )
nbetam = MAXVAL(nbeta(1:nsp))
nbetam = MAXVAL( upf(1:nsp)%nbeta )
if (lmaxkb > lmaxx) call errore(' pseudopotential_indexes ',' l > lmax ',lmaxkb)
lmaxq = 2*lmaxkb + 1
!
@ -117,7 +114,7 @@
! l of the beta functions but includes the l of the local potential
!
do is=1,nsp
nqlc(is) = MIN ( nqlc(is), lmaxq )
upf(is)%nqlc = MIN ( upf(is)%nqlc, lmaxq )
end do
if (nkb <= 0) call errore(' pseudopotential_indexes ',' not implemented ?',nkb)
@ -135,13 +132,13 @@
!
do is = 1, nsp
ind = 0
do iv = 1, nbeta(is)
lm = lll(iv,is)**2
do il = 1, 2*lll( iv, is ) + 1
do iv = 1, upf(is)%nbeta
lm = upf(is)%lll(iv)**2
do il = 1, 2* upf(is)%lll( iv ) + 1
lm = lm + 1
ind = ind + 1
nhtolm( ind, is ) = lm
nhtol( ind, is ) = lll( iv, is )
nhtol( ind, is ) = upf(is)%lll( iv )
indv( ind, is ) = iv
end do
end do
@ -240,7 +237,7 @@
use bhs, ONLY : rc1, rc2, wrc2, wrc1, rcl, al, bl, lloc
USE splines, ONLY : init_spline, allocate_spline, kill_spline, nullify_spline
USE pseudo_base, ONLY : formfn, formfa
USE uspp_param, only : vloc_at, oldvan
USE uspp_param, only : upf, oldvan
USE control_flags, only : tpre
use reciprocal_vectors, ONLY : g, gstart
USE cp_interfaces, ONLY : compute_xgtab, chkpstab
@ -284,17 +281,19 @@
if ( numeric(is) ) then
call formfn( vps_sp(is)%y, dvps_sp(is)%y, rgrid(is)%r, rgrid(is)%rab, vloc_at(:,is), &
zv(is), rcmax(is), xgtab, 1.0d0, tpiba2, rgrid(is)%mesh, &
call formfn( vps_sp(is)%y, dvps_sp(is)%y, rgrid(is)%r, &
rgrid(is)%rab, upf(is)%vloc(1), zv(is), &
rcmax(is), xgtab, 1.0d0, tpiba2, rgrid(is)%mesh, &
mmx, oldvan(is), tpre )
else
! bhs pseudopotentials
!
call formfa( vps_sp(is)%y, dvps_sp(is)%y, rc1(is), rc2(is), wrc1(is), wrc2(is), &
rcl(:,is,lloc(is)), al(:,is,lloc(is)), bl(:,is,lloc(is)), &
zv(is), rcmax(is), xgtab, 1.0d0, tpiba2, mmx, 2 , tpre )
call formfa( vps_sp(is)%y, dvps_sp(is)%y, rc1(is), rc2(is), &
wrc1(is), wrc2(is), rcl(:,is,lloc(is)), &
al(:,is,lloc(is)), bl(:,is,lloc(is)), zv(is), &
rcmax(is), xgtab, 1.0d0, tpiba2, mmx, 2 , tpre )
end if
@ -389,9 +388,8 @@
USE splines, ONLY : init_spline, allocate_spline, kill_spline, nullify_spline
USE pseudo_base, ONLY : nlin_base
USE pseudo_base, ONLY : nlin_stress_base
USE pseudo_types, ONLY : pseudo_ncpp, pseudo_upf
USE reciprocal_vectors, ONLY : g, gstart
USE uspp_param, ONLY : nbeta, nbetam
USE uspp_param, ONLY : upf, nbetam
USE read_pseudo_module_fpmd, ONLY : ap, nspnl
USE cp_interfaces, ONLY : compute_xgtab, chkpstab
USE pseudopotential, ONLY : wnl_sp, wnla_sp, xgtab
@ -436,7 +434,7 @@
CALL nullify_spline( wnla_sp( l, is ) )
END DO
DO l = 1, nbeta( is )
DO l = 1, upf(is)%nbeta
CALL allocate_spline( wnl_sp(l,is), mmx, xgmin, xgmax )
CALL allocate_spline( wnla_sp(l,is), mmx, xgmin, xgmax )
END DO
@ -445,17 +443,17 @@
!
CALL nlin_base(ap(is), xgtab(:), fintl)
!
DO l = 1, nbeta( is )
DO l = 1, upf(is)%nbeta
wnl_sp( l, is )%y = fintl(:,l)
END DO
!
CALL nlin_stress_base( ap(is), xgtab, fintl )
DO l = 1, nbeta( is )
DO l = 1, upf(is)%nbeta
wnla_sp( l, is )%y = fintl(:,l)
END DO
!
DO l = 1, nbeta( is )
DO l = 1, upf(is)%nbeta
CALL init_spline( wnl_sp( l, is ) )
CALL init_spline( wnla_sp( l, is ) )
END DO
@ -483,7 +481,7 @@
USE pseudo_base, ONLY: nlin_base
USE control_flags, ONLY: gamma_only
USE uspp, ONLY: dvan
USE uspp_param, ONLY: nh, nbeta
USE uspp_param, ONLY: upf, nh
USE constants, ONLY: pi
USE splines, ONLY: spline
USE read_pseudo_module_fpmd, ONLY: ap, nspnl
@ -523,7 +521,7 @@
!
IF( tpstab ) THEN
!
DO l = 1, nbeta( is )
DO l = 1, upf(is)%nbeta
!
IF( gstart == 2 ) THEN
wnl(1,l,is,ik) = wnl_sp( l, is )%y(1)
@ -563,7 +561,7 @@
USE cell_base, ONLY : tpiba
USE pseudo_base, ONLY : nlin_stress_base
USE splines, ONLY : spline
USE uspp_param, ONLY : nbeta
USE uspp_param, ONLY : upf
USE read_pseudo_module_fpmd, ONLY : ap, nspnl
USE reciprocal_vectors, ONLY : g, gstart
USE pseudopotential, ONLY : wnla_sp, tpstab
@ -583,7 +581,7 @@
!
IF ( tpstab ) THEN
!
DO l = 1, nbeta( is )
DO l = 1, upf(is)%nbeta
!
IF( gstart == 2 ) THEN
wnla(1,l,is) = wnla_sp( l, is )%y(1)
@ -617,7 +615,7 @@
!
USE kinds, ONLY : DP
USE ions_base, ONLY : nsp
USE uspp_param, ONLY : nh, kkbeta, betar, nhm, nbeta, oldvan
USE uspp_param, ONLY : upf, nh, nhm, oldvan
USE atom, ONLY : rgrid, numeric
USE uspp, ONLY : nhtol, indv
USE betax, only : refg, betagx, mmx, dbetagx
@ -635,11 +633,10 @@
!
ALLOCATE( betagx ( mmx, nhm, nsp ) )
IF ( tpre ) ALLOCATE( dbetagx( mmx, nhm, nsp ) )
!
do is = 1, nsp
!
nr = kkbeta( is )
nr = upf(is)%kkbeta
!
if ( tpre ) then
allocate( dfint( nr ) )
@ -667,7 +664,8 @@
! beta(ir)=r*beta(r)
!
do ir = 1, nr
fint(ir) = rgrid(is)%r(ir) * betar( ir, indv(iv,is), is ) * jl(ir)
fint(ir) = rgrid(is)%r(ir) * jl(ir) * &
upf(is)%beta( ir, indv(iv,is) )
end do
if (oldvan(is)) then
call herman_skillman_int(nr,fint,rgrid(is)%rab,betagx(il,iv,is))
@ -677,7 +675,8 @@
!
if(tpre) then
do ir = 1, nr
dfint(ir) = rgrid(is)%r(ir) * betar( ir, indv(iv,is), is ) * djl(ir)
dfint(ir) = rgrid(is)%r(ir) * djl(ir) * &
upf(is)%beta( ir, indv(iv,is) )
end do
if (oldvan(is)) then
call herman_skillman_int(nr,dfint,rgrid(is)%rab,dbetagx(il,iv,is))
@ -717,8 +716,7 @@
USE kinds, ONLY : DP
use io_global, only : stdout
USE ions_base, ONLY : nsp
USE uspp_param, ONLY : nh, kkbeta, betar, nhm, nbetam, nqlc, qqq, &
lmaxq, nbeta, oldvan
USE uspp_param, ONLY : upf, nh, nhm, nbetam, lmaxq, oldvan
USE atom, ONLY : rgrid, numeric
USE uspp, ONLY : indv
USE betax, only : refg, qradx, mmx, dqradx
@ -748,9 +746,9 @@
! as vanderbilts ppot-code prints them out
!
WRITE( stdout,*) ' nlinit nh(is), ngb, is, kkbeta, lmaxq = ', &
& nh(is), ngb, is, kkbeta(is), nqlc(is)
& nh(is), ngb, is, upf(is)%kkbeta, upf(is)%nqlc
!
nr = kkbeta(is)
nr = upf(is)%kkbeta
!
IF ( tpre ) THEN
ALLOCATE( djl ( nr ) )
@ -759,11 +757,11 @@
!
ALLOCATE( fint( nr ) )
ALLOCATE( jl ( nr ) )
ALLOCATE( qrl( nr, nbeta(is)*(nbeta(is)+1)/2, nqlc(is)) )
ALLOCATE( qrl( nr, upf(is)%nbeta*(upf(is)%nbeta+1)/2, upf(is)%nqlc) )
!
call fill_qrl ( is, qrl )
!
do l = 1, nqlc( is )
do l = 1, upf(is)%nqlc
!
do il = 1, mmx
!
@ -778,8 +776,8 @@
endif
!
!
do iv = 1, nbeta(is)
do jv = iv, nbeta(is)
do iv = 1, upf(is)%nbeta
do jv = iv, upf(is)%nbeta
ijv = jv * ( jv - 1 ) / 2 + iv
!
! note qrl(r)=r^2*q(r)
@ -827,8 +825,8 @@
WRITE( stdout,*)
WRITE( stdout,'(20x,a)') ' qqq '
!
do iv=1,nbeta(is)
WRITE( stdout,'(8f9.4)') (qqq(iv,jv,is),jv=1,nbeta(is))
do iv=1,upf(is)%nbeta
WRITE( stdout,'(8f9.4)') (upf(is)%qqq(iv,jv),jv=1,upf(is)%nbeta)
end do
WRITE( stdout,*)
!
@ -844,8 +842,7 @@
USE kinds, ONLY : DP
use io_global, only: stdout
USE ions_base, ONLY: nsp
USE uspp_param, ONLY: nh, kkbeta, betar, nhm, nbetam, nqlc, qqq, &
lmaxq, nbeta, oldvan
USE uspp_param, ONLY: upf, nh, nhm, nbetam, lmaxq, oldvan
use uspp_param, only: lmaxkb
USE atom, ONLY: rgrid, numeric
USE uspp, ONLY: indv
@ -889,9 +886,9 @@
! as vanderbilts ppot-code prints them out
!
WRITE( stdout,*) ' nlinit nh(is), ngb, is, kkbeta, lmaxq = ', &
& nh(is), ngb, is, kkbeta(is), nqlc(is)
& nh(is), ngb, is, upf(is)%kkbeta, upf(is)%nqlc
!
nr = kkbeta(is)
nr = upf(is)%kkbeta
!
IF ( tpre ) THEN
ALLOCATE( djl ( nr ) )
@ -900,12 +897,12 @@
!
ALLOCATE( fint( nr ) )
ALLOCATE( jl ( nr ) )
ALLOCATE( qrl( nr, nbeta(is)*(nbeta(is)+1)/2, nqlc(is)) )
ALLOCATE( qrl( nr, upf(is)%nbeta*(upf(is)%nbeta+1)/2, upf(is)%nqlc) )
!
call fill_qrl ( is, qrl )
! qrl = 0.0d0
!
do l = 1, nqlc( is )
do l = 1, upf(is)%nqlc
!
do il = 1, ngb
!
@ -920,8 +917,8 @@
endif
!
!
do iv = 1, nbeta(is)
do jv = iv, nbeta(is)
do iv = 1, upf(is)%nbeta
do jv = iv, upf(is)%nbeta
ijv = jv * ( jv - 1 ) / 2 + iv
!
! note qrl(r)=r^2*q(r)
@ -969,8 +966,8 @@
WRITE( stdout,*)
WRITE( stdout,'(20x,a)') ' qqq '
!
do iv=1,nbeta(is)
WRITE( stdout,'(8f9.4)') (qqq(iv,jv,is),jv=1,nbeta(is))
do iv=1, upf(is)%nbeta
WRITE( stdout,'(8f9.4)') (upf(is)%qqq(iv,jv),jv=1, upf(is)%nbeta)
end do
WRITE( stdout,*)
!
@ -991,11 +988,11 @@
!
c = fpi / omegab
!
do iv= 1,nbeta(is)
do jv = iv, nbeta(is)
do iv= 1, upf(is)%nbeta
do jv = iv, upf(is)%nbeta
ijv = jv*(jv-1)/2 + iv
do ig=1,ngb
do l=1,nqlc(is)
do l=1,upf(is)%nqlc
qradb(ig,ijv,l,is)= c*qradx(ig,ijv,l,is)
enddo
enddo
@ -1035,10 +1032,10 @@
!
do is=1,nvb
!
do iv= 1,nbeta(is)
do jv=iv,nbeta(is)
do iv= 1, upf(is)%nbeta
do jv=iv, upf(is)%nbeta
ijv = jv*(jv-1)/2 + iv
do l=1,nqlc(is)
do l=1,upf(is)%nqlc
do ig=1,ngb
dqradb(ig,ijv,l,is) = dqradx(ig,ijv,l,is)
enddo
@ -1146,7 +1143,7 @@
USE gvecw, only: ngw
USE ions_base, only: nsp
USE reciprocal_vectors, only: g, gx, gstart
USE uspp_param, only: lmaxq, nqlc, lmaxkb, kkbeta, nbeta, nh
USE uspp_param, only: upf, lmaxq, lmaxkb, nh
USE uspp, only: qq, nhtolm, beta
USE cell_base, only: ainv, omega, tpiba2, tpiba
USE betax, ONLY : refg, betagx, dbetagx
@ -1251,7 +1248,7 @@
use uspp, only: qq, nhtolm, beta
use constants, only: pi, fpi
use ions_base, only: nsp
use uspp_param, only: lmaxq, nqlc, lmaxkb, kkbeta, nbeta, nbetam, nh
use uspp_param, only: upf, lmaxq, lmaxkb, nbetam, nh
use qradb_mod, only: qradb
use qgb_mod, only: qgb
use gvecb, only: gb, gxb, ngb
@ -1286,16 +1283,16 @@
!
c = fpi / omegab
!
do iv= 1,nbeta(is)
do jv = iv, nbeta(is)
do iv= 1, upf(is)%nbeta
do jv = iv, upf(is)%nbeta
ijv = jv*(jv-1)/2 + iv
do l=1,nqlc(is)
do l=1, upf(is)%nqlc
qradb(1,ijv,l,is) = c * qradx(1,ijv,l,is)
end do
do ig=2,ngb
gg=gb(ig)*tpibab*tpibab/refg
jj=int(gg)+1
do l=1,nqlc(is)
do l=1,upf(is)%nqlc
if(jj.ge.mmx) then
qradb(ig,ijv,l,is)=0.d0
else
@ -1340,10 +1337,10 @@
!
do is=1,nvb
!
do iv= 1,nbeta(is)
do jv=iv,nbeta(is)
do iv= 1, upf(is)%nbeta
do jv=iv, upf(is)%nbeta
ijv = jv*(jv-1)/2 + iv
do l=1,nqlc(is)
do l=1,upf(is)%nqlc
dqradb(1,ijv,l,is) = dqradx(1,ijv,l,is)
do ig=2,ngb
gg=gb(ig)*tpibab*tpibab/refg
@ -1415,8 +1412,7 @@
USE io_global, only : stdout
USE gvecw, only : ngw
USE ions_base, only : nsp
USE uspp_param, only : lmaxq, nqlc, lmaxkb, kkbeta, nbeta, nh, &
betar, nhm, oldvan
USE uspp_param, only : upf, lmaxq, lmaxkb, nh, nhm, oldvan
USE uspp, only : qq, nhtolm, beta, nhtol, indv
USE cell_base, only : ainv, omega, tpiba2, tpiba
USE cdvan, ONLY : dbeta
@ -1444,7 +1440,7 @@
!
do is = 1, nsp
!
nr = kkbeta(is)
nr = upf(is)%kkbeta
!
if ( tpre ) then
allocate( dfint( nr ) )
@ -1472,7 +1468,8 @@
! beta(ir)=r*beta(r)
!
do ir = 1, nr
fint(ir) = rgrid(is)%r(ir) * betar( ir, indv(iv,is), is ) * jl(ir)
fint(ir) = rgrid(is)%r(ir) * jl(ir) * &
upf(is)%beta( ir, indv(iv,is) )
end do
if (oldvan(is)) then
call herman_skillman_int(nr,fint,rgrid(is)%rab,betagx(il,iv,is))
@ -1482,7 +1479,8 @@
!
if(tpre) then
do ir = 1, nr
dfint(ir) = rgrid(is)%r(ir) * betar( ir, indv(iv,is), is ) * djl(ir)
dfint(ir) = rgrid(is)%r(ir) * djl(ir) * &
upf(is)%beta( ir, indv(iv,is) )
end do
if (oldvan(is)) then
call herman_skillman_int(nr,dfint,rgrid(ir)%rab,dbetagx(il,iv,is))
@ -1576,7 +1574,7 @@
!
! fill l-components of Q(r) as in Vanderbilt's approach
!
USE uspp_param, ONLY: qfunc, nqf, qfcoef, rinner, lll, nbeta, kkbeta
USE uspp_param, ONLY: upf
USE atom, ONLY: rgrid
USE kinds, ONLY: DP
USE io_global, ONLY: stdout
@ -1593,14 +1591,14 @@
dim2 = SIZE( qrl, 2 )
dim3 = SIZE( qrl, 3 )
!
IF ( kkbeta(is) > dim1 ) &
IF ( upf(is)%kkbeta > dim1 ) &
CALL errore ('fill_qrl', 'bad 1st dimension for array qrl', 1)
!
qrl = 0.0d0
!
do iv = 1, nbeta(is)
do iv = 1, upf(is)%nbeta
!
do jv = iv, nbeta(is)
do jv = iv, upf(is)%nbeta
!
ijv = (jv-1)*jv/2 + iv
!
@ -1609,11 +1607,11 @@
! notice that L runs from 1 to Lmax+1
lmin = ABS (lll(jv,is) - lll(iv,is)) + 1
lmax = lll(jv,is) + lll(iv,is) + 1
lmin = ABS (upf(is)%lll(jv) - upf(is)%lll(iv)) + 1
lmax = upf(is)%lll(jv) + upf(is)%lll(iv) + 1
! WRITE( stdout, * ) 'QRL is, jv, iv = ', is, jv, iv
! WRITE( stdout, * ) 'QRL lll jv, iv = ', lll(jv,is), lll(iv,is)
! WRITE( stdout, * ) 'QRL lll jv, iv = ', upf(is)%lll(jv), upf(is)%lll(iv)
! WRITE( stdout, * ) 'QRL lmin, lmax = ', lmin, lmax
! WRITE( stdout, * ) '---------------- '
@ -1622,15 +1620,17 @@
CALL errore ('fill_qrl', 'bad 3rd dimension for array qrl', 3)
END IF
do l = lmin, lmax
do ir = 1, kkbeta(is)
if ( rgrid(is)%r(ir) >= rinner(l,is) ) then
qrl(ir,ijv,l)=qfunc(ir,ijv,is)
do ir = 1, upf(is)%kkbeta
if ( rgrid(is)%r(ir) >= upf(is)%rinner(l) ) then
! qrl(ir,ijv,l)=upf(is)%qfunc(ir,ijv) TEMP
qrl(ir,ijv,l)=upf(is)%qfunc(ir,iv,jv)
else
qrl(ir,ijv,l)=qfcoef(1,l,iv,jv,is)
do i = 2, nqf(is)
qrl(ir,ijv,l)=upf(is)%qfcoef(1,l,iv,jv)
do i = 2, upf(is)%nqf
qrl(ir,ijv,l)=qrl(ir,ijv,l) + &
qfcoef(i,l,iv,jv,is) * rgrid(is)%r(ir)**(2*i-2)
upf(is)%qfcoef(i,l,iv,jv)*rgrid(is)%r(ir)**(2*i-2)
end do
qrl(ir,ijv,l) = qrl(ir,ijv,l) * rgrid(is)%r(ir)**(l+1)
end if

16
CPV/read_pseudo.f90
View File

@ -136,13 +136,12 @@ SUBROUTINE check_types_order( )
USE ions_base, ONLY: nsp
IMPLICIT NONE
INTEGER :: is, il
LOGICAL :: tvanp
!
! With Vanderbilt, only UPF are allowed
!
IF( ANY( upf(1:nsp)%tvanp ) ) THEN
CALL errore( &
' check_types_order ', ' vanderbilt pseudo, not yet implemented in FPMD ', 1 )
CALL errore( ' check_types_order ', &
' vanderbilt pseudo, not yet implemented in FPMD ', 1 )
END IF
!
! non-local species must be ahead the local one,
@ -194,7 +193,7 @@ END FUNCTION calculate_dx
USE mp, ONLY: mp_bcast, mp_sum
USE io_global, ONLY: stdout, ionode, ionode_id
USE uspp, ONLY : okvan
USE uspp_param, ONLY : zp, tvanp, oldvan
USE uspp_param, ONLY : upf, oldvan
USE atom, ONLY: numeric, nlcc, oc, lchi, nchi
USE cvan, ONLY: nvb
use ions_base, only: zv, nsp
@ -364,8 +363,9 @@ END FUNCTION calculate_dx
! ... Zv = valence charge of the (pseudo-)atom, read from PP files,
! ... is set equal to Zp = pseudo-charge of the pseudopotential
! (should be moved out from here)
zv(is) = zp(is)
zv(is) = upf(is)%zp
CALL mp_sum( ierr )
IF( ierr /= 0 ) THEN
@ -388,7 +388,7 @@ END FUNCTION calculate_dx
! check on input ordering: US first, NC later
!
if(is > 1) then
if ( (.NOT. tvanp(is-1)) .AND. tvanp(is) ) then
if ( (.NOT. upf(is-1)%tvanp) .AND. upf(is)%tvanp ) then
call errore ('readpp', &
'ultrasoft PPs must precede norm-conserving',is)
endif
@ -396,7 +396,7 @@ END FUNCTION calculate_dx
!
! count u-s vanderbilt species
!
if (tvanp(is)) nvb=nvb+1
if (upf(is)%tvanp) nvb=nvb+1
!
END IF
@ -436,7 +436,7 @@ END FUNCTION calculate_dx
END IF
okvan = ( nvb > 0 )
!
RETURN
END SUBROUTINE readpp

2
Modules/make.depend
View File

@ -211,7 +211,6 @@ read_paw.o : atom.o
read_paw.o : constants.o
read_paw.o : functionals.o
read_paw.o : grid_paw_variables.o
read_paw.o : ions_base.o
read_paw.o : kind.o
read_paw.o : pseudo_types.o
read_paw.o : radial_grids.o
@ -255,6 +254,7 @@ upf_to_internal.o : uspp.o
uspp.o : constants.o
uspp.o : kind.o
uspp.o : parameters.o
uspp.o : pseudo_types.o
uspp.o : radial_grids.o
uspp.o : random_numbers.o
wannier.o : kind.o

1
Modules/parameters.f90
View File

@ -38,7 +38,6 @@ MODULE parameters
INTEGER, PARAMETER :: natx = 5000 ! maximum number of atoms
INTEGER, PARAMETER :: nbndxx = 10000 ! maximum number of electronic states
INTEGER, PARAMETER :: npkx = npk ! maximum number of K points
INTEGER, PARAMETER :: ncnsx = 101 ! maximum number of constraints
INTEGER, PARAMETER :: nspinx = 2 ! maximum number of spinors
INTEGER, PARAMETER :: nhclm = 4 ! maximum number NH chain length,

19
Modules/pseudo_types.f90
View File

@ -1,12 +1,9 @@
!
! Copyright (C) 2002-2003 PWSCF-FPMD-CP90 group
! Copyright (C) 2002-2007 Quantum-Espresso group
! This file is distributed under the terms of the
! GNU General Public License. See the file `License'
! in the root directory of the present distribution,
! or http://www.gnu.org/copyleft/gpl.txt .
!
! extracted from module "pseudo_types" of FPMD
!
MODULE pseudo_types
@ -17,9 +14,6 @@
USE kinds, ONLY: DP
USE parameters, ONLY: cp_lmax, lmaxx
use radial_grids, ONLY: ndmx, radial_grid_type
! USE ld1_parameters, ONLY: nwfsx
USE parameters, ONLY: nwfsx
! USE ld1_parameters, ONLY: nwfsx
IMPLICIT NONE
SAVE
@ -103,6 +97,11 @@ END TYPE paw_t
INTEGER :: mesh ! number of point in the radial mesh
INTEGER :: nwfc ! number of wavefunctions
INTEGER :: nbeta ! number of projectors
INTEGER :: kkbeta ! kkbeta=max(kbeta(:))
! kbeta<=mesh is the number of grid points for each beta function
! beta(r,nb) = 0 for r > r(kbeta(nb))
! kkbeta<=mesh is the largest of such number so that for all beta
! beta(r,nb) = 0 for r > r(kkbeta)
CHARACTER(LEN=2), POINTER :: els(:) ! els(nwfc)
CHARACTER(LEN=2), POINTER :: els_beta(:) ! els(nbeta)
INTEGER, POINTER :: lchi(:) ! lchi(nwfc)
@ -112,7 +111,7 @@ END TYPE paw_t
REAL(DP), POINTER :: rho_atc(:) ! rho_atc(mesh)
REAL(DP), POINTER :: vloc(:) ! vloc(mesh)
INTEGER, POINTER :: lll(:) ! lll(nbeta)
INTEGER, POINTER :: kkbeta(:) ! kkbeta(nbeta)
INTEGER, POINTER :: kbeta(:) ! kbeta(nbeta)
REAL(DP), POINTER :: beta(:,:) ! beta(mesh,nbeta)
INTEGER :: nd
REAL(DP), POINTER :: dion(:,:) ! dion(nbeta,nbeta)
@ -193,7 +192,7 @@ END TYPE paw_t
NULLIFY( upf%nn, upf%rcut)
NULLIFY( upf%els_beta)
NULLIFY( upf%rcutus, upf%epseu)
NULLIFY( upf%lll, upf%jjj, upf%kkbeta, upf%beta, upf%dion )
NULLIFY( upf%lll, upf%jjj, upf%kbeta, upf%beta, upf%dion )
NULLIFY( upf%rinner, upf%qqq, upf%qfunc, upf%qfcoef )
NULLIFY( upf%chi )
NULLIFY( upf%rho_at )
@ -229,7 +228,7 @@ END TYPE paw_t
IF( ASSOCIATED( upf%vloc ) ) DEALLOCATE( upf%vloc )
IF( ASSOCIATED( upf%lll ) ) DEALLOCATE( upf%lll )
IF( ASSOCIATED( upf%jjj ) ) DEALLOCATE( upf%jjj )
IF( ASSOCIATED( upf%kkbeta ) ) DEALLOCATE( upf%kkbeta )
IF( ASSOCIATED( upf%kbeta ) ) DEALLOCATE( upf%kbeta )
IF( ASSOCIATED( upf%beta ) ) DEALLOCATE( upf%beta )
IF( ASSOCIATED( upf%dion ) ) DEALLOCATE( upf%dion )
IF( ASSOCIATED( upf%rinner ) ) DEALLOCATE( upf%rinner )

57
Modules/read_paw.f90
View File

@ -333,14 +333,12 @@ subroutine set_pseudo_paw (is, pawset)
! PWSCF modules
!
USE radial_grids, ONLY: ndmx
USE atom, ONLY: rgrid, msh, &
chi, oc, nchi, lchi, jchi, rho_at, rho_atc, nlcc
USE atom, ONLY: rgrid, msh, chi, oc, nchi, lchi, jchi, rho_at, &
rho_atc, nlcc
! USE pseud, ONLY: lloc, lmax
USE uspp_param, ONLY: vloc_at, dion, betar, qqq, qfcoef, qfunc, nqf, nqlc, &
rinner, nbeta, kkbeta, lll, jjj, psd, tvanp, zp
USE uspp_param, ONLY: upf
USE funct, ONLY: set_dft_from_name, dft_is_meta, dft_is_hybrid
!
USE ions_base, ONLY: zv
! USE spin_orb, ONLY: lspinorb
USE pseudo_types
USE constants, ONLY: FPI
@ -375,9 +373,9 @@ subroutine set_pseudo_paw (is, pawset)
PRINT *, nrs, pawset%grid%r(nrs)
PRINT *, nrc, pawset%grid%r(nrc)
!
zp(is) = pawset%zval
psd (is)= pawset%symbol
tvanp(is)=.true.
upf(is)%zp = pawset%zval
upf(is)%psd = pawset%symbol
upf(is)%tvanp=.true.
tpawp(is)=.true.
nlcc(is) = pawset%nlcc
call set_dft_from_name( pawset%dft )
@ -415,30 +413,36 @@ subroutine set_pseudo_paw (is, pawset)
#endif
end do
!
nbeta(is)= pawset%nwfc
kkbeta(is)=0
upf(is)%nbeta= pawset%nwfc
allocate ( upf(is)%kbeta(pawset%nwfc) )
do nb=1,pawset%nwfc
kkbeta(is)=max(pawset%ikk(nb),kkbeta(is))
upf(is)%kbeta(nb)=pawset%ikk(nb)
end do
betar(1:pawset%grid%mesh, 1:pawset%nwfc, is) = pawset%proj(1:pawset%grid%mesh, 1:pawset%nwfc)
dion(1:pawset%nwfc, 1:pawset%nwfc, is) = pawset%dion(1:pawset%nwfc, 1:pawset%nwfc)
allocate (upf(is)%beta(1:pawset%grid%mesh, 1:pawset%nwfc))
upf(is)%beta(1:pawset%grid%mesh, 1:pawset%nwfc) = &
pawset%proj(1:pawset%grid%mesh, 1:pawset%nwfc)
allocate(upf(is)%dion(1:pawset%nwfc, 1:pawset%nwfc))
upf(is)%dion(1:pawset%nwfc, 1:pawset%nwfc) = pawset%dion(1:pawset%nwfc, 1:pawset%nwfc)
kdiff(1:pawset%nwfc, 1:pawset%nwfc, is) = pawset%kdiff(1:pawset%nwfc, 1:pawset%nwfc)
! HOPE!
! lmax(is) = pawset%lmax
nqlc(is) = 2*pawset%lmax+1
nqf (is) = 0 !! no rinner, all numeric
lll(1:pawset%nwfc,is) = pawset%l(1:pawset%nwfc)
rinner(1:nqlc(is),is) = 0._dp !! no rinner, all numeric
upf(is)%nqlc = 2*pawset%lmax+1
upf(is)%nqf = 0 !! no rinner, all numeric
allocate (upf(is)%lll(pawset%nwfc) )
upf(is)%lll(1:pawset%nwfc) = pawset%l(1:pawset%nwfc)
allocate (upf(is)%rinner(upf(is)%nqlc))
upf(is)%rinner(1:upf(is)%nqlc) = 0._dp !! no rinner, all numeric
!
! integral of augmentation charges vanishes for different values of l
!
allocate ( upf(is)%qqq(pawset%nwfc,pawset%nwfc))
do i = 1, pawset%nwfc
do j = 1, pawset%nwfc
if (pawset%l(i)==pawset%l(j)) then
qqq(i,j,is) = pawset%augmom(i,j,0) !!gf spherical approximation
upf(is)%qqq(i,j) = pawset%augmom(i,j,0) !!gf spherical approximation
else
qqq(i,j,is) = 0._dp
upf(is)%qqq(i,j) = 0._dp
end if
end do
end do
@ -455,10 +459,13 @@ subroutine set_pseudo_paw (is, pawset)
end do
! triangularize matrix of qfunc's
allocate ( upf(is)%qfunc(1:pawset%grid%mesh,pawset%nwfc,pawset%nwfc) )
do nb = 1, pawset%nwfc
do mb = nb, pawset%nwfc
ijv = mb * (mb-1) / 2 + nb
qfunc (1:pawset%grid%mesh, ijv, is) = pawset%augfun(1:pawset%grid%mesh,nb,mb,0)
!!! qfunc (1:pawset%grid%mesh, ijv, is) = pawset%augfun(1:pawset%grid%mesh,nb,mb,0)
upf(is)%qfunc (1:pawset%grid%mesh, nb,mb) = &
pawset%augfun(1:pawset%grid%mesh,nb,mb,0)
enddo
enddo
! augfun(1:pawset%grid%mesh,1:pawset%nwfc,1:pawset%nwfc,0:2*pawset%lmax,is) = &
@ -522,8 +529,9 @@ subroutine set_pseudo_paw (is, pawset)
!!$ jchi(1:pawset%nwfc, is) = pawset%jchi(1:pawset%nwfc)
!!$ jjj(1:pawset%nbeta, is) = pawset%jjj(1:pawset%nbeta)
!!$ else
jchi(1:pawset%nwfc, is) = 0._dp
jjj(1:pawset%nwfc, is) = 0._dp
jchi(1:pawset%nwfc,is) = 0._dp
allocate (upf(is)%jjj(1:pawset%nwfc))
upf(is)%jjj(1:pawset%nwfc) = 0._dp
!!$ endif
!
if ( pawset%nlcc) then
@ -548,7 +556,8 @@ subroutine set_pseudo_paw (is, pawset)
!!! answer (pltz): I don't, but it breaked dependencies (removed!)
! lloc(is) = 0
!!!
vloc_at(1:pawset%grid%mesh,is) = pawset%psloc(1:pawset%grid%mesh)
allocate (upf(is)%vloc(1:pawset%grid%mesh))
upf(is)%vloc(1:pawset%grid%mesh) = pawset%psloc(1:pawset%grid%mesh)
#if defined __DO_NOT_CUTOFF_PAW_FUNC
aevloc_at(1:pawset%grid%mesh,is) = pawset%aeloc(1:pawset%grid%mesh)
psvloc_at(1:pawset%grid%mesh,is) = pawset%psloc(1:pawset%grid%mesh)
@ -573,8 +582,6 @@ subroutine set_pseudo_paw (is, pawset)
!
5 msh (is) = 2 * ( (msh (is) + 1) / 2) - 1
zv(is) = zp(is) !!! maybe not needed: it is done in setup
end subroutine set_pseudo_paw
!=----------------------------------------------------------------------------=!

7
Modules/read_upf.f90
View File

@ -340,7 +340,7 @@ subroutine read_pseudo_nl (upf, iunps)
if ( upf%nbeta == 0) then
upf%nqf = 0
upf%nqlc= 0
ALLOCATE( upf%kkbeta( 1 ) )
ALLOCATE( upf%kbeta( 1 ) )
ALLOCATE( upf%lll( 1 ) )
ALLOCATE( upf%beta( 0:upf%mesh, 1 ) )
ALLOCATE( upf%dion( 1, 1 ) )
@ -353,7 +353,7 @@ subroutine read_pseudo_nl (upf, iunps)
ALLOCATE( upf%els_beta( 1 ) )
return
end if
ALLOCATE( upf%kkbeta( upf%nbeta ) )
ALLOCATE( upf%kbeta( upf%nbeta ) )
ALLOCATE( upf%lll( upf%nbeta ) )
ALLOCATE( upf%beta( 0:upf%mesh, upf%nbeta ) )
ALLOCATE( upf%dion( upf%nbeta, upf%nbeta ) )
@ -373,7 +373,8 @@ subroutine read_pseudo_nl (upf, iunps)
call scan_begin (iunps, "BETA", .false.)
read (iunps, *, err = 100, end = 100) idum, upf%lll(nb), dummy
read (iunps, *, err = 100, end = 100) ikk
upf%kkbeta(nb) = ikk
upf%kbeta(nb) = ikk
upf%kkbeta = MAX ( upf%kkbeta, upf%kbeta(nb) )
read (iunps, *, err = 100, end = 100) (upf%beta(ir,nb), ir=1,ikk)
read (iunps, *, err=200,iostat=ios) upf%rcut(nb), upf%rcutus(nb)

40
Modules/read_uspp.f90
View File

@ -124,7 +124,6 @@ CONTAINS
& lp, &! counter on Q angular momenta
& l, &! counter on angular momenta
& iv, jv, ijv, &! beta function counter
& kkbeta1, &! number of grid points for which betar.ne.0
& ir ! mesh points counter
!
character(len=20) :: title
@ -268,16 +267,14 @@ CONTAINS
! reads the number of beta functions
!
read( iunps, '(2i5)', err=100, iostat=ios ) &
upf%nbeta, kkbeta1
upf%nbeta, upf%kkbeta
!
! BEWARE: upf%kkbeta is an array (one per beta function)
!
ALLOCATE ( upf%kkbeta(upf%nbeta) )
upf%kkbeta(:) = kkbeta1
ALLOCATE ( upf%kbeta(upf%nbeta) )
upf%kbeta(:) = upf%kkbeta
!
if( upf%nbeta > nbrx .or. upf%nbeta <0 ) &
call errore( 'readvan','nbeta wrong or too large', is )
if( ANY (upf%kkbeta > upf%mesh) .or. ANY(upf%kkbeta < 0) ) &
if( upf%kkbeta > upf%mesh .or. upf%kkbeta < 0 ) &
call errore( 'readvan','kkbeta wrong or too large', is )
!
! Now reads the main Vanderbilt parameters
@ -290,8 +287,8 @@ CONTAINS
do iv=1,upf%nbeta
read( iunps, '(i5)',err=100, iostat=ios ) upf%lll(iv)
read( iunps, '(1p4e19.11)',err=100, iostat=ios ) &
eee(iv), ( upf%beta(ir,iv), ir=1,upf%kkbeta(iv) )
do ir=upf%kkbeta(iv)+1,upf%mesh
eee(iv), ( upf%beta(ir,iv), ir=1,upf%kkbeta )
do ir=upf%kkbeta+1,upf%mesh
upf%beta(ir,iv)=0.0_DP
enddo
if ( upf%lll(iv) > lmaxx .or. upf%lll(iv) < 0 ) &
@ -305,9 +302,9 @@ CONTAINS
ijv = jv * (jv-1) / 2 + iv
read( iunps, '(1p4e19.11)', err=100, iostat=ios ) &
upf%dion(iv,jv), ddd(iv,jv), upf%qqq(iv,jv), &
(upf%qfunc(ir,iv,jv),ir=1,upf%kkbeta(1)), &
(upf%qfunc(ir,iv,jv),ir=1,upf%kkbeta), &
((upf%qfcoef(i,lp,iv,jv),i=1,upf%nqf),lp=1,upf%nqlc)
do ir=upf%kkbeta(1)+1,upf%mesh
do ir=upf%kkbeta+1,upf%mesh
upf%qfunc(ir,iv,jv)=0.0_DP
enddo
!
@ -441,7 +438,7 @@ CONTAINS
end if
WRITE( stdout,1000)
1000 format(4x,'| new generation scheme:',32x,'|')
WRITE( stdout,1100) upf%nbeta, upf%kkbeta(1), rcloc
WRITE( stdout,1100) upf%nbeta, upf%kkbeta, rcloc
1100 format(4x,'| nbeta = ',i2,5x,'kkbeta =',i5,5x,'rcloc =',f10.4,4x,&
& '|'/4x,'| ibeta l epsilon rcut',25x,'|')
do iv = 1, upf%nbeta
@ -473,7 +470,7 @@ CONTAINS
!
allocate ( a(upf%nqf,upf%nqf), ainv(upf%nqf,upf%nqf) )
allocate ( b(upf%nqf), x(upf%nqf) )
ALLOCATE ( qrl(upf%kkbeta(1), upf%nqlc) )
ALLOCATE ( qrl(upf%kkbeta, upf%nqlc) )
!
do iv=1,upf%nbeta
do jv=iv,upf%nbeta
@ -491,7 +488,7 @@ CONTAINS
! read q_l(r) for all l
!
read(iunps,*, err=100) &
( (qrl(ir,l),ir=1,upf%kkbeta(1)), l=lmin,lmax)
( (qrl(ir,l),ir=1,upf%kkbeta), l=lmin,lmax)
!
!!! ijv = jv * (jv-1) / 2 + iv
!
@ -499,7 +496,7 @@ CONTAINS
!
! reconstruct rinner
!
do ir=upf%kkbeta(1),1,-1
do ir=upf%kkbeta,1,-1
if ( abs(qrl(ir,l)-upf%qfunc(ir,iv,jv)) > 1.0d-6) go to 10
end do
10 irinner = ir+1
@ -607,8 +604,7 @@ CONTAINS
pseudotype,&! the type of pseudopotential
ios, &! I/O control
ndum, &! dummy integer variable
l, &! counter on angular momentum
ikk ! the kkbeta for each beta
l ! counter on angular momentum
real(DP):: &
x, &! auxiliary variable
etotps, &! total energy of the pseudoatom
@ -695,15 +691,17 @@ CONTAINS
if ( upf%oc(nb) <= 0.0_DP) upf%oc(nb) = -1.0_DP
enddo
!
ALLOCATE ( upf%kkbeta(upf%nbeta) )
ALLOCATE ( upf%kbeta(upf%nbeta) )
ALLOCATE ( upf%dion(upf%nbeta,upf%nbeta), upf%qqq(upf%nbeta,upf%nbeta) )
ALLOCATE ( upf%beta(upf%mesh,upf%nbeta) )
ALLOCATE ( upf%qfunc(upf%mesh,upf%nbeta,upf%nbeta) )
upf%kkbeta = 0
do nb=1,upf%nbeta
read ( iunps, '(i6)',err=100, iostat=ios ) upf%kkbeta(nb)
read ( iunps, '(i6)',err=100, iostat=ios ) upf%kbeta(nb)
upf%kkbeta = MAX ( upf%kkbeta, upf%kbeta(nb) )
read ( iunps, '(1p4e19.11)',err=100, iostat=ios ) &
( upf%beta(ir,nb), ir=1,upf%kkbeta(nb))
do ir=upf%kkbeta(nb)+1,upf%mesh
( upf%beta(ir,nb), ir=1,upf%kbeta(nb))
do ir=upf%kbeta(nb)+1,upf%mesh
upf%beta(ir,nb)=0.0_DP
enddo
do mb=1,nb

5
Modules/upf_to_internal.f90
View File

@ -65,10 +65,7 @@ subroutine set_pseudo_upf (is, upf)
chi(1:upf%mesh, 1:upf%nwfc, is) = upf%chi(1:upf%mesh, 1:upf%nwfc)
!
nbeta(is)= upf%nbeta
kkbeta(is)=0
do nb=1,upf%nbeta
kkbeta(is)=max(upf%kkbeta(nb),kkbeta(is))
end do
kkbeta(is)=upf%kkbeta
betar(1:upf%mesh, 1:upf%nbeta, is) = upf%beta(1:upf%mesh, 1:upf%nbeta)
dion(1:upf%nbeta, 1:upf%nbeta, is) = upf%dion(1:upf%nbeta, 1:upf%nbeta)
!

3
Modules/uspp.f90
View File

@ -12,9 +12,12 @@ MODULE uspp_param
USE kinds, ONLY : DP
USE parameters, ONLY : lqmax, nbrx, npsx, nqfx
USE radial_grids, ONLY: ndmx
USE pseudo_types, ONLY: pseudo_upf
!
SAVE
!
TYPE (pseudo_upf), ALLOCATABLE, TARGET :: upf(:)
CHARACTER(LEN=2 ) :: psd(npsx) ! name of the pseudopotential
REAL(DP) :: &

4
PW/addusdens.f90
View File

@ -42,7 +42,7 @@ subroutine addusdens_g
USE noncollin_module, ONLY : noncolin
USE scf, ONLY : rho
USE uspp, ONLY : becsum, okvan
USE uspp_param, ONLY : lmaxq, tvanp, nh
USE uspp_param, ONLY : upf, lmaxq, nh
USE wvfct, ONLY : gamma_only
USE wavefunctions_module, ONLY : psic
!
@ -81,7 +81,7 @@ subroutine addusdens_g
qmod (ig) = sqrt (gg (ig) )
enddo
do nt = 1, ntyp
if (tvanp (nt) ) then
if ( upf(nt)%tvanp ) then
ijh = 0
do ih = 1, nh (nt)
do jh = ih, nh (nt)

6
PW/addusforce.f90
View File

@ -23,7 +23,7 @@ subroutine addusforce (forcenl)
USE lsda_mod, ONLY : nspin
USE scf, ONLY : vr, vltot
USE uspp, ONLY : becsum, okvan
USE uspp_param, ONLY : lmaxq, tvanp, nh, nhm
USE uspp_param, ONLY : upf, lmaxq, nh, nhm
USE wvfct, ONLY : gamma_only
!
implicit none
@ -77,13 +77,13 @@ subroutine addusforce (forcenl)
! I = sum_G i G_a exp(-iR.G) Q_nm v^*
!
do nt = 1, ntyp
if (tvanp (nt) ) then
if ( upf(nt)%tvanp ) then
ijh = 1
do ih = 1, nh (nt)
do jh = ih, nh (nt)
call qvan2 (ngm, ih, jh, nt, qmod, qgm, ylmk0)
do na = 1, nat
if (ityp (na) .eq.nt) then
if (ityp (na) == nt) then
!
! The product of potential, structure factor and iG
!

6
PW/addusstress.f90
View File

@ -21,7 +21,7 @@ subroutine addusstres (sigmanlc)
USE lsda_mod, ONLY : nspin
USE scf, ONLY : vr, vltot
USE uspp, ONLY : becsum, okvan
USE uspp_param, ONLY : lmaxq, tvanp, nh, nhm
USE uspp_param, ONLY : upf, lmaxq, nh, nhm
USE wvfct, ONLY : gamma_only
!
implicit none
@ -93,13 +93,13 @@ subroutine addusstres (sigmanlc)
do ipol = 1, 3
call dylmr2 (lmaxq * lmaxq, ngm, g, gg, dylmk0, ipol)
do nt = 1, ntyp
if (tvanp (nt) ) then
if ( upf(nt)%tvanp ) then
ijh = 1
do ih = 1, nh (nt)
do jh = ih, nh (nt)
call dqvan2 (ngm, ih, jh, nt, qmod, qgm, ylmk0, dylmk0, ipol)
do na = 1, nat
if (ityp (na) .eq.nt) then
if (ityp (na) == nt) then
!
do is = 1, nspin
do jpol = 1, ipol

12
PW/allocate_nlpot.f90
View File

@ -39,7 +39,7 @@ subroutine allocate_nlpot
nqxq, spline_ps
USE uspp, ONLY : indv, nhtol, nhtolm, qq, dvan, deeq, vkb, nkb, &
nkbus, nhtoj, becsum, qq_so, dvan_so, deeq_nc
USE uspp_param, ONLY : lmaxq, lmaxkb, lll, nbeta, nh, nhm, nbetam,tvanp
USE uspp_param, ONLY : upf, lmaxq, lmaxkb, nh, nhm, nbetam
USE spin_orb, ONLY : lspinorb, fcoef
!
implicit none
@ -64,16 +64,16 @@ subroutine allocate_nlpot
lmaxkb = - 1
do nt = 1, nsp
nh (nt) = 0
do nb = 1, nbeta (nt)
nh (nt) = nh (nt) + 2 * lll (nb, nt) + 1
lmaxkb = max (lmaxkb, lll (nb, nt) )
do nb = 1, upf(nt)%nbeta
nh (nt) = nh (nt) + 2 * upf(nt)%lll(nb) + 1
lmaxkb = max (lmaxkb, upf(nt)%lll(nb) )
enddo
enddo
!
! calculate the maximum number of beta functions
!
nhm = MAXVAL (nh (1:nsp))
nbetam = MAXVAL (nbeta (1:nsp))
nbetam = MAXVAL (upf(:)%nbeta)
!
! calculate the number of beta functions of the solid
!
@ -82,7 +82,7 @@ subroutine allocate_nlpot
do na = 1, nat
nt = ityp(na)
nkb = nkb + nh (nt)
if (tvanp(nt)) nkbus = nkbus + nh (nt)
if (upf(nt)%tvanp) nkbus = nkbus + nh (nt)
enddo
!
allocate (indv( nhm, nsp))

64
PW/average_pp.f90
View File

@ -12,7 +12,7 @@ SUBROUTINE average_pp ( ntyp )
USE kinds, ONLY : DP
USE atom, ONLY : chi, nchi, lchi, jchi, rgrid
USE spin_orb, ONLY : so
USE uspp_param, ONLY : betar, dion, jjj, lll, nbeta, tvanp
USE uspp_param, ONLY : upf
!
IMPLICIT NONE
!
@ -26,65 +26,67 @@ SUBROUTINE average_pp ( ntyp )
!
IF ( so(nt) ) THEN
!
IF ( tvanp(nt) ) &
IF ( upf(nt)%tvanp ) &
CALL errore( 'setup', 'US j-average not yet implemented', 1 )
!
nbe = 0
!
DO nb = 1, nbeta(nt)
DO nb = 1, upf(nt)%nbeta
!
nbe = nbe + 1
!
IF ( lll(nb,nt) /= 0 .AND. &
ABS( jjj(nb,nt) - lll(nb,nt) - 0.5D0 ) < 1.D-7 ) nbe = nbe - 1
IF ( upf(nt)%lll(nb) /= 0 .AND. &
ABS( upf(nt)%jjj(nb) - upf(nt)%lll(nb) - 0.5D0 ) < 1.D-7 ) &
nbe = nbe - 1
END DO
!
nbeta(nt) = nbe
upf(nt)%nbeta = nbe
!
nbe = 0
!
DO nb = 1, nbeta(nt)
DO nb = 1, upf(nt)%nbeta
!
nbe = nbe + 1
!
l = lll(nbe,nt)
l = upf(nt)%lll(nbe)
!
IF ( l /= 0 ) THEN
!
IF (ABS(jjj(nbe,nt)-lll(nbe,nt)+0.5d0).LT.1.d-7) THEN
IF (ABS(jjj(nbe+1,nt)-lll(nbe+1,nt)-0.5d0).GT.1.d-7) &
call errore('setup','wrong beta functions',1)
IF (ABS(upf(nt)%jjj(nbe)-upf(nt)%lll(nbe)+0.5d0) < 1.d-7) THEN
IF ( ABS( upf(nt)%jjj(nbe+1)-upf(nt)%lll(nbe+1)-0.5d0 ) &
> 1.d-7 ) call errore('setup','wrong beta functions',1)
ind=nbe+1
ind1=nbe
ELSE
IF (ABS(jjj(nbe+1,nt)-lll(nbe+1,nt)+0.5d0).GT.1.d-7) &
IF (ABS(upf(nt)%jjj(nbe+1)-upf(nt)%lll(nbe+1)+0.5d0) > 1.d-7) &
call errore('setup','wrong beta functions',1)
ind=nbe
ind1=nbe+1
ENDIF
!
vionl = ( ( l + 1.D0 ) * dion(ind,ind,nt) + &
l * dion(ind1,ind1,nt) ) / ( 2.D0 * l + 1.D0 )
vionl = ( ( l + 1.D0 ) * upf(nt)%dion(ind,ind) + &
l * upf(nt)%dion(ind1,ind1) ) / ( 2.D0 * l + 1.D0 )
!
betar(1:rgrid(nt)%mesh,nb,nt) = 1.D0 / ( 2.D0 * l + 1.D0 ) * &
( ( l + 1.D0 ) * SQRT( dion(ind,ind,nt) / vionl ) * &
betar(1:rgrid(nt)%mesh,ind,nt) + &
l * SQRT( dion(ind1,ind1,nt) / vionl ) * &
betar(1:rgrid(nt)%mesh,ind1,nt) )
upf(nt)%beta(1:rgrid(nt)%mesh,nb) = 1.D0 / ( 2.D0 * l + 1.D0 ) * &
( ( l + 1.D0 ) * SQRT( upf(nt)%dion(ind,ind) / vionl ) * &
upf(nt)%beta(1:rgrid(nt)%mesh,ind) + &
l * SQRT( upf(nt)%dion(ind1,ind1) / vionl ) * &
upf(nt)%beta(1:rgrid(nt)%mesh,ind1) )
!
dion(nb,nb,nt) = vionl
upf(nt)%dion(nb,nb) = vionl
!
nbe = nbe + 1
!
ELSE
!
betar(1:rgrid(nt)%mesh,nb,nt) = betar(1:rgrid(nt)%mesh,nbe,nt)
upf(nt)%beta(1:rgrid(nt)%mesh,nb) = &
upf(nt)%beta(1:rgrid(nt)%mesh,nbe)
!
dion(nb,nb,nt) = dion(nbe,nbe,nt)
upf(nt)%dion(nb,nb) = upf(nt)%dion(nbe,nbe)
!
END IF
!
lll(nb,nt)=lll(nbe,nt)
upf(nt)%lll(nb)=upf(nt)%lll(nbe)
!
END DO
!
@ -95,7 +97,8 @@ SUBROUTINE average_pp ( ntyp )
nbe = nbe + 1
!
IF ( lchi(nb,nt) /= 0 .AND. &
ABS(jchi(nb,nt)-lchi(nb,nt)-0.5D0 ) < 1.D-7 ) nbe = nbe - 1
ABS(jchi(nb,nt)-lchi(nb,nt)-0.5D0 ) < 1.D-7 ) &
nbe = nbe - 1
!
END DO
!
@ -111,19 +114,20 @@ SUBROUTINE average_pp ( ntyp )
!
IF ( l /= 0 ) THEN
!
IF (ABS(jchi(nbe,nt)-lchi(nbe,nt)+0.5d0).LT.1.d-7) THEN
IF (ABS(jchi(nbe+1,nt)-lchi(nbe+1,nt)-0.5d0).GT.1.d-7) &
call errore('setup','wrong chi functions',1)
IF (ABS(jchi(nbe,nt)-lchi(nbe,nt)+0.5d0) < 1.d-7) THEN
IF ( ABS(jchi(nbe+1,nt)-lchi(nbe+1,nt)-0.5d0) > &
1.d-7) call errore('setup','wrong chi functions',1)
ind=nbe+1
ind1=nbe
ELSE
IF (ABS(jchi(nbe+1,nt)-lchi(nbe+1,nt)+0.5d0).GT.1.d-7) &
call errore('setup','wrong chi functions',1)
IF ( ABS(jchi(nbe+1,nt)-lchi(nbe+1,nt)+0.5d0) > &
1.d-7) call errore('setup','wrong chi functions',1)
ind=nbe
ind1=nbe+1
END IF
!
chi(1:rgrid(nt)%mesh,nb,nt)=((l+1.D0) * chi(1:rgrid(nt)%mesh,ind,nt)+ &
chi(1:rgrid(nt)%mesh,nb,nt) = &
((l+1.D0) * chi(1:rgrid(nt)%mesh,ind,nt)+ &
l * chi(1:rgrid(nt)%mesh,ind1,nt)) / ( 2.D0 * l + 1.D0 )
!
nbe = nbe + 1

51
PW/bp_calc_btq.f90
View File

@ -18,63 +18,56 @@ SUBROUTINE calc_btq(ql,qr_k,idbes)
USE ions_base, ONLY : ntyp => nsp
USE cell_base, ONLY: omega
USE constants, ONLY: fpi
USE uspp_param, ONLY: lmaxq, qfunc, qfcoef, nqf, rinner, lll, &
nbeta, nbetam, kkbeta, tvanp
USE uspp_param, ONLY: upf, nbetam, lmaxq
!
IMPLICIT NONE
!
REAL(DP) :: ql, qr_k(nbetam,nbetam,lmaxq,ntyp)
INTEGER :: idbes
!
INTEGER :: msh_bp, i, np, l, ilmin, ilmax, iv, jv, ijv
INTEGER :: i, np, l, ilmin, ilmax, iv, jv, ijv
REAL(DP) :: qrk
REAL(DP), ALLOCATABLE :: jl(:), aux(:)
!
DO np=1,ntyp
IF (tvanp(np)) THEN
msh_bp=kkbeta(np)
ALLOCATE ( jl(msh_bp), aux(msh_bp) )
DO iv =1, nbeta(np)
DO jv =iv, nbeta(np)
!
IF ( upf(np)%tvanp ) THEN
!
ALLOCATE ( jl(upf(np)%kkbeta), aux(upf(np)%kkbeta) )
DO iv =1, upf(np)%nbeta
DO jv =iv, upf(np)%nbeta
ijv = jv * (jv-1) / 2 + iv
ilmin = iabs(lll(iv,np)-lll(jv,np))
ilmax = iabs(lll(iv,np)+lll(jv,np))
! only need to calculate for for lmin,lmin+2 ...lmax-2,lmax
ilmin = abs ( upf(np)%lll(iv) - upf(np)%lll(jv) )
ilmax = upf(np)%lll(iv) + upf(np)%lll(jv)
! only need to calculate for l=lmin,lmin+2 ...lmax-2,lmax
DO l = ilmin,ilmax,2
aux(:) = 0.0_DP
DO i = msh_bp,2,-1
IF (rgrid(np)%r(i) .LT. rinner(l+1,np)) GOTO 100
aux(i) = qfunc(i,ijv,np)
DO i = upf(np)%kkbeta,2,-1
IF (rgrid(np)%r(i) .LT. upf(np)%rinner(l+1)) GOTO 100
!!! aux(i) = qfunc(i,ijv,np) TEMP
aux(i) = upf(np)%qfunc(i,iv,jv)
ENDDO
100 CALL setqf(qfcoef(1,l+1,iv,jv,np),aux(1),rgrid(np)%r &
,nqf(np),l,i)
100 CALL setqf ( upf(np)%qfcoef(1,l+1,iv,jv), aux(1), &
rgrid(np)%r, upf(np)%nqf, l, i )
IF (idbes .EQ. 1) THEN
IF (idbes == 1) THEN
!
CALL sph_dbes( msh_bp, rgrid(np)%r, ql, l, jl )
!
! ... this is the old call
!
! CALL dbess( ql, l+1, msh_bp, r(1,np), jl )
CALL sph_dbes( upf(np)%kkbeta, rgrid(np)%r, ql, l, jl )
!
ELSE
!
CALL sph_bes( msh_bp, rgrid(np)%r, ql, l, jl )
!
! ... this is the old call
!
! CALL bess( ql, l+1, msh_bp, r(1,np), jl )
CALL sph_bes( upf(np)%kkbeta, rgrid(np)%r, ql, l, jl )
!
ENDIF
! jl is the Bessel function (or its derivative) calculated at ql
! now integrate qfunc*jl*r^2 = Bessel transform of qfunc
DO i=1, msh_bp
DO i=1, upf(np)%kkbeta
aux(i) = jl(i)*aux(i)
ENDDO
! if (tlog(np)) then
CALL radlg1(msh_bp,aux,rgrid(np)%rab,qrk)
CALL radlg1(upf(np)%kkbeta,aux,rgrid(np)%rab,qrk)
qr_k(iv,jv,l+1,np) = qrk*fpi/omega
qr_k(jv,iv,l+1,np) = qr_k(iv,jv,l+1,np)

4
PW/c_phase_field.f90
View File

@ -34,7 +34,7 @@ SUBROUTINE c_phase_field
USE gvect, ONLY : ngm, nr1, nr2, nr3, nrx1, nrx2, nrx3, &
ecutwfc, g, gcutm
USE uspp, ONLY : nkb, vkb, okvan
USE uspp_param, ONLY : lmaxq, nh, nhm, nbetam, tvanp
USE uspp_param, ONLY : upf, lmaxq, nh, nhm, nbetam
USE lsda_mod, ONLY : nspin
USE klist, ONLY : nelec, degauss, nks, xk, wk
USE wvfct, ONLY : npwx, npw, nbnd
@ -277,7 +277,7 @@ IF ((degauss > 0.01d0) .OR. (nbnd /= nelec/2)) &
! CALL setv(2*nhm*nhm*ntyp,0.d0,q_dk,1)
q_dk=(0.d0,0.d0)
DO np =1, ntyp
if(tvanp(np)) then
if( upf(nt)%tvanp ) then
DO iv = 1, nh(np)
DO jv = iv, nh(np)
call qvan3(iv,jv,np,pref,ylm_dk,qrad_dk)

55
PW/compute_qdipol.f90
View File

@ -15,8 +15,7 @@ SUBROUTINE compute_qdipol(dpqq)
USE atom, ONLY: rgrid
USE ions_base, ONLY: ntyp => nsp
USE uspp, only: nhtol, nhtolm, indv, nlx, ap
USE uspp_param, only: nbrx, nbeta, lll, kkbeta, qfunc, rinner, &
qfcoef, nqf, tvanp, nh, nhm
USE uspp_param, only: upf, nbetam, nh, nhm
implicit none
@ -26,51 +25,53 @@ SUBROUTINE compute_qdipol(dpqq)
integer :: nt, l, ir, nb, mb, ijv, ilast, ipol, ih, ivl, jh, jvl, lp, ndm
call start_clock('cmpt_qdipol')
ndm = MAXVAL (kkbeta(1:ntyp))
allocate (qrad2( nbrx , nbrx, ntyp))
ndm = MAXVAL ( upf(1:ntyp)%kkbeta )
allocate (qrad2( nbetam , nbetam, ntyp))
allocate (aux( ndm))
allocate (qtot( ndm, nbrx, nbrx))
allocate (qtot( ndm, nbetam, nbetam))
qrad2(:,:,:)=0.d0
dpqq=0.d0
do nt = 1, ntyp
if (tvanp (nt) ) then
if ( upf(nt)%tvanp ) then
l=1
!
! Only l=1 terms enter in the dipole of Q
!
do nb = 1, nbeta (nt)
do mb = nb, nbeta (nt)
do nb = 1, upf(nt)%nbeta
do mb = nb, upf(nt)%nbeta
ijv = mb * (mb-1) /2 + nb
if ((l.ge.abs(lll(nb,nt)-lll(mb,nt))) .and. &
(l.le.lll(nb,nt)+lll(mb,nt)) .and. &
(mod (l+lll(nb,nt)+lll(mb,nt),2) .eq.0) ) then
do ir = 1, kkbeta (nt)
if (rgrid(nt)%r(ir).ge.rinner(l+1, nt)) then
qtot(ir, nb, mb)=qfunc(ir,ijv,nt)
if ( ( l >= abs(upf(nt)%lll(nb) - upf(nt)%lll(mb)) ) .and. &
( l <= upf(nt)%lll(nb) + upf(nt)%lll(mb) ) .and. &
(mod (l+upf(nt)%lll(nb)+upf(nt)%lll(mb), 2) == 0) ) then
do ir = 1, upf(nt)%kkbeta
if (rgrid(nt)%r(ir) >= upf(nt)%rinner(l+1)) then
! qtot(ir, nb, mb)=qfunc(ir,ijv,nt) TEMP
qtot(ir, nb, mb)=upf(nt)%qfunc(ir,nb,mb)
else
ilast = ir
endif
enddo
if (rinner(l+1, nt).gt.0.d0) &
call setqf(qfcoef (1, l+1, nb, mb, nt), &
qtot(1,nb,mb), rgrid(nt)%r, nqf(nt),l,ilast)
if ( upf(nt)%rinner(l+1) > 0.0_dp) &
call setqf( upf(nt)%qfcoef (1, l+1, nb, mb), &
qtot(1,nb,mb), rgrid(nt)%r, upf(nt)%nqf, l, ilast)
endif
enddo
enddo
do nb=1, nbeta(nt)
do nb=1, upf(nt)%nbeta
!
! the Q are symmetric with respect to indices
!
do mb=nb, nbeta(nt)
if ( (l.ge.abs(lll(nb,nt)-lll(mb,nt) ) ) .and. &
(l.le.lll(nb,nt) + lll(mb,nt) ) .and. &
(mod(l+lll(nb,nt)+lll(mb,nt), 2).eq.0) ) then
do ir = 1, kkbeta (nt)
do mb=nb, upf(nt)%nbeta
if ( ( l >= abs(upf(nt)%lll(nb) - upf(nt)%lll(mb)) ) .and. &
( l <= upf(nt)%lll(nb) + upf(nt)%lll(mb) ) .and. &
(mod (l+upf(nt)%lll(nb)+upf(nt)%lll(mb), 2) == 0) ) then
do ir = 1, upf(nt)%kkbeta
aux(ir)=rgrid(nt)%r(ir)*qtot(ir, nb, mb)
enddo
call simpson (kkbeta(nt),aux,rgrid(nt)%rab,qrad2(nb,mb,nt))
call simpson ( upf(nt)%kkbeta, aux, rgrid(nt)%rab, &
qrad2(nb,mb,nt) )
endif
enddo
enddo
@ -87,7 +88,7 @@ SUBROUTINE compute_qdipol(dpqq)
fact=-fact
endif
do nt = 1,ntyp
if (tvanp(nt)) then
if ( upf(nt)%tvanp ) then
do ih = 1, nh(nt)
ivl = nhtolm(ih, nt)
mb = indv(ih, nt)
@ -96,8 +97,8 @@ SUBROUTINE compute_qdipol(dpqq)
nb=indv(jh,nt)
if (ivl > nlx) call errore('compute_qdipol',' ivl > nlx', ivl)
if (jvl > nlx) call errore('compute_qdipol',' jvl > nlx', jvl)
if (nb > nbrx) call errore('compute_qdipol',' nb > nbrx', nb)
if (mb > nbrx) call errore('compute_qdipol',' mb > nbrx', mb)
if (nb > nbetam) call errore('compute_qdipol',' nb > nbrx', nb)
if (mb > nbetam) call errore('compute_qdipol',' mb > nbrx', mb)
if (mb > nb) call errore('compute_qdipol',' mb > nb', 1)
dpqq(ih,jh,ipol,nt)=fact*ap(lp,ivl,jvl)*qrad2(mb,nb,nt)
dpqq(jh,ih,ipol,nt)=dpqq(ih,jh,ipol,nt)

4
PW/compute_qdipol_so.f90
View File

@ -18,7 +18,7 @@ SUBROUTINE compute_qdipol_so(dpqq,dpqq_so)
USE kinds, ONLY : DP
USE ions_base, ONLY : ntyp => nsp
USE lsda_mod, ONLY : nspin
USE uspp_param, ONLY : nh, tvanp, nhm
USE uspp_param, ONLY : upf, nh, nhm
USE spin_orb, ONLY : lspinorb, so, fcoef
!
IMPLICIT NONE
@ -34,7 +34,7 @@ SUBROUTINE compute_qdipol_so(dpqq,dpqq_so)
dpqq_so=(0.d0,0.d0)
DO ipol=1,3
DO nt = 1, ntyp
IF (tvanp (nt) ) THEN
IF ( upf(nt)%tvanp ) THEN
IF (so(nt)) THEN
DO ih=1,nh(nt)
DO jh=1,nh(nt)

6
PW/force_us.f90
View File

@ -21,7 +21,7 @@ SUBROUTINE force_us( forcenl )
USE klist, ONLY : nks, xk, ngk
USE gvect, ONLY : g
USE uspp, ONLY : nkb, vkb, qq, deeq, qq_so, deeq_nc
USE uspp_param, ONLY : nh, tvanp, newpseudo
USE uspp_param, ONLY : upf, nh, newpseudo
USE wvfct, ONLY : nbnd, npw, npwx, igk, wg, et
USE lsda_mod, ONLY : lsda, current_spin, isk
USE symme, ONLY : irt, s, nsym
@ -122,7 +122,7 @@ SUBROUTINE force_us( forcenl )
END DO
END DO
!
IF ( tvanp(nt) .OR. newpseudo(nt) ) THEN
IF ( upf(nt)%tvanp .OR. newpseudo(nt) ) THEN
!
! ... in US case there is a contribution for jh<>ih.
! ... We use here the symmetry in the interchange
@ -306,7 +306,7 @@ SUBROUTINE force_us( forcenl )
END IF
END DO
!
IF ( tvanp(nt) .OR. newpseudo(nt) ) THEN
IF ( upf(nt)%tvanp .OR. newpseudo(nt) ) THEN
!
! ... in US case there is a contribution for jh<>ih.
! ... We use here the symmetry in the interchange

40
PW/grid_paw_routines.f90
View File

@ -167,8 +167,7 @@ CONTAINS
USE lsda_mod, ONLY : nspin
USE us, ONLY : nqxq, dq, nqx, tab, qrad
USE uspp, ONLY : qq, qq_so
USE uspp_param, ONLY : lmaxq, betar, qfunc, qfcoef, rinner, nbeta, &
kkbeta, nqf, nqlc, lll, jjj, lmaxkb, nh, tvanp, nhm, tvanp
USE uspp_param, ONLY : upf, lmaxq, nh, nhm
USE spin_orb, ONLY : lspinorb, rot_ylm, fcoef
!
USE grid_paw_variables, ONLY: tpawp, pfunc, ptfunc, pp, ppt, prad, ptrad, &
@ -236,17 +235,17 @@ CONTAINS
CALL divide (nqxq, startq, lastq)
DO nt = 1, ntyp
IF (tpawp (nt) ) THEN
DO l = 0, nqlc (nt) - 1
DO l = 0, upf(nt)%nqlc - 1
!
! for each nb,mb,l we build the total P_l(|r|) function.
! l is the total (combined) angular momentum and
! the arrays have dimensions 1..l+1
!
DO nb = 1, nbeta (nt)
DO mb = nb, nbeta (nt)
IF ( (l >= ABS (lll(nb,nt) - lll(mb,nt) ) ) .AND. &
(l <= lll(nb,nt) + lll(mb,nt) ) .AND. &
(MOD (l + lll(nb,nt) + lll(mb,nt),2) == 0) ) THEN
DO nb = 1, upf(nt)%nbeta
DO mb = nb, upf(nt)%nbeta
IF ( (l >= ABS (upf(nt)%lll(nb) - upf(nt)%lll(mb) ) ) .AND. &
(l <= upf(nt)%lll(nb) + upf(nt)%lll(mb) ) .AND. &
(MOD (l + upf(nt)%lll(nb) + upf(nt)%lll(mb),2) == 0) ) THEN
qtot(1:msh(nt),nb,mb) = pfunc_(1:msh(nt),nb,mb,nt)
ENDIF
ENDDO ! mb
@ -260,13 +259,13 @@ CONTAINS
!
! and then we integrate with all the Q functions
!
DO nb = 1, nbeta (nt)
DO mb = nb, nbeta (nt)
DO nb = 1, upf(nt)%nbeta
DO mb = nb, upf(nt)%nbeta
! the P are symmetric with respect to indices
nmb = mb * (mb - 1) / 2 + nb
IF ( (l >= ABS( lll(nb,nt) - lll(mb,nt) ) ) .AND. &
(l <= lll(nb,nt) + lll(mb,nt) ) .AND. &
(MOD (l+lll(nb,nt) + lll(mb,nt), 2) == 0) ) THEN
IF ( (l >= ABS( upf(nt)%lll(nb) - upf(nt)%lll(mb) ) ) .AND. &
(l <= upf(nt)%lll(nb) + upf(nt)%lll(mb) ) .AND. &
(MOD (l+upf(nt)%lll(nb) + upf(nt)%lll(mb), 2) == 0) ) THEN
!!kk!! DO ir = 1, kkbeta (nt)
DO ir = 1, msh (nt)
aux1 (ir) = aux (ir) * qtot (ir, nb, mb)
@ -321,8 +320,8 @@ CONTAINS
! Compute the integrals of pfunc*r^2 (not in init_us_1)
DO nt = 1, ntyp
IF (tpawp(nt)) THEN
DO nb = 1, nbeta (nt)
DO mb = nb, nbeta (nt)
DO nb = 1, upf(nt)%nbeta
DO mb = nb, upf(nt)%nbeta
aux2(1:msh(nt)) = pfunc_(1:msh(nt), nb, mb, nt) * &
rgrid(nt)%r(1:msh(nt)) **2
! add augmentation charge if ps
@ -339,9 +338,9 @@ CONTAINS
pmultipole_ = 0.0_DP
DO nt = 1, ntyp
IF (tpawp(nt)) THEN
DO nb = 1, nbeta (nt)
DO mb = nb, nbeta (nt)
DO l = 0, nqlc(nt) - 1
DO nb = 1, upf(nt)%nbeta
DO mb = nb, upf(nt)%nbeta
DO l = 0, upf(nt)%nqlc - 1
aux2(1:msh(nt))= pfunc_(1:msh(nt), nb, mb, nt) * &
rgrid(nt)%r(1:msh(nt))**l
! add augmentation charge if ps
@ -468,7 +467,7 @@ CONTAINS
USE gvect, ONLY : nr1, nr2, nr3, nrx1, nrx2, nrx3, nrxx, &
ngm, nl, nlm, gg, g
USE lsda_mod, ONLY : nspin
USE uspp_param, ONLY : lmaxq, tvanp, nh, nhm
USE uspp_param, ONLY : lmaxq, nh, nhm
USE wvfct, ONLY : gamma_only
USE wavefunctions_module, ONLY : psic
!
@ -1294,8 +1293,7 @@ CONTAINS
g, gg, ngm, gstart, nl
USE lsda_mod, ONLY : nspin
USE scf, ONLY : vr, vltot
USE uspp, ONLY : deeq, dvan, deeq_nc, dvan_so, okvan
USE uspp_param, ONLY : lmaxq, nh, nhm, tvanp
USE uspp_param, ONLY : lmaxq, nh, nhm
USE wvfct, ONLY : gamma_only
USE wavefunctions_module, ONLY : psic
USE spin_orb, ONLY : lspinorb

10
PW/h_epsi_her_apply.f90
View File

@ -27,7 +27,7 @@ subroutine h_epsi_her_apply(lda, n,nbande, psi, hpsi)
USE scf, ONLY : vrs
USE gvect
USE uspp
USE uspp_param
USE uspp_param, ONLY: upf, nh, nhm, nbetam, lmaxq
USE bp
USE basis
USE klist
@ -68,7 +68,7 @@ subroutine h_epsi_her_apply(lda, n,nbande, psi, hpsi)
jkb_bp=0
DO nt=1,ntyp
DO na=1,nat
IF (ityp(na).eq.nt) THEN
IF (ityp(na)== nt) THEN
DO i=1, nh(nt)
jkb_bp=jkb_bp+1
nkbtona(jkb_bp) = na
@ -105,7 +105,7 @@ subroutine h_epsi_her_apply(lda, n,nbande, psi, hpsi)
jkb1 = jkb - nhjkb
DO j = 1,nhjkbm
pref = pref+CONJG(bec_evcel(jkb,mb))*becp0(jkb1+j,nb) &!bec_evcel is relative to ik
*qqq(nhjkb,j,np)
*upf(np)%qqq(nhjkb,j)
ENDDO
ENDDO
sca= sca + pref
@ -143,14 +143,14 @@ subroutine h_epsi_her_apply(lda, n,nbande, psi, hpsi)
ijkb0 = 0
do nt = 1, ntyp
do na = 1, nat
if (ityp (na) .eq.nt) then
if (ityp (na) == nt) then
do ibnd = 1, nbnd
do jh = 1, nh (nt)
jkb = ijkb0 + jh
do ih = 1, nh (nt)
ikb = ijkb0 + ih
ps (ikb, ibnd) = ps (ikb, ibnd) + &
qqq(ih,jh,ityp(na))* bec_evcel(jkb,ibnd)
upf(nt)%qqq(ih,jh)* bec_evcel(jkb,ibnd)
enddo
enddo
enddo

14
PW/h_epsi_her_set.f90
View File

@ -26,7 +26,7 @@ subroutine h_epsi_her_set
USE scf, ONLY : vrs
USE gvect
USE uspp
USE uspp_param
USE uspp_param, ONLY: upf, nh, nhm, nbetam, lmaxq
USE bp, ONLY : gdir,nppstr,efield,fact_hepsi,evcel,evcp=>evcelp,evcm=>evcelm
USE basis
USE klist
@ -269,7 +269,7 @@ subroutine h_epsi_her_set
! --- Form factor: 4 pi sum_LM c_ij^LM Y_LM(Omega) Q_ij^L(|r|) ---
q_dk=(0.d0,0.d0)
DO np =1, ntyp
if(tvanp(np)) then
if( upf(nt)%tvanp ) then
DO iv = 1, nh(np)
DO jv = iv, nh(np)
call qvan3(iv,jv,np,pref,ylm_dk,qrad_dk)
@ -287,7 +287,7 @@ subroutine h_epsi_her_set
q_dkp=(0.d0,0.d0)
DO np =1, ntyp
if(tvanp(np)) then
if( upf(nt)%tvanp ) then
DO iv = 1, nh(np)
DO jv = iv, nh(np)
call qvan3(iv,jv,np,pref,ylm_dk,qrad_dk)
@ -1091,14 +1091,14 @@ subroutine h_epsi_her_set
ijkb0 = 0
do nt = 1, ntyp
do na = 1, nat
if (ityp (na) .eq.nt) then
if (ityp (na) == nt) then
do ibnd = 1, nbnd
do jh = 1, nh (nt)
jkb = ijkb0 + jh
do ih = 1, nh (nt)
ikb = ijkb0 + ih
ps (ikb, ibnd) = ps (ikb, ibnd) + &
qqq(ih,jh,ityp(na))* becp1(jkb,ibnd)
upf(nt)%qqq(ih,jh)* becp1(jkb,ibnd)
enddo
enddo
enddo
@ -1114,14 +1114,14 @@ subroutine h_epsi_her_set
ijkb0 = 0
do nt = 1, ntyp
do na = 1, nat
if (ityp (na) .eq.nt) then
if (ityp (na) == nt) then
do ibnd = 1, nbnd
do jh = 1, nh (nt)
jkb = ijkb0 + jh
do ih = 1, nh (nt)
ikb = ijkb0 + ih
ps (ikb, ibnd) = ps (ikb, ibnd) + &
qqq(ih,jh,ityp(na))* becp1(jkb,ibnd)
upf(nt)%qqq(ih,jh)* becp1(jkb,ibnd)
enddo
enddo
enddo

127
PW/init_us_1.f90
View File

@ -1,5 +1,5 @@
!
! Copyright (C) 2001 PWSCF group
! Copyright (C) 2001-2007 Quantum-Espresso group
! This file is distributed under the terms of the
! GNU General Public License. See the file `License'
! in the root directory of the present distribution,
@ -39,9 +39,7 @@ subroutine init_us_1
USE splinelib
USE uspp, ONLY : nhtol, nhtoj, nhtolm, dvan, qq, indv, ap, aainit, &
qq_so, dvan_so, okvan
USE uspp_param, ONLY : lmaxq, dion, betar, qfunc, qfcoef, rinner, nbeta, &
kkbeta, nqf, nqlc, lll, jjj, lmaxkb, nh, tvanp, &
nbetam, nhm
USE uspp_param, ONLY : upf, lmaxq, nbetam, nh, nhm, lmaxkb
USE spin_orb, ONLY : lspinorb, so, rot_ylm, fcoef
!! NEW-AUG !!
USE grid_paw_variables, ONLY : really_do_paw, okpaw, tpawp, aug
@ -80,7 +78,7 @@ subroutine init_us_1
!
! Initialization of the variables
!
ndm = MAXVAL (kkbeta(1:ntyp))
ndm = MAXVAL ( upf(:)%kkbeta )
allocate (aux ( ndm))
allocate (aux1( ndm))
allocate (besr( ndm))
@ -94,8 +92,8 @@ subroutine init_us_1
! l of the beta functions but includes the l of the local potential
!
do nt=1,ntyp
nqlc(nt) = MIN ( nqlc(nt), lmaxq )
IF ( nqlc(nt) < 0 ) nqlc(nt) = 0
upf(nt)%nqlc = MIN ( upf(nt)%nqlc, lmaxq )
IF ( upf(nt)%nqlc < 0 ) upf(nt)%nqlc = 0
end do
prefr = fpi / omega
@ -132,17 +130,26 @@ subroutine init_us_1
!
do nt = 1, ntyp
ih = 1
do nb = 1, nbeta (nt)
l = lll (nb, nt)
j = jjj (nb, nt)
do nb = 1, upf(nt)%nbeta
l = upf(nt)%lll (nb)
do m = 1, 2 * l + 1
nhtol (ih, nt) = l
nhtolm(ih, nt) = l*l+m
nhtoj (ih, nt) = j
indv (ih, nt) = nb
ih = ih + 1
enddo
enddo
if ( so(nt) ) then
ih = 1
do nb = 1, upf(nt)%nbeta
l = upf(nt)%lll (nb)
j = upf(nt)%jjj (nb)
do m = 1, 2 * l + 1
nhtoj (ih, nt) = j
ih = ih + 1
enddo
enddo
endif
!
! From now on the only difference between KB and US pseudopotentials
! is in the presence of the q and Q functions.
@ -163,7 +170,7 @@ subroutine init_us_1
jk = nhtoj(kh, nt)
mk = nhtolm(kh, nt)-lk*lk
vk = indv (kh, nt)
if (li.eq.lk.and.abs(ji-jk).lt.1.d-7) then
if (li == lk .and. abs(ji-jk) < 1.d-7) then
do is1=1,2
do is2=1,2
coeff = (0.d0, 0.d0)
@ -190,7 +197,7 @@ subroutine init_us_1
do is1=1,2
do is2=1,2
ijs=ijs+1
dvan_so(ih,jh,ijs,nt) = dion(vi,vj,nt) * &
dvan_so(ih,jh,ijs,nt) = upf(nt)%dion(vi,vj) * &
fcoef(ih,jh,is1,is2,nt)
if (vi.ne.vj) fcoef(ih,jh,is1,is2,nt)=(0.d0,0.d0)
enddo
@ -205,10 +212,10 @@ subroutine init_us_1
ir = indv (ih, nt)
is = indv (jh, nt)
if (lspinorb) then
dvan_so (ih, jh, 1, nt) = dion (ir, is, nt)
dvan_so (ih, jh, 4, nt) = dion (ir, is, nt)
dvan_so (ih, jh, 1, nt) = upf(nt)%dion (ir, is)
dvan_so (ih, jh, 4, nt) = upf(nt)%dion (ir, is)
else
dvan (ih, jh, nt) = dion (ir, is, nt)
dvan (ih, jh, nt) = upf(nt)%dion (ir, is)
endif
endif
enddo
@ -224,35 +231,37 @@ subroutine init_us_1
! Q functions.
!
call divide (nqxq, startq, lastq)
!
do nt = 1, ntyp
!
if (tvanp (nt) ) then
do l = 0, nqlc (nt) - 1
if ( upf(nt)%tvanp ) then
do l = 0, upf(nt)%nqlc - 1
!
! first we build for each nb,mb,l the total Q(|r|) function
! note that l is the true (combined) angular momentum
! and that the arrays have dimensions 0..l (no more 1..l+1)
!
do nb = 1, nbeta (nt)
do mb = nb, nbeta (nt)
do nb = 1, upf(nt)%nbeta
do mb = nb, upf(nt)%nbeta
ijv = mb * (mb-1) / 2 + nb
paw:& ! in PAW formalism aug. charge is computed elsewhere
if (tpawp(nt)) then
qtot(1:kkbeta(nt),ijv) = aug(nt)%fun(1:kkbeta(nt),nb,mb,l)
qtot(1:upf(nt)%kkbeta,ijv) = aug(nt)%fun(1:upf(nt)%kkbeta,nb,mb,l)
else
if ( (l >= abs (lll (nb, nt) - lll (mb, nt) ) ) .and. &
(l <= lll (nb, nt) + lll (mb, nt) ) .and. &
(mod (l + lll (nb, nt) + lll (mb, nt), 2) == 0) ) then
do ir = 1, kkbeta (nt)
if (rgrid(nt)%r(ir) >= rinner (l + 1, nt) ) then
qtot (ir, ijv) = qfunc (ir, ijv, nt)
if ( ( l >= abs(upf(nt)%lll(nb) - upf(nt)%lll(mb)) ) .and. &
( l <= upf(nt)%lll(nb) + upf(nt)%lll(mb) ) .and. &
(mod (l+upf(nt)%lll(nb)+upf(nt)%lll(mb), 2) == 0) ) then
do ir = 1, upf(nt)%kkbeta
if (rgrid(nt)%r(ir) >=upf(nt)%rinner (l+1) ) then
! qtot (ir, ijv) = qfunc (ir, ijv, nt) TEMP
qtot (ir, ijv) = upf(nt)%qfunc(ir,nb,mb)
else
ilast = ir
endif
enddo
if (rinner (l + 1, nt) > 0.d0) &
call setqf(qfcoef (1, l+1, nb, mb, nt), &
qtot(1,ijv), rgrid(nt)%r, nqf(nt),l,ilast)
if ( upf(nt)%rinner (l+1) > 0.0_dp) &
call setqf(upf(nt)%qfcoef (1, l+1, nb, mb), &
qtot(1,ijv), rgrid(nt)%r, upf(nt)%nqf, &
l, ilast)
endif
endif paw
enddo
@ -262,24 +271,24 @@ subroutine init_us_1
!
do iq = startq, lastq
q = (iq - 1) * dq * tpiba
call sph_bes (kkbeta (nt), rgrid(nt)%r, q, l, aux)
call sph_bes ( upf(nt)%kkbeta, rgrid(nt)%r, q, l, aux)
!
! and then we integrate with all the Q functions
!
do nb = 1, nbeta (nt)
do nb = 1, upf(nt)%nbeta
!
! the Q are symmetric with respect to indices
!
do mb = nb, nbeta (nt)
do mb = nb, upf(nt)%nbeta
ijv = mb * (mb - 1) / 2 + nb
if ( (l >= abs (lll (nb, nt) - lll (mb, nt) ) ) .and. &
(l <= lll (nb, nt) + lll (mb, nt) ) .and. &
(mod (l + lll(nb, nt) + lll(mb, nt), 2) == 0) &
.or. tpawp(nt) ) then
do ir = 1, kkbeta (nt)
if ( ( l >= abs(upf(nt)%lll(nb) - upf(nt)%lll(mb)) ) .and. &
( l <= upf(nt)%lll(nb) + upf(nt)%lll(mb) ) .and. &
(mod (l+upf(nt)%lll(nb)+upf(nt)%lll(mb), 2) == 0) .or.&
tpawp(nt) ) then
do ir = 1, upf(nt)%kkbeta
aux1 (ir) = aux (ir) * qtot (ir, ijv)
enddo
call simpson (kkbeta(nt), aux1, rgrid(nt)%rab, &
call simpson ( upf(nt)%kkbeta, aux1, rgrid(nt)%rab, &
qrad(iq,ijv,l + 1, nt) )
endif
enddo
@ -305,7 +314,7 @@ subroutine init_us_1
#endif
call ylmr2 (lmaxq * lmaxq, 1, g, gg, ylmk0)
do nt = 1, ntyp
if (tvanp (nt) ) then
if ( upf(nt)%tvanp ) then
if (so(nt)) then
do ih=1,nh(nt)
do jh=1,nh(nt)
@ -361,39 +370,37 @@ subroutine init_us_1
call divide (nqx, startq, lastq)
tab (:,:,:) = 0.d0
do nt = 1, ntyp
do nb = 1, nbeta (nt)
l = lll (nb, nt)
do nb = 1, upf(nt)%nbeta
l = upf(nt)%lll (nb)
do iq = startq, lastq
qi = (iq - 1) * dq
call sph_bes (kkbeta (nt), rgrid(nt)%r, qi, l, besr)
do ir = 1, kkbeta (nt)
aux (ir) = betar (ir, nb, nt) * besr (ir) * rgrid(nt)%r(ir)
call sph_bes (upf(nt)%kkbeta, rgrid(nt)%r, qi, l, besr)
do ir = 1, upf(nt)%kkbeta
aux (ir) = upf(nt)%beta (ir, nb) * besr (ir) * rgrid(nt)%r(ir)
enddo
call simpson (kkbeta (nt), aux, rgrid(nt)%rab, vqint)
call simpson (upf(nt)%kkbeta, aux, rgrid(nt)%rab, vqint)
tab (iq, nb, nt) = vqint * pref
enddo
enddo
enddo
#ifdef __PARA
call reduce (nqx * nbetam * ntyp, tab)
#endif
! initialize spline interpolation
if (spline_ps) then
allocate( xdata(nqx) )
do iq = 1, nqx
xdata(iq) = (iq - 1) * dq
enddo
do nt = 1, ntyp
do nb = 1, nbeta (nt)
l = lll (nb, nt)
d1 = (tab(2,nb,nt) - tab(1,nb,nt)) / dq
call spline(xdata, tab(:,nb,nt), 0.d0, d1, tab_d2y(:,nb,nt))
allocate( xdata(nqx) )
do iq = 1, nqx
xdata(iq) = (iq - 1) * dq
enddo
enddo
deallocate(xdata)
do nt = 1, ntyp
do nb = 1, upf(nt)%nbeta
d1 = (tab(2,nb,nt) - tab(1,nb,nt)) / dq
call spline(xdata, tab(:,nb,nt), 0.d0, d1, tab_d2y(:,nb,nt))
enddo
enddo
deallocate(xdata)
endif
deallocate (ylmk0)

4
PW/init_us_2.f90
View File

@ -23,7 +23,7 @@ subroutine init_us_2 (npw_, igk_, q_, vkb_)
USE us, ONLY : nqx, dq, tab, tab_d2y, spline_ps
USE splinelib
USE uspp, ONLY : nkb, vkb, nhtol, nhtolm, indv
USE uspp_param, ONLY : lmaxkb, nbeta, nhm, nh
USE uspp_param, ONLY : upf, lmaxkb, nhm, nh
!
implicit none
!
@ -84,7 +84,7 @@ subroutine init_us_2 (npw_, igk_, q_, vkb_)
jkb = 0
do nt = 1, ntyp
! calculate beta in G-space using an interpolation table
do nb = 1, nbeta (nt)
do nb = 1, upf(nt)%nbeta
do ig = 1, npw_
if (spline_ps) then
vq(ig) = splint(xdata, tab(:,nb,nt), tab_d2y(:,nb,nt), qg(ig))

6
PW/init_vloc.f90
View File

@ -15,7 +15,7 @@ subroutine init_vloc()
! potential vloc(ig,it) for each type of atom
!
USE atom, ONLY : numeric, msh, rgrid
USE uspp_param, ONLY : vloc_at, zp
USE uspp_param, ONLY : upf
USE ions_base, ONLY : ntyp => nsp
USE cell_base, ONLY : omega, tpiba2
USE vlocal, ONLY : vloc
@ -35,8 +35,8 @@ subroutine init_vloc()
! compute V_loc(G) for a given type of atom
!
call vloc_of_g (lloc (nt), lmax (nt), numeric (nt), rgrid(nt)%mesh, &
msh (nt), rgrid(nt)%rab, rgrid(nt)%r, vloc_at (1, nt), cc (1, &
nt), alpc (1, nt), nlc (nt), nnl (nt), zp (nt), aps (1, 0, nt), &
msh (nt), rgrid(nt)%rab, rgrid(nt)%r, upf(nt)%vloc(1), cc (1, &
nt), alpc (1, nt), nlc (nt), nnl (nt), upf(nt)%zp, aps (1, 0, nt), &
alps (1, 0, nt), tpiba2, ngl, gl, omega, vloc (1, nt) )
enddo

2
PW/make.depend
View File

@ -783,7 +783,6 @@ paw.o : ../Modules/pseudo_types.o
paw.o : ../Modules/radial_grids.o
paw.o : ../Modules/read_upf.o
paw.o : ../Modules/splinelib.o
paw.o : ../Modules/uspp.o
paw_xc.o : ../Modules/cell_base.o
paw_xc.o : ../Modules/constants.o
paw_xc.o : ../Modules/functionals.o
@ -1011,7 +1010,6 @@ restart_in_ions.o : ../Modules/io_global.o
restart_in_ions.o : ../Modules/ions_base.o
restart_in_ions.o : ../Modules/kind.o
restart_in_ions.o : ../Modules/wavefunctions.o
restart_in_ions.o : noncol.o
restart_in_ions.o : pwcom.o
rgen.o : ../Modules/kind.o
rho2zeta.o : ../Modules/constants.o

4
PW/newd.f90
View File

@ -33,7 +33,7 @@ SUBROUTINE newd_g()
USE lsda_mod, ONLY : nspin
USE scf, ONLY : vr, vltot
USE uspp, ONLY : deeq, dvan, deeq_nc, dvan_so, okvan, indv
USE uspp_param, ONLY : lmaxq, nh, nhm, tvanp
USE uspp_param, ONLY : upf, lmaxq, nh, nhm
USE wvfct, ONLY : gamma_only
USE wavefunctions_module, ONLY : psic
USE spin_orb, ONLY : lspinorb, so, domag
@ -143,7 +143,7 @@ SUBROUTINE newd_g()
!
DO nt = 1, ntyp
!
IF ( tvanp(nt) ) THEN
IF ( upf(nt)%tvanp ) THEN
!
DO ih = 1, nh(nt)
!

1
PW/paw.f90
View File

@ -228,7 +228,6 @@ CONTAINS
use parameters, only : ntypx
USE io_global, ONLY : stdout
use splinelib
USE uspp_param, ONLY : vloc_at
IMPLICIT NONE

2
PW/rad_paw_trash.f90
View File

@ -33,7 +33,7 @@ SUBROUTINE integrate_pfunc
USE lsda_mod, ONLY : nspin
USE us, ONLY : nqxq, dq, nqx, tab, qrad
USE uspp
USE uspp_param
USE uspp_param
USE spin_orb, ONLY : lspinorb, rot_ylm, fcoef
!
USE grid_paw_variables, ONLY: tpawp, pfunc, ptfunc, pp, ppt, prad, ptrad, okpaw

10
PW/read_file.f90
View File

@ -36,7 +36,7 @@ SUBROUTINE read_file()
USE vlocal, ONLY : strf
USE io_files, ONLY : tmp_dir, prefix, iunpun, nwordwfc, iunwfc
USE buffers, ONLY : open_buffer, close_buffer
USE uspp_param, ONLY : nbeta, jjj, tvanp
USE uspp_param, ONLY : upf, nbeta, jjj
USE noncollin_module, ONLY : noncolin, npol
USE mp_global, ONLY : kunit
USE pw_restart, ONLY : pw_readfile
@ -149,17 +149,17 @@ SUBROUTINE read_file()
!
CALL readpp()
!
okvan = ANY ( tvanp(1:nsp) )
okvan = ANY ( upf(:)%tvanp )
!
! ... check for spin-orbit pseudopotentials
!
DO nt = 1, nsp
!
so(nt) = ( nbeta(nt) > 0 )
so(nt) = ( upf(nt)%nbeta > 0 )
!
DO nb = 1, nbeta(nt)
DO nb = 1, upf(nt)%nbeta
!
so(nt) = so(nt) .AND. ( ABS( jjj(nb,nt) ) > 1.D-7 )
so(nt) = so(nt) .AND. ( ABS( upf(nt)%jjj(nb) ) > 1.D-7 )
!
END DO
!

19
PW/read_ncpp.f90
View File

@ -146,19 +146,20 @@ subroutine read_ncpp (iunps, np, upf)
upf%r(ir) = exp (x) / upf%zmesh
upf%rab(ir) = upf%dx * upf%r(ir)
enddo
ALLOCATE ( upf%kkbeta(upf%nbeta) )
do ir = 1, upf%mesh
if ( upf%r(ir) > rcut) then
upf%kkbeta(:) = ir
upf%kkbeta = ir
go to 5
end if
end do
upf%kkbeta(:) = upf%mesh
upf%kkbeta = upf%mesh
!
! ... force kkbeta to be odd for simpson integration (obsolete?)
!
5 upf%kkbeta(:) = 2 * ( ( upf%kkbeta(:) + 1 ) / 2) - 1
5 upf%kkbeta = 2 * ( ( upf%kkbeta + 1 ) / 2) - 1
!
ALLOCATE ( upf%kbeta(upf%nbeta) )
upf%kbeta(:) = upf%kkbeta
ALLOCATE ( upf%vloc(upf%mesh) )
upf%vloc (:) = 0.d0
if (.not. numeric(np)) then
@ -169,7 +170,7 @@ subroutine read_ncpp (iunps, np, upf)
CALL bachel( alps(1,0,np), aps(1,0,np), 1, lmax(np) )
!
do i = 1, nlc (np)
do ir = 1, upf%kkbeta(1)
do ir = 1, upf%kkbeta
upf%vloc (ir) = upf%vloc (ir) - upf%zp * e2 * cc (i, np) * &
erf ( sqrt (alpc(i,np)) * upf%r(ir) ) / upf%r(ir)
end do
@ -217,14 +218,14 @@ subroutine read_ncpp (iunps, np, upf)
do l = 0, lmax (np)
if (l /= lloc (np) ) then
nb = nb + 1
! betar is used here as work space
do ir = 1, upf%kkbeta(1)
! upf%beta is used here as work space
do ir = 1, upf%kkbeta
upf%beta (ir, nb) = upf%chi(ir, l+1) **2 * vnl(ir, l)
end do
call simpson (upf%kkbeta (1), upf%beta (1, nb), upf%rab, vll )
call simpson (upf%kkbeta, upf%beta (1, nb), upf%rab, vll )
upf%dion (nb, nb) = 1.d0 / vll
! upf%beta stores projectors |beta(r)> = |V_nl(r)phi(r)>
do ir = 1, upf%kkbeta (1)
do ir = 1, upf%kkbeta
upf%beta (ir, nb) = vnl (ir, l) * upf%chi (ir, l + 1)
enddo
upf%lll (nb) = l

38
PW/read_pseudo.f90
View File

@ -20,14 +20,14 @@ subroutine readpp
USE upf_to_internal, ONLY : set_pseudo_upf
USE paw, ONLY : set_paw_upf
USE atom, ONLY : chi, nchi, oc, msh, numeric, rgrid
USE uspp_param, ONLY : zp, iver, tvanp, newpseudo
USE uspp_param, ONLY : iver, newpseudo
USE ions_base, ONLY : ntyp => nsp
USE funct, ONLY : get_iexch, get_icorr, get_igcx, get_igcc
USE io_files, ONLY : pseudo_dir, psfile
USE io_global, ONLY : stdout
USE ions_base, ONLY : zv
USE pseud, ONLY : lmax, lloc
USE uspp_param, ONLY : lll, nbeta
USE uspp_param, ONLY : upf
USE parameters, ONLY : nchix !PAW
USE grid_paw_variables, ONLY : tpawp
USE read_paw_module, ONLY : paw_io, allocate_pseudo_paw, deallocate_pseudo_paw
@ -36,7 +36,6 @@ subroutine readpp
!
real(DP), parameter :: rcut = 10.d0, eps = 1.0D-08
!
TYPE (pseudo_upf) :: upf
TYPE(paw_t) :: pawset
!
character(len=256) :: file_pseudo
@ -49,6 +48,7 @@ subroutine readpp
!
iunps = 4
l = len_trim (pseudo_dir)
ALLOCATE ( upf(ntyp) )
do nt = 1, ntyp
tpawp(nt) = .false.
!
@ -79,15 +79,13 @@ subroutine readpp
! read UPF pseudopotentials - the UPF format is detected via the
! presence of the keyword '<PP_HEADER>' at the beginning of the file
!
call read_pseudo_upf(iunps, upf, isupf)
call read_pseudo_upf(iunps, upf(nt), isupf)
!
if (isupf == 0) then
call set_pseudo_upf (nt, upf)
call set_paw_upf (nt, upf)
CALL deallocate_pseudo_upf( upf )
call set_pseudo_upf (nt, upf(nt)) ! TEMP
call set_paw_upf (nt, upf(nt))
! for compatibility with old formats
newpseudo (nt) = .true.
lmax(nt) = max ( lmax(nt), MAXVAL( lll( 1:nbeta(nt), nt) ) )
!
else
rewind (unit = iunps)
@ -107,14 +105,11 @@ subroutine readpp
newpseudo (nt) = ( pseudo_type (psfile (nt) ) == 2 )
!
IF ( newpseudo (nt) ) THEN
call readrrkj (iunps, nt, upf)
call readrrkj (iunps, nt, upf(nt))
ELSE
CALL readvan (iunps, nt, upf)
CALL readvan (iunps, nt, upf(nt))
ENDIF
CALL set_pseudo_upf (nt, upf) ! TEMP
CALL deallocate_pseudo_upf( upf )! TEMP
!
lmax(nt) = max ( lmax(nt), MAXVAL( lll( 1:nbeta(nt), nt) ) )
CALL set_pseudo_upf (nt, upf(nt)) ! TEMP
!
else if (pseudo_type (psfile (nt) ) ==3) then
!
@ -123,7 +118,6 @@ subroutine readpp
!tpaw(nt)=.true.
numeric (nt) = .true.
newpseudo (nt) = .true.
tvanp (nt) = .true.
open (unit = iunps, file = file_pseudo, status = 'old', &
form='formatted', iostat = ios)
call paw_io (pawset, iunps, "INP") !,ndmx,nchix,lmaxx)
@ -132,23 +126,19 @@ subroutine readpp
call deallocate_pseudo_paw (pawset)
!
else
tvanp (nt) = .false.
newpseudo (nt) = .false.
!
call read_ncpp (iunps, nt, upf)
call read_ncpp (iunps, nt, upf(nt))
!
CALL set_pseudo_upf (nt, upf) ! TEMP
CALL deallocate_pseudo_upf( upf )! TEMP
CALL set_pseudo_upf (nt, upf(nt)) ! TEMP
!
endif
! for compatibility with old formats - maybe obsolete?
lmax(nt) = max ( lmax(nt), &
MAXVAL( upf(nt)%lll( 1:upf(nt)%nbeta ) ) )
endif
close (iunps)
!
! ... Zv = valence charge of the (pseudo-)atom, read from PP files,
! ... is set equal to Zp = pseudo-charge of the pseudopotential
!
zv(nt) = zp(nt)
!
if (nt == 1) then
iexch_ = get_iexch()
icorr_ = get_icorr()

94
PW/realus.f90
View File

@ -1,5 +1,5 @@
!
! Copyright (C) 2001-2006 Quantum-ESPRESSO group
! Copyright (C) 2001-2007 Quantum-ESPRESSO group
! This file is distributed under the terms of the
! GNU General Public License. See the file `License'
! in the root directory of the present distribution,
@ -63,9 +63,7 @@ MODULE realus
USE cell_base, ONLY : at, bg, omega, alat
USE gvect, ONLY : nr1, nr2, nr3, nrx1, nrx2, nrx3, nrxx
USE uspp, ONLY : okvan, indv, nhtol, nhtolm, ap, nhtoj, lpx, lpl
USE uspp_param, ONLY : lmaxq, nh, nhm, tvanp, kkbeta, nbeta, &
qfunc, dion, lmaxkb, qfcoef, nqf, nqlc, &
lll, rinner
USE uspp_param, ONLY : upf, lmaxq, nh, nhm
USE atom, ONLY : rgrid
USE pfft, ONLY : npp
USE mp_global, ONLY : me_pool
@ -73,7 +71,7 @@ MODULE realus
!
IMPLICIT NONE
!
INTEGER :: qsdim, ia, it, mbia, iqs, iqsia
INTEGER :: qsdim, ia, mbia, iqs, iqsia
INTEGER :: indm, inbrx, idimension, &
ih, jh, ijh, lllnbnt, lllmbnt
INTEGER :: roughestimate, goodestimate, lamx2, l, nt
@ -107,13 +105,16 @@ MODULE realus
!
boxrad(:) = 0.D0
!
DO it = 1, nsp
DO inbrx = 1, nbeta(it)*(nbeta(it)+1)/2
DO indm = kkbeta(it), 1, -1
DO nt = 1, nsp
!DO inbrx = 1, upf(nt)%nbeta*(upf(nt)%nbeta+1)/2 TEMP
DO nb = 1, upf(nt)%nbeta
DO mb = nb, upf(nt)%nbeta
DO indm = upf(nt)%kkbeta, 1, -1
!
IF ( ABS( qfunc(indm,inbrx,it) ) > eps16 ) THEN
! IF ( ABS( qfunc(indm,inbrx,nt) ) > eps16 ) THEN TEMP
IF ( ABS( upf(nt)%qfunc(indm,nb,mb) ) > eps16 ) THEN
!
boxrad(it) = MAX( rgrid(it)%r(indm), boxrad(it) )
boxrad(nt) = MAX( rgrid(nt)%r(indm), boxrad(nt) )
!
CYCLE
!
@ -121,6 +122,7 @@ MODULE realus
!
END DO
END DO
END DO
END DO
!
boxrad(:) = boxrad(:) / alat
@ -174,9 +176,11 @@ MODULE realus
!
DO ia = 1, nat
!
IF ( .NOT. tvanp(ityp(ia)) ) CYCLE
nt = ityp(ia)
!
boxradsq_ia = boxrad(ityp(ia))**2
IF ( .NOT. upf(nt)%tvanp ) CYCLE
!
boxradsq_ia = boxrad(nt)**2
!
tau_ia(1) = tau(1,ia)
tau_ia(2) = tau(2,ia)
@ -258,7 +262,9 @@ MODULE realus
!
DO ia = 1, nat
!
IF ( .NOT. tvanp(ityp(ia)) ) CYCLE
nt = ityp(ia)
!
IF ( .NOT. upf(nt)%tvanp ) CYCLE
!
idimension = maxbox(ia)
!
@ -294,7 +300,7 @@ MODULE realus
mbia = maxbox(ia)
IF ( mbia == 0 ) CYCLE
nt = ityp(ia)
IF ( .NOT. tvanp(nt) ) CYCLE
IF ( .NOT. upf(nt)%tvanp ) CYCLE
DO ih = 1, nh(nt)
DO jh = ih, nh(nt)
qsdim = qsdim + mbia
@ -330,60 +336,63 @@ MODULE realus
!
nt = ityp(ia)
!
IF ( .NOT. tvanp(nt) ) CYCLE
IF ( .NOT. upf(nt)%tvanp ) CYCLE
!
ALLOCATE( qtot( kkbeta(nt), nbeta(nt), nbeta(nt) ) )
ALLOCATE( qtot( upf(nt)%kkbeta, upf(nt)%nbeta, upf(nt)%nbeta ) )
!
! ... variables used for spline interpolation
!
ALLOCATE( xsp( kkbeta(nt) ), ysp( kkbeta(nt) ), wsp( kkbeta(nt ) ) )
ALLOCATE( xsp( upf(nt)%kkbeta ), ysp( upf(nt)%kkbeta ), &
wsp( upf(nt)%kkbeta ) )
!
! ... the radii in x
!
xsp(:) = rgrid(nt)%r(1:kkbeta(nt))
xsp(:) = rgrid(nt)%r(1:upf(nt)%kkbeta)
!
DO l = 0, nqlc(nt) - 1
DO l = 0, upf(nt)%nqlc - 1
!
! ... first we build for each nb,mb,l the total Q(|r|) function
! ... note that l is the true (combined) angular momentum
! ... and that the arrays have dimensions 1..l+1
!
DO nb = 1, nbeta(nt)
DO mb = nb, nbeta(nt)
DO nb = 1, upf(nt)%nbeta
DO mb = nb, upf(nt)%nbeta
ijv = mb * (mb-1) /2 + nb
!
lllnbnt = lll(nb,nt)
lllmbnt = lll(mb,nt)
lllnbnt = upf(nt)%lll(nb)
lllmbnt = upf(nt)%lll(mb)
!
IF ( .NOT. ( l >= ABS( lllnbnt - lllmbnt ) .AND. &
l <= lllnbnt + lllmbnt .AND. &
MOD( l + lllnbnt + lllmbnt, 2 ) == 0 ) ) CYCLE
!
DO ir = 1, kkbeta(nt)
IF ( rgrid(nt)%r(ir) >= rinner(l+1,nt) ) THEN
qtot(ir,nb,mb) = qfunc(ir,ijv,nt) / rgrid(nt)%r(ir)**2
DO ir = 1, upf(nt)%kkbeta
IF ( rgrid(nt)%r(ir) >= upf(nt)%rinner(l+1) ) THEN
!qtot(ir,nb,mb) = qfunc(ir,ijv,nt) / rgrid(nt)%r(ir)**2
qtot(ir,nb,mb) = upf(nt)%qfunc(ir,nb,mb) / &
rgrid(nt)%r(ir)**2
ELSE
ilast = ir
END IF
END DO
!
IF ( rinner(l+1,nt) > 0.D0 ) &
CALL setqfcorr( qfcoef(1,l+1,nb,mb,nt), &
qtot(1,nb,mb), rgrid(nt)%r(1), nqf(nt), l, ilast )
IF ( upf(nt)%rinner(l+1) > 0.D0 ) &
CALL setqfcorr( upf(nt)%qfcoef(1,l+1,nb,mb), &
qtot(1,nb,mb), rgrid(nt)%r(1), upf(nt)%nqf, l, ilast )
!
! ... we save the values in y
!
ysp(:) = qtot(1:kkbeta(nt),nb,mb)
ysp(:) = qtot(1:upf(nt)%kkbeta,nb,mb)
!
! ... compute the first derivative in first point
!
CALL setqfcorrptfirst( qfcoef(1,l+1,nb,mb,nt), &
first, rgrid(nt)%r(1), nqf(nt), l )
CALL setqfcorrptfirst( upf(nt)%qfcoef(1,l+1,nb,mb), &
first, rgrid(nt)%r(1), upf(nt)%nqf, l )
!
! ... compute the second derivative in second point
!
CALL setqfcorrptsecond( qfcoef(1,l+1,nb,mb,nt), &
second, rgrid(nt)%r(1), nqf(nt), l )
CALL setqfcorrptsecond( upf(nt)%qfcoef(1,l+1,nb,mb), &
second, rgrid(nt)%r(1), upf(nt)%nqf, l )
!
! ... call spline
!
@ -391,13 +400,13 @@ MODULE realus
!
DO ir = 1, maxbox(ia)
!
IF ( boxdist(ir,ia) < rinner(l+1,nt) ) THEN
IF ( boxdist(ir,ia) < upf(nt)%rinner(l+1) ) THEN
!
! ... if in the inner radius just compute the
! ... polynomial
!
CALL setqfcorrpt( qfcoef(1,l+1,nb,mb,nt), &
qtot_int, boxdist(ir,ia), nqf(nt), l )
CALL setqfcorrpt( upf(nt)%qfcoef(1,l+1,nb,mb), &
qtot_int, boxdist(ir,ia), upf(nt)%nqf, l )
!
ELSE
!
@ -461,9 +470,8 @@ MODULE realus
USE gvect, ONLY : nr1, nr2, nr3, nrxx
USE lsda_mod, ONLY : nspin
USE scf, ONLY : vr, vltot
USE uspp, ONLY : okvan
USE uspp, ONLY : deeq, deeq_nc, dvan, dvan_so
USE uspp_param, ONLY : nh, nhm, tvanp
USE uspp, ONLY : okvan, deeq, deeq_nc, dvan, dvan_so
USE uspp_param, ONLY : upf, nh, nhm
USE noncollin_module, ONLY : noncolin
USE spin_orb, ONLY : so, domag, lspinorb
!
@ -539,7 +547,7 @@ MODULE realus
!
nt = ityp(ia)
!
IF ( .NOT. tvanp(nt) ) CYCLE
IF ( .NOT. upf(nt)%tvanp ) CYCLE
!
nhnt = nh(nt)
!
@ -871,7 +879,7 @@ MODULE realus
USE klist, ONLY : nelec
USE gvect, ONLY : nr1, nr2, nr3
USE uspp, ONLY : okvan, becsum
USE uspp_param, ONLY : tvanp, nh
USE uspp_param, ONLY : upf, nh
USE noncollin_module, ONLY : noncolin
USE spin_orb, ONLY : domag
!
@ -901,7 +909,7 @@ MODULE realus
!
nt = ityp(ia)
!
IF ( .NOT. tvanp(nt) ) CYCLE
IF ( .NOT. upf(nt)%tvanp ) CYCLE
!
nhnt = nh(nt)
!

6
PW/s_psi.f90
View File

@ -30,7 +30,7 @@ SUBROUTINE s_psi( lda, n, m, psi, spsi )
USE kinds, ONLY : DP
USE wvfct, ONLY : gamma_only
USE uspp, ONLY : vkb, nkb, qq, okvan
USE uspp_param, ONLY : nh, tvanp
USE uspp_param, ONLY : upf, nh
USE wvfct, ONLY : igk, g2kin
USE gsmooth, ONLY : nls, nr1s, nr2s, nr3s, nrx1s, nrx2s, nrx3s, nrxxs
USE ldaU, ONLY : lda_plus_u
@ -93,7 +93,7 @@ SUBROUTINE s_psi( lda, n, m, psi, spsi )
!
ijkb0 = 0
DO nt = 1, nsp
IF ( tvanp (nt) ) THEN
IF ( upf(nt)%tvanp ) THEN
DO na = 1, nat
IF ( ityp(na) == nt ) THEN
DO ibnd = 1, m
@ -166,7 +166,7 @@ SUBROUTINE s_psi( lda, n, m, psi, spsi )
!
ijkb0 = 0
DO nt = 1, nsp
IF ( tvanp(nt) ) THEN
IF ( upf(nt)%tvanp ) THEN
DO na = 1, nat
IF ( ityp(na) == nt ) THEN
DO ibnd = 1, m

4
PW/s_psi_nc.f90
View File

@ -24,7 +24,7 @@ subroutine s_psi_nc (lda, n, m, psi, spsi )
! spsi S*psi
!
USE ions_base, ONLY: nat, ityp, ntyp => nsp
USE uspp_param, ONLY: nh, tvanp
USE uspp_param, ONLY: upf, nh
USE uspp, ONLY: nkb, vkb, qq, qq_so, okvan
use wvfct, ONLY: igk, g2kin
use gsmooth, ONLY: nls, nr1s, nr2s, nr3s, nrx1s, nrx2s, nrx3s, nrxxs
@ -61,7 +61,7 @@ subroutine s_psi_nc (lda, n, m, psi, spsi )
!
ijkb0 = 0
do nt = 1, ntyp
if (tvanp (nt) ) then
if ( upf(nt)%tvanp ) then
do na = 1, nat
if (ityp (na) .eq.nt) then
do ih = 1,nh(nt)

34
PW/setup.f90
View File

@ -68,7 +68,7 @@ SUBROUTINE setup()
USE control_flags, ONLY : para_diago_dim ! debug
USE relax, ONLY : starting_diag_threshold
USE cellmd, ONLY : calc
USE uspp_param, ONLY : zp, psd, nbeta, jjj, tvanp
USE uspp_param, ONLY : upf
USE uspp, ONLY : okvan
USE ldaU, ONLY : d1, d2, d3, lda_plus_u, Hubbard_U, &
Hubbard_l, Hubbard_alpha, Hubbard_lmax
@ -98,12 +98,8 @@ SUBROUTINE setup()
LOGICAL, EXTERNAL :: lchk_tauxk ! tests that atomic coordinates do not overlap
!
!
IF (dft_is_meta()) THEN
DO nt=1,ntyp
IF ( tvanp(nt) ) &
CALL errore( 'setup', 'US and Meta-GGA not yet implemented', 1 )
END DO
END IF
IF ( dft_is_meta() .AND. ANY ( upf(:)%tvanp ) ) &
CALL errore( 'setup', 'US and Meta-GGA not yet implemented', 1 )
ALLOCATE( m_loc( 3, nat ) )
!
@ -112,7 +108,7 @@ SUBROUTINE setup()
!
! ... Compute the ionic charge for each atom type
!
zv(1:ntyp) = zp(1:ntyp)
zv(1:ntyp) = upf(1:ntyp)%zp
!
#if defined (__PGI)
ionic_charge = - tot_charge
@ -406,18 +402,12 @@ SUBROUTINE setup()
! ... if this is not a spin-orbit calculation, all spin-orbit pseudopotentials
! ... are transformed into standard pseudopotentials
!
IF ( lspinorb .AND. ALL ( .NOT. upf(:)%has_so ) ) &
CALL infomsg ('setup','At least one non s.o. pseudo')
!
DO nt = 1, ntyp
!
so(nt) = ( nbeta(nt) > 0 )
!
IF ( lspinorb .AND. ALL ( ABS( jjj(1:nbeta(nt),nt) ) < 1.D-7 ) ) &
CALL infomsg ('setup','At least one non s.o. pseudo')
!
DO nb = 1, nbeta(nt)
!
so(nt) = so(nt) .AND. ( ABS( jjj(nb,nt) ) > 1.D-7 )
!
END DO
so(nt) = upf(nt)%has_so
!
END DO
!
@ -700,7 +690,7 @@ SUBROUTINE setup()
!
! ... okvan = .TRUE. : at least one pseudopotential is US
!
okvan = ANY( tvanp(1:ntyp) )
okvan = ANY( upf(:)%tvanp )
okpaw = ANY( tpawp(1:ntyp) )
!
! ... Needed for LDA+U
@ -715,12 +705,12 @@ SUBROUTINE setup()
!
IF ( Hubbard_U(nt) /= 0.D0 .OR. Hubbard_alpha(nt) /= 0.D0 ) THEN
!
Hubbard_l(nt) = set_Hubbard_l( psd(nt) )
Hubbard_l(nt) = set_Hubbard_l( upf(nt)%psd )
!
Hubbard_lmax = MAX( Hubbard_lmax, Hubbard_l(nt) )
!
WRITE( UNIT = stdout, &
FMT = * ) ' HUBBARD L FOR TYPE ',psd(nt),' IS ', Hubbard_l(nt)
WRITE( UNIT = stdout, FMT = * ) &
' HUBBARD L FOR TYPE ',upf(nt)%psd,' IS ', Hubbard_l(nt)
!
END IF
!

6
PW/stres_loc.f90
View File

@ -23,7 +23,7 @@ subroutine stres_loc (sigmaloc)
USE vlocal, ONLY : strf, vloc
USE wvfct, ONLY : gamma_only
USE wavefunctions_module, ONLY : psic
USE uspp_param, ONLY : vloc_at, zp
USE uspp_param, ONLY : upf
!
implicit none
!
@ -67,8 +67,8 @@ subroutine stres_loc (sigmaloc)
do nt = 1, ntyp
! dvloc contains dV_loc(G)/dG
call dvloc_of_g (lloc (nt), lmax (nt), numeric (nt), rgrid(nt)%mesh, &
msh (nt), rgrid(nt)%rab, rgrid(nt)%r, vloc_at (1, nt), &
cc (1, nt), alpc (1, nt), nlc (nt), nnl (nt), zp (nt), &
msh (nt), rgrid(nt)%rab, rgrid(nt)%r, upf(nt)%vloc(1), &
cc (1, nt), alpc (1, nt), nlc (nt), nnl (nt), upf(nt)%zp, &
aps (1, 0, nt), alps (1, 0, nt), tpiba2, ngl, gl, omega, dvloc)
! no G=0 contribution
do ng = 1, ngm

14
PW/stres_us.f90
View File

@ -19,7 +19,7 @@ SUBROUTINE stres_us( ik, gk, sigmanlc )
USE klist, ONLY : nks, xk
USE lsda_mod, ONLY : current_spin, lsda, isk
USE wvfct, ONLY : gamma_only, npw, npwx, nbnd, igk, wg, et
USE uspp_param, ONLY : lmaxkb, nh, tvanp, newpseudo
USE uspp_param, ONLY : upf, lmaxkb, nh, newpseudo
USE uspp, ONLY : nkb, vkb, qq, deeq, deeq_nc, qq_so
USE wavefunctions_module, ONLY : evc
USE spin_orb, ONLY : lspinorb
@ -111,7 +111,7 @@ SUBROUTINE stres_us( ik, gk, sigmanlc )
et(ibnd,ik) * qq(ih,ih,np)
evps = evps + fac * ps * ABS( becp(ikb,ibnd) )**2
!
IF ( tvanp(np) .OR. newpseudo(np) ) THEN
IF ( upf(np)%tvanp .OR. newpseudo(np) ) THEN
!
! ... only in the US case there is a contribution
! ... for jh<>ih
@ -149,7 +149,7 @@ SUBROUTINE stres_us( ik, gk, sigmanlc )
IF ( ityp(na) == np ) THEN
DO ih = 1, nh(np)
ikb = ijkb0 + ih
IF ( .NOT. ( tvanp(np) .OR. newpseudo(np) ) ) THEN
IF ( .NOT. ( upf(np)%tvanp .OR. newpseudo(np) ) ) THEN
ps = becp(ikb,ibnd) * &
( deeq(ih,ih,na,current_spin) - &
et(ibnd,ik) * qq(ih,ih,np) )
@ -202,7 +202,7 @@ SUBROUTINE stres_us( ik, gk, sigmanlc )
IF ( ityp(na) == np ) THEN
DO ih = 1, nh(np)
ikb = ijkb0 + ih
IF ( .NOT. ( tvanp(np) .OR. newpseudo(np) ) ) THEN
IF ( .NOT. ( upf(np)%tvanp .OR. newpseudo(np) ) ) THEN
ps = becp(ikb,ibnd) * &
( deeq(ih,ih,na,current_spin) - &
et(ibnd,ik) * qq(ih,ih,np ) )
@ -342,7 +342,7 @@ SUBROUTINE stres_us( ik, gk, sigmanlc )
evps = evps + fac * ps * ABS( becp(ikb,ibnd) )**2
END IF
!
IF ( tvanp(np) .OR. newpseudo(np) ) THEN
IF ( upf(np)%tvanp .OR. newpseudo(np) ) THEN
!
! ... only in the US case there is a contribution
! ... for jh<>ih
@ -408,7 +408,7 @@ SUBROUTINE stres_us( ik, gk, sigmanlc )
IF ( ityp(na) == np ) THEN
DO ih = 1, nh(np)
ikb = ijkb0 + ih
IF ( .NOT. ( tvanp(np) .OR. newpseudo(np) ) ) THEN
IF ( .NOT. ( upf(np)%tvanp .OR. newpseudo(np) ) ) THEN
IF (noncolin) THEN
if (lspinorb) call errore('stres_us','wrong case',1)
ijs=0
@ -518,7 +518,7 @@ SUBROUTINE stres_us( ik, gk, sigmanlc )
IF ( ityp(na) == np ) THEN
DO ih = 1, nh(np)
ikb = ijkb0 + ih
IF ( .NOT. ( tvanp(np) .OR. newpseudo(np) ) ) THEN
IF ( .NOT. ( upf(np)%tvanp .OR. newpseudo(np) ) ) THEN
IF (noncolin) THEN
ijs=0
ps_nc = (0.D0,0.D0)

8
PW/sum_band.f90
View File

@ -32,7 +32,7 @@ SUBROUTINE sum_band()
USE io_files, ONLY : iunwfc, nwordwfc, iunigk
USE buffers, ONLY : get_buffer
USE uspp, ONLY : nkb, vkb, becsum, nhtol, nhtoj, indv, okvan
USE uspp_param, ONLY : nh, tvanp, nhm
USE uspp_param, ONLY : upf, nh, nhm
USE wavefunctions_module, ONLY : evc, psic, psic_nc
USE noncollin_module, ONLY : noncolin, npol
USE spin_orb, ONLY : lspinorb, domag, so, fcoef
@ -343,7 +343,7 @@ SUBROUTINE sum_band()
!
DO np = 1, ntyp
!
IF ( tvanp(np) ) THEN
IF ( upf(np)%tvanp ) THEN
!
DO na = 1, nat
!
@ -572,7 +572,7 @@ SUBROUTINE sum_band()
!
DO np = 1, ntyp
!
IF ( tvanp(np) ) THEN
IF ( upf(np)%tvanp ) THEN
!
DO na = 1, nat
!
@ -665,7 +665,7 @@ SUBROUTINE sum_band()
IF (noncolin.and.okvan) THEN
DO np = 1, ntyp
IF ( tvanp(np) ) THEN
IF ( upf(np)%tvanp ) THEN
DO na = 1, nat
IF (ityp(na)==np) THEN
IF (so(np)) THEN

38
PW/summary.f90
View File

@ -22,7 +22,6 @@ SUBROUTINE summary()
USE constants, ONLY : amconv
USE cell_base, ONLY : alat, ibrav, omega, at, bg, celldm
USE ions_base, ONLY : nat, atm, zv, tau, ntyp => nsp, ityp
USE uspp_param, ONLY : psd
USE char, ONLY : title, sname
USE cellmd, ONLY : calc, cmass
USE ions_base, ONLY : amass
@ -48,6 +47,7 @@ SUBROUTINE summary()
USE funct, ONLY : write_dft_name
USE bp, ONLY : lelfield, gdir, nppstr, efield, nberrycyc
USE fixed_occ, ONLY : f_inp, tfixed_occ
USE uspp_param, ONLY : upf
USE wvfct, ONLY : nbnd
USE lsda_mod, ONLY : nspin
USE mp_global, ONLY : intra_pool_comm
@ -181,10 +181,10 @@ SUBROUTINE summary()
DO nt = 1, ntyp
IF (calc.EQ.' ') THEN
WRITE( stdout, '(5x,a6,6x,f10.2,2x,f10.5,5x,5 (a2,"(",f5.2,")"))') &
atm(nt), zv(nt), amass(nt), psd(nt), xp
atm(nt), zv(nt), amass(nt), upf(nt)%psd, xp
ELSE
WRITE( stdout, '(5x,a6,6x,f10.2,2x,f10.5,5x,5 (a2,"(",f5.2,")"))') &
atm(nt), zv(nt), amass(nt)/amconv, psd(nt), xp
atm(nt), zv(nt), amass(nt)/amconv, upf(nt)%psd, xp
END IF
ENDDO
@ -404,11 +404,10 @@ SUBROUTINE print_ps_info
USE io_global, ONLY : stdout
USE io_files, ONLY : psfile
USE ions_base, ONLY : ntyp => nsp
USE atom, ONLY : nlcc, rgrid
USE atom, ONLY : rgrid, nlcc
USE pseud, ONLY : alps, alpc, cc, aps, nlc, nnl, lmax, lloc, &
a_nlcc, b_nlcc, alpha_nlcc
USE uspp_param, ONLY : zp, nqf, rinner, nqlc, nbeta, iver, lll, &
psd, tvanp
USE uspp_param, ONLY : upf, iver
USE grid_paw_variables, ONLY: tpawp
!
INTEGER :: nt
@ -416,7 +415,7 @@ SUBROUTINE print_ps_info
!
DO nt = 1, ntyp
!
IF ( tvanp(nt) ) THEN
IF ( upf(nt)%tvanp ) THEN
ps='Ultrasoft'
ELSE IF ( tpawp (nt) ) THEN
ps="Projector augmented-wave"
@ -424,34 +423,34 @@ SUBROUTINE print_ps_info
ps='Norm-conserving'
END IF
!
IF ( nlcc (nt) ) ps = TRIM(ps) // ' + core correction'
IF ( nlcc(nt) ) ps = TRIM(ps) // ' + core correction'
!
WRITE( stdout, '(/5x,"PseudoPot. #",i2," for ",a2," read from file ",a)')&
nt, psd(nt), TRIM (psfile(nt))
nt, upf(nt)%psd, TRIM (psfile(nt))
!
WRITE( stdout, '( 5x,"Pseudo is ",a,", Zval =",f5.1)') &
TRIM (ps), zp (nt)
TRIM (ps), upf(nt)%zp
!
IF ( iver (1, nt) > 0 ) &
WRITE( stdout, '(5x,"Generated by v. ", 3i3, " of US pseudo code")') &
(iver (i, nt) , i = 1, 3)
WRITE( stdout, '(5x,"Using radial grid of ", i4, " points, ", &
&i2," beta functions with: ")') rgrid(nt)%mesh, nbeta (nt)
DO ib = 1, nbeta (nt)
&i2," beta functions with: ")') rgrid(nt)%mesh, upf(nt)%nbeta
DO ib = 1, upf(nt)%nbeta
IF (ib<10) THEN
WRITE( stdout, '(15x," l(",i1,") = ",i3)') ib, lll (ib, nt)
WRITE( stdout, '(15x," l(",i1,") = ",i3)') ib, upf(nt)%lll(ib)
ELSE
WRITE( stdout, '(14x," l(",i2,") = ",i3)') ib, lll (ib, nt)
WRITE( stdout, '(14x," l(",i2,") = ",i3)') ib, upf(nt)%lll(ib)
ENDIF
END DO
IF ( tvanp(nt) ) THEN
IF (nqf(nt)==0) THEN
IF ( upf(nt)%tvanp ) THEN
IF (upf(nt)%nqf==0) THEN
WRITE( stdout, '(5x,"Q(r) pseudized with 0 coefficients ",/)')
ELSE
WRITE( stdout, '(5x,"Q(r) pseudized with ", &
& i2," coefficients, rinner = ",3f8.3,/ &
& 52x,3f8.3,/ &
& 52x,3f8.3)') nqf(nt), (rinner(i,nt), i=1,nqlc(nt) )
& 52x,3f8.3,/ 52x,3f8.3)') &
& upf(nt)%nqf, (upf(nt)%rinner(i), i=1,upf(nt)%nqlc)
END IF
!
ELSE
@ -473,7 +472,8 @@ SUBROUTINE print_ps_info
WRITE( stdout, '(5x,"a(i) =",4x,3g13.5)') (aps (i, l, nt) , i = 1,3)
WRITE( stdout, '(5x,"a(i+3)=",4x,3g13.5)') (aps (i, l, nt) , i= 4, 6)
ENDDO
IF ( nlcc(nt) ) WRITE( stdout, 200) a_nlcc(nt), b_nlcc(nt), alpha_nlcc(nt)
IF ( nlcc(nt) ) &
WRITE( stdout, 200) a_nlcc(nt), b_nlcc(nt), alpha_nlcc(nt)
200 FORMAT(/5x,'nonlinear core correction: ', &
& 'rho(r) = ( a + b r^2) exp(-alpha r^2)', &
& /,5x,'a =',4x,g11.5, &

20
PW/tabd.f90
View File

@ -15,7 +15,7 @@ subroutine tabd (nt, occ_loc)
! (PPs usually are built on non physical configurations)
!
USE kinds, ONLY: DP
USE uspp_param, ONLY: psd
USE uspp_param, ONLY: upf
implicit none
real(DP) :: occ_loc
! output: the total number of d electrons
@ -24,29 +24,29 @@ subroutine tabd (nt, occ_loc)
!
! TRANSITION METALS
!
if (psd (nt) .eq.'Mn') then
if (upf(nt)%psd .eq.'Mn') then
occ_loc = 5.d0
elseif (psd (nt) .eq.'Fe') then
elseif (upf(nt)%psd .eq.'Fe') then
occ_loc = 6.d0
elseif (psd (nt) .eq.'Co') then
elseif (upf(nt)%psd .eq.'Co') then
occ_loc = 7.d0
elseif (psd (nt) .eq.'Ni') then
elseif (upf(nt)%psd .eq.'Ni') then
occ_loc = 8.d0
elseif (psd (nt) .eq.'Cu') then
elseif (upf(nt)%psd .eq.'Cu') then
occ_loc = 10.d0
!
! RARE EARTHS
!
elseif (psd (nt) .eq.'Ce') then
elseif (upf(nt)%psd .eq.'Ce') then
occ_loc = 2.d0
!
! OTHER ELEMENTS
!
elseif (psd (nt) .eq.'C') then
elseif (upf(nt)%psd .eq.'C') then
occ_loc = 2.d0
elseif (psd (nt) .eq.'O') then
elseif (upf(nt)%psd .eq.'O') then
occ_loc = 4.d0
elseif (psd (nt) .eq.'H') then
elseif (upf(nt)%psd .eq.'H') then
occ_loc = 1.d0
else
occ_loc = 0.d0

4
PW/usnldiag.f90
View File

@ -18,7 +18,7 @@ subroutine usnldiag (h_diag, s_diag)
USE wvfct, ONLY: npw, npwx
USE lsda_mod, ONLY: current_spin
USE uspp, ONLY: deeq, vkb, qq, qq_so, deeq_nc
USE uspp_param, ONLY: nh, tvanp, newpseudo
USE uspp_param, ONLY: upf, nh, newpseudo
USE spin_orb, ONLY: lspinorb
USE noncollin_module, ONLY: noncolin, npol
!
@ -69,7 +69,7 @@ subroutine usnldiag (h_diag, s_diag)
s_diag (ig,ipol) = s_diag (ig,ipol) + ps2(ipol) * ar
enddo
enddo
if (tvanp (nt) .or.newpseudo (nt) ) then
if ( upf(nt)%tvanp .or.newpseudo (nt) ) then
do jh = 1, nh (nt)
if (jh.ne.ih) then
jkb = ijkb0 + jh