mirror of https://gitlab.com/scemama/eplf
149 lines
2.5 KiB
Fortran
149 lines
2.5 KiB
Fortran
|
|
!recursive double precision function Boys(x,n) result(res)
|
|
! implicit none
|
|
! include 'constants.F'
|
|
!
|
|
! real, intent(in) :: x
|
|
! integer, intent(in) :: n
|
|
!
|
|
! ASSERT (x > 0.)
|
|
! if (n == 0) then
|
|
! res = sqrt(pi/(4.*x))*erf(sqrt(x))
|
|
! else
|
|
! res = (dble(2*n-1) * Boys(x,(n-1)) - exp(-x) )/(2.*x)
|
|
! endif
|
|
!
|
|
!end function
|
|
|
|
|
|
double precision function binom(k,n)
|
|
implicit none
|
|
integer, intent(in) :: k,n
|
|
double precision, save :: memo(0:20,0:20)
|
|
!DEC$ ATTRIBUTES ALIGN : 32 :: memo
|
|
integer, save :: ifirst
|
|
double precision :: fact
|
|
if (ifirst /= 0) then
|
|
|
|
if (n<20) then
|
|
binom = memo(k,n)
|
|
return
|
|
else
|
|
binom=fact(n)/(fact(k)*fact(n-k))
|
|
return
|
|
endif
|
|
|
|
else
|
|
|
|
ifirst = 1
|
|
integer :: i,j
|
|
do j=0,20
|
|
do i=0,j
|
|
memo(i,j) = fact(j)/(fact(i)*fact(j-i))
|
|
enddo
|
|
enddo
|
|
|
|
if (n<20) then
|
|
binom = memo(k,n)
|
|
return
|
|
else
|
|
binom=fact(n)/(fact(k)*fact(n-k))
|
|
endif
|
|
endif
|
|
end
|
|
|
|
double precision function fact2(n)
|
|
implicit none
|
|
integer :: n
|
|
double precision, save :: memo(100)
|
|
integer, save :: memomax = 1
|
|
|
|
ASSERT (mod(n,2) /= 0)
|
|
memo(1) = 1.d0
|
|
if (n<2) then
|
|
fact2 = 1.d0
|
|
else if (n<=memomax) then
|
|
fact2 = memo(n)
|
|
else
|
|
|
|
integer :: i
|
|
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 *= float(i)
|
|
enddo
|
|
|
|
endif
|
|
|
|
end function
|
|
|
|
|
|
double precision function fact(n)
|
|
implicit none
|
|
integer :: n
|
|
double precision, save :: memo(0:100)
|
|
integer, save :: memomax = 1
|
|
|
|
memo(0) = 1.d0
|
|
memo(1) = 1.d0
|
|
|
|
if (n<=memomax) then
|
|
fact = memo(n)
|
|
else
|
|
integer :: i
|
|
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 *= float(i)
|
|
enddo
|
|
endif
|
|
end function
|
|
|
|
|
|
double precision function rintgauss(n)
|
|
implicit none
|
|
include 'constants.F'
|
|
|
|
integer :: n
|
|
|
|
rintgauss = dsqrt_pi
|
|
if ( n == 0 ) then
|
|
return
|
|
else if ( n == 1 ) then
|
|
rintgauss = 0.
|
|
else if ( mod(n,2) == 1) then
|
|
rintgauss = 0.
|
|
else
|
|
double precision :: fact2
|
|
rintgauss = rintgauss/(2.**(n/2))
|
|
rintgauss = rintgauss * fact2(n-1)
|
|
endif
|
|
end function
|
|
|
|
double precision function goverlap(gamA,gamB,nA)
|
|
implicit none
|
|
|
|
real :: gamA, gamB
|
|
integer :: nA(3)
|
|
|
|
double precision :: gamtot
|
|
gamtot = gamA+gamB
|
|
|
|
goverlap=1.0
|
|
|
|
integer :: l
|
|
double precision :: rintgauss
|
|
do l=1,3
|
|
goverlap = goverlap * rintgauss(nA(l)+nA(l))/ (gamtot**(0.5+float(nA(l))))
|
|
enddo
|
|
|
|
end function
|
|
|