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-2021 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_(dict.lookup<scalar>("alpha")),
46  Tsss_(dict.lookup<scalar>("Tsss")),
47  Ts_(dict.lookup<scalar>("Ts")),
48  Tss_(dict.lookup<scalar>("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  const scalar logPr = log10(max(Pr, small));
61 
62  const scalar Fcent =
63  (1 - alpha_)*exp(-T/Tsss_)
64  + alpha_*exp(-T/Ts_)
65  + exp(-Tss_/T);
66 
67  const scalar logFcent = log10(max(Fcent, small));
68 
69  const scalar c = -0.4 - 0.67*logFcent;
70  static const scalar d = 0.14;
71  const scalar n = 0.75 - 1.27*logFcent;
72 
73  const scalar x1 = n - d*(logPr + c);
74  const scalar x2 = (logPr + c)/x1;
75  const scalar x3 = 1 + sqr(x2);
76  const scalar x4 = logFcent/x3;
77 
78  return pow(10, x4);
79 }
80 
81 
82 inline Foam::scalar Foam::TroeFallOffFunction::ddT
83 (
84  const scalar T,
85  const scalar Pr,
86  const scalar F
87 ) const
88 {
89  static const scalar logTen = log(scalar(10));
90 
91  const scalar logPr = log10(max(Pr, small));
92 
93  const scalar Fcent =
94  (1 - alpha_)*exp(-T/Tsss_)
95  + alpha_*exp(-T/Ts_)
96  + exp(-Tss_/T);
97  const scalar dFcentdT =
98  - (1 - alpha_)/Tsss_*exp(-T/Tsss_)
99  - alpha_/Ts_*exp(-T/Ts_)
100  + Tss_/sqr(T)*exp(-Tss_/T);
101 
102  const scalar logFcent = log10(max(Fcent, small));
103  const scalar dlogFcentdT = Fcent >= small ? dFcentdT/Fcent/logTen : 0;
104 
105  const scalar c = -0.4 - 0.67*logFcent;
106  const scalar dcdT = -0.67*dlogFcentdT;
107  static const scalar d = 0.14;
108  const scalar n = 0.75 - 1.27*logFcent;
109  const scalar dndT = - 1.27*dlogFcentdT;
110 
111  const scalar x1 = n - d*(logPr + c);
112  const scalar dx1dT = dndT - d*dcdT;
113  const scalar x2 = (logPr + c)/x1;
114  const scalar dx2dT = (dcdT - x2*dx1dT)/x1;
115  const scalar x3 = 1 + sqr(x2);
116  const scalar dx3dT = 2*x2*dx2dT;
117  const scalar x4 = logFcent/x3;
118  const scalar dx4dT = (dlogFcentdT - x4*dx3dT)/x3;
119 
120  return logTen*F*dx4dT;
121 }
122 
123 
124 inline Foam::scalar Foam::TroeFallOffFunction::ddPr
125 (
126  const scalar T,
127  const scalar Pr,
128  const scalar F
129 ) const
130 {
131  static const scalar logTen = log(scalar(10));
132 
133  const scalar logPr = log10(max(Pr, small));
134  const scalar dlogPrdPr = Pr >= small ? 1/(logTen*Pr) : 0;
135 
136  const scalar Fcent =
137  (1 - alpha_)*exp(-T/Tsss_)
138  + alpha_*exp(-T/Ts_)
139  + exp(-Tss_/T);
140 
141  const scalar logFcent = log10(max(Fcent, small));
142 
143  const scalar c = -0.4 - 0.67*logFcent;
144  static const scalar d = 0.14;
145  const scalar n = 0.75 - 1.27*logFcent;
146 
147  const scalar x1 = n - d*(logPr + c);
148  const scalar dx1dPr = -d*dlogPrdPr;
149  const scalar x2 = (logPr + c)/x1;
150  const scalar dx2dPr = (dlogPrdPr - x2*dx1dPr)/x1;
151  const scalar x3 = 1 + sqr(x2);
152  const scalar dx3dPr = 2*x2*dx2dPr;
153  const scalar x4 = logFcent/x3;
154  const scalar dx4dPr = -x4*dx3dPr/x3;
155 
156  return logTen*F*dx4dPr;
157 }
158 
159 
161 {
162  writeEntry(os, "alpha", alpha_);
163  writeEntry(os, "Tsss", Tsss_);
164  writeEntry(os, "Ts", Ts_);
165  writeEntry(os, "Tss", Tss_);
166 }
167 
168 
169 // * * * * * * * * * * * * * * * Ostream Operator * * * * * * * * * * * * * //
170 
171 inline Foam::Ostream& Foam::operator<<
172 (
173  Foam::Ostream& os,
174  const Foam::TroeFallOffFunction& tfof
175 )
176 {
177  tfof.write(os);
178  return os;
179 }
180 
181 
182 // ************************************************************************* //
dictionary dict
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
virtual Ostream & write(const char)=0
Write character.
scalar ddPr(const scalar T, const scalar Pr, const scalar F) const
dimensionedScalar log(const dimensionedScalar &ds)
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
dimensionedSymmTensor sqr(const dimensionedVector &dv)
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
The Troe fall-off function.
void write(Ostream &os) const
Write to stream.
const dimensionedScalar c
Speed of light in a vacuum.
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:54
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 ddT(const scalar T, const scalar Pr, const scalar F) const
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
label n
dimensionedScalar log10(const dimensionedScalar &ds)