derivylm Subroutine

public subroutine derivylm(v1, v2, v3, lmax, Rabs, Ylm, dYdth, dYdfi, d2Ydth2, d2Ydfi2, d2Ydthdfi)

Calculate the 1st and 2nd derivatives of real spherical harmonics with respect to , .

Use recursion relations for the assoc. Legendre functions P[l,m] to generate the derivatives. These are (taken from Abramowitz and Stegun, Handbook of Mathematical Functions, chapt. 8.):

P[l,m+1] = (x^2-1)^(-1/2) ( (l-m)xP[l,m] - (l+m)*P[l-1,m] ) (8.5.1)

(x^2-1)dP[l,m]/dx = (l+m)(l-m+1)(x^2-1)^(1/2) P[l,m-1] - mx*P[l,m] (8.5.2)

(x^2-1)dP[l,m]/dx = lxP[l,m] - (l+m)P[l-1,m] (8.5.4)

where x=cos(th), (x^2-1)^(1/2) = -sin(th), d/dth = -sin(th) d/dx.

Adding (8.5.2)+(8.5.4) and using (8.5.1) results in:

dPl,m / dth = (1/2) * ( -(l+m)(l-m+1)P[l,m-1] + P[l,m+1] ) (A)

It is implied that P[l,m]=0 if m>l or m<-l. Also, the term (x^2-1)^(1/2) is ambiguous for real x, 0<x<1; here it is interpreted as (x^2-1)^(1/2)=-sin(th), but (x^2-1)^(-1/2)=+1/sin(th) (in 8.5.1), otherwise the result (A) (which is cross-checked and correct) does not follow.

For the 2nd derivative apply (A) twice. Result:

ddPl,m/dth/dth = (1/4) * ( (l+m)(l-m+1)(l+m-1)(l-m+2) * P[l,m-2] - ( (l-m)(l+m+1)+(l+m)(l-m+1) ) P[l,m]
+ P[l,m+2] ) (B)

The fi-derivatives act on cos(fi),sin(fi) and are trivial.

For the associated Legendre functions use the recursion formulas:

(l-m+1)P[l+1,m] = (2l+1)cos(th)P[l,m] - (l+m)P[l-1,m] (8.5.3)

P[l+1,m] = P[l-1,m] - (2l+1)sin(th)*P[l,m-1] (8.5.5)

( with x=cos(th) ).

Recursion algorithm for the calculation of P[l,m] and calculation of Ylm taken over from subr. ymy of KKR program (implemented there by M. Weinert, B. Drittler).

For m<0, use P[l,-m] = P[l,m] (l-m)!/(l+m)! (C)

Taking into account the lm-prefactors of the spherical harmonics, we construct and use the functions

Q[l,m] = sqrt((2l+1)/(4pi)) * sqrt((l-m)!/(l+m)!) * P[l,m]

whence (A) and (B) become

dQ[l,m]/dth = (1/2) * ( -sqrt((l+m)(l-m+1))Q[l,m-1] + sqrt((l+m+1)(l-m))Q[l,m+1] ) (A1)

ddQ[l,m]/dth/dth = (1/4) * ( sqrt((l+m)(l+m-1)(l-m+1)(l-m+2)) * Q[l,m-2] + ((l-m)(l+m+1)+(l+m)(l-m+1)) * Q[l,m] + sqrt((l-m)(l-m-1)(l+m+1)(l+m+2)) * Q[l,m+2] ) (A2)

Note on sign convension:

For the needs of GGA PW91 as implemented here, ylm and derivatives come with a different sign convention compared to the usual in the program: . Thus some signs change.

Arguments

Type IntentOptional Attributes Name
real(kind=8) :: v1
real(kind=8) :: v2
real(kind=8) :: v3
integer :: lmax
real(kind=8) :: Rabs
real(kind=8) :: Ylm(*)
real(kind=8) :: dYdth(*)
real(kind=8) :: dYdfi(*)
real(kind=8) :: d2Ydth2(*)
real(kind=8) :: d2Ydfi2(*)
real(kind=8) :: d2Ydthdfi(*)