Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

Introduce BulkRate class template #1187

Merged
merged 10 commits into from
Feb 9, 2022
Merged

Introduce BulkRate class template #1187

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
10 commits
Select commit Hold shift + click to select a range
880c24e
[Kinetics] Introduce Arrhenius3
ischoegl Feb 6, 2022
bde21ed
Replace ArrheniusBase by Arrhenius3 as appropriate
ischoegl Feb 6, 2022
342981e
[Kinetics] Re-introduce BlowersMasel class
ischoegl Feb 6, 2022
3d521aa
[Kinetics] Introduce BulkRate<ReactionType,DataType> template
ischoegl Feb 6, 2022
7a68a11
[Kinetics] Define BlowersMaselRate as BulkRate
ischoegl Feb 6, 2022
3441bc2
[Kinetics] Define TwoTempPlasmaRate as BulkRate
ischoegl Feb 6, 2022
62979a0
Eliminate C++ TwoTempPlasmaReaction class specialization
ischoegl Feb 7, 2022
d610a2c
[UnitTests] Update tests and introduce pytest features
ischoegl Feb 7, 2022
49ff30b
[Kinetics] Eliminate multiple inheritance
ischoegl Feb 7, 2022
e0095a4
[Kinetics] Eliminate chained function calls
ischoegl Feb 8, 2022
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
227 changes: 97 additions & 130 deletions include/cantera/kinetics/Arrhenius.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -30,13 +30,8 @@ class AnyMap;
* This base class provides a minimally functional interface that allows for parameter
* access from derived classes as well as classes that use Arrhenius-type expressions
* internally, for example FalloffRate and PlogRate.
*
* @todo supersedes Arrhenius2 and will replace Arrhenius(2) after Cantera 2.6,
* The class should be renamed to Arrhenius after removal of Arrhenius2. The new
* behavior can be forced in self-compiled Cantera installations by defining
* CT_NO_LEGACY_REACTIONS_26 via the 'no_legacy_reactions' option in SCons.
*/
class ArrheniusBase
class ArrheniusBase : public ReactionRate
{
public:
//! Default constructor.
Expand Down Expand Up @@ -72,23 +67,7 @@ class ArrheniusBase
void getRateParameters(AnyMap& node) const;

//! Check rate expression
void checkRate(const std::string& equation, const AnyMap& node);

//! Evaluate reaction rate
/*!
* @internal Non-virtual method that should not be overloaded
*/
double evalRate(double logT, double recipT) const {
return m_A * std::exp(m_b * logT - m_Ea_R * recipT);
}

//! Evaluate natural logarithm of the rate constant.
/*!
* @internal Non-virtual method that should not be overloaded
*/
double evalLog(double logT, double recipT) const {
return m_logA + m_b * logT - m_Ea_R * recipT;
}
virtual void check(const std::string& equation, const AnyMap& node) override;

//! Return the pre-exponential factor *A* (in m, kmol, s to powers depending
//! on the reaction order)
Expand All @@ -109,6 +88,12 @@ class ArrheniusBase
return m_Ea_R * GasConstant;
}

//! Return the activation energy divided by the gas constant (i.e. the
//! activation temperature) [K]
double activationEnergy_R() const {
return m_Ea_R;
}

