Merge branch 'development' into prob_params_cleaning

AMReX-Astro · Jan 6, 2024 · feb54ef · feb54ef
2 parents d8d0832 + 8e7784c
commit feb54ef
Show file tree

Hide file tree

Showing 2 changed files with 28 additions and 9 deletions.
diff --git a/Exec/science/massive_star/problem_initialize_state_data.H b/Exec/science/massive_star/problem_initialize_state_data.H
@@ -7,6 +7,7 @@
 #ifdef NSE_TABLE
 #include <nse_table_type.H>
 #include <nse_table_check.H>
+#include <nse_eos.H>
 #endif
 
 AMREX_GPU_HOST_DEVICE AMREX_INLINE
@@ -128,8 +129,33 @@ void problem_initialize_state_data (int i, int j, int k,
     }
 
     if (problem::interpolate_pres == 1) {
-        eos(eos_input_rp, eos_state);
-        state(i,j,k,UTEMP) = eos_state.T;
+
+	// we need to get T from P, rho, but consistent with NSE (if
+	// we are in NSE)
+	if (nse_check) {
+	    // this will also change the composition, since the new T
+	    // alters the NSE
+	    nse_T_abar_from_p(eos_state.rho, eos_state.p, eos_state.aux[AuxZero::iye],
+			      eos_state.T, eos_state.aux[AuxZero::iabar]);
+	    state(i,j,k,UTEMP) = eos_state.T;
+	    state(i,j,k,UFX+AuxZero::iabar) = eos_state.aux[AuxZero::iabar];
+
+	    // now call the EOS with the new T and abar to get e
+	    eos(eos_input_rt, eos_state);
+
+	    // finally, get the updated B/A
+	    nse_table_t nse_state;
+	    nse_state.T = state(i,j,k,UTEMP);
+	    nse_state.rho = state(i,j,k,URHO);
+	    nse_state.Ye = state(i,j,k,UFX+AuxZero::iye);
+	    nse_interp(nse_state);
+
+	    state(i,j,k,UFX+AuxZero::ibea) = nse_state.bea;
+
+	} else {
+	    eos(eos_input_rp, eos_state);
+	    state(i,j,k,UTEMP) = eos_state.T;
+	}
     } else {
         eos(eos_input_rt, eos_state);
     }

diff --git a/Source/driver/Castro.cpp b/Source/driver/Castro.cpp
@@ -324,13 +324,6 @@ Castro::read_params ()
     }
 #endif
 
-#ifdef REACTIONS
-#ifdef SIMPLIFIED_SDC
-    if (jacobian == 1) {
-      amrex::Abort("Simplified SDC requires the numerical Jacobian now (jacobian = 2)");
-    }
-#endif
-#endif
     // sanity checks
 
     if (grown_factor < 1) {