diff --git a/src/nuclide.cpp b/src/nuclide.cpp
index b64a1681ab..ee1a35111f 100644
--- a/src/nuclide.cpp
+++ b/src/nuclide.cpp
@@ -893,6 +893,8 @@ void check_data_version(hid_t file_id)
 //==============================================================================
 
 extern "C" void extend_nuclides();
+extern "C" void nuclide_from_hdf5(hid_t group, const Nuclide* ptr,
+  const double* temps, int n, int n_nuclide);
 
 extern "C" int openmc_load_nuclide(const char* name)
 {
@@ -917,6 +919,11 @@ extern "C" int openmc_load_nuclide(const char* name)
       int i_nuclide = data::nuclides.size();
       data::nuclides.push_back(std::make_unique<Nuclide>(
         group, temperature, i_nuclide));
+
+      // Read from Fortran too
+      nuclide_from_hdf5(group, data::nuclides.back().get(),
+        &temperature.front(), temperature.size(), i_nuclide + 1);
+
       close_group(group);
       file_close(file_id);
 
diff --git a/tests/regression_tests/deplete/example_geometry.py b/tests/regression_tests/deplete/example_geometry.py
index ca10c1f725..e54c087f67 100644
--- a/tests/regression_tests/deplete/example_geometry.py
+++ b/tests/regression_tests/deplete/example_geometry.py
@@ -47,10 +47,8 @@ def generate_initial_number_density():
     fuel_dict['U238'] = 2.21371e22
     fuel_dict['O16'] = 4.62954e22
     fuel_dict['O17'] = 1.127684e20
-    fuel_dict['I135'] = 1.0e10
     fuel_dict['Xe135'] = 1.0e10
     fuel_dict['Xe136'] = 1.0e10
-    fuel_dict['Cs135'] = 1.0e10
     fuel_dict['Gd156'] = 1.0e10
     fuel_dict['Gd157'] = 1.0e10
     # fuel_dict['O18'] = 9.51352e19 # Does not exist in ENDF71, merged into 17
@@ -62,10 +60,8 @@ def generate_initial_number_density():
     fuel_gd_dict['Gd156'] = 3.95517E+10
     fuel_gd_dict['Gd157'] = 1.08156e20
     fuel_gd_dict['O16'] = 4.64035e22
-    fuel_dict['I135'] = 1.0e10
     fuel_dict['Xe136'] = 1.0e10
     fuel_dict['Xe135'] = 1.0e10
-    fuel_dict['Cs135'] = 1.0e10
     # There are a whole bunch of 1e-10 stuff here.
 
     # Concentration to be used for cladding
diff --git a/tests/regression_tests/deplete/test_reference.h5 b/tests/regression_tests/deplete/test_reference.h5
index 8750924b11..2f1ef375a9 100644
Binary files a/tests/regression_tests/deplete/test_reference.h5 and b/tests/regression_tests/deplete/test_reference.h5 differ
diff --git a/tests/unit_tests/test_capi.py b/tests/unit_tests/test_capi.py
index e69d15a270..eb9a577a59 100644
--- a/tests/unit_tests/test_capi.py
+++ b/tests/unit_tests/test_capi.py
@@ -390,7 +390,11 @@ def test_restart(capi_init):
 
 
 def test_load_nuclide(capi_init):
+    # load multiple nuclides
     openmc.capi.load_nuclide('H3')
+    assert 'H3' in openmc.capi.nuclides
     openmc.capi.load_nuclide('Pu239')
+    assert 'Pu239' in openmc.capi.nuclides
+    # load non-existent nuclide
     with pytest.raises(exc.DataError):
         openmc.capi.load_nuclide('Pu3')
diff --git a/tests/unit_tests/test_deplete_resultslist.py b/tests/unit_tests/test_deplete_resultslist.py
index 994ee531af..5a066be6e4 100644
--- a/tests/unit_tests/test_deplete_resultslist.py
+++ b/tests/unit_tests/test_deplete_resultslist.py
@@ -20,8 +20,7 @@ def test_get_atoms(res):
     t, n = res.get_atoms("1", "Xe135")
 
     t_ref = [0.0, 1296000.0, 2592000.0, 3888000.0]
-    n_ref = [6.674732823364922e+08, 4.022377919003244e+14,
-             3.469325382312551e+14, 3.639463751297369e+14]
+    n_ref = [6.67473282e+08, 3.76986925e+14, 3.68587383e+14, 3.91338675e+14]
 
     np.testing.assert_allclose(t, t_ref)
     np.testing.assert_allclose(n, n_ref)
@@ -32,13 +31,11 @@ def test_get_reaction_rate(res):
     t, r = res.get_reaction_rate("1", "Xe135", "(n,gamma)")
 
     t_ref = [0.0, 1296000.0, 2592000.0, 3888000.0]
-    n_ref = [6.674732823364922e+08, 4.022377919003244e+14,
-             3.469325382312551e+14, 3.639463751297369e+14]
-    xs_ref = np.array([3.272630849911638e-05, 2.663333206501429e-05,
-                       3.378999816491878e-05, 3.277013860196171e-05])
+    n_ref = [6.67473282e+08, 3.76986925e+14, 3.68587383e+14, 3.91338675e+14]
+    xs_ref = [3.32282266e-05, 2.76207120e-05, 4.10986677e-05, 3.72453665e-05]
 
     np.testing.assert_allclose(t, t_ref)
-    np.testing.assert_allclose(r, n_ref * xs_ref)
+    np.testing.assert_allclose(r, np.array(n_ref) * xs_ref)
 
 
 def test_get_eigenvalue(res):
@@ -46,8 +43,7 @@ def test_get_eigenvalue(res):
     t, k = res.get_eigenvalue()
 
     t_ref = [0.0, 1296000.0, 2592000.0, 3888000.0]
-    k_ref = [1.2162584685809341, 1.1021840248590875, 1.223040353281542,
-             1.2360979537312164]
+    k_ref = [1.16984322, 1.19097427, 1.03012572, 1.20045627]
 
     np.testing.assert_allclose(t, t_ref)
     np.testing.assert_allclose(k, k_ref)