Working on PR updates, have to fix a few more bugs

in regards to equations, handling volumes, etc.

Author: anthony-walker

doc/sphinx/develop/index.md

@@ -1,6 +1,7 @@
 # Develop
 
 (sec-compiling)=
+
 ## Compiling Cantera from Source
 
 If you're interested in contributing new features to Cantera, or you want to try the
@@ -12,9 +13,9 @@ Cantera](compiling/configure-build) on your computer.
 
 The following additional references may also be useful:
 
-- [](compiling/dependencies.md)
-- [](compiling/config-options)
-- [](compiling/special-cases)
+-   [](compiling/dependencies.md)
+-   [](compiling/config-options)
+-   [](compiling/special-cases)
 
 ```{toctree}
 :caption: Compiling Cantera from Source
@@ -35,7 +36,7 @@ compiling/special-cases
 This section is a work in progress.
 ```
 
-- [](reactor-integration)
+-   [](reactor-integration)
 
 ```{toctree}
 :caption: How Cantera Works
@@ -47,13 +48,13 @@ reactor-integration
 
 ## Adding New Features to Cantera
 
-- [](CONTRIBUTING)
-- [](style-guidelines)
-- [](vscode-tips)
-- [](writing-tests)
-- [](running-tests)
-- [](writing-examples)
-- [](doc-formatting)
+-   [](CONTRIBUTING)
+-   [](style-guidelines)
+-   [](vscode-tips)
+-   [](writing-tests)
+-   [](running-tests)
+-   [](writing-examples)
+-   [](doc-formatting)
 
 ```{toctree}
 :caption: Adding New Features to Cantera

include/cantera/zeroD/FlowDevice.h

@@ -141,9 +141,10 @@ class FlowDevice
     //! @warning This function is an experimental part of the %Cantera API and may be
     //! changed
     //! or removed without notice.
-    //! @since New in %Cantera 3.0.
+    //! @since New in %Cantera 3.1.
     //!
-    virtual void buildReactorJacobian(ReactorBase* r, vector<Eigen::Triplet<double>>& jacVector) {
+    virtual void buildReactorJacobian(ReactorBase* r,
+        vector<Eigen::Triplet<double>>& jacVector) {
         throw NotImplementedError(type() + "::buildReactorJacobian");
     }
 
@@ -155,29 +156,12 @@ class FlowDevice
     //! @warning This function is an experimental part of the %Cantera API and may be
     //! changed
     //! or removed without notice.
-    //! @since New in %Cantera 3.0.
+    //! @since New in %Cantera 3.1.
     //!
     virtual void buildNetworkJacobian(vector<Eigen::Triplet<double>>& jacVector) {
-        if (!m_jac_calculated) {
-            throw NotImplementedError(type() + "::buildNetworkJacobian");
-        }
+        throw NotImplementedError(type() + "::buildNetworkJacobian");
     }
 
-    //! Specify the jacobian terms have been calculated and should not be recalculated.
-    //! @warning This function is an experimental part of the %Cantera API and may be
-    //! changed
-    //! or removed without notice.
-    //! @since New in %Cantera 3.0.
-    //!
-    void jacobianCalculated() { m_jac_calculated = true; };
-
-    //! Specify that jacobian terms have not been calculated and should be recalculated.
-    //! @warning This function is an experimental part of the %Cantera API and may be changed
-    //! or removed without notice.
-    //! @since New in %Cantera 3.0.
-    //!
-    void jacobianNotCalculated() { m_jac_calculated = false; };
-
 protected:
     string m_name;  //!< Flow device name.
     bool m_defaultNameSet = false;  //!< `true` if default name has been previously set.

include/cantera/zeroD/IdealGasConstPressureMoleReactor.h

@@ -47,9 +47,7 @@ class IdealGasConstPressureMoleReactor : public ConstPressureMoleReactor
 
     bool preconditionerSupported() const override { return true; };
 
-    double moleDerivative(size_t index) override;
-
-    double moleRadiationDerivative(size_t index) override;
+    double temperature_ddni(size_t index) override;
 
     size_t speciesOffset() const override { return m_sidx; };
 

include/cantera/zeroD/IdealGasMoleReactor.h

@@ -43,9 +43,7 @@ class IdealGasMoleReactor : public MoleReactor
 
     bool preconditionerSupported() const override {return true;};
 
-    double moleDerivative(size_t index) override;
-
-    double moleRadiationDerivative(size_t index) override;
+    double temperature_ddni(size_t index) override;
 
     size_t speciesOffset() const override { return m_sidx; };
 

include/cantera/zeroD/Reactor.h

@@ -70,12 +70,6 @@ class Reactor : public ReactorBase
         m_chem = cflag;
     }
 
