Compare commits

...

1 Commits

Author SHA1 Message Date
ysyapa 0f3bcdc95e Interop with COSMO (same results) 2024-04-16 19:41:55 +00:00
1 changed files with 3 additions and 11 deletions

View File

@ -1,5 +1,4 @@
using SparseArrays using SparseArrays
using NuclearToolkit
using SpecialFunctions using SpecialFunctions
include("helper.jl") include("helper.jl")
@ -84,15 +83,6 @@ function sp_V_matrix(V_l, ns, ls; mask=trues(length(ns),length(ns)), dtype=Float
end end
function Moshinsky_transform(Es, n1s, l1s, n2s, l2s, Λ) function Moshinsky_transform(Es, n1s, l1s, n2s, l2s, Λ)
E_max = maximum(Es)
j_max = 2 * E_max + 1
l_max = j_max
to = 0 # unused
dtri = NuclearToolkit.prep_dtri(l_max + 1);
dcgm0 = NuclearToolkit.prep_dcgm0(l_max);
d6j = nothing # NuclearToolkit.prep_d6j_int(E_max, j_max, to);
mat = spzeros(length(Es), length(Es)) mat = spzeros(length(Es), length(Es))
s = hcat(Es, n1s, l1s, n2s, l2s) s = hcat(Es, n1s, l1s, n2s, l2s)
for idx in CartesianIndices(mat) for idx in CartesianIndices(mat)
@ -100,9 +90,11 @@ function Moshinsky_transform(Es, n1s, l1s, n2s, l2s, Λ)
(Elhs, N, L, n, l) = s[i, :] (Elhs, N, L, n, l) = s[i, :]
(Erhs, n1, l1, n2, l2) = s[j, :] (Erhs, n1, l1, n2, l2) = s[j, :]
if Elhs == Erhs && triangle_ineq(L, l, Λ) && triangle_ineq(l1, l2, Λ) if Elhs == Erhs && triangle_ineq(L, l, Λ) && triangle_ineq(l1, l2, Λ)
mat[i, j] = NuclearToolkit.HObracket_d6j(N, L, n, l, n1, l1, n2, l2, Λ, 1.0, dtri, dcgm0, d6j) mat[i, j] = get_bracket(N, L, n, l, n1, l1, n2, l2, Λ, 1.0)
end end
end end
return mat return mat
end end
get_bracket(N, L, n, l, n1, l1, n2, l2, Λ, d) = @ccall "../cosmo/bin/shared.so"._Z24MoshinskyBracketsWrapperiiiiiiiiid(n::Cint,l::Cint,N::Cint,L::Cint,n1::Cint,l1::Cint,n2::Cint,l2::Cint,Λ::Cint,d::Cdouble)::Cdouble