squeeze some more space out of VODE's memory footprint

integration/VODE/vode_dvhin.H

@@ -13,7 +13,8 @@
 
 template <typename BurnT, typename DvodeT>
 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
-void dvhin (BurnT& state, DvodeT& vstate, amrex::Real& H0, int& NITER, int& IER)
+int
+dvhin (BurnT& state, DvodeT& vstate, amrex::Real& H0, int& NITER)
 {
     // This routine computes the step size, H0, to be attempted on the
     // first step, when the user has not supplied a value for this.
@@ -34,8 +35,7 @@ void dvhin (BurnT& state, DvodeT& vstate, amrex::Real& H0, int& NITER, int& IER)
 
     if (TDIST < 2.0_rt * TROUND) {
         // Error return for vstate.tout - vstate.t too small.
-        IER = -1;
-        return;
+        return -1;
     }
 
     // Set a lower bound on h based on the roundoff level in vstate.t and vstate.tout.
@@ -132,7 +132,7 @@ void dvhin (BurnT& state, DvodeT& vstate, amrex::Real& H0, int& NITER, int& IER)
     // apply sign
     H0 = std::copysign(H0, vstate.tout - vstate.t);
     NITER = iter;
-    IER = 0;
+    return 0;
 
 }
 

integration/VODE/vode_dvjust.H

@@ -77,7 +77,7 @@ void dvjust (int IORD, BurnT& state, DvodeT& vstate)
                 ALPH0 -= 1.0_rt / (j + 1);
                 ALPH1 += 1.0_rt / XI;
                 for (int iback = 1; iback <= j+1; ++iback) {
-                    int i = (j + 4) - iback;
+                    const int i = (j + 4) - iback;
                     vstate.el(i) = vstate.el(i) * XIOLD + vstate.el(i-1);
                 }
                 XIOLD = XI;

integration/VODE/vode_dvode.H

@@ -27,7 +27,7 @@ int dvode (BurnT& state, DvodeT& vstate)
 
     // Local variables
     amrex::Real H0{}, S{};
-    int IER{}, NITER{};
+    int NITER{};
 
     // Flag determining if we were successful.
 
@@ -84,7 +84,7 @@ int dvode (BurnT& state, DvodeT& vstate)
 
     // Call DVHIN to set initial step size H0 to be attempted.
     H0 = 0.0_rt;
-    dvhin(state, vstate, H0, NITER, IER);
+    auto IER = dvhin(state, vstate, H0, NITER);
     vstate.n_rhs += NITER;
 
     if (IER != 0) {
@@ -189,7 +189,7 @@ int dvode (BurnT& state, DvodeT& vstate)
 
        }
 
-       int kflag = dvstep(state, vstate);
+       const int kflag = dvstep(state, vstate);
 
 
        // Branch on KFLAG. KFLAG can be 0, -1, or -2.

integration/VODE/vode_type.H

@@ -24,11 +24,11 @@ constexpr amrex::Real vode_decrease_change_factor = 0.25_rt;
 
 // For the backward differentiation formula (BDF) integration
 // the maximum order should be no greater than 5.
-constexpr int VODE_MAXORD = 5;
-constexpr int VODE_LMAX = VODE_MAXORD + 1;
+constexpr std::uint8_t VODE_MAXORD = 5;
+constexpr std::uint8_t VODE_LMAX = VODE_MAXORD + 1;
 
 // How many timesteps should pass before refreshing the Jacobian
-constexpr int max_steps_between_jacobian_evals = 50;
+constexpr std::uint8_t max_steps_between_jacobian_evals = 50;
 
 // Type dvode_t contains the integration solution and control variables
 template<int int_neqs>
@@ -86,31 +86,31 @@ struct dvode_t
     //                   (recoverable error)
     //            2 means convergence failure with current Jacobian or
     //                   singular matrix (unrecoverable error)
-    short ICF;
+    std::uint8_t ICF;
 
     // IPUP   = Saved flag to signal updating of Newton matrix
-    short IPUP;
+    bool IPUP;
 
     // JCUR   = Output flag from DVJAC showing Jacobian status:
     //            JCUR = 0 means J is not current
     //            JCUR = 1 means J is current
-    short JCUR;
+    bool JCUR;
 
     // L      = Integer variable, NQ + 1, current order plus one
-    short L;
+    std::uint8_t L;
 
     // NEWH   = Saved integer to flag change of H
-    short NEWH;
+    bool NEWH;
 
     // NEWQ   = The method order to be used on the next step
-    short NEWQ;
+    std::uint8_t NEWQ;
 
     // NQ     = Integer variable, the current integration method order
-    short NQ;
+    std::uint8_t NQ;
 
     // NQWAIT = A counter controlling the frequency of order changes.
     //          An order change is about to be considered if NQWAIT = 1.
-    short NQWAIT;
+    std::uint8_t NQWAIT;
 
     // NSLJ   = The number of steps taken as of the last Jacobian update
     int NSLJ;
@@ -119,7 +119,7 @@ struct dvode_t
     int NSLP;
 
     // jacobian_type = the type of Jacobian to use (1 = analytic, 2 = numerical)
-    short jacobian_type;
+    std::uint8_t jacobian_type;
 
     // EL     = Real array of integration coefficients.  See DVSET
     amrex::Array1D<amrex::Real, 1, VODE_LMAX> el;
@@ -155,7 +155,7 @@ struct dvode_t
 
     amrex::Array1D<amrex::Real, 1, int_neqs> ewt, savf;
 
-    amrex::Array1D<short, 1, int_neqs> pivot;
+    IArray1D pivot;
 
     // Array of size NEQ used for the accumulated corrections on each
     // step, scaled in the output to represent the estimated local

integration/integrator.H

@@ -13,7 +13,7 @@ void integrator_wrapper (BurnT& state, amrex::Real dt)
 {
 
     if constexpr (enable_retry) {
-        burn_t old_state{state};
+        const burn_t old_state{state};
 
         actual_integrator(state, dt);
 

integration/integrator_data.H

@@ -1,12 +1,18 @@
 #ifndef INTEGRATOR_DATA_H
 #define INTEGRATOR_DATA_H
 
+#include <cstdint>
+
 #include <burn_type.H>
 
 // Define the size of the ODE system that VODE will integrate
 
 constexpr int INT_NEQS = NumSpec + 1;
 
+// the datatype we will use for integer quantities that hold
+// indexing into our INT_NEQS
+using smallint_t = short;  //std::conditional_t<(INT_NEQS < 255), std::uint8_t, short>;
+
 // We will use this parameter to determine if a given species
 // abundance is unreasonably small or large (each X must satisfy
 // -failure_tolerance <= X <= 1.0 + failure_tolerance).
@@ -50,7 +56,7 @@ constexpr int integrator_neqs ()
     return int_neqs;
 }
 
-using IArray1D = amrex::Array1D<short, 1, INT_NEQS>;
+using IArray1D = amrex::Array1D<smallint_t, 1, INT_NEQS>;
 using RArray1D = amrex::Array1D<amrex::Real, 1, INT_NEQS>;
 using RArray2D = ArrayUtil::MathArray2D<jac_t, 1, INT_NEQS, 1, INT_NEQS>;