35 inline Foam::scalar Foam::fluxLimitedLangmuirHinshelwoodReactionRate::Av
62 additionalAdsorbableSpecieNames_
66 "additionalAdsorbableSpecies",
71 a_(
dict.lookupOrDefault(
"a", scalar(1))),
73 Ta_(
dict.lookup(
"Ta")),
77 AvUniform_(limited_ ?
dict.lookup(
"Av")[0].isNumber() : true),
78 Av_(limited_ && AvUniform_ ?
dict.lookup<scalar>(
"Av") : 0),
79 AvName_(limited_ && !AvUniform_ ?
dict.lookup(
"Av") :
word::null),
97 nReactants_ = lhs.
size();
99 ra_.
setSize(nReactants_ + additionalAdsorbableSpecieNames_.
size());
105 ra_[i] = lhs[i].index;
106 nu_[i] = lhs[i].stoichCoeff;
107 exp_[i] = lhs[i].exponent;
110 forAll(additionalAdsorbableSpecieNames_, i)
112 ra_[nReactants_ + i] = st[additionalAdsorbableSpecieNames_[i]];
121 1 + nReactants_ + additionalAdsorbableSpecieNames_.
size(),
132 1 + nReactants_ + additionalAdsorbableSpecieNames_.
size(),
137 if (!
dict.found(
"m"))
142 s_ =
dict.lookupOrDefault(
"s",
scalarList(nReactants_, scalar(0)));
149 const scalar nCoeffs =
150 1 + nReactants_ + additionalAdsorbableSpecieNames_.
size();
152 if (A_.
size() != nCoeffs)
155 <<
"Number of A coefficients != " << nCoeffs
159 if (Ta_.
size() != nCoeffs)
162 <<
"Number of Ta coefficients != " << nCoeffs
166 if (beta_.
size() != nCoeffs)
169 <<
"Number of beta coefficients != " << nCoeffs
173 if (m_.
size() != nCoeffs)
176 <<
"Number of m coefficients != " << nCoeffs
180 if (s_.
size() != nReactants_)
183 <<
"Number of s coefficients != "
187 if (W_.
size() != nReactants_)
190 <<
"Number of W coefficients != "
219 Foam::fluxLimitedLangmuirHinshelwoodReactionRate::operator()
230 const label ip1 = i + 1;
233 A_[ip1]*
pow(
T, beta_[ip1])*
exp(-Ta_[ip1]/
T)
234 *
pow(
c[ra_[i]], m_[ip1]);
239 const scalar TaByT0 = Ta_[0]/
T;
240 const scalar k0 = A_[0]*
pow(
T, beta_[0])*
exp(-TaByT0);
242 scalar r = k0/
max(
pow(a_ + sumkc, m_[0]), rootVSmall);
249 rc *=
pow(
c[ra_[i]], exp_[i]);
265 return this->Av(li)*r;
278 scalar sumBetaKc = 0;
281 const label ip1 = i + 1;
284 A_[ip1]*
pow(
T, beta_[ip1])*
exp(-Ta_[ip1]/
T)
285 *
pow(
c[ra_[i]], m_[ip1]);
288 sumBetaKc += (beta_[ip1] + Ta_[ip1]/
T)*kc;
291 const scalar TaByT0 = Ta_[0]/
T;
292 const scalar k0 = A_[0]*
pow(
T, beta_[0])*
exp(-TaByT0);
295 ((beta_[0] + TaByT0)*k0*(a_ + sumkc) - m_[0]*k0*sumBetaKc)
296 /(
max(
pow(a_ + sumkc, m_[0] + 1), rootVSmall)*
T);
303 rc *=
pow(
c[ra_[i]], exp_[i]);
306 scalar r = k0/
pow(a_ + sumkc, m_[0]);
327 (s_[l]*
c[ra_[l]]/(nu_[l]*rc))
333 return this->Av(li)*ddT;
364 "additionalAdsorbableSpecies",
365 additionalAdsorbableSpecieNames_
#define forAll(list, i)
Loop across all elements in list.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Input from memory buffer stream.
void size(const label)
Override size to be inconsistent with allocated storage.
void setSize(const label)
Reset size of List.
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....
Langmuir-Hinshelwood reaction rate for gaseous reactions on surfaces including the optional flux limi...
void preEvaluate() const
Pre-evaluation hook.
void postEvaluate() const
Post-evaluation hook.
void write(Ostream &os) const
Write to stream.
fluxLimitedLangmuirHinshelwoodReactionRate(const speciesTable &species, const objectRegistry &ob, const dictionary &dict)
Construct from dictionary.
void ddc(const scalar p, const scalar T, const scalarField &c, const label li, scalarField &ddc) const
scalar ddT(const scalar p, const scalar T, const scalarField &c, const label li) const
A wordList with hashed indices for faster lookup by name.
Registry of regIOobjects.
static void setLRhs(Istream &, const speciesTable &, List< specieCoeffs > &lhs, List< specieCoeffs > &rhs)
Construct the left- and right-hand-side reaction coefficients.
A class for managing temporary objects.
A class for handling words, derived from string.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
const scalar twoPi(2 *pi)
const dimensionedScalar RR
Universal gas constant: default SI units: [J/kmol/K].
Thermodynamic scalar constants.
const dimensionedScalar c
Speed of light in a vacuum.
errorManipArg< error, int > exit(error &err, const int errNo=1)
dimensionedScalar exp(const dimensionedScalar &ds)
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
List< scalar > scalarList
A List of scalars.
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
dimensionedScalar sqrt(const dimensionedScalar &ds)
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)