mirror of
https://gitlab.com/scemama/eplf
synced 2024-06-26 23:22:20 +02:00
99 lines
1.6 KiB
FortranFixed
99 lines
1.6 KiB
FortranFixed
|
|
||
|
|
double precision function Boys(x,n)
|
||
|
|
implicit none
|
||
|
|
include 'constants.F'
|
||
|
|
|
||
|
|
real, intent(in) :: x
|
||
|
|
integer, intent(in) :: n
|
||
|
|
|
||
|
|
integer :: k
|
||
|
|
|
||
|
|
real, parameter :: thr = 6.
|
||
|
|
integer ,parameter :: Nmax = 20
|
||
|
|
|
||
|
|
double precision :: fact,fact2
|
||
|
|
|
||
|
|
if (n == 0) then
|
||
|
|
if (x > thr) then
|
||
|
|
Boys = 0.5d0*sqrt(pi/x)
|
||
|
|
else
|
||
|
|
Boys = 1.d0/dble(2*n+1)
|
||
|
|
do k=1,Nmax
|
||
|
|
Boys = Boys + (-x)**k/dble(fact(k)*(2*n+2*k+1))
|
||
|
|
enddo
|
||
|
|
|
||
|
|
endif
|
||
|
|
|
||
|
|
else
|
||
|
|
if (x > thr) then
|
||
|
|
Boys = fact2(2*n-1)*0.5d0**(n+1)*sqrt(pi/x**(2*n+1))
|
||
|
|
else
|
||
|
|
Boys = 1.d0/dble(2*n+1)
|
||
|
|
do k=1,Nmax
|
||
|
|
Boys = Boys + (-x)**k/dble(fact(k)*(2*n+2*k+1))
|
||
|
|
enddo
|
||
|
|
endif
|
||
|
|
endif
|
||
|
|
|
||
|
|
end function
|
||
|
|
|
||
|
|
double precision function fact2(n)
|
||
|
|
implicit none
|
||
|
|
integer :: n
|
||
|
|
double precision, save :: memo(1:100)
|
||
|
|
integer, save :: memomax = 1
|
||
|
|
|
||
|
|
ASSERT (mod(n,2) /= 0)
|
||
|
|
if (n<=memomax) then
|
||
|
|
if (n<3) then
|
||
|
|
fact2 = 1.d0
|
||
|
|
else
|
||
|
|
fact2 = memo(n)
|
||
|
|
endif
|
||
|
|
return
|
||
|
|
endif
|
||
|
|
|
||
|
|
integer :: i
|
||
|
|
memo(1) = 1.d0
|
||
|
|
do i=memomax+2,min(n,99),2
|
||
|
|
memo(i) = memo(i-2)* float(i)
|
||
|
|
enddo
|
||
|
|
memomax = min(n,99)
|
||
|
|
fact2 = memo(memomax)
|
||
|
|
|
||
|
|
do i=101,n,2
|
||
|
|
fact2 = fact2*float(i)
|
||
|
|
enddo
|
||
|
|
|
||
|
|
end function
|
||
|
|
|
||
|
|
|
||
|
|
double precision function fact(n)
|
||
|
|
implicit none
|
||
|
|
integer :: n
|
||
|
|
double precision, save :: memo(1:100)
|
||
|
|
integer, save :: memomax = 1
|
||
|
|
|
||
|
|
if (n<=memomax) then
|
||
|
|
if (n<2) then
|
||
|
|
fact = 1.d0
|
||
|
|
else
|
||
|
|
fact = memo(n)
|
||
|
|
endif
|
||
|
|
return
|
||
|
|
endif
|
||
|
|
|
||
|
|
integer :: i
|
||
|
|
memo(1) = 1.d0
|
||
|
|
do i=memomax+1,min(n,100)
|
||
|
|
memo(i) = memo(i-1)*float(i)
|
||
|
|
enddo
|
||
|
|
memomax = min(n,100)
|
||
|
|
fact = memo(memomax)
|
||
|
|
do i=101,n
|
||
|
|
fact = fact*float(i)
|
||
|
|
enddo
|
||
|
|
end function
|
||
|
|
|
||
|
|
|