Merge branch 'development' into ROCK4

.github/workflows/test_rhs.yml

@@ -47,3 +47,122 @@ jobs:
         run: |
           cd unit_test/test_rhs
           diff test.out ci-benchmarks/ignition_simple.out
+
+      - name: Compile, test_rhs (VODE, iso7)
+        run: |
+          cd unit_test/test_rhs
+          make realclean
+          make NETWORK_DIR=iso7 -j 4
+
+      - name: Run test_rhs (VODE, iso7)
+        run: |
+          cd unit_test/test_rhs
+          ./main3d.gnu.ex inputs_iso7
+          ../../external/amrex/Tools/Plotfile/fextrema.gnu.ex react_iso7_test_rhs.VODE > test.out
+
+      - name: Compare to stored output (VODE, iso7)
+        run: |
+          cd unit_test/test_rhs
+          diff test.out ci-benchmarks/iso7.out
+
+      - name: Compile, test_rhs (VODE, aprox13)
+        run: |
+          cd unit_test/test_rhs
+          make realclean
+          make NETWORK_DIR=aprox13 -j 4
+
+      - name: Run test_rhs (VODE, aprox13)
+        run: |
+          cd unit_test/test_rhs
+          ./main3d.gnu.ex inputs_aprox13
+          ../../external/amrex/Tools/Plotfile/fextrema.gnu.ex react_aprox13_test_rhs.VODE > test.out
+
+      - name: Compare to stored output (VODE, aprox13)
+        run: |
+          cd unit_test/test_rhs
+          diff test.out ci-benchmarks/aprox13.out
+
+      - name: Compile, test_rhs (VODE, aprox19)
+        run: |
+          cd unit_test/test_rhs
+          make realclean
+          make NETWORK_DIR=aprox19 -j 4
+
+      - name: Run test_rhs (VODE, aprox19)
+        run: |
+          cd unit_test/test_rhs
+          ./main3d.gnu.ex inputs_aprox19
+          ../../external/amrex/Tools/Plotfile/fextrema.gnu.ex react_aprox19_test_rhs.VODE > test.out
+
+      - name: Compare to stored output (VODE, aprox19)
+        run: |
+          cd unit_test/test_rhs
+          diff test.out ci-benchmarks/aprox19.out
+
+      - name: Compile, test_rhs (VODE, aprox21)
+        run: |
+          cd unit_test/test_rhs
+          make realclean
+          make NETWORK_DIR=aprox21 -j 4
+
+      - name: Run test_rhs (VODE, aprox21)
+        run: |
+          cd unit_test/test_rhs
+          ./main3d.gnu.ex inputs_aprox21
+          ../../external/amrex/Tools/Plotfile/fextrema.gnu.ex react_aprox21_test_rhs.VODE > test.out
+
+      - name: Compare to stored output (VODE, aprox21)
+        run: |
+          cd unit_test/test_rhs
+          diff test.out ci-benchmarks/aprox21.out
+
+      - name: Compile, test_rhs (VODE, rprox)
+        run: |
+          cd unit_test/test_rhs
+          make realclean
+          make NETWORK_DIR=rprox -j 4
+
+      - name: Run test_rhs (VODE, rprox)
+        run: |
+          cd unit_test/test_rhs
+          ./main3d.gnu.ex inputs_rprox
+          ../../external/amrex/Tools/Plotfile/fextrema.gnu.ex react_rprox_test_rhs.VODE > test.out
+
+      - name: Compare to stored output (VODE, rprox)
+        run: |
+          cd unit_test/test_rhs
+          diff test.out ci-benchmarks/rprox.out
+
+      - name: Compile, test_rhs (VODE, powerlaw)
+        run: |
+          cd unit_test/test_rhs
+          make realclean
+          make NETWORK_DIR=powerlaw -j 4
+
+      - name: Run test_rhs (VODE, powerlaw)
+        run: |
+          cd unit_test/test_rhs
+          ./main3d.gnu.ex inputs_powerlaw
+          ../../external/amrex/Tools/Plotfile/fextrema.gnu.ex react_powerlaw_test_rhs.VODE > test.out
+
+      - name: Compare to stored output (VODE, powerlaw)
+        run: |
+          cd unit_test/test_rhs
+          diff test.out ci-benchmarks/powerlaw.out
+
+      - name: Compile, test_rhs (VODE, triple_alpha_plus_cago)
+        run: |
+          cd unit_test/test_rhs
+          make realclean
+          make NETWORK_DIR=triple_alpha_plus_cago -j 4
+
+      - name: Run test_rhs (VODE, triple_alpha_plus_cago)
+        run: |
+          cd unit_test/test_rhs
+          ./main3d.gnu.ex inputs_triple_alpha_plus_cago
+          ../../external/amrex/Tools/Plotfile/fextrema.gnu.ex react_triple_alpha_plus_cago_test_rhs.VODE > test.out
+
+      - name: Compare to stored output (VODE, triple_alpha_plus_cago)
+        run: |
+          cd unit_test/test_rhs
+          diff test.out ci-benchmarks/triple_alpha_plus_cago.out

