1 2! Copyright (C) 2008 F. Bultmark, F. Cricchio, L. Nordstrom and J. K. Dewhurst. 3! This file is distributed under the terms of the GNU General Public License. 4! See the file COPYING for license details. 5 6!BOP 7! !ROUTINE: genfdu 8! !INTERFACE: 9subroutine genfdu(i,u,j,f) 10! !USES: 11use moddftu 12use modmpi 13! !INPUT/OUTPUT PARAMETERS: 14! i : DFT+U entry (in,integer) 15! u : parameter U (inout,real) 16! j : parameter J (inout,real) 17! f : Slater parameters (inout,real) 18! !DESCRIPTION: 19! Calculate the Slater parameters for DFT+$U$ calculation with different 20! approaches, see {\it Phys. Rev. B} {\bf 80}, 035121 (2009). The relations 21! among Slater and Racah parameters are from E. U. Condon and G. H. Shortley, 22! {\it The Theory of Atomic Spectra}, The University Press, Cambridge (1935). 23! 24! !REVISION HISTORY: 25! Created July 2008 (Francesco Cricchio) 26!EOP 27!BOC 28implicit none 29! arguments 30integer, intent(in) :: i 31real(8), intent(inout) :: u,j 32real(8), intent(inout) :: f(0:2*lmaxdm) 33! local variables 34integer is,l,k,q 35real(8) r1,r2 36real(8) lambda,ufix 37! automatic arrays 38real(8) :: e(0:lmaxdm) 39real(8) :: a(3,3),v1(3),v2(3) 40! external functions 41real(8), external :: fyukawa,fyukawa0 42is=idftu(1,i) 43l=idftu(2,i) 44! load input parameters to calculate Slater integrals 45u=ujdu(1,i) 46j=ujdu(2,i) 47f(:)=fdu(:,i) 48e(:)=edu(:,i) 49lambda=lambdadu(i) 50if (inpdftu.lt.4) then 51! F(0) = U for any l-shell 52 if (inpdftu.eq.1) f(0)=u 53 select case(l) 54 case(0) 55! s electrons only f(0)=u 56 if (inpdftu.eq.3) then 57 f(0)=e(0) 58 u=f(0) 59 end if 60 case(1) 61! p electrons 62 if (inpdftu.eq.1) then 63! F(2) = 5.0 * J 64 f(2)=5.d0*j 65 else if (inpdftu.eq.3) then 66! F(0) = E(0) + J= E(0) + 5/3 * E(1) 67 f(0)=e(0)+(5.d0/3.d0)*e(1) 68! F(2) = 5 * J = 25/3 * E1, Eq. 101 69 f(2)=(25.d0/3.d0)*e(1) 70 end if 71 case(2) 72! d electrons 73 if (inpdftu.eq.1) then 74! r1 = F(4)/F(2), see PRB 52, R5467 (1995) 75 r1=0.625d0 76 f(2)=(14.d0*j)/(1.d0+r1) 77 f(4)=f(2)*r1 78 else if (inpdftu.eq.3) then 79! copy Racah parameters 80 v1(1:3)=e(0:2) 81! transformation matrix from Racah to Slater parameters 82! obtained from inversion of Eq. 110-112, LN Notes 29-12-08 83 a(1,1)=1.d0 84 a(1,2)=1.4d0 85 a(1,3)=0.d0 86 a(2,1)=0.d0 87 a(2,2)=0.1428571428571428d0 88 a(2,3)=1.285714285714286d0 89 a(3,1)=0.d0 90 a(3,2)=2.8571428571428571d-2 91 a(3,3)=-0.1428571428571428d0 92! multiply transformation matrix by Racah parameters 93 call r3mv(a,v1,v2) 94! Slater parameters, Eq. 104-105, LN Notes 29-12-08 95 f(0)=v2(1) 96 f(2)=49.d0*v2(2) 97 f(4)=441.d0*v2(3) 98 end if 99 case(3) 100! f electrons 101 if (inpdftu.eq.1) then 102! r2 = F(6)/F(2), r1 = F(4)/F(2), see PRB 50, 16861 (1994) 103 r1=451.d0/675.d0 104 r2=1001.d0/2025.d0 105 f(2)=6435.d0*j/(286.d0+195.d0*r1+250.d0*r2) 106 f(4)=f(2)*r1 107 f(6)=f(2)*r2 108 else if (inpdftu.eq.3) then 109! F(0) = E(0) + 9/7 * E(1) , Eq. 119, LN Notes 29-12-08 110 f(0)=e(0)+(9.d0/7.d0)*e(1) 111! copy Racah parameters 112 v1(1:3)=e(1:3) 113! transformation matrix from Racah to Slater parameters 114! obtained from inversion of Eq. 120-122, LN Notes 29-12-08 115 a(1,1)=2.3809523809523808d-2 116 a(1,2)=3.404761904761904d0 117 a(1,3)=0.2619047619047619d0 118 a(2,1)=1.2987012987012984d-2 119 a(2,2)=-1.688311688311688d0 120 a(2,3)=5.1948051948051951d-2 121 a(3,1)=2.1645021645021645d-3 122 a(3,2)=7.5757575757575760d-2 123 a(3,3)=-1.5151515151515152d-2 124! multiply transformation matrix by Racah parameters 125 call r3mv(a,v1,v2) 126! Slater parameters, Eq. 115-117, LN Notes 29-12-08 127 f(2)=225.d0*v2(1) 128 f(4)=1089.d0*v2(2) 129 f(6)=(184041.d0/25.d0)*v2(3) 130 end if 131 case default 132 write(*,*) 133 write(*,'("Error(genfdu): invalid l : ",I8)') l 134 write(*,*) 135 stop 136 end select 137else if (inpdftu.ge.4) then 138! define energies for Slater parameters 139 call energyfdu 140! write energies for Slater parameters to a file 141 if (mp_mpi) call writeefdu 142! calculate radial functions for Slater parameters 143 call genfdufr 144 if (inpdftu.eq.5) then 145 ufix=udufix(i) 146! calculate the lambda corresponding to udufix 147! lambdadu0 is in/out and is initialized to 0 in readinput 148 call findlambdadu(is,l,ufix,lambdadu0(i),lambda) 149 end if 150 do q=0,l 151 k=2*q 152 if (lambda.lt.1.d-2) then 153! unscreened Slater parameters 154 f(k)=fyukawa0(is,l,k) 155 else 156! screened Slater parameters 157 f(k)=fyukawa(is,l,k,lambda) 158 end if 159 end do 160end if 161! calculate U and J from Slater integrals 162if (inpdftu.ne.1) then 163 u=f(0) 164 select case(l) 165 case(0) 166 j=0.d0 167 case(1) 168! J = 1/5 * F(2) 169 j=(1.d0/5.d0)*f(2) 170 case(2) 171! J = 1/14 * ( F(2) + F(4) ), Eq. 106, LN Notes 29-12-08 172 j=(1.d0/14.d0)*(f(2)+f(4)) 173 case(3) 174! J= 2/45 * F(2) + 1/33 * F(4) + 50/1287 * F(6), Eq. 118, LN Notes 29-12-08 175 j=(2.d0/45.d0)*f(2)+(1.d0/33.d0)*f(4)+(50.d0/1287.d0)*f(6) 176 case default 177 write(*,*) 178 write(*,'("Error(genfdu): invalid l : ",I8)') l 179 write(*,*) 180 stop 181 end select 182end if 183! save calculated parameters 184! (except Racah parameters that are provided only as input) 185ujdu(1,i)=u 186ujdu(2,i)=j 187fdu(:,i)=f(:) 188lambdadu(i)=lambda 189end subroutine 190!EOC 191 192