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
98 const scalar k0 = k0_(p, T, c, li);
99 const scalar kInf = kInf_(p, T, c, li);
100 const scalar Pr = k0*thirdBodyEfficiencies_.M(c)/kInf;
102 return k0*(1/(1 + Pr))*F_(T, Pr);
106 template<
class ReactionRate,
class ChemicallyActivationFunction>
110 ChemicallyActivationFunction
119 const scalar k0 = k0_(p, T, c, li);
120 const scalar kInf = kInf_(p, T, c, li);
121 const scalar Pr = k0*thirdBodyEfficiencies_.M(c)/kInf;
123 return (1/(1 + Pr))*F_(T, Pr)*k0_.ddT(p, T, c, li);
127 template<
class ReactionRate,
class ChemicallyActivationFunction>
131 ChemicallyActivationFunction
141 const scalar
M = thirdBodyEfficiencies_.M(c);
145 const scalar k0 = k0_(p, T, c, li);
146 const scalar kInf = kInf_(p, T, c, li);
148 const scalar Pr = k0*M/kInf;
149 const scalar F = F_(T, Pr);
153 const scalar dPrdci = -beta_[i].second()*k0/kInf;
154 const scalar dFdci = F_.ddc(Pr, F, dPrdci, T);
155 dcidc[i] = (-dPrdci/(1 + Pr) + dFdci/F);
168 template<
class ReactionRate,
class ChemicallyActivationFunction>
172 ChemicallyActivationFunction
181 const scalar
M = thirdBodyEfficiencies_.M(c);
185 const scalar k0 = k0_(p, T, c, li);
186 const scalar kInf = kInf_(p, T, c, li);
188 const scalar Pr = k0*thirdBodyEfficiencies_.M(c)/kInf;
189 const scalar F = F_(T, Pr);
191 Pr*(k0_.ddT(p, T, c, li)/k0 - kInf_.ddT(p, T, c, li)/kInf - 1/
T);
192 const scalar dFdT = F_.ddT(Pr, F, dPrdT, T);
194 return (-dPrdT/(1 + Pr) + dFdT/F);
203 template<
class ReactionRate,
class ChemicallyActivationFunction>
207 ChemicallyActivationFunction
213 thirdBodyEfficiencies_.write(os);
217 template<
class ReactionRate,
class ChemicallyActivationFunction>
222 <ReactionRate, ChemicallyActivationFunction>& carr
#define forAll(list, i)
Loop across all elements in list.
virtual Ostream & write(const char)=0
Write character.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
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.
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.
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.