# Simple Gaussian benchmark, with units
# Expected 2-body bound state: E₂ ≈ -1.00 MeV

include("../Hamiltonian.jl")
mode = cpu_tensor
T = Float32 # single-precision mode

V_0 = -38.4671 # MeV
r_0 = 2 # fm

V_benchmark(r2) = V_0 * exp(−r2 / r_0^2)

m_N = 939.565 # MeV/c²
ħ = 197.327 # MeV·fm/c

μ = (m_N / 2) / ħ^2

ϕ = 0
d = 3
n = 2
N = 80
L = 40 # fm
n_imag = 1

s = system{T}(d, n, N, L, μ)
H = Hamiltonian{T}(s, V_benchmark, ϕ, n_imag, mode)
@time evals, _, info = eig(H, 3, resonances = false)

sort!(evals)
display(evals[1:3])

# Result = -1.001821 MeV