-    //! Returns `true` if changes in the reactor composition due to chemical reactions
-    //! are enabled.
-    bool chemistryEnabled() const {
-        return m_chem;
-    }
-
     void setEnergy(int eflag=1) override {
         if (eflag > 0) {
             m_energy = true;
@@ -84,11 +78,6 @@ class Reactor : public ReactorBase
         }
     }
 
-    //! Returns `true` if solution of the energy equation is enabled.
-    bool energyEnabled() const {
-        return m_energy;
-    }
-
     //! Number of equations (state variables) for this reactor
     size_t neq() {
         if (!m_nv) {
@@ -195,15 +184,7 @@ class Reactor : public ReactorBase
     //!
     //! @warning  This method is an experimental part of the %Cantera
     //! API and may be changed or removed without notice.
-    virtual Eigen::SparseMatrix<double> jacobian() {
-        m_jac_trips.clear();
-        // Add before, during, after evals
-        buildJacobian(m_jac_trips);
-        // construct jacobian from vector
-        Eigen::SparseMatrix<double> jac(m_nv, m_nv);
-        jac.setFromTriplets(m_jac_trips.begin(), m_jac_trips.end());
-        return jac;
-    }
+    virtual Eigen::SparseMatrix<double> jacobian();
 
     //! Calculate the Jacobian of a specific Reactor specialization.
     //! @param jacVector vector where jacobian triplets are added
@@ -321,8 +302,6 @@ class Reactor : public ReactorBase
 
     vector<double> m_wdot; //!< Species net molar production rates
     vector<double> m_uk; //!< Species molar internal energies
-    bool m_chem = false;
-    bool m_energy = true;
     size_t m_nv = 0;
     size_t m_nv_surf; //!!< Number of variables associated with reactor surfaces
 

include/cantera/zeroD/ReactorBase.h

@@ -266,26 +266,15 @@ class ReactorBase
     //! Set the ReactorNet that this reactor belongs to.
     void setNetwork(ReactorNet* net);
 
-    //! Calculate the derivative of T with respect to the ith species in the heat
-    //! transfer equation based on the reactor specific equation of state.
+    //! Calculate the derivative of T with respect to the ith species in the energy
+    //! conservation equation based on the reactor specific equation of state.
     //! @param index index of the species the derivative is with respect too
-    //! @warning This function is an experimental part of the %Cantera API and may be changed
-    //! or removed without notice.
-    //! @since New in %Cantera 3.0.
-    //!
-    virtual double moleDerivative(size_t index) {
-        throw NotImplementedError("Reactor::moleDerivative");
-    }
-
-    //! Calculate the derivative of T with respect to the ith species in the heat
-    //! transfer radiation equation based on the reactor specific equation of state.
-    //! @param index index of the species the derivative is with respect too
-    //! @warning This function is an experimental part of the %Cantera API and may be changed
-    //! or removed without notice.
-    //! @since New in %Cantera 3.0.
+    //! @warning This function is an experimental part of the %Cantera API and may
+    //! be changed or removed without notice.
+    //! @since New in %Cantera 3.1.
     //!
-    virtual double moleRadiationDerivative(size_t index) {
-        throw NotImplementedError("Reactor::moleRadiationDerivative");
+    virtual double temperature_ddni(size_t index) {
+        throw NotImplementedError("Reactor::temperature_ddni");
     }
 
     //! Return the index associated with energy of the system
@@ -294,6 +283,17 @@ class ReactorBase
     //! Return the offset between species and state variables
     virtual size_t speciesOffset() const { return m_sidx; };
 
+    //! Returns `true` if solution of the energy equation is enabled.
+    virtual bool energyEnabled() const {
+        return m_energy;
+    }
+
+    //! Returns `true` if changes in the reactor composition due to chemical reactions
+    //! are enabled.
+    bool chemistryEnabled() const {
+        return m_chem;
+    }
+
 protected:
     //! Specify the mixture contained in the reactor. Note that a pointer to
     //! this substance is stored, and as the integration proceeds, the state of
@@ -338,6 +338,12 @@ class ReactorBase
 
     //! Composite thermo/kinetics/transport handler
     shared_ptr<Solution> m_solution;
+
+    //! A bool that enables the energy equation
+    bool m_energy = true;
+
+    //! A bool that enables the chemical kinetics equations
+    bool m_chem = false;
 };
 }
 

include/cantera/zeroD/ReactorDelegator.h

@@ -74,7 +74,8 @@ class ReactorDelegator : public Delegator, public R, public ReactorAccessor
         install("updateState", m_updateState,
             [this](std::array<size_t, 1> sizes, double* y) { R::updateState(y); });
         install("updateSurfaceState", m_updateSurfaceState,
-            [this](std::array<size_t, 1> sizes, double* y) { R::updateSurfaceState(y); });
+            [this](std::array<size_t, 1> sizes, double* y)
+            { R::updateSurfaceState(y); });
         install("getSurfaceInitialConditions", m_getSurfaceInitialConditions,
             [this](std::array<size_t, 1> sizes, double* y) {
                 R::getSurfaceInitialConditions(y);
@@ -89,8 +90,8 @@ class ReactorDelegator : public Delegator, public R, public ReactorAccessor
         );
         install("evalWalls", m_evalWalls, [this](double t) { R::evalWalls(t); });
         install("evalSurfaces", m_evalSurfaces,
-            [this](std::array<size_t, 3> sizes, double* LHS, double* RHS, double* sdot) {
-                R::evalSurfaces(LHS, RHS, sdot);
+            [this](std::array<size_t, 3> sizes, double* LHS, double* RHS, double* sd) {
+                R::evalSurfaces(LHS, RHS, sd);
             }
         );
         install("componentName", m_componentName,

include/cantera/zeroD/ReactorNet.h

@@ -7,6 +7,8 @@
 #define CT_REACTORNET_H
 
 #include "Reactor.h"
+#include "Wall.h"
+#include "FlowDevice.h"
 #include "cantera/numerics/FuncEval.h"
 
 
@@ -239,16 +241,7 @@ class ReactorNet : public FuncEval
     //! Return the index corresponding to the start of the reactor specific state
     //! vector in the reactor with index *reactor* in the global state vector for the
     //! reactor network.
-    size_t globalStartIndex(Reactor* curr_reactor) {
-        for (size_t i = 0; i < m_reactors.size(); i++) {
-            if (curr_reactor == m_reactors[i]) {
-                return m_start[i];
-            }
-        }
-        throw CanteraError("ReactorNet::globalStartIndex: ",
-                curr_reactor->name(), " not found in network.");
-        return npos;
-    }
+    size_t globalStartIndex(ReactorBase* curr_reactor);
 
     //! Return the name of the i-th component of the global state vector. The
     //! name returned includes both the name of the reactor and the specific
@@ -417,6 +410,10 @@ class ReactorNet : public FuncEval
     //! derivative settings
     bool m_jac_skip_walls = false;
     bool m_jac_skip_flow_devices = false;
+    //! set to store walls for Jacobian calculation
+    set<WallBase*> m_walls;
+    //! set to store flow devices for Jacobian calculation
+    set<FlowDevice*> m_flow_devices;
 };
 }
 

include/cantera/zeroD/Wall.h

@@ -108,47 +108,27 @@ class WallBase
     //! Build the Jacobian terms specific to the flow device for the given connected
     //! reactor.
     //! @param r a pointer to the calling reactor
-    //! @param jacVector a vector of triplets to be added to the jacobian for the
-    //! reactor
+    //! @param jacVector a vector of triplets to be added to the reactor Jacobian
     //! @warning This function is an experimental part of the %Cantera API and may be
-    //! changed
-    //! or removed without notice.
-    //! @since New in %Cantera 3.0.
+    //! changed or removed without notice.
+    //! @since New in %Cantera 3.1.
     //!
-    virtual void buildReactorJacobian(ReactorBase* r, vector<Eigen::Triplet<double>>& jacVector) {
+    virtual void buildReactorJacobian(ReactorBase* r,
+        vector<Eigen::Triplet<double>>& jacVector) {
         throw NotImplementedError("WallBase::buildReactorJacobian");
     }
 
-    //! Build the Jacobian terms specific to the flow device for the network. These
-    //! terms
-    //! will be adjusted to the networks indexing system outside of the reactor.
-    //! @param jacVector a vector of triplets to be added to the jacobian for the
-    //! reactor
+    //! Build the Jacobian cross-reactor terms specific to the flow device for the
+    //! network.
+    //! @param jacVector a vector of triplets to be added to the network Jacobian
     //! @warning This function is an experimental part of the %Cantera API and may be
-    //! changed
-    //! or removed without notice.
-    //! @since New in %Cantera 3.0.
+    //! changed or removed without notice.
+    //! @since New in %Cantera 3.1.
     //!
     virtual void buildNetworkJacobian(vector<Eigen::Triplet<double>>& jacVector) {
         throw NotImplementedError("WallBase::buildNetworkJacobian");
     }
 
-    //! Specify the jacobian terms have been calculated and should not be recalculated.
-    //! @warning This function is an experimental part of the %Cantera API and may be
-    //! changed
-    //! or removed without notice.
-    //! @since New in %Cantera 3.0.
-    //!
-    void jacobianCalculated() { m_jac_calculated = true; };
-
-    //! Specify that jacobian terms have not been calculated and should be recalculated.
-    //! @warning This function is an experimental part of the %Cantera API and may be
-    //! changed
-    //! or removed without notice.
-    //! @since New in %Cantera 3.0.
-    //!
-    void jacobianNotCalculated() { m_jac_calculated = false; };
-
 protected:
     string m_name;  //!< Wall name.
     bool m_defaultNameSet = false;  //!< `true` if default name has been previously set.
@@ -160,10 +140,6 @@ class WallBase
     double m_time = 0.0;
 
     double m_area = 1.0;
-
-    //! a variable to switch on and off so calculations are not doubled by the calling
-    //! reactor or network
-    bool m_jac_calculated = false;
 };
 
 //! Represents a wall between between two ReactorBase objects.
@@ -268,9 +244,11 @@ class Wall : public WallBase
         return m_k;
     }
 
-    virtual void buildReactorJacobian(ReactorBase* r, vector<Eigen::Triplet<double>>& jacVector) override;
+    void buildReactorJacobian(ReactorBase* r,
+        vector<Eigen::Triplet<double>>& jacVector) override;
 
-    virtual void buildNetworkJacobian(vector<Eigen::Triplet<double>>& jacVector) override;
+    void buildNetworkJacobian(vector<Eigen::Triplet<double>>& jacVector)
+        override;
 
 protected:
 

interfaces/cython/cantera/delegator.pyx

@@ -319,6 +319,7 @@ cdef int assign_delegates(obj, CxxDelegator* delegator) except -1:
 
         if when is None:
             continue
+
         cxx_name = stringify(options[0])
         callback = options[1].replace(' ', '')