TABLE OF CONTENTS
ABINIT/xcwign [ Functions ]
NAME
xcwign
FUNCTION
Returns exc, vxc, and eventually d(vxc)/d($\rho$) from input $\rho$. Wigner exchange and correlation (xc)--see e.g. David Pines, Elementary Excitations in Solids, p. 94, NY 1964. Expression is exc=-(0.44)/(rs+7.8)-efac/rs (hartree), efac below. rs = $(3/(4\pi))^{1/3}* \rho (r)^{-1/3}$.
COPYRIGHT
Copyright (C) 1998-2018 ABINIT group (DCA, XG, GMR) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt . For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
INPUTS
npt=number of real space points on which density is provided order=gives the maximal derivative of Exc computed. rspts(npt)=corresponding Wigner-Seitz radii, precomputed
OUTPUT
exc(npt)=exchange-correlation energy density (hartree) vxc(npt)=xc potential (d($\rho$*exc)/d($\rho$)) (hartree) if(order>1) dvxc(npt)=derivative d(vxc)/d($\rho$) (hartree*bohr^3)
PARENTS
drivexc
CHILDREN
SOURCE
37 #if defined HAVE_CONFIG_H 38 #include "config.h" 39 #endif 40 41 #include "abi_common.h" 42 43 44 subroutine xcwign(exc,npt,order,rspts,vxc,& !Mandatory arguments 45 & dvxc) !Optional arguments 46 47 use defs_basis 48 use m_errors 49 use m_profiling_abi 50 51 !This section has been created automatically by the script Abilint (TD). 52 !Do not modify the following lines by hand. 53 #undef ABI_FUNC 54 #define ABI_FUNC 'xcwign' 55 !End of the abilint section 56 57 implicit none 58 59 !Arguments ------------------------------------ 60 !scalars 61 integer,intent(in) :: npt,order 62 !arrays 63 real(dp),intent(in) :: rspts(npt) 64 real(dp),intent(out) :: exc(npt),vxc(npt) 65 real(dp),intent(out),optional :: dvxc(npt) 66 67 !Local variables------------------------------- 68 !c1 and c2 are the Wigner parameters in hartree and bohr resp. 69 !scalars 70 integer :: ipt 71 real(dp),parameter :: c1=0.44_dp,c2=7.8_dp,c4_3=4.0_dp/3.0_dp 72 real(dp),parameter :: c8_27=8.0_dp/27.0_dp 73 real(dp) :: dfac,efac,rs,rsc2m1,rsm1,vfac,vxcnum 74 character(len=500) :: message 75 76 ! ************************************************************************* 77 78 !Checks the values of order 79 if(order<0 .or. order>2)then 80 write(message, '(a,a,a,i0)' )& 81 & 'With Wigner xc functional, the only',ch10,& 82 & 'allowed values for order are 0, 1 or 2, while it is found to be',order 83 MSG_BUG(message) 84 end if 85 86 !Checks the compatibility between the order and the presence of the optional arguments 87 if(order <= 1 .and. present(dvxc))then 88 write(message, '(a,a,a,i3)' )& 89 & 'The order chosen does not need the presence',ch10,& 90 & 'of the vector dvxc, that is needed only with order=2 , while we have',order 91 MSG_BUG(message) 92 end if 93 94 !Compute vfac=(3/(2*Pi))^(2/3) 95 vfac=(1.5_dp/pi)**(2.0_dp/3.0_dp) 96 !Compute efac=(3/4)*vfac 97 efac=0.75_dp*vfac 98 !Compute dfac=(4*Pi/9)*vfac 99 dfac=(4.0_dp*pi/9.0_dp)*vfac 100 101 !separate cases with respect to order 102 if (order==2) then 103 104 ! Loop over grid points 105 do ipt=1,npt 106 rs=rspts(ipt) 107 rsm1=1.0_dp/rs 108 rsc2m1=1.0_dp/(rs+c2) 109 ! compute energy density (hartree) 110 exc(ipt)=-c1*rsc2m1-efac*rsm1 111 vxcnum=-(c4_3*rs+c2)*c1 112 ! compute potential (hartree) 113 vxc(ipt)=vxcnum*rsc2m1**2-vfac*rsm1 114 ! compute d(vxc)/d(rho) (hartree*bohr^3) 115 dvxc(ipt)=-(c8_27*pi)*(c1*rs**4)*(rs+rs+c2)*rsc2m1**3-dfac*rs**2 116 end do 117 else 118 119 ! Loop over grid points 120 do ipt=1,npt 121 rs=rspts(ipt) 122 rsm1=1.0_dp/rs 123 rsc2m1=1.0_dp/(rs+c2) 124 ! compute energy density (hartree) 125 exc(ipt)=-c1*rsc2m1-efac*rsm1 126 vxcnum=-(c4_3*rs+c2)*c1 127 ! compute potential (hartree) 128 vxc(ipt)=vxcnum*rsc2m1**2-vfac*rsm1 129 end do 130 131 end if 132 ! 133 end subroutine xcwign