EOS/helmholtz/actual_eos.H

@@ -1486,6 +1486,8 @@ template <typename I>
 AMREX_GPU_HOST_DEVICE AMREX_INLINE
 bool is_input_valid (I input)
 {
+  amrex::ignore_unused(input);
+
   static_assert(std::is_same<I, eos_input_t>::value, "input must be an eos_input_t");
 
   bool valid = true;

integration/QSS/actual_integrator.H

@@ -6,7 +6,6 @@
 #include <rhs.H>
 #include <actual_rhs.H>
 #include <burn_type.H>
-#include <rate_type.H>
 #include <eos_type.H>
 #include <eos.H>
 #include <extern_parameters.H>

integration/RKC/rkc.H

@@ -125,8 +125,8 @@ int rkclow (BurnT& state, RkcT& rstate)
 
     // note: the original Fortran code only had a single rtol, so
     // I needed to generalize this
-    int mmax = std::round(std::sqrt(std::max(rstate.rtol_spec, rstate.rtol_enuc) /
-                                    (10.0_rt * UROUND)));
+    int mmax = static_cast<int>(std::round(std::sqrt(std::max(rstate.rtol_spec,
+                                                              rstate.rtol_enuc) / (10.0_rt * UROUND))));
     mmax = std::max(mmax, 2);
     bool newspc = true;
     bool jacatt = false;

integration/utils/rkc_util.H

@@ -3,19 +3,19 @@
 
 template <typename BurnT, typename IntT>
 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
-Real rkc_init_dt (BurnT& state, IntT& rstate, const Real max_dt, const Real sprad)
+Real rkc_init_dt (BurnT& state, IntT& rstate, const Real max_timestep, const Real sprad)
 {
 
     // estimate the initial timestep
     //
     // this requires that yn hold the initial state and fn the corresponding RHS
     // it also requires temporary space yjm1 and yjm2
     // finally, we require the initial spectral radius, sprad, and the maximum
-    // timestep (tout - tstart), max_dt.
+    // timestep (tout - tstart), max_timestep.
 
-    Real hmin = 10.0_rt * UROUND * std::max(std::abs(rstate.t), max_dt);
+    Real hmin = 10.0_rt * UROUND * std::max(std::abs(rstate.t), max_timestep);
 
-    Real absh = max_dt;
+    Real absh = max_timestep;
 
     if (sprad * absh > 1.0_rt) {
         absh = 1.0_rt / sprad;
@@ -53,10 +53,10 @@ Real rkc_init_dt (BurnT& state, IntT& rstate, const Real max_dt, const Real spra
     }
     est = absh * std::sqrt(est / INT_NEQS);
 
-    if (0.1_rt * absh < max_dt * std::sqrt(est)) {
+    if (0.1_rt * absh < max_timestep * std::sqrt(est)) {
         absh = std::max(0.1_rt * absh / std::sqrt(est), hmin);
     } else {
-        absh = max_dt;
+        absh = max_timestep;
     }
 
     return absh;
@@ -65,7 +65,7 @@ Real rkc_init_dt (BurnT& state, IntT& rstate, const Real max_dt, const Real spra
 
 template <typename BurnT, typename IntT>
 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
-int rkcrho (BurnT& state, IntT& rstate, const Real max_dt, Real& sprad)
+int rkcrho (BurnT& state, IntT& rstate, const Real max_timestep, Real& sprad)
 {
 
     // RKCRHO attempts to compute a close upper bound, SPRAD, on
@@ -85,7 +85,7 @@ int rkcrho (BurnT& state, IntT& rstate, const Real max_dt, Real& sprad)
     // rstate.fn to the RHS called with yn for nsteps == 0
     // rstate.sprad to the previous eigenvector for nsteps > 0
     //
-    // max_dt is the maximum timestep, typically tout - tbegin
+    // max_timestep is the maximum timestep, typically tout - tbegin
     //
     // we use rstate.yjm1 and rstate.yjm2 as internal storage
 
@@ -94,7 +94,7 @@ int rkcrho (BurnT& state, IntT& rstate, const Real max_dt, Real& sprad)
     // sprad smaller than small = 1/hmax are not
     // interesting because they do not constrain the step size.
 
-    Real small = 1.0_rt / max_dt;
+    Real small = 1.0_rt / max_timestep;
 
     // The initial slope is used as guess when nsteps = 0 and
     // thereafter the last computed eigenvector.  Some care

interfaces/eos_type.H

@@ -375,10 +375,7 @@ std::ostream& print_state (std::ostream& o, T const& eos_state)
     return o;
 }
 
-template <typename T, typename U = typename std::enable_if<std::is_same<T, eos_t>::value ||
-                                                           std::is_same<T, eos_re_t>::value ||
-                                                           std::is_same<T, eos_rep_t>::value
-                                                           >::type>
+template <typename T, typename U = typename std::enable_if<std::is_base_of<eos_base_t, T>::value>::type>
 inline
 std::ostream& operator<< (std::ostream& o, T const& eos_state)
 {

networks/general_null/actual_rhs.H

@@ -4,7 +4,6 @@
 #include <extern_parameters.H>
 #include <actual_network.H>
 #include <jacobian_utilities.H>
-#include <rate_type.H>
 
 using namespace amrex;
 using namespace network_rp;

networks/ignition_chamulak/actual_rhs.H

@@ -11,7 +11,6 @@
 #include <microphysics_math.H>
 #include <actual_network.H>
 #include <burn_type.H>
-#include <rate_type.H>
 #include <screen.H>
 #include <jacobian_utilities.H>
 #include <fundamental_constants.H>

networks/ignition_simple/Make.package

@@ -4,8 +4,8 @@ ifeq ($(USE_REACT),TRUE)
   CEXE_sources += actual_network_data.cpp
   CEXE_headers += actual_network.H
 
-  CEXE_sources += actual_rhs_data.cpp
   CEXE_headers += actual_rhs.H
 
+  USE_RATES     = TRUE
   USE_SCREENING = TRUE
 endif

networks/ignition_simple/actual_network.H

@@ -1,13 +1,17 @@
 #ifndef actual_network_H
 #define actual_network_H
 
+#define NEW_NETWORK_IMPLEMENTATION
+
 #include <AMReX_REAL.H>
 #include <AMReX_Array.H>
 
+using namespace amrex;
+
 #include <fundamental_constants.H>
 #include <network_properties.H>
-
-using namespace amrex;
+#include <rhs_type.H>
+#include <aprox_rates.H>
 
 void actual_network_init();
 
@@ -21,8 +25,107 @@ namespace network
 
 namespace Rates
 {
-    const int NumRates = 1;
-    const int NumGroups = 2;
+    enum NetworkRates
+    {
+        C12_C12_to_Mg24 = 1,
+        NumRates = C12_C12_to_Mg24
+    };
 }
 
+namespace RHS
+{
+    AMREX_GPU_HOST_DEVICE AMREX_INLINE
+    constexpr rhs_t rhs_data (int rate)
+    {
+        using namespace Species;
+        using namespace Rates;
+
+        rhs_t data{};
+
+        data.species_A = -1;
+        data.species_B = -1;
+        data.species_C = -1;
+        data.species_D = -1;
+        data.species_E = -1;
+        data.species_F = -1;
+
+        data.number_A = 0;
+        data.number_B = 0;
+        data.number_C = 0;
+        data.number_D = 0;
+        data.number_E = 0;
+        data.number_F = 0;
+
+        data.exponent_A = 0;
+        data.exponent_B = 0;
+        data.exponent_C = 0;
+        data.exponent_D = 0;
+        data.exponent_E = 0;
+        data.exponent_F = 0;
+
+        data.forward_branching_ratio = 1.0_rt;
+        data.reverse_branching_ratio = 1.0_rt;
+
+        data.apply_identical_particle_factor = 1;
+
+        data.rate_can_be_tabulated = 1;
+
+        data.screen_forward_reaction = 1;
+        data.screen_reverse_reaction = 1;
+
+        data.additional_reaction_1 = -1;
+        data.additional_reaction_2 = -1;
+        data.additional_reaction_3 = -1;
+
+        switch (rate) {
+
+        case C12_C12_to_Mg24:
+            data.species_A = C12;
+            data.species_D = Mg24;
+
+            data.number_A = 2;
+            data.number_D = 1;
+
+            data.exponent_A = 2;
+            data.exponent_D = 1;
+
+            // This network doesn't allow the reverse reaction
+            data.reverse_branching_ratio = 0.0_rt;
+            break;
+        }
+
+        return data;
+    }
+
+    template<int rate>
+    AMREX_GPU_HOST_DEVICE AMREX_INLINE
+    void evaluate_analytical_rate (const rhs_state_t& state, rate_t& rates)
+    {
+        using namespace Species;
+        using namespace Rates;
+
+        if constexpr (rate == C12_C12_to_Mg24) {
+            rate_c12c12(state.tf, 1.0_rt, rates.fr, rates.frdt, rates.rr, rates.rrdt);
+        }
+    }
+
+    template<int rate>
+    AMREX_GPU_HOST_DEVICE AMREX_INLINE
+    void postprocess_rate (const rhs_state_t& state, rate_t& rates,
+                           rate_t& rates1, [[maybe_unused]] rate_t& rates2, [[maybe_unused]] rate_t& rates3)
+    {
+        using namespace Species;
+        using namespace Rates;
+
+        // Nothing to do for this network.
+    }
+
+    template<int spec>
+    AMREX_GPU_HOST_DEVICE AMREX_INLINE
+    Real ener_gener_rate (Real const& dydt)
+    {
+        return dydt * network::mion(spec) * C::Legacy::enuc_conv2;
+    }
+} // namespace RHS
+
 #endif