Enforce chabrier1998 or null screening with NSE_NET (#1398)

Make.Microphysics_extern

@@ -145,6 +145,13 @@ ifeq ($(USE_NSE_TABLE), TRUE)
   endif
 endif
 
+# Only allow null chabrier1998 screening with NSE_NET
+ifeq ($(USE_NSE_NET),TRUE)
+  ifneq ($(or $(filter-out null chabrier1998,$(SCREEN_METHOD)),),)
+    $(error USE_NSE_NET can only be used with null or chabrier1998 screening)
+  endif
+endif
+
 clean::
 	@if [ -L helm_table.dat ]; then rm -f helm_table.dat; fi
 	@if [ -L reaclib_rate_metadata.dat ]; then rm -f reaclib_rate_metadata.dat; fi

nse_solver/make_table/GNUmakefile

@@ -24,6 +24,8 @@ EOS_DIR     := helmholtz
 # This sets the network directory
 NETWORK_DIR := aprox21
 
+SCREEN_METHOD := chabrier1998
+
 CONDUCTIVITY_DIR := stellar
 
 INTEGRATOR_DIR =  VODE

nse_solver/nse_solver.H

@@ -89,10 +89,10 @@ void apply_nse_exponent(T& nse_state) {
   Real T_in = nse_state.T_fixed > 0.0_rt ? nse_state.T_fixed : nse_state.T;
   Real exponent;
 
-  // if we use a non-null screening method.
+  // if we use chabrier1998 screening
   // Get the required terms to calculate coulomb correction term, u_c
 
-#if SCREEN_METHOD != 0
+#if SCREEN_METHOD == 4
 
   // Find n_e for original state;
   const amrex::Real n_e = nse_state.rho * nse_state.y_e / C::m_u;
@@ -112,8 +112,8 @@ void apply_nse_exponent(T& nse_state) {
 #endif
       // term for calculating u_c
 
-      // if use non-null screening, calculate the coulomb correction term
-#if SCREEN_METHOD != 0
+      // if use chabrier1998 screening, calculate the coulomb correction term
+#if SCREEN_METHOD == 4
       gamma = std::pow(zion[n], 5.0_rt/3.0_rt) * Gamma_e;
 
       // chemical potential for coulomb correction

unit_test/test_nse/GNUmakefile

@@ -26,6 +26,8 @@ EOS_DIR     := helmholtz
 # This sets the network directory
 NETWORK_DIR := aprox21
 
+SCREEN_METHOD := chabrier1998
+
 CONDUCTIVITY_DIR := stellar
 
 INTEGRATOR_DIR =  VODE