SCOPELaminarFlameSpeed.C
Go to the documentation of this file.
1/*---------------------------------------------------------------------------*\
2 ========= |
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4 \\ / O peration |
5 \\ / A nd | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2011-2017 OpenFOAM Foundation
9-------------------------------------------------------------------------------
10License
11 This file is part of OpenFOAM.
12
13 OpenFOAM is free software: you can redistribute it and/or modify it
14 under the terms of the GNU General Public License as published by
15 the Free Software Foundation, either version 3 of the License, or
16 (at your option) any later version.
17
18 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21 for more details.
22
23 You should have received a copy of the GNU General Public License
24 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25
26\*---------------------------------------------------------------------------*/
27
28#include "IFstream.H"
31
32// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
33
34namespace Foam
35{
36namespace laminarFlameSpeedModels
37{
38 defineTypeNameAndDebug(SCOPE, 0);
39
41 (
42 laminarFlameSpeed,
43 SCOPE,
44 dictionary
45 );
46}
47}
48
49
50// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
51
53(
54 const dictionary& polyDict
55)
56:
57 FixedList<scalar, 7>(polyDict.lookup("coefficients")),
58 ll(polyDict.get<scalar>("lowerLimit")),
59 ul(polyDict.get<scalar>("upperLimit")),
60 llv(polyPhi(ll, *this)),
61 ulv(polyPhi(ul, *this)),
62 lu(0)
63{}
64
65
67(
68 const dictionary& dict,
69 const psiuReactionThermo& ct
70)
71:
72 laminarFlameSpeed(dict, ct),
73
74 coeffsDict_
75 (
76 dictionary
77 (
78 IFstream
79 (
80 dict.get<fileName>("fuelFile")
81 )()
82 ).optionalSubDict(typeName + "Coeffs")
83 ),
84 LFL_
85 (
86 coeffsDict_.getCompat<scalar>
87 (
88 "lowerFlammabilityLimit",
89 {{"lowerFlamabilityLimit", 1712}}
90 )
91 ),
92 UFL_
93 (
94 coeffsDict_.getCompat<scalar>
95 (
96 "upperFlammabilityLimit",
97 {{"upperFlamabilityLimit", 1712}}
98 )
99 ),
100 SuPolyL_(coeffsDict_.subDict("lowerSuPolynomial")),
101 SuPolyU_(coeffsDict_.subDict("upperSuPolynomial")),
102 Texp_(coeffsDict_.get<scalar>("Texp")),
103 pexp_(coeffsDict_.get<scalar>("pexp")),
104 MaPolyL_(coeffsDict_.subDict("lowerMaPolynomial")),
105 MaPolyU_(coeffsDict_.subDict("upperMaPolynomial"))
106{
107 SuPolyL_.ll = max(SuPolyL_.ll, LFL_) + SMALL;
108 SuPolyU_.ul = min(SuPolyU_.ul, UFL_) - SMALL;
109
110 SuPolyL_.lu = 0.5*(SuPolyL_.ul + SuPolyU_.ll);
111 SuPolyU_.lu = SuPolyL_.lu - SMALL;
112
113 MaPolyL_.lu = 0.5*(MaPolyL_.ul + MaPolyU_.ll);
114 MaPolyU_.lu = MaPolyL_.lu - SMALL;
115
116 if (debug)
117 {
118 Info<< "phi Su (T = Tref, p = pref)" << endl;
119 label n = 200;
120 for (int i=0; i<n; i++)
121 {
122 scalar phi = (2.0*i)/n;
123 Info<< phi << token::TAB << SuRef(phi) << endl;
124 }
125 }
126}
127
128
129// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
130
132{}
133
134
135// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
136
137inline Foam::scalar Foam::laminarFlameSpeedModels::SCOPE::polyPhi
138(
139 scalar phi,
140 const polynomial& a
141)
142{
143 scalar x = phi - 1.0;
144
145 return
146 a[0]
147 *(
148 scalar(1)
149 + x*(a[1] + x*(a[2] + x*(a[3] + x*(a[4] + x*(a[5] + x*a[6])))))
150 );
151}
152
153
154inline Foam::scalar Foam::laminarFlameSpeedModels::SCOPE::SuRef
155(
156 scalar phi
157) const
158{
159 if (phi < LFL_ || phi > UFL_)
160 {
161 // Return 0 beyond the flammability limits
162 return scalar(0);
163 }
164 else if (phi < SuPolyL_.ll)
165 {
166 // Use linear interpolation between the low end of the
167 // lower polynomial and the lower flammability limit
168 return SuPolyL_.llv*(phi - LFL_)/(SuPolyL_.ll - LFL_);
169 }
170 else if (phi > SuPolyU_.ul)
171 {
172 // Use linear interpolation between the upper end of the
173 // upper polynomial and the upper flammability limit
174 return SuPolyU_.ulv*(UFL_ - phi)/(UFL_ - SuPolyU_.ul);
175 }
176 else if (phi < SuPolyL_.lu)
177 {
178 // Evaluate the lower polynomial
179 return polyPhi(phi, SuPolyL_);
180 }
181 else if (phi > SuPolyU_.lu)
182 {
183 // Evaluate the upper polynomial
184 return polyPhi(phi, SuPolyU_);
185 }
186 else
187 {
189 << "phi = " << phi
190 << " cannot be handled by SCOPE function with the "
191 "given coefficients"
192 << exit(FatalError);
193
194 return scalar(0);
195 }
196}
197
198
200(
201 scalar phi
202) const
203{
204 if (phi < MaPolyL_.ll)
205 {
206 // Beyond the lower limit assume Ma is constant
207 return MaPolyL_.llv;
208 }
209 else if (phi > MaPolyU_.ul)
210 {
211 // Beyond the upper limit assume Ma is constant
212 return MaPolyU_.ulv;
213 }
214 else if (phi < SuPolyL_.lu)
215 {
216 // Evaluate the lower polynomial
217 return polyPhi(phi, MaPolyL_);
218 }
219 else if (phi > SuPolyU_.lu)
220 {
221 // Evaluate the upper polynomial
222 return polyPhi(phi, MaPolyU_);
223 }
224 else
225 {
227 << "phi = " << phi
228 << " cannot be handled by SCOPE function with the "
229 "given coefficients"
230 << exit(FatalError);
231
232 return scalar(0);
233 }
234}
235
236
237inline Foam::scalar Foam::laminarFlameSpeedModels::SCOPE::Su0pTphi
238(
239 scalar p,
240 scalar Tu,
241 scalar phi
242) const
243{
244 static const scalar Tref = 300.0;
245 static const scalar pRef = 1.013e5;
246
247 return SuRef(phi)*pow((Tu/Tref), Texp_)*pow((p/pRef), pexp_);
248}
249
250
251Foam::tmp<Foam::volScalarField> Foam::laminarFlameSpeedModels::SCOPE::Su0pTphi
252(
253 const volScalarField& p,
254 const volScalarField& Tu,
255 scalar phi
256) const
257{
258 tmp<volScalarField> tSu0
259 (
261 (
262 IOobject
263 (
264 "Su0",
265 p.time().timeName(),
266 p.db(),
269 ),
270 p.mesh(),
272 )
273 );
274
275 volScalarField& Su0 = tSu0.ref();
276
277 forAll(Su0, celli)
278 {
279 Su0[celli] = Su0pTphi(p[celli], Tu[celli], phi);
280 }
281
282 volScalarField::Boundary& Su0Bf = Su0.boundaryFieldRef();
283
284 forAll(Su0Bf, patchi)
285 {
286 scalarField& Su0p = Su0Bf[patchi];
287 const scalarField& pp = p.boundaryField()[patchi];
288 const scalarField& Tup = Tu.boundaryField()[patchi];
289
290 forAll(Su0p, facei)
291 {
292 Su0p[facei] = Su0pTphi(pp[facei], Tup[facei], phi);
293 }
294 }
295
296 return tSu0;
297}
298
299
300Foam::tmp<Foam::volScalarField> Foam::laminarFlameSpeedModels::SCOPE::Su0pTphi
301(
302 const volScalarField& p,
303 const volScalarField& Tu,
304 const volScalarField& phi
305) const
306{
307 tmp<volScalarField> tSu0
308 (
310 (
311 IOobject
312 (
313 "Su0",
314 p.time().timeName(),
315 p.db(),
318 ),
319 p.mesh(),
321 )
322 );
323
324 volScalarField& Su0 = tSu0.ref();
325
326 forAll(Su0, celli)
327 {
328 Su0[celli] = Su0pTphi(p[celli], Tu[celli], phi[celli]);
329 }
330
331 volScalarField::Boundary& Su0Bf = Su0.boundaryFieldRef();
332
333 forAll(Su0Bf, patchi)
334 {
335 scalarField& Su0p = Su0Bf[patchi];
336 const scalarField& pp = p.boundaryField()[patchi];
337 const scalarField& Tup = Tu.boundaryField()[patchi];
338 const scalarField& phip = phi.boundaryField()[patchi];
339
340 forAll(Su0p, facei)
341 {
342 Su0p[facei] =
343 Su0pTphi
344 (
345 pp[facei],
346 Tup[facei],
347 phip[facei]
348 );
349 }
350 }
351
352 return tSu0;
353}
354
355
357(
358 const volScalarField& phi
359) const
360{
361 tmp<volScalarField> tMa
362 (
364 (
365 IOobject
366 (
367 "Ma",
368 phi.time().timeName(),
369 phi.db(),
372 ),
373 phi.mesh(),
375 )
376 );
377
378 volScalarField& ma = tMa.ref();
379
380 forAll(ma, celli)
381 {
382 ma[celli] = Ma(phi[celli]);
383 }
384
385 volScalarField::Boundary& maBf = ma.boundaryFieldRef();
386
387 forAll(maBf, patchi)
388 {
389 scalarField& map = maBf[patchi];
390 const scalarField& phip = phi.boundaryField()[patchi];
391
392 forAll(map, facei)
393 {
394 map[facei] = Ma(phip[facei]);
395 }
396 }
397
398 return tMa;
399}
400
401
404{
405 if (psiuReactionThermo_.composition().contains("ft"))
406 {
407 const volScalarField& ft = psiuReactionThermo_.composition().Y("ft");
408
409 return Ma
410 (
412 (
413 "stoichiometricAirFuelMassRatio", dimless, psiuReactionThermo_
414 )*ft/(scalar(1) - ft)
415 );
416 }
417 else
418 {
419 const fvMesh& mesh = psiuReactionThermo_.p().mesh();
420
421 return tmp<volScalarField>
422 (
424 (
425 IOobject
426 (
427 "Ma",
428 mesh.time().timeName(),
429 mesh,
432 ),
433 mesh,
434 dimensionedScalar("Ma", dimless, Ma(equivalenceRatio_))
435 )
436 );
437 }
438}
439
440
443{
444 if (psiuReactionThermo_.composition().contains("ft"))
445 {
446 const volScalarField& ft = psiuReactionThermo_.composition().Y("ft");
447
448 return Su0pTphi
449 (
450 psiuReactionThermo_.p(),
451 psiuReactionThermo_.Tu(),
453 (
454 "stoichiometricAirFuelMassRatio", dimless, psiuReactionThermo_
455 )*ft/(scalar(1) - ft)
456 );
457 }
458 else
459 {
460 return Su0pTphi
461 (
462 psiuReactionThermo_.p(),
463 psiuReactionThermo_.Tu(),
464 equivalenceRatio_
465 );
466 }
467}
468
469
470// ************************************************************************* //
label n
Y[inertIndex] max(0.0)
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
surfaceScalarField & phi
Internal & ref(const bool updateAccessTime=true)
Return a reference to the dimensioned internal field.
GeometricBoundaryField< scalar, fvPatchField, volMesh > Boundary
Type of boundary fields.
Laminar flame speed obtained from the SCOPE correlation.
tmp< volScalarField > operator()() const
Return the laminar flame speed [m/s].
tmp< volScalarField > Ma() const
Return the Markstein number.
Polynomial equation for the saturation vapour temperature in terms of the vapour pressure (in Pa).
Definition: polynomial.H:66
A class for managing temporary objects.
Definition: tmp.H:65
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition: className.H:121
volScalarField & p
dynamicFvMesh & mesh
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
List< ReturnType > get(const UPtrList< T > &list, const AccessOp &aop)
Namespace for OpenFOAM.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
const dimensionSet dimless
Dimensionless.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:82
const dimensionSet dimVelocity
messageStream Info
Information stream (stdout output on master, null elsewhere)
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:33
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
error FatalError
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
dictionary dict
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333