yapanuwan
/
BergEC-jl
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Numerical accuracy

This commit is contained in:
Nuwan Yapa 2024-09-23 14:37:30 -04:00
parent cc774b05da
commit e5477bfeb1
2 changed files with 3 additions and 8 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

9
math.jl
View File

@ -22,13 +22,8 @@ gauss_int(a, n) = double_factorial(n - 1) / (2 * a)^((n + 1)/2) * (iseven(n) ? s
"Gives ∫dp p u' u where u' and u' may have different l (Ref: worked out in Mathematica)" "Gives ∫dp p u' u where u' and u' may have different l (Ref: worked out in Mathematica)"
function integral(np, lp, n, l, μω) function integral(np, lp, n, l, μω)
s = 0.0 s = [(-1)^(m + mp) * gauss_int(1, 2 * m + 2 * mp + l + lp + 3) * N_lnk(l, n, m) * N_lnk(lp, np, mp) for (m, mp) in Iterators.product(0:n, 0:np)]
for m in 0 : n return 2 * sqrt(μω) * better_sum(s)
for mp in 0 : np
s += (-1)^(m + mp) * gauss_int(1, 2 * m + 2 * mp + l + lp + 3) / (factorial(n - m) * gamma(l + 3/2 + m) * factorial(m) * factorial(np - mp) * gamma(lp + 3/2 + mp) * factorial(mp))
end
end
return 2 * sqrt(μω * factorial(n) * factorial(np) * gamma(n + l + 3/2) * gamma(np + lp + 3/2)) * s
end end
"Gives <n' l'|| p ||n l> for the HO basis" "Gives <n' l'|| p ||n l> for the HO basis"

2
test/racahs_reduction_formula.jl
View File

@ -9,7 +9,7 @@ function explicit_integral(np, lp, n, l, μω)
return integral return integral
end end
n_list = 0 : 5 n_list = 0 : 10
l_list = 0 : 5 l_list = 0 : 5
μω_list = [1, 2.66] μω_list = [1, 2.66]