28 template<
class ReactionRate,
class ChemicallyActivationFunction>
32 ChemicallyActivationFunction
33 >::ChemicallyActivatedReactionRate
35 const ReactionRate& k0,
36 const ReactionRate& kInf,
37 const ChemicallyActivationFunction& F,
44 thirdBodyEfficiencies_(tbes)
53 template<
class ReactionRate,
class ChemicallyActivationFunction>
57 ChemicallyActivationFunction
67 thirdBodyEfficiencies_(species, dict)
69 forAll(thirdBodyEfficiencies_, i)
76 thirdBodyEfficiencies_[i]
85 template<
class ReactionRate,
class ChemicallyActivationFunction>
89 ChemicallyActivationFunction
97 const scalar k0 = k0_(p, T, c);
98 const scalar kInf = kInf_(p, T, c);
99 const scalar
Pr = k0*thirdBodyEfficiencies_.M(c)/kInf;
101 return k0*(1/(1 +
Pr))*F_(T, Pr);
105 template<
class ReactionRate,
class ChemicallyActivationFunction>
109 ChemicallyActivationFunction
117 const scalar k0 = k0_(p, T, c);
118 const scalar kInf = kInf_(p, T, c);
119 const scalar
Pr = k0*thirdBodyEfficiencies_.M(c)/kInf;
121 return (1/(1 + Pr))*F_(T, Pr)*k0_.ddT(p, T, c);
125 template<
class ReactionRate,
class ChemicallyActivationFunction>
129 ChemicallyActivationFunction
138 const scalar
M = thirdBodyEfficiencies_.M(c);
142 const scalar k0 = k0_(p, T, c);
143 const scalar kInf = kInf_(p, T, c);
145 const scalar
Pr = k0*M/kInf;
146 const scalar F = F_(T, Pr);
150 const scalar dPrdci = -beta_[i].second()*k0/kInf;
151 const scalar dFdci = F_.ddc(Pr, F, dPrdci, T);
152 dcidc[i] = (-dPrdci/(1 +
Pr) + dFdci/F);
165 template<
class ReactionRate,
class ChemicallyActivationFunction>
169 ChemicallyActivationFunction
177 const scalar
M = thirdBodyEfficiencies_.M(c);
181 const scalar k0 = k0_(p, T, c);
182 const scalar kInf = kInf_(p, T, c);
184 const scalar
Pr = k0*thirdBodyEfficiencies_.M(c)/kInf;
185 const scalar F = F_(T, Pr);
187 Pr*(k0_.ddT(p, T, c)/k0 - kInf_.ddT(p, T, c)/kInf - 1/
T);
188 const scalar dFdT = F_.ddT(Pr, F, dPrdT, T);
190 return (-dPrdT/(1 + Pr) + dFdT/F);
199 template<
class ReactionRate,
class ChemicallyActivationFunction>
203 ChemicallyActivationFunction
209 thirdBodyEfficiencies_.write(os);
213 template<
class ReactionRate,
class ChemicallyActivationFunction>
218 <ReactionRate, ChemicallyActivationFunction>& carr
dimensionedScalar Pr("Pr", dimless, laminarTransport)
#define forAll(list, i)
Loop across all elements in list.
A list of keyword definitions, which are a keyword followed by any number of values (e...
A 2-tuple for storing two objects of different types.
General class for handling chemically-activated bimolecular reactions.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
A wordList with hashed indices for faster lookup by name.
const dimensionedScalar c
Speed of light in a vacuum.
virtual Ostream & write(const token &)=0
Write next token to stream.