Extensible rate: unit handling, serialization, and validation

doc/sphinx/cython/utilities.rst

@@ -6,8 +6,13 @@ Utilities
 .. contents::
    :local:
 
-YAML Output
------------
+YAML Input/Output
+-----------------
+
+AnyMap
+^^^^^^
+
+.. autoclass:: AnyMap
 
 YamlWriter
 ^^^^^^^^^^

doc/sphinx/yaml/general.rst

@@ -40,6 +40,8 @@ following structure::
     - equation: A + C <=> 2 D
       # Additional fields come after this
 
+.. _sec-yaml-units:
+
 Units
 -----
 

include/cantera/base/AnyMap.h

@@ -621,6 +621,9 @@ class AnyMap : public AnyBase
     //! Return the default units that should be used to convert stored values
     const UnitSystem& units() const { return *m_units; }
 
+    //! @copydoc units()
+    shared_ptr<UnitSystem> unitsShared() const { return m_units; }
+
     //! Use the supplied UnitSystem to set the default units, and recursively
     //! process overrides from nodes named `units`.
     /*!

include/cantera/base/Delegator.h

@@ -7,6 +7,7 @@
 #define CT_DELEGATOR_H
 
 #include "cantera/base/global.h"
+#include "cantera/base/Units.h"
 #include "cantera/base/ctexceptions.h"
 #include "cantera/base/ExtensionManager.h"
 #include <array>
@@ -145,6 +146,17 @@ class Delegator
         *m_funcs_v_d[name] = makeDelegate(func, when, *m_funcs_v_d[name]);
     }
 
+    //! set delegates for member functions with the signature `void(AnyMap&)`
+    void setDelegate(const string& name, const function<void(AnyMap&)>& func,
+                     const string& when)
+    {
+        if (!m_funcs_v_AMr.count(name)) {
+            throw NotImplementedError("Delegator::setDelegate",
+                "for function named '{}' with signature 'void(AnyMap&)'.", name);
+        }
+        *m_funcs_v_AMr[name] = makeDelegate(func, when, *m_funcs_v_AMr[name]);
+    }
+
     //! set delegates for member functions with the signature
     //! `void(AnyMap&, UnitStack&)`
     void setDelegate(const std::string& name,
@@ -160,6 +172,19 @@ class Delegator
         *m_funcs_v_cAMr_cUSr[name] = makeDelegate(func, when, *m_funcs_v_cAMr_cUSr[name]);
     }
 
+    //! set delegates for member functions with the signature
+    //! `void(const string&, void*)`
+    void setDelegate(const string& name,
+                     const function<void(const string&, void*)>& func,
+                     const string& when)
+    {
+        if (!m_funcs_v_csr_vp.count(name)) {
+            throw NotImplementedError("Delegator::setDelegate",
+                "for function named '{}' with signature 'void(const string&, void*)'.");
+        }
+        *m_funcs_v_csr_vp[name] = makeDelegate(func, when, *m_funcs_v_csr_vp[name]);
+    }
+
     //! Set delegates for member functions with the signature `void(double*)`
     void setDelegate(const std::string& name,
                      const std::function<void(std::array<size_t, 1>, double*)>& func,
@@ -254,8 +279,21 @@ class Delegator
         *m_funcs_sz_csr[name] = makeDelegate(name, func, when, m_base_sz_csr[name]);
     }
 
-    void holdExternalHandle(const shared_ptr<ExternalHandle>& handle) {
-        m_handles.push_back(handle);
+    //! Store a handle to a wrapper for the delegate from an external language interface
+    void holdExternalHandle(const string& name,
+                            const shared_ptr<ExternalHandle>& handle) {
+        m_handles[name] = handle;
+    }
+
+    //! Get the handle for a wrapper for the delegate from the external language
+    //! interface specified by *name*.
+    //! Returns a null pointer if the requested handle does not exist.
+    shared_ptr<ExternalHandle> getExternalHandle(const string& name) const {
+        if (m_handles.count(name)) {
+            return m_handles.at(name);
+        } else {
+            return shared_ptr<ExternalHandle>();
+        }
     }
 
 protected:
@@ -283,6 +321,14 @@ class Delegator
         m_funcs_v_d[name] = &target;
     }
 
+    //! Install a function with the signature `void(AnyMap&)` as being delegatable
+    void install(const string& name, function<void(AnyMap&)>& target,
+                 const function<void(AnyMap&)>& func)
+    {
+        target = func;
+        m_funcs_v_AMr[name] = &target;
+    }
+
     //! Install a function with the signature `void(const AnyMap&, const UnitStack&)`
     //! as being delegatable
     void install(const std::string& name,
@@ -293,7 +339,14 @@ class Delegator
         m_funcs_v_cAMr_cUSr[name] = &target;
     }
 
-
+    //! Install a function with the signature `void(const string&, void*) as being
+    //! delegatable
+    void install(const string& name, function<void(const string&, void*)>& target,
+                 const function<void(const string&, void*)>& func)
+    {
+        target = func;
+        m_funcs_v_csr_vp[name] = &target;
+    }
 
     //! Install a function with the signature `void(double*)` as being delegatable
     void install(const std::string& name,
@@ -470,8 +523,10 @@ class Delegator
     std::map<std::string, std::function<void()>*> m_funcs_v;
     std::map<std::string, std::function<void(bool)>*> m_funcs_v_b;
     std::map<std::string, std::function<void(double)>*> m_funcs_v_d;
+    map<string, function<void(AnyMap&)>*> m_funcs_v_AMr;
     std::map<std::string,
         std::function<void(const AnyMap&, const UnitStack&)>*> m_funcs_v_cAMr_cUSr;
+    map<string, function<void(const string&, void*)>*> m_funcs_v_csr_vp;
     std::map<std::string,
         std::function<void(std::array<size_t, 1>, double*)>*> m_funcs_v_dp;
     std::map<std::string,
@@ -496,8 +551,10 @@ class Delegator
         std::function<size_t(const std::string&)>*> m_funcs_sz_csr;
     //! @}
 
-    //! Cleanup functions to be called from the destructor
-    std::list<shared_ptr<ExternalHandle>> m_handles;
+    //! Handles to wrappers for the delegated object in external language interfaces.
+    //! Used to provide access to these wrappers as well as managing cleanup functions
+    //! called from the destructor of the derived ExternalHandle class.
+    map<string, shared_ptr<ExternalHandle>> m_handles;
 
     //! Name of the class in the extension language
     std::string m_delegatorName;

include/cantera/base/ExtensionManager.h

@@ -88,14 +88,22 @@ class ExtensionManager
 
     //! Register a function that can be used to create wrappers for ReactionData objects
     //! in an external language and link them to the corresponding C++ object
-    static void registerReactionDataLinker(const std::string& rateName,
-        std::function<void(ReactionDataDelegator&)> link);
+    //!
+    //! @param rateName  The name of the reaction rate type
+    //! @param wrapperName  The name used for Solution wrappers to be used with this
+    //!     object, corresponding to a type registered with registerSolutionLinker().
+    static void registerReactionDataLinker(const string& rateName,
+        const string& wrapperName, function<void(ReactionDataDelegator&)> link);
 
     //! Register a function that can be used to create wrappers for Solution objects in
     //! an external language and link it to the corresponding C++ objects
     static void registerSolutionLinker(const std::string& wrapperName,
         std::function<shared_ptr<ExternalHandle>(shared_ptr<Solution>)> link);
 
+    //! Get the Solution wrapper type corresponding to the specified user-defined
+    //! reaction rate type.
+    static string getSolutionWrapperType(const string& userType);
+
 protected:
     //! Functions for wrapping and linking ReactionData objects
     static std::map<std::string,
@@ -105,6 +113,8 @@ class ExtensionManager
     static std::map<std::string,
         std::function<shared_ptr<ExternalHandle>(shared_ptr<Solution>)>> s_Solution_linkers;
 
+    //! Mapping from user-defined rate types to Solution wrapper types
+    static map<string, string> s_userTypeToWrapperType;
 };
 
 }

include/cantera/base/Solution.h

@@ -137,7 +137,7 @@ class Solution : public std::enable_shared_from_this<Solution>
 
     AnyMap m_header;  //!< Additional input fields; usually from a YAML input file
 
-    //! Wrappers for this Kinetics object in extension languages, for evaluation
+    //! Wrappers for this Solution object in extension languages, for evaluation
     //! of user-defined reaction rates
     std::map<std::string, shared_ptr<ExternalHandle>> m_externalHandles;
 

include/cantera/base/Units.h

@@ -106,6 +106,8 @@ struct UnitStack
     UnitStack(std::initializer_list<std::pair<Units, double>> units)
         : stack(units) {}
 
+    UnitStack() = default;
+
     //! Size of UnitStack
     size_t size() const { return stack.size(); }
 
@@ -220,6 +222,13 @@ class UnitSystem
     double convert(const AnyValue& val, const std::string& dest) const;
     double convert(const AnyValue& val, const Units& dest) const;
 
+    //! Convert a generic AnyValue node representing a reaction rate coefficient to the
+    //! units specified in `dest`. Works like `convert(AnyValue&, Units&)` but with
+    //! special handling for the case where the destination units are undefined.
+    //!
+    //! @since New in Cantera 3.0
+    double convertRateCoeff(const AnyValue& val, const Units& dest) const;
+
     //! Convert an array of AnyValue nodes to the units specified in `dest`. For
     //! each node, if the value is a double, convert it using the default units,
     //! and if it is a string, treat it as a value with the given dimensions.

include/cantera/kinetics/Arrhenius.h

@@ -121,24 +121,18 @@ class ArrheniusBase : public ReactionRate
         return m_Ea_R * GasConstant;
     }
 
-    // Return units of the reaction rate expression
-    const Units& rateUnits() const {
-        return m_rate_units;
-    }
-
     //! Return reaction order associated with the reaction rate
     double order() const {
         return m_order;
     }
 
     //! Set units of the reaction rate expression
-    void setRateUnits(const UnitStack& rate_units) {
+    void setRateUnits(const UnitStack& rate_units) override {
+        ReactionRate::setRateUnits(rate_units);
         if (rate_units.size() > 1) {
-            m_rate_units = rate_units.product();
-            m_order = 1 - m_rate_units.dimension("quantity");
+            m_order = 1 - rate_units.product().dimension("quantity");
         } else {
             m_order = NAN;
-            m_rate_units = rate_units.standardUnits();
         }
     }
 
@@ -164,7 +158,6 @@ class ArrheniusBase : public ReactionRate
     std::string m_b_str = "b"; //!< The string for temperature exponent
     std::string m_Ea_str = "Ea"; //!< The string for activation energy
     std::string m_E4_str = ""; //!< The string for an optional 4th parameter
-    Units m_rate_units{0.}; //!< Reaction rate units
 };
 
 //! Arrhenius reaction rate type depends only on temperature

include/cantera/kinetics/ChebyshevRate.h

@@ -232,8 +232,6 @@ class ChebyshevRate final : public ReactionRate
 
     Array2D m_coeffs; //!<< coefficient array
     vector_fp dotProd_; //!< dot product of coeffs with the reduced pressure polynomial
-
-    Units m_rate_units = Units(0.); //!< Reaction rate units
 };
 
 }

include/cantera/kinetics/InterfaceRate.h

@@ -456,12 +456,10 @@ class StickingRate : public RateType, public StickingCoverage
 
     //! Constructor based on AnyMap content
     StickingRate(const AnyMap& node, const UnitStack& rate_units) {
-        // sticking coefficients are dimensionless
-        setParameters(node, Units(1.0));
+        setParameters(node, rate_units);
     }
     explicit StickingRate(const AnyMap& node) {
-        // sticking coefficients are dimensionless
-        setParameters(node, Units(1.0));
+        setParameters(node, {});
     }
 
     unique_ptr<MultiRateBase> newMultiRate() const override {
@@ -473,10 +471,16 @@ class StickingRate : public RateType, public StickingCoverage
         return "sticking-" + RateType::type();
     }
 
+    void setRateUnits(const UnitStack& rate_units) override {
+        // Sticking coefficients are dimensionless
+        RateType::m_conversion_units = Units(1.0);
+    }
+
     virtual void setParameters(
         const AnyMap& node, const UnitStack& rate_units) override
     {
         InterfaceRateBase::setParameters(node);
+        setRateUnits(rate_units);
         RateType::m_negativeA_ok = node.getBool("negative-A", false);
         setStickingParameters(node);
         if (!node.hasKey("sticking-coefficient")) {

include/cantera/kinetics/KineticsFactory.h

@@ -58,10 +58,12 @@ shared_ptr<Kinetics> newKinetics(const string& model);
  *     to the Kinetics object.
  * @param rootNode   The root node of the file containing the phase definition,
  *     which will be treated as the default source for reactions
+ * @param soln       The Solution object that this Kinetics object is being added to.
  */
 shared_ptr<Kinetics> newKinetics(const vector<shared_ptr<ThermoPhase>>& phases,
                                  const AnyMap& phaseNode,
-                                 const AnyMap& rootNode=AnyMap());
+                                 const AnyMap& rootNode=AnyMap(),
+                                 shared_ptr<Solution> soln={});
 
 //! @copydoc newKinetics(const vector<shared_ptr<ThermoPhase>>&, const AnyMap&, const AnyMap&)
 //! @deprecated  To be removed after Cantera 3.0;

