qe_erf to erf in CPV

This commit is contained in:
fabrizio22 2021-02-18 18:03:40 +01:00
parent a33f94f6c6
commit 07ac37218d
4 changed files with 16 additions and 24 deletions

View File

@ -167,7 +167,7 @@ SUBROUTINE EFERMI(NEL,NBANDS,DEL,NKPTS,OCC,EF,EIGVAL, &
REAL(kind=DP) :: weight(nkpts), sort(nbands*nkpts)
REAL(kind=DP), EXTERNAL :: qe_erfc,FERMID,DELTHM,POSHM,POSHM2, EFERMI_SPLINE
REAL(kind=DP), EXTERNAL :: FERMID,DELTHM,POSHM,POSHM2, EFERMI_SPLINE
INTEGER, PARAMETER :: JMAX =300
REAL(kind=DP), PARAMETER :: XACC=1.0D-17
@ -387,7 +387,7 @@ SUBROUTINE EFERMI(NEL,NBANDS,DEL,NKPTS,OCC,EF,EIGVAL, &
DO J = 1,NBANDS
X = (XE2 - EIGVAL(J,ISPPT))/DEL
IF(ISMEAR.EQ.1) THEN
Z1 = Z1 + WEIGHT(ISPPT)*( 2.d0 - qe_erfc(X) )/fspin
Z1 = Z1 + WEIGHT(ISPPT)*( 2.d0 - erfc(X) )/fspin
ELSEIF(ISMEAR.EQ.2) THEN
Z1 = Z1 + WEIGHT(ISPPT)*FERMID(-X)/fspin
ELSEIF(ISMEAR.EQ.3) THEN
@ -413,7 +413,7 @@ SUBROUTINE EFERMI(NEL,NBANDS,DEL,NKPTS,OCC,EF,EIGVAL, &
DO J = 1,NBANDS
X = (XE1 - EIGVAL(J,ISPPT))/DEL
IF(ISMEAR.EQ.1) THEN
Z1 = Z1 + WEIGHT(ISPPT)*( 2.d0 - qe_erfc(X) )/fspin
Z1 = Z1 + WEIGHT(ISPPT)*( 2.d0 - erfc(X) )/fspin
ELSEIF(ISMEAR.EQ.2) THEN
Z1 = Z1 + WEIGHT(ISPPT)*FERMID(-X)/fspin
ELSEIF(ISMEAR.EQ.3) THEN
@ -454,7 +454,7 @@ SUBROUTINE EFERMI(NEL,NBANDS,DEL,NKPTS,OCC,EF,EIGVAL, &
DO J2 = 1,NBANDS
X = (XMID - EIGVAL(J2,ISPPT))/DEL
IF(ISMEAR.EQ.1) THEN
Z1 = Z1 + WEIGHT(ISPPT)*( 2.d0 - qe_erfc(X) )/fspin
Z1 = Z1 + WEIGHT(ISPPT)*( 2.d0 - erfc(X) )/fspin
ELSEIF(ISMEAR.EQ.2) THEN
Z1 = Z1 + WEIGHT(ISPPT)*FERMID(-X)/fspin
ELSEIF(ISMEAR.EQ.3) THEN
@ -488,7 +488,7 @@ SUBROUTINE EFERMI(NEL,NBANDS,DEL,NKPTS,OCC,EF,EIGVAL, &
DO J = 1,NBANDS
X = ( EF-EIGVAL(J,ISPPT))/DEL
IF(ISMEAR.EQ.1) THEN
OCC(J,ISPPT) = 2.d0 - qe_erfc(X)
OCC(J,ISPPT) = 2.d0 - erfc(X)
ELSEIF(ISMEAR.EQ.2) THEN
OCC(J,ISPPT) = FERMID(-X)
ELSEIF(ISMEAR.EQ.3) THEN
@ -530,7 +530,7 @@ SUBROUTINE EFERMI(NEL,NBANDS,DEL,NKPTS,OCC,EF,EIGVAL, &
& *(2.0d0*x*x-1.d0)*exp(-x*x)/(2.0d0*sqrt(pi))
ELSEIF(ISMEAR.EQ.4) THEN
x=abs(x)
zeta=eesh*abs(x)*exp(-(x+sq2i)**2)+piesqq*qe_erfc(x+sq2i)
zeta=eesh*abs(x)*exp(-(x+sq2i)**2)+piesqq*erfc(x+sq2i)
delcor=delcor-del*WEIGHT(ISPPT)*zeta
ELSEIF(ISMEAR.EQ.5) THEN
a=-0.5634d0
@ -629,7 +629,6 @@ FUNCTION delthm(xx)
REAL(kind=DP) :: delthm
REAL(kind=DP), INTENT(in) :: xx
REAL(kind=DP), EXTERNAL :: qe_erfc
REAL(kind=DP) :: pi
@ -639,7 +638,7 @@ FUNCTION delthm(xx)
ELSEIF(XX .LT. -10.D0) THEN
DELTHM=0.D0
ELSE
DELTHM=(2.D0-qe_erfc(XX))+XX*EXP(-XX*XX)/SQRT(PI)
DELTHM=(2.D0-erfc(XX))+XX*EXP(-XX*XX)/SQRT(PI)
ENDIF
!
RETURN
@ -680,7 +679,7 @@ FUNCTION poshm(x)
REAL(kind=DP) :: poshm
REAL(kind=DP), INTENT(in) :: x
REAL(kind=DP), EXTERNAL :: qe_erfc
REAL(kind=DP), EXTERNAL :: erfc
REAL(kind=DP) :: pi,a
@ -692,7 +691,7 @@ FUNCTION poshm(x)
ELSEIF(X .LT. -10.D0) THEN
POSHM=0.D0
ELSE
POSHM=(2.D0-qe_erfc(X))+(-2.d0*a*x*x+2*x+a)*EXP(-X*X)/SQRT(PI)/2.d0
POSHM=(2.D0-erfc(X))+(-2.d0*a*x*x+2*x+a)*EXP(-X*X)/SQRT(PI)/2.d0
ENDIF
!
RETURN
@ -710,7 +709,6 @@ FUNCTION poshm2(x)
REAL(kind=DP) :: poshm2
REAL(kind=DP), INTENT(in) :: x
REAL(kind=DP), EXTERNAL :: qe_erfc
REAL(kind=DP) :: pi
@ -720,7 +718,7 @@ FUNCTION poshm2(x)
ELSEIF(X .LT. -10.D0) THEN
POSHM2=0.D0
ELSE
POSHM2=(2.D0-qe_erfc(X-1.d0/sqrt(2.d0)))+ &
POSHM2=(2.D0-erfc(X-1.d0/sqrt(2.d0)))+ &
& sqrt(2.d0)*exp(-x*x+sqrt(2.d0)*x-0.5d0)/sqrt(pi)
ENDIF
!

View File

@ -45,8 +45,6 @@
REAL(DP), INTENT(IN) :: omega, hmat( 3, 3 )
COMPLEX(DP) :: screen_coul( dfftp%ngm )
REAL(DP), EXTERNAL :: qe_erf
! ... Locals
!
REAL(DP), ALLOCATABLE :: grr(:,:)
@ -87,7 +85,7 @@
IF( rmod < gsmall ) THEN
grr( ir, 1 ) = fact * 2.0d0 * rc / SQRT( pi )
ELSE
grr( ir, 1 ) = fact * qe_erf( rc * rmod ) / rmod
grr( ir, 1 ) = fact * erf( rc * rmod ) / rmod
END IF
END DO
END DO
@ -427,8 +425,6 @@
LOGICAL, INTENT(IN) :: TSTRESS
REAL(DP), INTENT(in) :: hmat( 3, 3 )
REAL(DP), EXTERNAL :: qe_erfc
INTEGER, EXTERNAL :: ldim_block, gind_block
@ -525,7 +521,7 @@
ELSE
ESRTZERO=1.D0
END IF
ADDESR = ZV2_KJ * qe_erfc(ARG) / RLM
ADDESR = ZV2_KJ * erfc(ARG) / RLM
ESR = ESR + ESRTZERO*ADDESR
ADDPRE = FACT_PRE * EXP(-ARG*ARG)
REPAND = ESRTZERO*(ADDESR + ADDPRE)/ERRE2

View File

@ -212,7 +212,6 @@
integer :: ig, ir, irmax
real(DP), allocatable:: f(:),vscr(:), figl(:)
real(DP), allocatable:: df(:), dfigl(:)
real(DP), external :: qe_erf, qe_erfc
!
allocate( vscr(mesh), figl(ngs) )
if (tpre) then
@ -227,7 +226,7 @@
end do
!
do ir = 1, irmax
vscr(ir) = 0.5d0 * r(ir) * vloc_at(ir) + zv * qe_erf( r(ir) / rcmax )
vscr(ir) = 0.5d0 * r(ir) * vloc_at(ir) + zv * erf( r(ir) / rcmax )
end do
do ir = irmax + 1, mesh
vscr(ir)=0.0d0
@ -250,7 +249,7 @@
allocate( df(mesh) )
end if
DO ir = 1, irmax
f(ir) = fpi * ( zv * qe_erfc( r(ir)/rcmax ) ) * r(ir)
f(ir) = fpi * ( zv * erfc( r(ir)/rcmax ) ) * r(ir)
END DO
DO ir = irmax + 1, mesh
f(ir)=0.0d0

View File

@ -42,7 +42,6 @@ function wgauss (x, n)
integer :: i, ni
! counter on the n indices
! counter on 2n
real(DP) :: gauss_freq
real(DP), parameter :: maxarg = 200.d0
! maximum value for the argument of the exponential
real(DP), parameter :: c = 0.7071067811865475_DP
@ -70,8 +69,8 @@ function wgauss (x, n)
endif
! Methfessel-Paxton
gauss_freq = 0.5_DP * ERFC( - x * c)
wgauss = gauss_freq (x * sqrt (2.0d0) )
!gauss_freq(x) = 0.5_DP * ERFC( - x * c)
wgauss = 0.5_DP * ERFC( - x * sqrt (2.0d0) * c)
if (n.eq.0) return
hd = 0.d0
arg = min (maxarg, x**2)