NonEquilibriumReversibleReaction.C
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-2018 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 
30 template
31 <
32  template<class> class ReactionType,
33  class ReactionThermo,
34  class ReactionRate
35 >
37 <
38  ReactionType,
39  ReactionThermo,
40  ReactionRate
41 >::
42 NonEquilibriumReversibleReaction
43 (
44  const ReactionType<ReactionThermo>& reaction,
45  const ReactionRate& forwardReactionRate,
46  const ReactionRate& reverseReactionRate
47 )
48 :
49  ReactionType<ReactionThermo>(reaction),
50  fk_(forwardReactionRate),
51  rk_(reverseReactionRate)
52 {}
53 
54 
55 template
56 <
57  template<class> class ReactionType,
58  class ReactionThermo,
59  class ReactionRate
60 >
62 <
63  ReactionType,
64  ReactionThermo,
65  ReactionRate
66 >::
68 (
69  const speciesTable& species,
70  const HashPtrTable<ReactionThermo>& thermoDatabase,
71  const dictionary& dict
72 )
73 :
74  ReactionType<ReactionThermo>(species, thermoDatabase, dict),
75  fk_(species, dict.subDict("forward")),
76  rk_(species, dict.subDict("reverse"))
77 {}
78 
79 
80 template
81 <
82  template<class> class ReactionType,
83  class ReactionThermo,
84  class ReactionRate
85 >
87 <
88  ReactionType,
89  ReactionThermo,
90  ReactionRate
91 >::
93 (
95  <
96  ReactionType,
97  ReactionThermo,
98  ReactionRate
99  >& nerr,
100  const speciesTable& species
101 )
102 :
103  ReactionType<ReactionThermo>(nerr, species),
104  fk_(nerr.fk_),
105  rk_(nerr.rk_)
106 {}
107 
108 
109 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
110 
111 template
112 <
113  template<class> class ReactionType,
114  class ReactionThermo,
115  class ReactionRate
116 >
117 Foam::scalar
119 <
120  ReactionType,
121  ReactionThermo,
122  ReactionRate
123 >::kf
124 (
125  const scalar p,
126  const scalar T,
127  const scalarField& c
128 ) const
129 {
130  return fk_(p, T, c);
131 }
132 
133 
134 template
135 <
136  template<class> class ReactionType,
137  class ReactionThermo,
138  class ReactionRate
139 >
140 Foam::scalar
142 <
143  ReactionType,
144  ReactionThermo,
145  ReactionRate
146 >::kr
147 (
148  const scalar,
149  const scalar p,
150  const scalar T,
151  const scalarField& c
152 ) const
153 {
154  return rk_(p, T, c);
155 }
156 
157 
158 template
159 <
160  template<class> class ReactionType,
161  class ReactionThermo,
162  class ReactionRate
163 >
164 Foam::scalar
166 <
167  ReactionType,
168  ReactionThermo,
169  ReactionRate
170 >::kr
171 (
172  const scalar p,
173  const scalar T,
174  const scalarField& c
175 ) const
176 {
177  return rk_(p, T, c);
178 }
179 
180 
181 template
182 <
183  template<class> class ReactionType,
184  class ReactionThermo,
185  class ReactionRate
186 >
187 Foam::scalar
189 <
190  ReactionType,
191  ReactionThermo,
192  ReactionRate
193 >::dkfdT
194 (
195  const scalar p,
196  const scalar T,
197  const scalarField& c
198 ) const
199 {
200  return fk_.ddT(p, T, c);
201 }
202 
203 
204 template
205 <
206  template<class> class ReactionType,
207  class ReactionThermo,
208  class ReactionRate
209 >
210 Foam::scalar
212 <
213  ReactionType,
214  ReactionThermo,
215  ReactionRate
216 >::dkrdT
217 (
218  const scalar p,
219  const scalar T,
220  const scalarField& c,
221  const scalar dkfdT,
222  const scalar kr
223 ) const
224 {
225  return rk_.ddT(p, T, c);
226 }
227 
228 
229 template
230 <
231  template<class> class ReactionType,
232  class ReactionThermo,
233  class ReactionRate
234 >
237 <
238  ReactionType,
239  ReactionThermo,
240  ReactionRate
241 >::beta() const
242 {
243  return fk_.beta();
244 }
245 
246 
247 template
248 <
249  template<class> class ReactionType,
250  class ReactionThermo,
251  class ReactionRate
252 >
254 <
255  ReactionType,
256  ReactionThermo,
257  ReactionRate
258 >::dcidc
259 (
260  const scalar p,
261  const scalar T,
262  const scalarField& c,
263  scalarField& dcidc
264 ) const
265 {
266  fk_.dcidc(p, T, c, dcidc);
267 }
268 
269 
270 template
271 <
272  template<class> class ReactionType,
273  class ReactionThermo,
274  class ReactionRate
275 >
277 <
278  ReactionType,
279  ReactionThermo,
280  ReactionRate
281 >::dcidT
282 (
283  const scalar p,
284  const scalar T,
285  const scalarField& c
286 ) const
287 {
288  return fk_.dcidT(p, T, c);
289 }
290 
291 
292 template
293 <
294  template<class> class ReactionType,
295  class ReactionThermo,
296  class ReactionRate
297 >
299 <
300  ReactionType,
301  ReactionThermo,
302  ReactionRate
303 >::write
304 (
305  Ostream& os
306 ) const
307 {
309 
310  os << indent << "forward" << nl;
311  os << indent << token::BEGIN_BLOCK << nl;
312  os << incrIndent;
313  fk_.write(os);
314  os << decrIndent;
315  os << indent << token::END_BLOCK << nl;
316 
317  os << indent << "reverse" << nl;
318  os << indent << token::BEGIN_BLOCK << nl;
319  os << incrIndent;
320  rk_.write(os);
321  os << decrIndent;
322  os << indent << token::END_BLOCK << nl;
323 }
324 
325 
326 // ************************************************************************* //
virtual scalar dcidT(const scalar p, const scalar T, const scalarField &c) const
Temperature derivative of the pressure dependent term.
dictionary dict
Ostream & indent(Ostream &os)
Indent stream.
Definition: Ostream.H:226
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
virtual void dcidc(const scalar p, const scalar T, const scalarField &c, scalarField &dcidc) const
Species concentration derivative of the pressure dependent term.
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: HashTable.H:60
A HashTable specialization for hashing pointers.
Definition: HashPtrTable.H:50
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
Definition: dictionary.C:692
void write(Ostream &, const label, const dictionary &)
Write with dictionary lookup.
Simple extension of Reaction to handle reversible reactions using equilibrium thermodynamics.
CombustionModel< rhoReactionThermo > & reaction
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
static const char nl
Definition: Ostream.H:265
Ostream & decrIndent(Ostream &os)
Decrement the indent level.
Definition: Ostream.H:240
A wordList with hashed indices for faster lookup by name.
dimensionedScalar beta("beta", dimless/dimTemperature, laminarTransport)
virtual Ostream & write(const token &)=0
Write next token to stream.
Ostream & incrIndent(Ostream &os)
Increment the indent level.
Definition: Ostream.H:233