Commit 32c41c5d authored by Matthias Redies's avatar Matthias Redies

xcwgn moved

parent bf734f31
......@@ -22,11 +22,11 @@ xc-pot/vxcwb91.f90
xc-pot/xcbh.f90
xc-pot/xcpz.f90
xc-pot/xcvwn.f90
xc-pot/xcwgn.f90
)
set(fleur_F77 ${fleur_F77}
xc-pot/xch91.F
xc-pot/xcwgn.f
xc-pot/xcxal.f
)
......@@ -49,11 +49,11 @@ xc-pot/vxcwb91.f90
xc-pot/xcbh.f90
xc-pot/xcpz.f90
xc-pot/xcvwn.f90
xc-pot/xcwgn.f90
)
set(inpgen_F77 ${inpgen_F77}
xc-pot/xch91.F
xc-pot/xcwgn.f
xc-pot/xcxal.f
)
......
......@@ -4,77 +4,75 @@
! of the MIT license as expressed in the LICENSE file in more detail.
!--------------------------------------------------------------------------------
MODULE m_xcwgn
MODULE m_xcwgn
!-----------------------------------------------------------------------
! Called in case of icorr=1 : non-spinpolarized exchange-correlation
! Called in case of icorr=1 : non-spinpolarized exchange-correlation
! from wigners interpolation formula.
!
! krla=1: Relativistic correction of exchange energy and potential
! krla=1: Relativistic correction of exchange energy and potential
! related to Dirac kinetic energy, according to:
! A.H. MacDonald and S.H. Vosko, J. Phys. C12, 2977 (1979)
!
! be careful: calculation in rydberg!
!
! vxc calculates the XC-potential and
! exc calculates the XC-energy
!
! based on a subroutine by S. Bluegel; r.pentcheva 22.01.96
!-----------------------------------------------------------------------
USE m_constants, ONLY : pi_const
USE m_relcor
IMPLICIT NONE
USE m_constants, ONLY : pi_const
USE m_relcor
IMPLICIT NONE
REAL, PARAMETER, PRIVATE :: one = 1.0 , three = 3.0 , four = 4.0
REAL, PARAMETER, PRIVATE :: thrd = one/three , d_15 = 1.e-15
REAL, PARAMETER, PRIVATE :: cex = 0.91633058742 ! 3/2 * ( 3/(2*pi) )^(2/3)
REAL, PARAMETER, PRIVATE :: one = 1.0 , three = 3.0 , four = 4.0
REAL, PARAMETER, PRIVATE :: thrd = one/three , d_15 = 1.e-15
REAL, PARAMETER, PRIVATE :: cex = 0.91633058742 ! 3/2 * ( 3/(2*pi) )^(2/3)
CONTAINS
CONTAINS
!************************************************************************
SUBROUTINE vxcwgn(
> krla,jspins,
> mgrid,ngrid,rh,
< vx,vxc)
SUBROUTINE vxcwgn( &
krla,jspins, &
mgrid,ngrid,rh, &
vx,vxc)
!************************************************************************
!
! .. Scalar Arguments ..
INTEGER, INTENT (IN) :: jspins
INTEGER, INTENT (IN) :: krla ! run mode parameters
INTEGER, INTENT (IN) :: mgrid,ngrid ! mesh points
!
! .. Array Arguments ..
REAL, INTENT (IN) :: rh(mgrid,jspins) ! charge density
REAL, INTENT (OUT) :: vx(mgrid,jspins) ! x potential
REAL, INTENT (OUT) :: vxc(mgrid,jspins) ! xc potential
!
! .. Local Scalars ..
REAL fothrd,vxp,vcp
REAL :: fothrd,vxp,vcp
REAL :: rs, rho, thfpi, exp, ecp
INTEGER :: ir
!
! .. Local Arrays ..
REAL, ALLOCATABLE :: psi(:) ! relativistic exchange potential corr.
!
!-----s Intrinsic Functions
INTRINSIC alog,max
thfpi = three / ( four * pi_const )
fothrd = four/three
!
!-----> evaluate relativistic corrections for exchange
!
ALLOCATE ( psi(ngrid) )
CALL relcor(
> mgrid,ngrid,jspins,krla,.true.,rh,
< psi)
CALL relcor( &
mgrid,ngrid,jspins,krla, .TRUE. ,rh, &
psi)
DO ir = 1,ngrid
IF (jspins.EQ.1) THEN
rho = max(d_15,rh(ir,1))
IF (jspins == 1) THEN
rho = max(d_15,rh(ir,1))
ELSE
rho = max(d_15,rh(ir,1)) + max(d_15,rh(ir,jspins))
rho = max(d_15,rh(ir,1)) + max(d_15,rh(ir,jspins))
ENDIF
rs = (thfpi/rho)**thrd
exp = - cex / rs ! exchange energy = -0.9163306/rs
......@@ -84,13 +82,13 @@
vcp = fothrd * ecp + 2.288 / ( rs + 7.8 ) ** 2
vxc(ir,1) = vcp + vxp * psi(ir)
vx(ir,1) = vxp * psi(ir)
ENDDO
IF ( jspins .EQ. 2 ) THEN ! spinpolarized calculation
IF ( jspins == 2 ) THEN ! spinpolarized calculation
DO ir = 1,ngrid
vxc(ir,jspins) = vxc(ir,1)
vx(ir,jspins) = vx(ir,1)
ENDDO
ENDIF
......@@ -98,42 +96,42 @@
DEALLOCATE (psi)
RETURN
END SUBROUTINE vxcwgn
C***********************************************************************
SUBROUTINE excwgn(
> iofile,krla,jspins,
> mgrid,ngrid,rh,
< exc)
C***********************************************************************
!
END SUBROUTINE vxcwgn
!***********************************************************************
SUBROUTINE excwgn( &
iofile,krla,jspins, &
mgrid,ngrid,rh, &
exc)
!***********************************************************************
! .. Scalar Arguments ..
INTEGER, INTENT (IN) :: jspins
INTEGER, INTENT (IN) :: krla ! run mode parameters
INTEGER, INTENT (IN) :: iofile ! file number for read and write
INTEGER, INTENT (IN) :: mgrid,ngrid ! mesh points
!
! .. Array Arguments ..
REAL, INTENT (IN) :: rh(mgrid,jspins) ! charge density
REAL, INTENT (OUT) :: exc(mgrid) ! xc energy
!
! .. Local Arrays ..
REAL, ALLOCATABLE :: phi(:) ! relativistic exchange energy correct.
REAL :: rs, rho, thfpi, exp, ecp
INTEGER :: ir
!
thfpi = three / ( four * pi_const )
ALLOCATE ( phi(ngrid) )
CALL relcor(
> mgrid,ngrid,jspins,krla,.false.,rh,
< phi)
CALL relcor( &
mgrid,ngrid,jspins,krla, .FALSE. ,rh, &
phi)
DO ir = 1,ngrid
IF (jspins.EQ.1) THEN
rho = max(d_15,rh(ir,1))
IF (jspins == 1) THEN
rho = max(d_15,rh(ir,1))
ELSE
rho = max(d_15,rh(ir,1)) + max(d_15,rh(ir,jspins))
rho = max(d_15,rh(ir,1)) + max(d_15,rh(ir,jspins))
ENDIF
rs = (thfpi/rho)**thrd
exp = - cex / rs ! exchange energy = -0.9163306/rs
......@@ -141,13 +139,13 @@ C***********************************************************************
exc(ir) = ecp + exp * phi(ir)
ENDDO
IF (jspins.EQ.2) THEN
WRITE (iofile,'('' WARNING: Wigner correlation !'',
+ ''only applicable for non-spinpolarized calculations'')')
IF (jspins == 2) THEN
WRITE (iofile,'('' WARNING: Wigner correlation !'', &
''only applicable for non-spinpolarized calculations'')')
ENDIF
DEALLOCATE (phi)
RETURN
END SUBROUTINE excwgn
END MODULE m_xcwgn
END SUBROUTINE excwgn
END MODULE m_xcwgn
Markdown is supported
0%
or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment