mirror of
https://github.com/QuantumPackage/qp2.git
synced 2024-11-19 20:42:36 +01:00
51 lines
1.1 KiB
Fortran
51 lines
1.1 KiB
Fortran
|
double precision function plgndr(l,m,x)
|
||
|
|
integer, intent(in) :: l,m
|
||
|
|
double precision, intent(in) :: x
|
||
|
|
BEGIN_DOC
|
||
|
|
! associated Legenre polynom P_l,m(x). Used for the Y_lm(theta,phi)
|
||
|
|
! Taken from https://iate.oac.uncor.edu/~mario/materia/nr/numrec/f6-8.pdf
|
||
|
|
END_DOC
|
||
|
|
integer :: i,ll
|
||
|
|
double precision :: fact,pll,pmm,pmmp1,somx2
|
||
|
|
if(m.lt.0.or.m.gt.l.or.dabs(x).gt.1.d0)then
|
||
|
|
print*,'bad arguments in plgndr'
|
||
|
|
pause
|
||
|
|
endif
|
||
|
|
pmm=1.d0
|
||
|
|
if(m.gt.0) then
|
||
|
|
somx2=dsqrt((1.d0-x)*(1.d0+x))
|
||
|
|
fact=1.d0
|
||
|
|
do i=1,m
|
||
|
|
pmm=-pmm*fact*somx2
|
||
|
|
fact=fact+2.d0
|
||
|
|
enddo
|
||
|
|
endif ! m > 0
|
||
|
|
if(l.eq.m) then
|
||
|
|
plgndr=pmm
|
||
|
|
else
|
||
|
|
pmmp1=x*(2*m+1)*pmm ! Compute P_m+1^m
|
||
|
|
if(l.eq.m+1) then
|
||
|
|
plgndr=pmmp1
|
||
|
|
else ! Compute P_l^m, l> m+1
|
||
|
|
do ll=m+2,l
|
||
|
|
pll=(x*dble(2*ll-1)*pmmp1-dble(ll+m-1)*pmm)/(ll-m)
|
||
|
|
pmm=pmmp1
|
||
|
|
pmmp1=pll
|
||
|
|
enddo
|
||
|
|
plgndr=pll
|
||
|
|
endif ! l.eq.m+1
|
||
|
|
endif ! l.eq.m
|
||
|
|
return
|
||
|
|
end
|
||
|
|
|
||
|
|
double precision function ortho_assoc_gaus_pol(l1,m1,l2)
|
||
|
|
implicit none
|
||
|
|
integer, intent(in) :: l1,m1,l2
|
||
|
|
double precision :: fact
|
||
|
|
if(l1.ne.l2)then
|
||
|
|
ortho_assoc_gaus_pol= 0.d0
|
||
|
|
else
|
||
|
|
ortho_assoc_gaus_pol = 2.d0*fact(l1+m1) / (dble(2*l1+1)*fact(l1-m1))
|
||
|
|
endif
|
||
|
|
end
|