Individual Reactor Initialization - Segmentation Fault Fix.

include/cantera/zeroD/Reactor.h

@@ -193,7 +193,8 @@ class Reactor : public ReactorBase
     //!     `true` for reactors where the pressure is a dependent property,
     //!     calculated from the state, and `false` when the pressure is constant
     //!     or an independent variable.
-    virtual void updateConnected(bool updatePressure);
+    //! @param t0 initialization time for the reactor if not obtained from the network
+    virtual void updateConnected(bool updatePressure, double t0=0.0);
 
     //! Get initial conditions for SurfPhase objects attached to this reactor
     virtual void getSurfaceInitialConditions(double* y);

src/zeroD/Reactor.cpp

@@ -88,7 +88,7 @@ void Reactor::initialize(doublereal t0)
     m_thermo->restoreState(m_state);
     m_sdot.resize(m_nsp, 0.0);
     m_wdot.resize(m_nsp, 0.0);
-    updateConnected(true);
+    updateConnected(true, t0);
 
     for (size_t n = 0; n < m_wall.size(); n++) {
         WallBase* W = m_wall[n];
@@ -177,7 +177,7 @@ void Reactor::updateSurfaceState(double* y)
     }
 }
 
-void Reactor::updateConnected(bool updatePressure) {
+void Reactor::updateConnected(bool updatePressure, double t0) {
     // save parameters needed by other connected reactors
     m_enthalpy = m_thermo->enthalpy_mass();
     if (updatePressure) {
@@ -187,7 +187,7 @@ void Reactor::updateConnected(bool updatePressure) {
     m_thermo->saveState(m_state);
 
     // Update the mass flow rate of connected flow devices
-    double time = m_net->time();
+    double time = (m_net != nullptr) ? m_net->time() : t0;
     for (size_t i = 0; i < m_outlet.size(); i++) {
         m_outlet[i]->updateMassFlowRate(time);
     }

test/SConscript

@@ -290,6 +290,7 @@ addTestProgram('thermo', 'thermo')
 addTestProgram('equil', 'equil')
 addTestProgram('kinetics', 'kinetics')
 addTestProgram('transport', 'transport')
+addTestProgram('zeroD', 'zeroD')
 
 python_subtests = ['']
 test_root = '#interfaces/cython/cantera/test'

test/zeroD/test_zeroD.cpp

@@ -0,0 +1,69 @@
+#include <cstdio>
+#include <string>
+#include "gtest/gtest.h"
+#include "cantera/thermo.h"
+#include "cantera/kinetics.h"
+#include "cantera/zerodim.h"
+
+using namespace Cantera;
+
+// This test ensures that prior reactor initialization of a reactor does
+// not affect later integration within a network. This example was
+// adapted from test_reactor.py::test_equilibrium_HP.
+TEST(ZeroDim, test_individual_reactor_initialization)
+{
+    // Initial conditions
+    double T0 = 1100.0;
+    double P0 = 1013250;
+    std::string X0 = "H2:1.0, O2:0.5, AR:8.0";
+    // Reactor solution, phase, and kinetics objects
+    std::shared_ptr<Solution> sol1 = newSolution("h2o2.yaml");
+    std::shared_ptr<ThermoPhase> gas1 = sol1->thermo();
+    std::shared_ptr<Kinetics> kin1 = sol1->kinetics();
+    gas1->setState_TPX(T0, P0, X0);
+    // Set up reactor object
+    Reactor reactor1;
+    reactor1.setKineticsMgr(*kin1);
+    reactor1.setThermoMgr(*gas1);
+    // Initialize reactor at arbitrary value prior to integration to ensure no impact
+    reactor1.initialize(0.1);
+    // Setup reactor network and integrate
+    ReactorNet network;
+    network.addReactor(reactor1);
+    network.advance(1.0);
+    // Secondary gas for comparison
+    std::shared_ptr<Solution> sol2 = newSolution("h2o2.yaml");
+    std::shared_ptr<ThermoPhase> gas2 = sol2->thermo();
+    std::shared_ptr<Kinetics> kin2 = sol2->kinetics();
+    gas2->setState_TPX(T0, P0, X0);
+    gas2->equilibrate("UV");
+    // Secondary reactor for comparison
+    Reactor reactor2;
+    reactor2.setKineticsMgr(*kin2);
+    reactor2.setThermoMgr(*gas2);
+    reactor2.initialize(0);
+    // Get state of reactors
+    std::vector<double> state1 (reactor1.neq());
+    std::vector<double> state2 (reactor2.neq());
+    reactor1.getState(state1.data());
+    reactor2.getState(state2.data());
+    // Compare the reactors.
+    EXPECT_EQ(reactor1.neq(), reactor2.neq());
+    double tol = 1e-14;
+    EXPECT_NEAR(state1[0], state2[0], tol);
+    EXPECT_NEAR(state1[1], state2[1], tol);
+    for(size_t i = 3; i < reactor1.neq(); i++)
+    {
+        EXPECT_NEAR(state1[i], state2[i], tol);
+    }
+}
+
+int main(int argc, char** argv)
+{
+    printf("Running main() from test_zeroD.cpp\n");
+    testing::InitGoogleTest(&argc, argv);
+    Cantera::make_deprecation_warnings_fatal();
+    int result = RUN_ALL_TESTS();
+    Cantera::appdelete();
+    return result;
+}