Cantera  3.1.0a1
Loading...
Searching...
No Matches
BlowersMaselRate.h
Go to the documentation of this file.
1//! @file BlowersMaselRate.h Header for Blowers-Masel reaction rates
2
3// This file is part of Cantera. See License.txt in the top-level directory or
4// at https://cantera.org/license.txt for license and copyright information.
5
6#ifndef CT_BLOWERSMASELRATE_H
7#define CT_BLOWERSMASELRATE_H
8
9#include "Arrhenius.h"
10
11namespace Cantera
12{
13
14//! Data container holding shared data specific to BlowersMaselRate
15/**
16 * The data container `BlowersMaselData` holds precalculated data common to
17 * all `BlowersMaselRate` objects.
18 */
20{
21 BlowersMaselData() = default;
22
23 void update(double T) override;
24 bool update(const ThermoPhase& phase, const Kinetics& kin) override;
26
27 void resize(size_t nSpecies, size_t nReactions, size_t nPhases) override {
28 partialMolarEnthalpies.resize(nSpecies, 0.);
29 ready = true;
30 }
31
32 bool ready = false; //!< boolean indicating whether vectors are accessible
33 double density = NAN; //!< used to determine if updates are needed
34 vector<double> partialMolarEnthalpies; //!< partial molar enthalpies
35
36protected:
37 int m_state_mf_number = -1; //!< integer that is incremented when composition changes
38};
39
40
41//! Blowers Masel reaction rate type depends on the enthalpy of reaction
42/**
43 * The Blowers Masel approximation @cite blowers2000 adjusts the activation energy
44 * based on enthalpy change of a reaction:
45 *
46 * @f{eqnarray*}{
47 * E_a &=& 0\; &\text{if }\Delta H < -4E_0 \\
48 * E_a &=& \Delta H\; &\text{if }\Delta H > 4E_0 \\
49 * E_a &=& \frac{(w + \Delta H / 2)(V_P - 2w +
50 * \Delta H)^2}{(V_P^2 - 4w^2 + (\Delta H)^2)}\; &\text{otherwise}
51 * @f}
52 * where
53 * @f[
54 * V_P = \frac{2w (w + E_0)}{w - E_0},
55 * @f]
56 * @f$ w @f$ is the average bond dissociation energy of the bond breaking
57 * and that being formed in the reaction. Since the expression is
58 * very insensitive to @f$ w @f$ for @f$ w >= 2 E_0 @f$, @f$ w @f$
59 * can be approximated to an arbitrary high value like 1000 kJ/mol.
60 *
61 * After the activation energy is determined by Blowers-Masel approximation,
62 * it can be plugged into Arrhenius function to calculate the rate constant.
63 * @f[
64 * k_f = A T^b \exp (-E_a/RT)
65 * @f]
66 *
67 * @ingroup arrheniusGroup
68 */
70{
71public:
72 //! Default constructor.
74
75 //! Constructor.
76 /*!
77 * @param A Pre-exponential factor. The unit system is (kmol, m, s); actual units
78 * depend on the reaction order and the dimensionality (surface or bulk).
79 * @param b Temperature exponent (non-dimensional)
80 * @param Ea0 Intrinsic activation energy in energy units [J/kmol]
81 * @param w Average bond dissociation energy of the bond being formed and
82 * broken in the reaction, in energy units [J/kmol]
83 */
84 BlowersMaselRate(double A, double b, double Ea0, double w);
85
86 explicit BlowersMaselRate(const AnyMap& node,
87 const UnitStack& rate_units={});
88
89 unique_ptr<MultiRateBase> newMultiRate() const override {
90 return make_unique<MultiRate<BlowersMaselRate, BlowersMaselData>>();
91 }
92
93 const string type() const override {
94 return "Blowers-Masel";
95 }
96
97 void setContext(const Reaction& rxn, const Kinetics& kin) override;
98
99 //! Evaluate reaction rate
100 double evalRate(double logT, double recipT) const {
102 return m_A * std::exp(m_b * logT - Ea_R * recipT);
103 }
104
105 //! Update information specific to reaction
106 void updateFromStruct(const BlowersMaselData& shared_data) {
107 if (shared_data.ready) {
108 m_deltaH_R = 0.;
109 for (const auto& [k, multiplier] : m_stoich_coeffs) {
110 m_deltaH_R += shared_data.partialMolarEnthalpies[k] * multiplier;
111 }
113 }
114 }
115
116 //! Evaluate reaction rate
117 /*!
118 * @param shared_data data shared by all reactions of a given type
119 */
120 double evalFromStruct(const BlowersMaselData& shared_data) const {
122 return m_A * std::exp(m_b * shared_data.logT - Ea_R * shared_data.recipT);
123 }
124
125 //! Evaluate derivative of reaction rate with respect to temperature
126 //! divided by reaction rate
127 /*!
128 * This method does not consider potential changes due to a changed reaction
129 * enthalpy. A corresponding warning is raised.
130 * @param shared_data data shared by all reactions of a given type
131 */
132 double ddTScaledFromStruct(const BlowersMaselData& shared_data) const;
133
134protected:
135 //! Return the effective activation energy (a function of the delta H of reaction)
136 //! divided by the gas constant (that is, the activation temperature) [K]
137 double effectiveActivationEnergy_R(double deltaH_R) const {
138 if (deltaH_R < -4 * m_Ea_R) {
139 return 0.;
140 }
141 if (deltaH_R > 4 * m_Ea_R) {
142 return deltaH_R;
143 }
144 // m_E4_R is the bond dissociation energy "w" (in temperature units)
145 double vp = 2 * m_E4_R * ((m_E4_R + m_Ea_R) / (m_E4_R - m_Ea_R)); // in Kelvin
146 double vp_2w_dH = (vp - 2 * m_E4_R + deltaH_R); // (Vp - 2 w + dH)
147 return (m_E4_R + deltaH_R / 2) * (vp_2w_dH * vp_2w_dH) /
148 (vp * vp - 4 * m_E4_R * m_E4_R + deltaH_R * deltaH_R); // in Kelvin
149 }
150
151public:
152 double activationEnergy() const override {
154 }
155
156 //! Return the bond dissociation energy *w* [J/kmol]
157 double bondEnergy() const {
158 return m_E4_R * GasConstant;
159 }
160
161 //! Return current enthalpy change of reaction [J/kmol]
162 double deltaH() const {
163 return m_deltaH_R * GasConstant;
164 }
165
166 //! Set current enthalpy change of reaction [J/kmol]
167 /*!
168 * @note used for testing purposes only; this quantity is not an
169 * independent variable and will be overwritten during an update of the state.
170 *
171 * @warning This method is an experimental part of the %Cantera API and
172 * may be changed or removed without notice.
173 */
174 void setDeltaH(double deltaH) {
176 }
177
178protected:
179 //! Pairs of species indices and multipliers to calculate enthalpy change
180 vector<pair<size_t, double>> m_stoich_coeffs;
181
182 double m_deltaH_R = 0.0; //!< enthalpy change of reaction (in temperature units)
183};
184
185}
186
187#endif
Header for reaction rates that involve Arrhenius-type kinetics.
A map of string keys to values whose type can vary at runtime.
Definition AnyMap.h:427
Base class for Arrhenius-type Parameterizations.
Definition Arrhenius.h:44
double m_E4_R
Optional 4th energy parameter (in temperature units)
Definition Arrhenius.h:150
double m_A
Pre-exponential factor.
Definition Arrhenius.h:147
double m_b
Temperature exponent.
Definition Arrhenius.h:148
double m_Ea_R
Activation energy (in temperature units)
Definition Arrhenius.h:149
Blowers Masel reaction rate type depends on the enthalpy of reaction.
double evalFromStruct(const BlowersMaselData &shared_data) const
Evaluate reaction rate.
void setContext(const Reaction &rxn, const Kinetics &kin) override
Set context of reaction rate evaluation.
unique_ptr< MultiRateBase > newMultiRate() const override
Create a rate evaluator for reactions of a particular derived type.
double evalRate(double logT, double recipT) const
Evaluate reaction rate.
double deltaH() const
Return current enthalpy change of reaction [J/kmol].
vector< pair< size_t, double > > m_stoich_coeffs
Pairs of species indices and multipliers to calculate enthalpy change.
double bondEnergy() const
Return the bond dissociation energy w [J/kmol].
BlowersMaselRate()
Default constructor.
double effectiveActivationEnergy_R(double deltaH_R) const
Return the effective activation energy (a function of the delta H of reaction) divided by the gas con...
double m_deltaH_R
enthalpy change of reaction (in temperature units)
double ddTScaledFromStruct(const BlowersMaselData &shared_data) const
Evaluate derivative of reaction rate with respect to temperature divided by reaction rate.
void updateFromStruct(const BlowersMaselData &shared_data)
Update information specific to reaction.
void setDeltaH(double deltaH)
Set current enthalpy change of reaction [J/kmol].
const string type() const override
String identifying reaction rate specialization.
double activationEnergy() const override
Return the activation energy Ea [J/kmol] The value corresponds to the constant specified by input par...
Public interface for kinetics managers.
Definition Kinetics.h:125
Abstract base class which stores data about a reaction and its rate parameterization so that it can b...
Definition Reaction.h:25
Base class for a phase with thermodynamic properties.
const double GasConstant
Universal Gas Constant [J/kmol/K].
Definition ct_defs.h:120
Namespace for the Cantera kernel.
Definition AnyMap.cpp:564
Data container holding shared data specific to BlowersMaselRate.
int m_state_mf_number
integer that is incremented when composition changes
vector< double > partialMolarEnthalpies
partial molar enthalpies
void update(double T) override
Update data container based on temperature T
void resize(size_t nSpecies, size_t nReactions, size_t nPhases) override
Update number of species, reactions and phases.
bool ready
boolean indicating whether vectors are accessible
double density
used to determine if updates are needed
Data container holding shared data used for ReactionRate calculation.
double recipT
inverse of temperature
virtual void update(double T)
Update data container based on temperature T
double logT
logarithm of temperature
Unit aggregation utility.
Definition Units.h:105