From 5df79c53745fde5d6c3340a2979b1429cd5892c1 Mon Sep 17 00:00:00 2001 From: Henrik Rydberg Date: Sat, 8 Oct 2011 20:30:28 +0200 Subject: Initial import of htcd system 1.0 Signed-off-by: Henrik Rydberg --- src/labat/ldauc.f | 32 ++++++++++++++++++++++++++++++++ 1 file changed, 32 insertions(+) create mode 100644 src/labat/ldauc.f (limited to 'src/labat/ldauc.f') diff --git a/src/labat/ldauc.f b/src/labat/ldauc.f new file mode 100644 index 0000000..4285c9a --- /dev/null +++ b/src/labat/ldauc.f @@ -0,0 +1,32 @@ + subroutine ldauc(rs,zet,ec,ecrs,eczet,vcup,vcdn) +c uniform-gas correlation of perdew and wang 1991 +c calculates the local correlation potential within the lda approx. +c input rs : seitz radius +c input zet : relative spin polarization +c output ec : correlation energy +c output ecrs : derivative of ec with respect to rs +c output eczet : derivative of ec with respect to zet +c output vcup, vcdn : up- and down-spin potentials + implicit double precision (a-h,o-z) + data gam,fzz/0.5198421d0,1.709921d0/ + data thrd,thrd4/0.333333333333d0,1.333333333333d0/ + f = ((1.d0+zet)**thrd4+(1.d0-zet)**thrd4-2.d0)/gam + call gcor(0.0310907d0,0.21370d0,7.5957d0,3.5876d0,1.6382d0, + 1 0.49294d0,1.00d0,rs,eu,eurs) + call gcor(0.01554535d0,0.20548d0,14.1189d0,6.1977d0,3.3662d0, + 1 0.62517d0,1.00d0,rs,ep,eprs) + call gcor(0.0168869d0,0.11125d0,10.357d0,3.6231d0,0.88026d0, + 1 0.49671d0,1.00d0,rs,alfm,alfrsm) +c alfm is minus the spin stiffness alfc + z4 = zet**4 + ec = eu*(1.d0-f*z4)+ep*f*z4-alfm*f*(1.d0-z4)/fzz +c energy done. now the potential: + ecrs = eurs*(1.d0-f*z4)+eprs*f*z4-alfrsm*f*(1.d0-z4)/fzz + fz = thrd4*((1.d0+zet)**thrd-(1.d0-zet)**thrd)/gam + eczet = 4.d0*(zet**3)*f*(ep-eu+alfm/fzz)+fz*(z4*ep-z4*eu + 1 -(1.d0-z4)*alfm/fzz) + comm = ec -rs*ecrs/3.d0-zet*eczet + vcup = comm + eczet + vcdn = comm - eczet + return + end -- cgit v1.2.3