// Return units of the reaction rate expression
const Units& rateUnits() const {
return m_rate_units;
Expand Down Expand Up @@ -165,45 +150,36 @@ class ArrheniusBase
* \f]
*
* @ingroup arrheniusGroup
*
* @todo supersedes Arrhenius2 and will replace Arrhenius(2) after Cantera 2.6,
* The class should be renamed to Arrhenius after removal of Arrhenius2. The new
* behavior can be forced in self-compiled Cantera installations by defining
* CT_NO_LEGACY_REACTIONS_26 via the 'no_legacy_reactions' option in SCons.
*/
class ArrheniusRate final : public ArrheniusBase, public ReactionRate
class Arrhenius3 : public ArrheniusBase
{
public:
ArrheniusRate() = default;
using ArrheniusBase::ArrheniusBase; // inherit constructors

//! Constructor based on AnyMap content
ArrheniusRate(const AnyMap& node, const UnitStack& rate_units={}) {
setParameters(node, rate_units);
}

unique_ptr<MultiRateBase> newMultiRate() const override {
return unique_ptr<MultiRateBase>(new MultiRate<ArrheniusRate, ArrheniusData>);
}

//! Identifier of reaction rate type
virtual const std::string type() const override {
return "Arrhenius";
}

//! Perform object setup based on AnyMap node information
/*!
* @param node AnyMap containing rate information
* @param rate_units Unit definitions specific to rate information
*/
virtual void setParameters(const AnyMap& node, const UnitStack& rate_units) override;

virtual void getParameters(AnyMap& node) const override;
//! Evaluate reaction rate
double evalRate(double logT, double recipT) const {
return m_A * std::exp(m_b * logT - m_Ea_R * recipT);
}

void check(const std::string& equation, const AnyMap& node) override {
checkRate(equation, node);
//! Evaluate natural logarithm of the rate constant.
double evalLog(double logT, double recipT) const {
return m_logA + m_b * logT - m_Ea_R * recipT;
}

//! Evaluate reaction rate
/*!
* @param shared_data data shared by all reactions of a given type
*/
double evalFromStruct(const ReactionData& shared_data) const {
double evalFromStruct(const ArrheniusData& shared_data) const {
return m_A * std::exp(m_b * shared_data.logT - m_Ea_R * shared_data.recipT);
}

Expand All @@ -212,15 +188,9 @@ class ArrheniusRate final : public ArrheniusBase, public ReactionRate
/*!
* @param shared_data data shared by all reactions of a given type
*/
virtual double ddTScaledFromStruct(const ReactionData& shared_data) const {
double ddTScaledFromStruct(const ArrheniusData& shared_data) const {
return (m_Ea_R * shared_data.recipT + m_b) * shared_data.recipT;
}

//! Return the activation energy divided by the gas constant (i.e. the
//! activation temperature) [K]
double activationEnergy_R() const {
return m_Ea_R;
}
};


Expand All @@ -240,12 +210,13 @@ class ArrheniusRate final : public ArrheniusBase, public ReactionRate
* Kinetic scheme of the non-equilibrium discharge in nitrogen-oxygen mixtures.
* Plasma Sources Science and Technology, 1(3), 207.
* doi: 10.1088/0963-0252/1/3/011
*
* @ingroup arrheniusGroup
*/
class TwoTempPlasmaRate final : public ArrheniusBase, public ReactionRate
class TwoTempPlasma : public ArrheniusBase
{
public:
TwoTempPlasmaRate();
TwoTempPlasma();

//! Constructor.
/*!
Expand All @@ -255,55 +226,34 @@ class TwoTempPlasmaRate final : public ArrheniusBase, public ReactionRate
* @param Ea Activation energy in energy units [J/kmol]
* @param EE Activation electron energy in energy units [J/kmol]
*/
TwoTempPlasmaRate(double A, double b, double Ea=0.0, double EE=0.0);

unique_ptr<MultiRateBase> newMultiRate() const override {
return unique_ptr<MultiRateBase>(
new MultiRate<TwoTempPlasmaRate, TwoTempPlasmaData>);
}

//! Constructor based on AnyMap content
TwoTempPlasmaRate(const AnyMap& node, const UnitStack& rate_units={})
: TwoTempPlasmaRate()
{
setParameters(node, rate_units);
}
TwoTempPlasma(double A, double b, double Ea=0.0, double EE=0.0);

//! Identifier of reaction rate type
virtual const std::string type() const override {
return "two-temperature-plasma";
}

//! Perform object setup based on AnyMap node information
/*!
* @param node AnyMap containing rate information
* @param rate_units Unit definitions specific to rate information
*/
virtual void setParameters(const AnyMap& node, const UnitStack& rate_units) override;

virtual void getParameters(AnyMap& node) const override;

void check(const std::string& equation, const AnyMap& node) override {
checkRate(equation, node);
}
virtual void setContext(const Reaction& rxn, const Kinetics& kin) override;
ischoegl marked this conversation as resolved.
Show resolved Hide resolved

//! Evaluate reaction rate
/*!
* @param shared_data data shared by all reactions of a given type
*/
double evalFromStruct(const TwoTempPlasmaData& shared_data) {
double evalFromStruct(const TwoTempPlasmaData& shared_data) const {
// m_E4_R is the electron activation (in temperature units)
return m_A * std::exp(m_b * shared_data.logTe -
m_Ea_R * shared_data.recipT +
m_E4_R * (shared_data.electronTemp - shared_data.temperature)
* shared_data.recipTe * shared_data.recipT);
}

//! Return the activation energy divided by the gas constant (i.e. the
//! activation temperature) [K]
double activationEnergy_R() const {
return m_Ea_R;
}
//! Evaluate derivative of reaction rate with respect to temperature
//! divided by reaction rate
/*!
* This method does not consider changes of electron temperature.
* A corresponding warning is raised.
* @param shared_data data shared by all reactions of a given type
*/
double ddTScaledFromStruct(const TwoTempPlasmaData& shared_data) const;

//! Return the electron activation energy *Ea* [J/kmol]
double activationElectronEnergy() const {
Expand Down Expand Up @@ -347,11 +297,11 @@ class TwoTempPlasmaRate final : public ArrheniusBase, public ReactionRate
*
* @ingroup arrheniusGroup
*/
class BlowersMaselRate final : public ArrheniusBase, public ReactionRate
class BlowersMasel : public ArrheniusBase
{
public:
//! Default constructor.
BlowersMaselRate();
BlowersMasel();

//! Constructor.
/*!
Expand All @@ -362,45 +312,15 @@ class BlowersMaselRate final : public ArrheniusBase, public ReactionRate
* @param w Average bond dissociation energy of the bond being formed and
* broken in the reaction, in energy units [J/kmol]
*/
BlowersMaselRate(double A, double b, double Ea0, double w);

unique_ptr<MultiRateBase> newMultiRate() const override {
return unique_ptr<MultiRateBase>(
new MultiRate<BlowersMaselRate, BlowersMaselData>);
}

//! Constructor based on AnyMap content
BlowersMaselRate(const AnyMap& node, const UnitStack& rate_units={})
: BlowersMaselRate()
{
setParameters(node, rate_units);
}
BlowersMasel(double A, double b, double Ea0, double w);

//! Identifier of reaction rate type
virtual const std::string type() const override {
return "Blowers-Masel";
}

//! Perform object setup based on AnyMap node information
/*!
* @param node AnyMap containing rate information
* @param rate_units Unit definitions specific to rate information
*/
virtual void setParameters(
const AnyMap& node, const UnitStack& rate_units) override;

virtual void getParameters(AnyMap& node) const override;

void check(const std::string& equation, const AnyMap& node) override {
checkRate(equation, node);
}

virtual void setContext(const Reaction& rxn, const Kinetics& kin) override;

//! Update information specific to reaction
/*!
* @param shared_data data shared by all reactions of a given type
*/
void updateFromStruct(const BlowersMaselData& shared_data) {
if (shared_data.ready) {
m_deltaH_R = 0.;
Expand All @@ -417,24 +337,22 @@ class BlowersMaselRate final : public ArrheniusBase, public ReactionRate
/*!
* @param shared_data data shared by all reactions of a given type
*/
double evalFromStruct(const BlowersMaselData& shared_data) {
double Ea_R = activationEnergy_R(m_deltaH_R);
double evalFromStruct(const BlowersMaselData& shared_data) const {
double Ea_R = effectiveActivationEnergy_R(m_deltaH_R);
return m_A * std::exp(m_b * shared_data.logT - Ea_R * shared_data.recipT);
}

//! Evaluate derivative of reaction rate with respect to temperature
//! divided by reaction rate
/*!
* This method is used to override the numerical derivative, which does not
* consider potential changes due to a changed reaction enthalpy. A corresponding
* warning is raised.
* This method does not consider potential changes due to a changed reaction
* enthalpy. A corresponding warning is raised.
* @param shared_data data shared by all reactions of a given type
*/
virtual double ddTScaledFromStruct(const BlowersMaselData& shared_data) const;
double ddTScaledFromStruct(const BlowersMaselData& shared_data) const;

//! Return the effective activation energy (a function of the delta H of reaction)
//! divided by the gas constant (i.e. the activation temperature) [K]
double activationEnergy_R(double deltaH_R) const {
double effectiveActivationEnergy_R(double deltaH_R) const {
if (deltaH_R < -4 * m_Ea_R) {
return 0.;
}
Expand All @@ -453,7 +371,7 @@ class BlowersMaselRate final : public ArrheniusBase, public ReactionRate
* @param deltaH Enthalpy change of reaction [J/kmol]
*/
double effectiveActivationEnergy(double deltaH) const {
return activationEnergy_R(deltaH / GasConstant) * GasConstant;
return effectiveActivationEnergy_R(deltaH / GasConstant) * GasConstant;
}

//! Return the bond dissociation energy *w* [J/kmol]
Expand All @@ -468,7 +386,56 @@ class BlowersMaselRate final : public ArrheniusBase, public ReactionRate
double m_deltaH_R; //!< enthalpy change of reaction (in temperature units)
};

}

//! A class template for bulk phase reaction rate specifications
template <class RateType, class DataType>
class BulkRate : public RateType
{
public:
BulkRate() = default;
using RateType::RateType; // inherit constructors

//! Constructor based on AnyMap content
BulkRate(const AnyMap& node, const UnitStack& rate_units={}) {
setParameters(node, rate_units);
}

unique_ptr<MultiRateBase> newMultiRate() const override {
return unique_ptr<MultiRateBase>(
new MultiRate<BulkRate<RateType, DataType>, DataType>);
}

virtual void setParameters(
const AnyMap& node, const UnitStack& rate_units) override
{
RateType::m_negativeA_ok = node.getBool("negative-A", false);
if (!node.hasKey("rate-constant")) {
RateType::setRateParameters(AnyValue(), node.units(), rate_units);
return;
}
RateType::setRateParameters(node["rate-constant"], node.units(), rate_units);
}

virtual void getParameters(AnyMap& node) const override {
if (RateType::m_negativeA_ok) {
node["negative-A"] = true;
}
AnyMap rateNode;
RateType::getRateParameters(rateNode);
if (!rateNode.empty()) {
// RateType object is configured
node["rate-constant"] = std::move(rateNode);
}
if (RateType::type() != "Arrhenius") {
node["type"] = RateType::type();
}
}
};

typedef BulkRate<Arrhenius3, ArrheniusData> ArrheniusRate;
typedef BulkRate<TwoTempPlasma, TwoTempPlasmaData> TwoTempPlasmaRate;
typedef BulkRate<BlowersMasel, BlowersMaselData> BlowersMaselRate;

}

#endif
Loading