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-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 
26 #include "ReversibleReaction.H"
27 
28 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
29 
30 template<class ThermoType, class ReactionRate>
32 (
34  const ReactionRate& k
35 )
36 :
37  Reaction<ThermoType>(reaction),
38  k_(k)
39 {}
40 
41 
42 template<class ThermoType, class ReactionRate>
44 (
45  const speciesTable& species,
46  const PtrList<ThermoType>& speciesThermo,
47  const dictionary& dict
48 )
49 :
50  Reaction<ThermoType>(species, speciesThermo, dict),
51  k_(species, this->kfDims(), dict)
52 {}
53 
54 
55 template<class ThermoType, class ReactionRate>
57 (
58  const speciesTable& species,
59  const PtrList<ThermoType>& speciesThermo,
60  const objectRegistry& ob,
61  const dictionary& dict
62 )
63 :
64  Reaction<ThermoType>(species, speciesThermo, dict),
65  k_(species, ob, this->kfDims(), dict)
66 {}
67 
68 
69 template<class ThermoType, class ReactionRate>
71 (
73  const speciesTable& species
74 )
75 :
76  Reaction<ThermoType>(rr, species),
77  k_(rr.k_)
78 {}
79 
80 
81 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
82 
83 template<class ThermoType, class ReactionRate>
85 {
86  k_.preEvaluate();
87 }
88 
89 
90 template<class ThermoType, class ReactionRate>
92 {
93  k_.postEvaluate();
94 }
95 
96 
97 template<class ThermoType, class ReactionRate>
99 (
100  const scalar p,
101  const scalar T,
102  const scalarField& c,
103  const label li
104 ) const
105 {
106  return k_(p, T, c, li);
107 }
108 
109 
110 template<class ThermoType, class ReactionRate>
112 (
113  const scalar kfwd,
114  const scalar p,
115  const scalar T,
116  const scalarField& c,
117  const label li
118 ) const
119 {
120  return kfwd/max(this->Kc(p, T), rootSmall);
121 }
122 
123 
124 template<class ThermoType, class ReactionRate>
126 (
127  const scalar p,
128  const scalar T,
129  const scalarField& c,
130  const label li
131 ) const
132 {
133  return kr(kf(p, T, c, li), p, T, c, li);
134 }
135 
136 
137 template<class ThermoType, class ReactionRate>
139 (
140  const scalar p,
141  const scalar T,
142  const scalarField& c,
143  const label li
144 ) const
145 {
146  return k_.ddT(p, T, c, li);
147 }
148 
149 
150 template<class ThermoType, class ReactionRate>
152 (
153  const scalar p,
154  const scalar T,
155  const scalarField& c,
156  const label li,
157  const scalar dkfdT,
158  const scalar kr
159 ) const
160 {
161  const scalar Kc = max(this->Kc(p, T), rootSmall);
162 
163  return dkfdT/Kc - (Kc > rootSmall ? kr*this->dKcdTbyKc(p, T) : 0);
164 }
165 
166 
167 template<class ThermoType, class ReactionRate>
168 bool
170 {
171  return k_.hasDdc();
172 }
173 
174 
175 template<class ThermoType, class ReactionRate>
177 (
178  const scalar p,
179  const scalar T,
180  const scalarField& c,
181  const label li,
182  scalarField& dkfdc
183 ) const
184 {
185  k_.ddc(p, T, c, li, dkfdc);
186 }
187 
188 
189 template<class ThermoType, class ReactionRate>
191 (
192  const scalar p,
193  const scalar T,
194  const scalarField& c,
195  const label li,
196  const scalarField& dkfdc,
197  const scalar kr,
198  scalarField& dkrdc
199 ) const
200 {
201  const scalar Kc = max(this->Kc(p, T), rootSmall);
202 
203  dkrdc = dkfdc/Kc;
204 }
205 
206 
207 template<class ThermoType, class ReactionRate>
209 (
210  Ostream& os
211 ) const
212 {
214  k_.write(os);
215 }
216 
217 
218 // ************************************************************************* //
label k
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:57
Simple extension of ThermoType to handle reaction kinetics in addition to the equilibrium thermodynam...
Definition: Reaction.H:72
virtual void write(Ostream &) const
Write.
Definition: Reaction.C:550
Extension of Reaction to handle reversible reactions using equilibrium thermodynamics.
virtual scalar kr(const scalar kfwd, const scalar p, const scalar T, const scalarField &c, const label li) const
Reverse rate constant from the given forward rate constant.
virtual scalar dkfdT(const scalar p, const scalar T, const scalarField &c, const label li) const
Temperature derivative of forward rate.
virtual void write(Ostream &) const
Write.
virtual void preEvaluate() const
Pre-evaluation hook.
virtual void postEvaluate() const
Post-evaluation hook.
virtual scalar dkrdT(const scalar p, const scalar T, const scalarField &c, const label li, const scalar dkfdT, const scalar kr) const
Temperature derivative of backward rate.
void dkfdc(const scalar p, const scalar T, const scalarField &c, const label li, scalarField &dkfdc) const
Concentration derivative of forward rate.
virtual bool hasDkdc() const
Does this reaction have concentration-dependent rate constants?
virtual scalar kf(const scalar p, const scalar T, const scalarField &c, const label li) const
Forward rate constant.
void dkrdc(const scalar p, const scalar T, const scalarField &c, const label li, const scalarField &dkfdc, const scalar kr, scalarField &dkrdc) const
Concentration derivative of reverse rate.
ReversibleReaction(const Reaction< ThermoType > &reaction, const ReactionRate &k)
Construct from components.
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
A wordList with hashed indices for faster lookup by name.
Registry of regIOobjects.
Reaction base-class holding the specie names and coefficients.
Definition: reaction.H:57
const dimensionedScalar c
Speed of light in a vacuum.
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
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
dictionary dict
volScalarField & p