diff --git a/include/cantera/zeroD/Reactor.h b/include/cantera/zeroD/Reactor.h
index 94dde02c5..6c31d5de2 100644
--- a/include/cantera/zeroD/Reactor.h
+++ b/include/cantera/zeroD/Reactor.h
@@ -115,6 +115,8 @@ public:
     virtual void evalEqs(doublereal t, doublereal* y,
                          doublereal* ydot, doublereal* params);
 
+    virtual void syncState();
+
     //! Set the state of the reactor to correspond to the state vector *y*.
     virtual void updateState(doublereal* y);
 
diff --git a/include/cantera/zeroD/ReactorBase.h b/include/cantera/zeroD/ReactorBase.h
index 9fbe77c9b..650127015 100644
--- a/include/cantera/zeroD/ReactorBase.h
+++ b/include/cantera/zeroD/ReactorBase.h
@@ -130,6 +130,11 @@ public:
         m_thermo->restoreState(m_state);
     }
 
+    //! Set the state of the reactor to correspond to the state of the
+    //! associated ThermoPhase object. This is the inverse of restoreState().
+    //! Calling this will trigger integrator reinitialization.
+    virtual void syncState();
+
     //! return a reference to the contents.
     thermo_t& contents() {
         return *m_thermo;
diff --git a/include/cantera/zeroD/ReactorNet.h b/include/cantera/zeroD/ReactorNet.h
index 40186d87f..5a3d8b86e 100644
--- a/include/cantera/zeroD/ReactorNet.h
+++ b/include/cantera/zeroD/ReactorNet.h
@@ -36,7 +36,7 @@ public:
      */
     void setInitialTime(doublereal time) {
         m_time = time;
-        m_init = false;
+        m_integrator_init = false;
     }
 
     //! Set the maximum time step.
@@ -205,6 +205,18 @@ public:
         return m_paramNames.at(p);
     }
 
