diff --git a/nucleons.jl b/nucleons.jl
index cc82506..b57c66c 100644
--- a/nucleons.jl
+++ b/nucleons.jl
@@ -120,7 +120,7 @@ l_κ(κ::Int) = abs(κ) - (κ < 0) # since true = 1 and false = 0
 "Calculate scalar and vector densities of a nucleon species on [0,r_max] divided into (divs+1) points and returns them as vectors (ρ_s, ρ_v) where
     the arrays κs, Es, occs tabulate the energies and occupation numbers corresponding to each κ,
     the other parameters are defined above"
-function calculateNucleonDensity(κs, Es, occs, p, Φ0, W0, B0, A0, r_max, divs, E_min=0, E_max=(p ? M_p : M_n))
+function calculateNucleonDensity(κs, Es, occs, p, Φ0, W0, B0, A0, r_max, divs)
     ρr2 = zeros(2, divs + 1) # ρ×r² values
     
     for (κ, E, occ) in zip(κs, Es, occs)
@@ -141,8 +141,8 @@ end
 
 "Solve the Dirac equation and calculate scalar and vector densities of a nucleon species where
     the other parameters are defined above"
-function solveNucleonDensity(N, p, Φ0, W0, B0, A0, r_max, divs, E_min=0, E_max=(p ? M_p : M_n))
+function solveNucleonDensity(N, p, Φ0, W0, B0, A0, r_max, divs, E_min=800, E_max=939)
     κs, Es = findAllOrbitals(p, Φ0, W0, B0, A0, r_max, E_min, E_max)
     occs = fillNucleons(N, κs, Es)
-    return calculateNucleonDensity(κs, Es, occs, p, Φ0, W0, B0, A0, r_max, divs, E_min, E_max)
+    return calculateNucleonDensity(κs, Es, occs, p, Φ0, W0, B0, A0, r_max, divs)
 end