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