diff --git a/src/lib.rs b/src/lib.rs
index 45121f4..76db817 100644
--- a/src/lib.rs
+++ b/src/lib.rs
@@ -11,18 +11,17 @@ fn neutron_energy_mean_and_std_dev(
 ) -> PyResult<(f64, f64)> {
     // values from Ballabio paper
     let (a_1, a_2, a_3, a_4, mean) = match reaction {
-        Some("D+D=n+He3") => (4.69515, -0.040729, 0.47, 0.81844, 2.4495e6),
-        Some("D+T=n+a") => (5.30509, 2.4736e-3, 1.84, 1.3818, 14.021e6),
+        Some("D+D=n+He3") => (4.69515, -0.040729, 0.47, 0.81844, 2.4486858678216934e6),
+        Some("D+T=n+a") => (5.30509, 0.0024736, 1.84, 1.3818, 14028394.744466662),
         _ => return Err(PyErr::new::<pyo3::exceptions::PyValueError, _>("reaction must be either 'D+D=n+He3' or 'D+T=n+a'")),
     };
 
     let ion_temperature_kev: f64 = scale_temperature_units_to_kev(ion_temperature, temperature_units); // Ballabio equation accepts KeV units
 
-    // units of mean_delta are in ev
-    let mean_delta = a_1 * ion_temperature_kev.powf(2.0 / 3.0) / (1.0 + a_2 * ion_temperature_kev.powf(a_3)) + a_4 * ion_temperature_kev;
+    // units of mean_delta are in put into ev with the 1000 multiplication
+    let mean_delta = 1000.0 *( a_1 * ion_temperature_kev.powf(0.66666666) / (1.0 + a_2 * ion_temperature_kev.powf(a_3)) + a_4 * ion_temperature_kev);
 
     let mean_adjusted = mean + mean_delta;  
-    print!("{}", mean);
 
     let mean_scaled =  scale_energy_in_kev_to_requested_units(mean_adjusted/1e3, neutron_energy_units);
 
diff --git a/tests/test_average_neutron_energy.py b/tests/test_average_neutron_energy.py
index b405a97..d8557c4 100644
--- a/tests/test_average_neutron_energy.py
+++ b/tests/test_average_neutron_energy.py
@@ -22,10 +22,10 @@ def test_mean_energy_with_nesst():
         temperature_units='eV',
         neutron_energy_units='eV'
     )
-    # assert nesst_dt_mean == fnu_dt_mean
-    assert nesst_dd_mean == fnu_dd_mean
-    assert nest_dt_std_dev == fnu_dt_std_dev
-    assert nest_dd_std_dev == fnu_dd_std_dev
+    assert nesst_dt_mean == approx(fnu_dt_mean, rel=1e-6)
+    assert nesst_dd_mean == approx(fnu_dd_mean, rel=1e-6)
+    assert nest_dt_std_dev == approx(fnu_dt_std_dev, rel=1e-6)
+    assert nest_dd_std_dev == approx(fnu_dd_std_dev, rel=1e-6)
 
 def test_mean_energy():
     mean, std_dev = neutron_energy_mean_and_std_dev(