Changed default E_min and E_max

This commit is contained in:
Nuwan Yapa 2025-02-18 18:22:58 -05:00
parent 459d35e08e
commit f49207fc42
5 changed files with 8 additions and 17 deletions

View File

@ -87,10 +87,10 @@ function boundaryValueFunc(κ, p, s::system; dtype=Float64, algo=Tsit5())
return func return func
end end
"Find all bound energies between E_min (=0) and E_max (=mass) where "Find all bound energies between E_min (=850) and E_max (=938) where
tol_digits determines the precision for root finding and the threshold for identifying duplicates, tol_digits determines the precision for root finding and the threshold for identifying duplicates,
the other parameters are the same from dirac!(...)." the other parameters are the same from dirac!(...)."
function findEs(κ, p, s::system, E_min=0, E_max=(p ? M_p : M_n), tol_digits=5) function findEs(κ, p, s::system, E_min=850, E_max=938, tol_digits=5)
func = boundaryValueFunc(κ, p, s) func = boundaryValueFunc(κ, p, s)
Es = find_zeros(func, (E_min, E_max); xatol=1/10^tol_digits) Es = find_zeros(func, (E_min, E_max); xatol=1/10^tol_digits)
return unique(E -> round(E; digits=tol_digits), Es) return unique(E -> round(E; digits=tol_digits), Es)
@ -98,7 +98,7 @@ end
"Find all orbitals and return two lists containing κ values and corresponding energies for a single species where "Find all orbitals and return two lists containing κ values and corresponding energies for a single species where
the other parameters are defined above" the other parameters are defined above"
function findAllOrbitals(p, s::system, E_min=0, E_max=(p ? M_p : M_n)) function findAllOrbitals(p, s::system, E_min=850, E_max=938)
κs = Int[] κs = Int[]
Es = Float64[] Es = Float64[]
# start from κ=-1 and go both up and down # start from κ=-1 and go both up and down
@ -160,7 +160,7 @@ end
"Solve the Dirac equation and calculate scalar and vector densities of a nucleon species where "Solve the Dirac equation and calculate scalar and vector densities of a nucleon species where
the other parameters are defined above" the other parameters are defined above"
function solveNucleonDensity(p, s::system, E_min=800, E_max=939) function solveNucleonDensity(p, s::system, E_min=850, E_max=938)
κs, Es = findAllOrbitals(p, s, E_min, E_max) κs, Es = findAllOrbitals(p, s, E_min, E_max)
occs = fillNucleons(Z_or_N(s, p), κs, Es) occs = fillNucleons(Z_or_N(s, p), κs, Es)
return calculateNucleonDensity(κs, Es, occs, p, s) return calculateNucleonDensity(κs, Es, occs, p, s)

View File

@ -20,10 +20,7 @@ s.W0 = Vs
s.B0 = Rs s.B0 = Rs
s.A0 = As s.A0 = As
E_min = 860 κs, Es = findAllOrbitals(p, s)
E_max = 939
κs, Es = findAllOrbitals(p, s, E_min, E_max)
Ebinds = (p ? M_p : M_n) .- Es Ebinds = (p ? M_p : M_n) .- Es
occ = fillNucleons(Z_or_N(s, p), κs, Es) occ = fillNucleons(Z_or_N(s, p), κs, Es)

View File

@ -22,7 +22,7 @@ s.B0 = Rs
s.A0 = As s.A0 = As
E_min = 850 E_min = 850
E_max = 939 E_max = 938
boundEs = findEs(κ, p, s, E_min, E_max) boundEs = findEs(κ, p, s, E_min, E_max)

View File

@ -22,10 +22,7 @@ s.W0 = Vs
s.B0 = Rs s.B0 = Rs
s.A0 = As s.A0 = As
E_min = 800 κs, Es = findAllOrbitals(p, s)
E_max = 939
κs, Es = findAllOrbitals(p, s, E_min, E_max)
Ebinds = (p ? M_p : M_n) .- Es Ebinds = (p ? M_p : M_n) .- Es
occ = fillNucleons(Z_or_N(s, p), κs, Es) occ = fillNucleons(Z_or_N(s, p), κs, Es)

View File

@ -21,10 +21,7 @@ s.W0 = Vs
s.B0 = Rs s.B0 = Rs
s.A0 = As s.A0 = As
E_min = 800 groundE = findEs(κ, p, s) |> minimum
E_max = 939
groundE = findEs(κ, p, s, E_min, E_max) |> minimum
println("ground state E = $groundE") println("ground state E = $groundE")
wf = solveNucleonWf(κ, p, groundE, s) wf = solveNucleonWf(κ, p, groundE, s)