26 #define CHEMICALLY_ACTIVATION_FUNCTION_JACOBIAN false
30 template<
class ReactionRate,
class ChemicallyActivationFunction>
34 ChemicallyActivationFunction
37 const ReactionRate& k0,
38 const ReactionRate& kInf,
39 const ChemicallyActivationFunction&
F,
46 thirdBodyEfficiencies_(tbes)
50 template<
class ReactionRate,
class ChemicallyActivationFunction>
54 ChemicallyActivationFunction
61 k0_(species,
dict.subDict(
"k0")),
62 kInf_(species,
dict.subDict(
"kInf")),
63 F_(
dict.subDict(
"F")),
64 thirdBodyEfficiencies_(species,
dict.subDict(
"thirdBodyEfficiencies"))
70 template<
class ReactionRate,
class ChemicallyActivationFunction>
74 ChemicallyActivationFunction
75 >::preEvaluate()
const
82 template<
class ReactionRate,
class ChemicallyActivationFunction>
86 ChemicallyActivationFunction
87 >::postEvaluate()
const
94 template<
class ReactionRate,
class ChemicallyActivationFunction>
98 ChemicallyActivationFunction
107 const scalar k0 = k0_(
p,
T,
c, li);
108 const scalar kInf = kInf_(
p,
T,
c, li);
109 const scalar
M = thirdBodyEfficiencies_.M(
c);
110 const scalar Pr = k0/kInf*
M;
111 const scalar
F = F_(
T, Pr);
113 return k0/(1 + Pr)*
F;
117 template<
class ReactionRate,
class ChemicallyActivationFunction>
121 ChemicallyActivationFunction
130 const scalar k0 = k0_(
p,
T,
c, li);
131 const scalar kInf = kInf_(
p,
T,
c, li);
132 const scalar
M = thirdBodyEfficiencies_.M(
c);
133 const scalar Pr = k0/kInf*
M;
134 const scalar
F = F_(
T, Pr);
136 const scalar dk0dT = k0_.ddT(
p,
T,
c, li);
138 #if CHEMICALLY_ACTIVATION_FUNCTION_JACOBIAN
140 const scalar dkInfdT = kInf_.ddT(
p,
T,
c, li);
141 const scalar dPrdT = (
M*dk0dT - Pr*dkInfdT)/kInf;
142 const scalar dFdT = F_.ddT(
T, Pr,
F) + F_.ddPr(
T, Pr,
F)*dPrdT;
146 - k0*dPrdT/
sqr(1 + Pr)*
F
151 return dk0dT/(1 + Pr)*
F;
157 template<
class ReactionRate,
class ChemicallyActivationFunction>
161 ChemicallyActivationFunction
168 template<
class ReactionRate,
class ChemicallyActivationFunction>
172 ChemicallyActivationFunction
182 const scalar k0 = k0_(
p,
T,
c, li);
183 const scalar kInf = kInf_(
p,
T,
c, li);
184 const scalar
M = thirdBodyEfficiencies_.M(
c);
185 const scalar Pr = k0/kInf*
M;
186 const scalar
F = F_(
T, Pr);
188 #if CHEMICALLY_ACTIVATION_FUNCTION_JACOBIAN
191 k0_.ddc(
p,
T,
c, li, dk0dc);
193 kInf_.ddc(
p,
T,
c, li, dkInfdc);
196 scalarField dPrdc((
M*dk0dc - Pr*dkInfdc + k0*dMdc)/kInf);
197 const scalar dFdPr = F_.ddPr(
T, Pr,
F);
201 - k0*dPrdc/
sqr(1 + Pr)*
F
202 + k0/(1 + Pr)*dFdPr*dPrdc;
206 k0_.ddc(
p,
T,
c, li, ddc);
214 template<
class ReactionRate,
class ChemicallyActivationFunction>
218 ChemicallyActivationFunction
242 os <<
indent <<
"thirdBodyEfficiencies" <<
nl;
245 thirdBodyEfficiencies_.
write(os);
251 template<
class ReactionRate,
class ChemicallyActivationFunction>
256 <ReactionRate, ChemicallyActivationFunction>& carr
General class for handling chemically-activated bimolecular reactions.
void preEvaluate() const
Pre-evaluation hook.
void postEvaluate() const
Post-evaluation hook.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
virtual Ostream & write(const char)=0
Write character.
A list of keyword definitions, which are a keyword followed by any number of values (e....
A wordList with hashed indices for faster lookup by name.
A class for managing temporary objects.
const dimensionedScalar F
Faraday constant: default SI units: [C/mol].
const dimensionedScalar c
Speed of light in a vacuum.
void write(std::ostream &os, const bool binary, List< floatScalar > &fField)
Write floats ascii or binary.
Ostream & decrIndent(Ostream &os)
Decrement the indent level.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Ostream & incrIndent(Ostream &os)
Increment the indent level.
Ostream & indent(Ostream &os)
Indent stream.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)