+    //! Reinitialize the integrator. Used to solve a new problem (different
+    //! initial conditions) but with the same configuration of the reactor
+    //! network. Can be called manually, or automatically after calling
+    //! setInitialTime or modifying a reactor's contents.
+    void reinitialize();
+
+    //! Called to trigger integrator reinitialization before further
+    //! integration.
+    void setNeedsReinit() {
+        m_integrator_init = false;
+    }
+
 protected:
     void connect(size_t i, size_t j) {
         m_connect[j*m_reactors.size() + i] = 1;
@@ -225,6 +237,7 @@ protected:
     Integrator* m_integ;
     doublereal m_time;
     bool m_init;
+    bool m_integrator_init; //! True if integrator initialization is current
     size_t m_nv;
 
     //! m_start[n] is the starting point in the state vector for reactor n
diff --git a/interfaces/cython/cantera/_cantera.pxd b/interfaces/cython/cantera/_cantera.pxd
index 00a68cbf9..b7cdf9ce8 100644
--- a/interfaces/cython/cantera/_cantera.pxd
+++ b/interfaces/cython/cantera/_cantera.pxd
@@ -251,6 +251,7 @@ cdef extern from "cantera/zeroD/ReactorBase.h" namespace "Cantera":
         CxxReactorBase()
         void setThermoMgr(CxxThermoPhase&) except +
         void restoreState() except +
+        void syncState() except +
         double volume()
         string name()
         void setName(string)
@@ -327,6 +328,7 @@ cdef extern from "cantera/zeroD/ReactorNet.h":
         void addReactor(CxxReactor&)
         void advance(double) except +
         double step(double) except +
+        void reinitialize() except +
         double time()
         void setInitialTime(double)
         void setTolerances(double, double)
diff --git a/interfaces/cython/cantera/reactor.pyx b/interfaces/cython/cantera/reactor.pyx
index cf4324a47..437674f15 100644
--- a/interfaces/cython/cantera/reactor.pyx
+++ b/interfaces/cython/cantera/reactor.pyx
@@ -44,6 +44,14 @@ cdef class ReactorBase:
         def __set__(self, name):
             self.rbase.setName(stringify(name))
 
+    def syncState(self):
+        """
+        Set the state of the Reactor to match that of the associated
+        `ThermoPhase` object. After calling syncState(), call
+        ReactorNet.reinitialize() before further integration.
+        """
+        self.rbase.syncState()
+
     property thermo:
         """The `ThermoPhase` object representing the reactor's contents."""
         def __get__(self):
@@ -748,6 +756,14 @@ cdef class ReactorNet:
         """
         return self.net.step(t)
 
+    def reinitialize(self):
+        """
+        Reinitialize the integrator after making changing to the state of the
+        system. Changes to Reactor contents will automatically trigger
+        reinitialization.
+        """
+        self.net.reinitialize()
+
     property time:
         """The current time [s]."""
         def __get__(self):
diff --git a/interfaces/cython/cantera/test/test_reactor.py b/interfaces/cython/cantera/test/test_reactor.py
index 3a91e4b6c..86b7b2c7c 100644
--- a/interfaces/cython/cantera/test/test_reactor.py
+++ b/interfaces/cython/cantera/test/test_reactor.py
@@ -474,6 +474,28 @@ class TestReactor(utilities.CanteraTest):
         self.assertNear(p1a, p1b)
         self.assertNear(p2a, p2b)
 
+    def test_reinitialize(self):
+        self.make_reactors(T1=300, T2=1000, independent=False)
+        self.add_wall(U=200, A=1.0)
+        self.net.advance(1.0)
+        T1a = self.r1.T
+        T2a = self.r2.T
+
+        self.r1.thermo.TD = 300, None
+        self.r1.syncState()
+
+        self.r2.thermo.TD = 1000, None
+        self.r2.syncState()
+
+        self.assertNear(self.r1.T, 300)
+        self.assertNear(self.r2.T, 1000)
+        self.net.advance(2.0)
+        T1b = self.r1.T
+        T2b = self.r2.T
+
+        self.assertNear(T1a, T1b)
+        self.assertNear(T2a, T2b)
+
     def test_unpicklable(self):
         self.make_reactors()
         import pickle
diff --git a/src/zeroD/Reactor.cpp b/src/zeroD/Reactor.cpp
index ef1a4357f..3af7ffc4f 100644
--- a/src/zeroD/Reactor.cpp
+++ b/src/zeroD/Reactor.cpp
@@ -116,6 +116,12 @@ size_t Reactor::nSensParams()
     return m_nsens;
 }
 
+void Reactor::syncState()
+{
+    ReactorBase::syncState();
+    m_mass = m_thermo->density() * m_vol;
+}
+
 void Reactor::updateState(doublereal* y)
 {
     for (size_t i = 0; i < m_nv; i++) {
diff --git a/src/zeroD/ReactorBase.cpp b/src/zeroD/ReactorBase.cpp
index 47dc0f06b..9c60e43c7 100644
--- a/src/zeroD/ReactorBase.cpp
+++ b/src/zeroD/ReactorBase.cpp
@@ -7,6 +7,7 @@
 #include "cantera/zeroD/ReactorBase.h"
 #include "cantera/zeroD/FlowDevice.h"
 #include "cantera/zeroD/Wall.h"
+#include "cantera/zeroD/ReactorNet.h"
 
 using namespace std;
 namespace Cantera
@@ -39,6 +40,17 @@ void ReactorBase::setThermoMgr(thermo_t& thermo)
     m_pressure = m_thermo->pressure();
 }
 
+void ReactorBase::syncState()
+{
+    m_thermo->saveState(m_state);
+    m_enthalpy = m_thermo->enthalpy_mass();
+    m_intEnergy = m_thermo->intEnergy_mass();
+    m_pressure = m_thermo->pressure();
+    if (m_net) {
+        m_net->setNeedsReinit();
+    }
+}
+
 void ReactorBase::addInlet(FlowDevice& inlet)
 {
     m_inlet.push_back(&inlet);
diff --git a/src/zeroD/ReactorNet.cpp b/src/zeroD/ReactorNet.cpp
index 49c31da1f..769dd2b86 100644
--- a/src/zeroD/ReactorNet.cpp
+++ b/src/zeroD/ReactorNet.cpp
@@ -12,7 +12,7 @@ namespace Cantera
 {
 
 ReactorNet::ReactorNet() : Cantera::FuncEval(),
-    m_integ(0), m_time(0.0), m_init(false),
+    m_integ(0), m_time(0.0), m_init(false), m_integrator_init(false),
     m_nv(0), m_rtol(1.0e-9), m_rtolsens(1.0e-4),
     m_atols(1.0e-15), m_atolsens(1.0e-4),
     m_maxstep(-1.0), m_maxErrTestFails(0),
@@ -132,9 +132,21 @@ void ReactorNet::initialize()
         writelog(buf);
     }
     m_integ->initialize(m_time, *this);
+    m_integrator_init = true;
     m_init = true;
 }
 
+void ReactorNet::reinitialize()
+{
+    if (m_init) {
+        writelog("Re-initializing reactor network.\n", m_verbose);
+        m_integ->reinitialize(m_time, *this);
+        m_integrator_init = true;
+    } else {
+        initialize();
+    }
+}
+
 void ReactorNet::advance(doublereal time)
 {
     if (!m_init) {
@@ -142,6 +154,8 @@ void ReactorNet::advance(doublereal time)
             m_maxstep = time - m_time;
         }
         initialize();
+    } else if (!m_integrator_init) {
+        reinitialize();
     }
     m_integ->integrate(time);
     m_time = time;
@@ -155,6 +169,8 @@ double ReactorNet::step(doublereal time)
             m_maxstep = time - m_time;
         }
         initialize();
+    } else if (!m_integrator_init) {
+        reinitialize();
     }
     m_time = m_integ->step(time);
     updateState(m_integ->solution());