1
2! Copyright (C) 2008 F. Bultmark, F. Cricchio and L. Nordstrom.
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: fyukawa0
8! !INTERFACE:
9real(8) function fyukawa0(is,l1,k)
10! !USES:
11use modmain
12use moddftu
13! !INPUT/OUTPUT PARAMETERS:
14!   is : species type (in,integer)
15!   l  : an angular momentum (in,integer)
16!   k  : order of Slater parameter (in,integer)
17! !DESCRIPTION:
18!   Calculates the Slater parameters in the unscreened case. See {\it Phys. Rev.
19!   B} {\bf 52}, 1421 (1995) and {\it Phys. Rev. B} {\bf 80}, 035121 (2009).
20!
21! !REVISION HISTORY:
22!   Created April 2008 (LN)
23!   Modified and tested July 2008 (LN and FC)
24!EOP
25!BOC
26implicit none
27! arguments
28integer, intent(in) :: is
29integer, intent(in) :: l1
30integer, intent(in) :: k
31! local variables
32integer, parameter :: nstart=1
33integer ir,nr,ias
34integer ir1,ir2,nr1,nr2
35real(8) r2,x
36! automatic arrays
37real(8) clow(nrmtmax),chigh(nrmtmax)
38real(8) blow(nrmtmax),bhigh(nrmtmax),fint(nrmtmax)
39real(8) bint(nrmtmax),cint(nrmtmax)
40! allocatable arrays
41real(8), allocatable :: a(:,:),b(:,:)
42ias=idxas(1,is)
43nr=nrmt(is)
44allocate(a(0:k,nr),b(0:k,nr))
45a(:,:)=0.d0
46b(:,:)=0.d0
47! calculate unscreened Slater parameters
48do ir=1,nr
49  r2=rlmt(ir,2,is)
50  bint(ir)=fdufr(ir,l1,ias)*fdufr(ir,l1,ias)*r2
51  x=rsp(ir,is)**k
52  a(k,ir)=x
53  b(k,ir)=1.d0/(x*rsp(ir,is))
54end do
55do ir=nstart,nr
56  nr1=ir-nstart+1
57  nr2=nr-ir+1
58  do ir1=1,nr1
59    ir2=ir1+nstart-1
60    blow(ir1)=bint(ir2)*a(k,ir2)
61  end do
62  call fderiv(-1,nr1,rsp(nstart,is),blow,clow)
63  do ir1=1,nr2
64    ir2=ir1+ir-1
65    bhigh(ir1)=bint(ir2)*b(k,ir2)
66  end do
67  call fderiv(-1,nr2,rsp(ir,is),bhigh,chigh)
68  cint(ir-nstart+1)=bint(ir)*(b(k,ir)*clow(nr1)+a(k,ir)*chigh(nr2))
69end do
70nr1=nr-nstart+1
71call fderiv(-1,nr1,rsp(nstart,is),cint,fint)
72fyukawa0=fint(nr1)
73deallocate(a,b)
74end function
75!EOC
76
77