All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
FallOffReactionRateI.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-2019 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 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
27 
28 template<class ReactionRate, class FallOffFunction>
31 (
32  const ReactionRate& k0,
33  const ReactionRate& kInf,
34  const FallOffFunction& F,
35  const thirdBodyEfficiencies& tbes
36 )
37 :
38  k0_(k0),
39  kInf_(kInf),
40  F_(F),
41  thirdBodyEfficiencies_(tbes)
42 {
43  forAll(tbes, i)
44  {
45  beta_.append(Tuple2<label, scalar>(i, tbes[i]));
46  }
47 }
48 
49 
50 template<class ReactionRate, class FallOffFunction>
53 (
54  const speciesTable& species,
55  const dictionary& dict
56 )
57 :
58  k0_(species, dict.subDict("k0")),
59  kInf_(species, dict.subDict("kInf")),
60  F_(dict.subDict("F")),
61  thirdBodyEfficiencies_(species, dict.subDict("thirdBodyEfficiencies"))
62 {
63  forAll(thirdBodyEfficiencies_, i)
64  {
65  beta_.append
66  (
68  (
69  i,
70  thirdBodyEfficiencies_[i]
71  )
72  );
73  }
74 }
75 
76 
77 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
78 
79 template<class ReactionRate, class FallOffFunction>
80 inline Foam::scalar
82 (
83  const scalar p,
84  const scalar T,
85  const scalarField& c,
86  const label li
87 ) const
88 {
89  const scalar k0 = k0_(p, T, c, li);
90  const scalar kInf = kInf_(p, T, c, li);
91  const scalar Pr = k0*thirdBodyEfficiencies_.M(c)/kInf;
92 
93  return kInf*(Pr/(1 + Pr))*F_(T, Pr);
94 }
95 
96 
97 template<class ReactionRate, class FallOffFunction>
98 inline Foam::scalar
100 (
101  const scalar p,
102  const scalar T,
103  const scalarField& c,
104  const label li
105 ) const
106 {
107  const scalar k0 = k0_(p, T, c, li);
108  const scalar kInf = kInf_(p, T, c, li);
109  const scalar Pr = k0*thirdBodyEfficiencies_.M(c)/kInf;
110 
111  return (Pr/(1 + Pr))*F_(T, Pr)*kInf_.ddT(p, T, c, li);
112 }
113 
114 
115 template<class ReactionRate, class FallOffFunction>
117 (
118  const scalar p,
119  const scalar T,
120  const scalarField& c,
121  const label li,
122  scalarField& dcidc
123 ) const
124 {
125  const scalar M = thirdBodyEfficiencies_.M(c);
126 
127  if (M > small)
128  {
129  const scalar k0 = k0_(p, T, c, li);
130  const scalar kInf = kInf_(p, T, c, li);
131  const scalar Pr = k0*M/kInf;
132  const scalar F = F_(T, Pr);
133 
134  forAll(beta_, i)
135  {
136  const scalar dPrdci = -beta_[i].second()*k0/kInf;
137  const scalar dFdci = F_.ddc(Pr, F, dPrdci, T);
138  dcidc[i] = (dPrdci/(Pr*(1 + Pr)) + dFdci/F);
139  }
140  }
141  else
142  {
143  forAll(beta_, i)
144  {
145  dcidc[i] = 0;
146  }
147  }
148 }
149 
150 
151 template<class ReactionRate, class FallOffFunction>
152 inline Foam::scalar
154 (
155  const scalar p,
156  const scalar T,
157  const scalarField& c,
158  const label li
159 ) const
160 {
161  const scalar M = thirdBodyEfficiencies_.M(c);
162 
163  if (M > small)
164  {
165  const scalar k0 = k0_(p, T, c, li);
166  const scalar kInf = kInf_(p, T, c, li);
167 
168  const scalar Pr = k0*thirdBodyEfficiencies_.M(c)/kInf;
169  const scalar F = F_(T, Pr);
170  const scalar dPrdT =
171  Pr*(k0_.ddT(p, T, c, li)/k0 - kInf_.ddT(p, T, c, li)/kInf - 1/T);
172  const scalar dFdT = F_.ddT(Pr, F, dPrdT, T);
173 
174  return (dPrdT/(Pr*(1 + Pr)) + dFdT/F);
175  }
176  else
177  {
178  return 0;
179  }
180 }
181 
182 
183 template<class ReactionRate, class FallOffFunction>
185 (
186  Ostream& os
187 ) const
188 {
189  os << indent << "k0" << nl;
190  os << indent << token::BEGIN_BLOCK << nl;
191  os << incrIndent;
192  k0_.write(os);
193  os << decrIndent;
194  os << indent << token::END_BLOCK << nl;
195 
196  os << indent << "kInf" << nl;
197  os << indent << token::BEGIN_BLOCK << nl;
198  os << incrIndent;
199  kInf_.write(os);
200  os << decrIndent;
201  os << indent << token::END_BLOCK << nl;
202 
203  os << indent << "F" << nl;
204  os << indent << token::BEGIN_BLOCK << nl;
205  os << incrIndent;
206  F_.write(os);
207  os << decrIndent;
208  os << indent << token::END_BLOCK << nl;
209 
210  os << indent << "thirdBodyEfficiencies" << nl;
211  os << indent << token::BEGIN_BLOCK << nl;
212  os << incrIndent;
213  thirdBodyEfficiencies_.write(os);
214  os << decrIndent;
215  os << indent << token::END_BLOCK << nl;
216 }
217 
218 
219 template<class ReactionRate, class FallOffFunction>
220 inline Foam::Ostream& Foam::operator<<
221 (
222  Ostream& os,
224 )
225 {
226  forr.write(os);
227  return os;
228 }
229 
230 
231 // ************************************************************************* //
const dimensionedScalar & F
Faraday constant: default SI units: [C/mol].
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
virtual Ostream & write(const char)=0
Write character.
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
Ostream & indent(Ostream &os)
Indent stream.
Definition: Ostream.H:221
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
A 2-tuple for storing two objects of different types.
Definition: HashTable.H:65
const dimensionedScalar & c
Speed of light in a vacuum.
General class for handling unimolecular/recombination fall-off reactions.
void dcidc(const scalar p, const scalar T, const scalarField &c, const label li, scalarField &dcidc) const
Species concentration derivative of the pressure dependent term.
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
Definition: dictionary.C:934
scalar dcidT(const scalar p, const scalar T, const scalarField &c, const label li) const
Temperature derivative of the pressure dependent term.
scalar ddT(const scalar p, const scalar T, const scalarField &c, const label li) const
void write(Ostream &os) const
Write to stream.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
static const char nl
Definition: Ostream.H:260
Ostream & decrIndent(Ostream &os)
Decrement the indent level.
Definition: Ostream.H:235
const volScalarField & T
FallOffReactionRate(const ReactionRate &k0, const ReactionRate &kInf, const FallOffFunction &F, const thirdBodyEfficiencies &tbes)
Construct from components.
friend Ostream & operator(Ostream &, const FallOffReactionRate< ReactionRate, FallOffFunction > &)
A wordList with hashed indices for faster lookup by name.
volScalarField & p
Third body efficiencies.
Ostream & incrIndent(Ostream &os)
Increment the indent level.
Definition: Ostream.H:228
#define M(I)