ReversibleReaction.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 
26 #include "ReversibleReaction.H"
27 
28 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
29 
30 template
31 <
32  template<class> class ReactionType,
33  class ReactionThermo,
34  class ReactionRate
35 >
38 (
39  const ReactionType<ReactionThermo>& reaction,
40  const ReactionRate& k
41 )
42 :
43  ReactionType<ReactionThermo>(reaction),
44  k_(k)
45 {}
46 
47 
48 template
49 <
50  template<class> class ReactionType,
51  class ReactionThermo,
52  class ReactionRate
53 >
56 (
57  const speciesTable& species,
58  const HashPtrTable<ReactionThermo>& thermoDatabase,
59  const dictionary& dict
60 )
61 :
62  ReactionType<ReactionThermo>(species, thermoDatabase, dict),
63  k_(species, dict)
64 {}
65 
66 
67 template
68 <
69  template<class> class ReactionType,
70  class ReactionThermo,
71  class ReactionRate
72 >
75 (
77  const speciesTable& species
78 )
79 :
80  ReactionType<ReactionThermo>(rr, species),
81  k_(rr.k_)
82 {}
83 
84 
85 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
86 
87 template
88 <
89  template<class> class ReactionType,
90  class ReactionThermo,
91  class ReactionRate
92 >
93 Foam::scalar Foam::ReversibleReaction
94 <
95  ReactionType,
96  ReactionThermo,
97  ReactionRate
98 >::kf
99 (
100  const scalar p,
101  const scalar T,
102  const scalarField& c
103 ) const
104 {
105  return k_(p, T, c);
106 }
107 
108 
109 template
110 <
111  template<class> class ReactionType,
112  class ReactionThermo,
113  class ReactionRate
114 >
115 Foam::scalar Foam::ReversibleReaction
116 <
117  ReactionType,
118  ReactionThermo,
119  ReactionRate
120 >::kr
121 (
122  const scalar kfwd,
123  const scalar p,
124  const scalar T,
125  const scalarField& c
126 ) const
127 {
128  return kfwd/max(this->Kc(p, T), rootSmall);
129 }
130 
131 
132 template
133 <
134  template<class> class ReactionType,
135  class ReactionThermo,
136  class ReactionRate
137 >
138 Foam::scalar Foam::ReversibleReaction
139 <
140  ReactionType,
141  ReactionThermo,
142  ReactionRate
143 >::kr
144 (
145  const scalar p,
146  const scalar T,
147  const scalarField& c
148 ) const
149 {
150  return kr(kf(p, T, c), p, T, c);
151 }
152 
153 
154 template
155 <
156  template<class> class ReactionType,
157  class ReactionThermo,
158  class ReactionRate
159 >
160 Foam::scalar Foam::ReversibleReaction
161 <
162  ReactionType,
163  ReactionThermo,
164  ReactionRate
165 >::dkfdT
166 (
167  const scalar p,
168  const scalar T,
169  const scalarField& c
170 ) const
171 {
172  return k_.ddT(p, T, c);
173 }
174 
175 
176 template
177 <
178  template<class> class ReactionType,
179  class ReactionThermo,
180  class ReactionRate
181 >
182 Foam::scalar Foam::ReversibleReaction
183 <
184  ReactionType,
185  ReactionThermo,
186  ReactionRate
187 >::dkrdT
188 (
189  const scalar p,
190  const scalar T,
191  const scalarField& c,
192  const scalar dkfdT,
193  const scalar kr
194 ) const
195 {
196  scalar Kc = max(this->Kc(p, T), rootSmall);
197 
198  return dkfdT/Kc - kr*this->dKcdTbyKc(p, T);
199 }
200 
201 
202 template
203 <
204  template<class> class ReactionType,
205  class ReactionThermo,
206  class ReactionRate
207 >
210 <
211  ReactionType,
212  ReactionThermo,
213  ReactionRate
214 >::beta() const
215 {
216  return k_.beta();
217 }
218 
219 
220 template
221 <
222  template<class> class ReactionType,
223  class ReactionThermo,
224  class ReactionRate
225 >
227 <
228  ReactionType,
229  ReactionThermo,
230  ReactionRate
231 >::dcidc
232 (
233  const scalar p,
234  const scalar T,
235  const scalarField& c,
236  scalarField& dcidc
237 ) const
238 {
239  k_.dcidc(p, T, c, dcidc);
240 }
241 
242 
243 template
244 <
245  template<class> class ReactionType,
246  class ReactionThermo,
247  class ReactionRate
248 >
249 Foam::scalar Foam::ReversibleReaction
250 <
251  ReactionType,
252  ReactionThermo,
253  ReactionRate
254 >::dcidT
255 (
256  const scalar p,
257  const scalar T,
258  const scalarField& c
259 ) const
260 {
261  return k_.dcidT(p, T, c);
262 }
263 
264 
265 template
266 <
267  template<class> class ReactionType,
268  class ReactionThermo,
269  class ReactionRate
270 >
272 <
273  ReactionType,
274  ReactionThermo,
275  ReactionRate
276 >::write
277 (
278  Ostream& os
279 ) const
280 {
282  k_.write(os);
283 }
284 
285 
286 // ************************************************************************* //
dictionary dict
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
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
ReversibleReaction(const ReactionType< ReactionThermo > &reaction, const ReactionRate &k)
Construct from components.
virtual scalar dcidT(const scalar p, const scalar T, const scalarField &c) const
Temperature derivative of the pressure dependent term.
A HashTable specialization for hashing pointers.
Definition: HashPtrTable.H:50
void write(Ostream &, const label, const dictionary &)
Write with dictionary lookup.
CombustionModel< rhoReactionThermo > & reaction
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
A wordList with hashed indices for faster lookup by name.
dimensionedScalar beta("beta", dimless/dimTemperature, laminarTransport)
virtual void dcidc(const scalar p, const scalar T, const scalarField &c, scalarField &dcidc) const
Species concentration derivative of the pressure dependent term.
Simple extension of Reaction to handle reversible reactions using equilibrium thermodynamics.