include/cantera/kinetics/ReactionRate.h

@@ -53,6 +53,7 @@ class ReactionRate
         : m_input(other.m_input)
         , m_rate_index(other.m_rate_index)
         , m_valid(other.m_valid)
+        , m_conversion_units(other.m_conversion_units)
     {}
 
     ReactionRate& operator=(const ReactionRate& other) {
@@ -62,6 +63,7 @@ class ReactionRate
         m_input = other.m_input;
         m_rate_index = other.m_rate_index;
         m_valid = other.m_valid;
+        m_conversion_units = other.m_conversion_units;
         return *this;
     }
 
@@ -94,6 +96,7 @@ class ReactionRate
     //! @param node  AnyMap object containing reaction rate specification
     //! @param units  unit definitions specific to rate information
     virtual void setParameters(const AnyMap& node, const UnitStack& units) {
+        setRateUnits(units);
         m_input = node;
     }
 
@@ -102,10 +105,35 @@ class ReactionRate
     //! handled by the getParameters() method.
     AnyMap parameters() const {
         AnyMap out;
+        out["type"] = type();
         getParameters(out);
         return out;
     }
 
+    //! Get the units for converting the leading term in the reaction rate expression.
+    //!
+    //! These units are often the same as the units of the rate expression itself, but
+    //! not always; sticking coefficients are a notable exception.
+    //! @since New in Cantera 3.0
+    const Units& conversionUnits() const {
+        return m_conversion_units;
+    }
+
+    //! Set the units of the reaction rate expression
+    //!
+    //! Used to determine the units that should be used for converting terms in the
+    //! reaction rate expression, which often have the same units (for example, the
+    //! Arrhenius pre-exponential) but may also be different (for example, sticking
+    //! coefficients).
+    //! @since New in Cantera 3.0
+    virtual void setRateUnits(const UnitStack& rate_units) {
+        if (rate_units.size() > 1) {
+            m_conversion_units = rate_units.product();
+        } else {
+            m_conversion_units = rate_units.standardUnits();
+        }
+    }
+
     //! Check basic syntax and settings of reaction rate expression
     virtual void check(const std::string& equation) {}
 
@@ -224,6 +252,9 @@ class ReactionRate
     //! Flag indicating composition dependent rate
     bool m_composition_dependent_rate = false;
 
+    //! Units of the leading term in the reaction rate expression
+    Units m_conversion_units{0.};
+
 private:
     //! Return an object that be used to evaluate the rate by converting general input
     //! such as temperature and pressure into the `DataType` struct that is particular