LangmuirHinshelwoodReactionRateI.H
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2024 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
27 
28 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
29 
31 (
32  const speciesTable& st,
33  const dimensionSet& dims,
34  const dictionary& dict
35 )
36 :
37  reactantNames_(dict.lookup("reactants")),
38  a_(dict.lookupOrDefault<scalar>("a", 1)),
39  A_(dict.lookup("A")),
40  Ta_(dict.lookup("Ta")),
41  beta_
42  (
43  dict.lookupOrDefault<FixedList<scalar, 3>>
44  (
45  "beta",
46  FixedList<scalar, 3>({0, 0, 0})
47  )
48  ),
49  m_
50  (
52  (
53  "m",
54  FixedList<scalar, 3>({2, 1, 1})
55  )
56  )
57 {
58  forAll(reactantNames_, i)
59  {
60  r_[i] = st[reactantNames_[i]];
61  }
62 }
63 
64 
65 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
66 
68 {}
69 
70 
72 {}
73 
74 
75 inline Foam::scalar Foam::LangmuirHinshelwoodReactionRate::operator()
76 (
77  const scalar p,
78  const scalar T,
79  const scalarField& c,
80  const label li
81 ) const
82 {
83  const scalar c1m = pow(c[r_[0]], m_[1]);
84  const scalar c2m = pow(c[r_[1]], m_[2]);
85 
86  const scalar k0 = A_[0]*pow(T, beta_[0])*exp(-Ta_[0]/T);
87  const scalar k1 = A_[1]*pow(T, beta_[1])*exp(-Ta_[1]/T);
88  const scalar k2 = A_[2]*pow(T, beta_[2])*exp(-Ta_[2]/T);
89 
90  const scalar b = a_ + k1*c1m + k2*c2m;
91 
92  return k0/pow(b, m_[0]);
93 }
94 
95 
97 (
98  const scalar p,
99  const scalar T,
100  const scalarField& c,
101  const label li
102 ) const
103 {
104  const scalar c1m = pow(c[r_[0]], m_[1]);
105  const scalar c2m = pow(c[r_[1]], m_[2]);
106 
107  const scalar k0 = A_[0]*pow(T, beta_[0])*exp(-Ta_[0]/T);
108  const scalar k1 = A_[1]*pow(T, beta_[1])*exp(-Ta_[1]/T);
109  const scalar k2 = A_[2]*pow(T, beta_[2])*exp(-Ta_[2]/T);
110 
111  const scalar dk0dT = k0/T*(beta_[0] + Ta_[0]/T);
112  const scalar dk1dT = k1/T*(beta_[1] + Ta_[1]/T);
113  const scalar dk2dT = k2/T*(beta_[2] + Ta_[2]/T);
114 
115  const scalar b = a_ + k1*c1m + k2*c2m;
116  const scalar dbdT = dk1dT*c1m + dk2dT*c2m;
117 
118  return (dk0dT - k0*m_[0]*dbdT/b)/pow(b, m_[0]);
119 }
120 
121 
123 {
124  return true;
125 }
126 
127 
129 (
130  const scalar p,
131  const scalar T,
132  const scalarField& c,
133  const label li,
134  scalarField& ddc
135 ) const
136 {
137  const scalar c1m = pow(c[r_[0]], m_[1]);
138  const scalar c2m = pow(c[r_[1]], m_[2]);
139 
140  const scalar dc1mdc1 = m_[1]*c1m/c[r_[0]];
141  const scalar dc2mdc2 = m_[2]*c2m/c[r_[1]];
142 
143  const scalar k0 = A_[0]*pow(T, beta_[0])*exp(-Ta_[0]/T);
144  const scalar k1 = A_[1]*pow(T, beta_[1])*exp(-Ta_[1]/T);
145  const scalar k2 = A_[2]*pow(T, beta_[2])*exp(-Ta_[2]/T);
146 
147  const scalar b = a_ + k1*c1m + k2*c2m;
148  const scalar dbdc1 = k1*dc1mdc1;
149  const scalar dbdc2 = k2*dc2mdc2;
150 
151  const scalar k = k0/pow(b, m_[0]);
152  const scalar dkdb = - k*m_[0]/b;
153 
154  ddc = 0;
155  ddc[r_[0]] = dkdb*dbdc1;
156  ddc[r_[1]] = dkdb*dbdc2;
157 }
158 
159 
161 {
162  writeEntry(os, "reactants", reactantNames_);
163  writeEntry(os, "a", a_);
164  writeEntry(os, "A", A_);
165  writeEntry(os, "Ta", Ta_);
166  writeEntry(os, "beta", beta_);
167  writeEntry(os, "m", m_);
168 }
169 
170 
171 inline Foam::Ostream& Foam::operator<<
172 (
173  Ostream& os,
175 )
176 {
177  lhrr.write(os);
178  return os;
179 }
180 
181 
182 // ************************************************************************* //
label k
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
A 1D vector of objects of type <T> with a fixed size <Size>.
Definition: FixedList.H:78
Langmuir-Hinshelwood reaction rate for gaseous reactions on surfaces.
LangmuirHinshelwoodReactionRate(const speciesTable &species, const dimensionSet &dims, const dictionary &dict)
Construct from dictionary.
void write(Ostream &os) const
Write to stream.
void ddc(const scalar p, const scalar T, const scalarField &c, const label li, scalarField &ddc) const
The derivative of the rate w.r.t. concentration.
scalar ddT(const scalar p, const scalar T, const scalarField &c, const label li) const
The derivative of the rate w.r.t. temperature.
bool hasDdc() const
Is the rate a function of concentration?
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:57
virtual Ostream & write(const char)=0
Write character.
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T, if not found return the given default.
Dimension set for the base types.
Definition: dimensionSet.H:125
A wordList with hashed indices for faster lookup by name.
volScalarField & b
Definition: createFields.H:25
const dimensionedScalar c
Speed of light in a vacuum.
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.
Definition: label.H:59
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
dictionary dict
volScalarField & p