TroeFallOffFunctionI.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 
29 (
30  const scalar alpha,
31  const scalar Tsss,
32  const scalar Ts,
33  const scalar Tss
34 )
35 :
36  alpha_(alpha),
37  Tsss_(Tsss),
38  Ts_(Ts),
39  Tss_(Tss)
40 {}
41 
42 
44 :
45  alpha_(readScalar(dict.lookup("alpha"))),
46  Tsss_(readScalar(dict.lookup("Tsss"))),
47  Ts_(readScalar(dict.lookup("Ts"))),
48  Tss_(readScalar(dict.lookup("Tss")))
49 {}
50 
51 
52 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
53 
54 inline Foam::scalar Foam::TroeFallOffFunction::operator()
55 (
56  const scalar T,
57  const scalar Pr
58 ) const
59 {
60  scalar logFcent = log10
61  (
62  max
63  (
64  (1 - alpha_)*exp(-T/Tsss_) + alpha_*exp(-T/Ts_) + exp(-Tss_/T),
65  small
66  )
67  );
68 
69  scalar c = -0.4 - 0.67*logFcent;
70  static const scalar d = 0.14;
71  scalar n = 0.75 - 1.27*logFcent;
72 
73  scalar logPr = log10(max(Pr, small));
74  return pow(10, logFcent/(1 + sqr((logPr + c)/(n - d*(logPr + c)))));
75 }
76 
77 
78 inline Foam::scalar Foam::TroeFallOffFunction::ddT
79 (
80  const scalar Pr,
81  const scalar F,
82  const scalar dPrdT,
83  const scalar T
84 ) const
85 {
86  scalar logPr = log10(max(Pr, small));
87  scalar logTen = log(10.0);
88  scalar Fcent =
89  (
90  max
91  (
92  (1 - alpha_)*exp(-T/Tsss_) + alpha_*exp(-T/Ts_) + exp(-Tss_/T),
93  small
94  )
95  );
96  scalar logFcent = log10(Fcent);
97 
98  scalar dFcentdT =
99  (
100  (alpha_ - 1)*exp(-T/Tsss_)/Tsss_
101  - alpha_*exp(-T/Ts_)/Ts_
102  + Tss_*exp(-Tss_/T)/sqr(T)
103  );
104 
105  scalar d = 0.14;
106  scalar dlogFcentdT = dFcentdT/Fcent/logTen;
107  scalar c = -0.4 - 0.67*logFcent;
108  scalar dcdT = -0.67*dlogFcentdT;
109  scalar n = 0.75 - 1.27*logFcent;
110  scalar dndT = -1.27*dlogFcentdT;
111 
112  scalar dlogPrdT = dPrdT/Pr/logTen;
113 
114  scalar dParentdT =
115  2.0*(logPr + c)/sqr(n - d*(logPr + c))
116  *(
117  (dlogPrdT + dcdT)
118  - (logPr + c)*(dndT - d*(dlogPrdT + dcdT))/(n - d*(logPr + c))
119  );
120 
121  return
122  (
123  F*logTen
124  *(
125  dlogFcentdT/(1.0 + sqr((logPr + c)/(n - d*(logPr + c))))
126  - logFcent*dParentdT/sqr(1.0 + sqr((logPr + c)/(n - d*(logPr + c))))
127  )
128  );
129 }
130 
131 
132 inline Foam::scalar Foam::TroeFallOffFunction::ddc
133 (
134  const scalar Pr,
135  const scalar F,
136  const scalar dPrdc,
137  const scalar T
138 ) const
139 {
140  scalar logPr = log10(max(Pr, small));
141  scalar logTen = log(10.0);
142  scalar logFcent = log10
143  (
144  max
145  (
146  (1 - alpha_)*exp(-T/Tsss_) + alpha_*exp(-T/Ts_) + exp(-Tss_/T),
147  small
148  )
149  );
150 
151  scalar dlogPrdc = dPrdc/Pr/logTen;
152  scalar d = 0.14;
153  scalar c = -0.4 - 0.67*logFcent;
154  scalar n = 0.75 - 1.27*logFcent;
155 
156  scalar dParentdc =
157  2.0*(logPr + c)/sqr(n - d*(logPr + c))
158  *(
159  (dlogPrdc)
160  - (logPr + c)*(-d*(dlogPrdc))/(n - d*(logPr + c))
161  );
162 
163  return
164  (
165  F*logTen
166  *(
167  - logFcent*dParentdc/sqr(1.0 + sqr((logPr + c)/(n - d*(logPr + c))))
168  )
169  );
170 }
171 
172 
174 {
175  writeEntry(os, "alpha", alpha_);
176  writeEntry(os, "Tsss", Tsss_);
177  writeEntry(os, "Ts", Ts_);
178  writeEntry(os, "Tss", Tss_);
179 }
180 
181 
182 // * * * * * * * * * * * * * * * Ostream Operator * * * * * * * * * * * * * //
183 
184 inline Foam::Ostream& Foam::operator<<
185 (
186  Foam::Ostream& os,
187  const Foam::TroeFallOffFunction& tfof
188 )
189 {
190  tfof.write(os);
191  return os;
192 }
193 
194 
195 // ************************************************************************* //
#define readScalar
Definition: doubleScalar.C:38
dictionary dict
dimensionedScalar log(const dimensionedScalar &ds)
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
The Troe fall-off function.
void write(Ostream &os) const
Write to stream.
stressControl lookup("compactNormalStress") >> compactNormalStress
dimensionedScalar exp(const dimensionedScalar &ds)
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
TroeFallOffFunction(const scalar alpha, const scalar Tsss, const scalar Ts, const scalar Tss)
Construct from components.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
scalar ddc(const scalar Pr, const scalar F, const scalar dPrdc, const scalar T) const
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
const dimensionedScalar c
Speed of light in a vacuum.
label n
virtual Ostream & write(const token &)=0
Write next token to stream.
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
scalar ddT(const scalar Pr, const scalar F, const scalar dPrdT, const scalar T) const
dimensionedScalar log10(const dimensionedScalar &ds)