GuldersEGR.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-2018 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 "GuldersEGR.H"
28 
29 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
30 
31 namespace Foam
32 {
33 namespace laminarFlameSpeedModels
34 {
35  defineTypeNameAndDebug(GuldersEGR, 0);
36 
38  (
39  laminarFlameSpeed,
40  GuldersEGR,
41  dictionary
42  );
43 }
44 }
45 
46 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
47 
48 Foam::laminarFlameSpeedModels::GuldersEGR::GuldersEGR
49 (
50  const dictionary& dict,
51  const psiuReactionThermo& ct
52 )
53 :
54  laminarFlameSpeed(dict, ct),
55 
56  coeffsDict_(dict.optionalSubDict(typeName + "Coeffs").subDict(fuel_)),
57  W_(readScalar(coeffsDict_.lookup("W"))),
58  eta_(readScalar(coeffsDict_.lookup("eta"))),
59  xi_(readScalar(coeffsDict_.lookup("xi"))),
60  f_(readScalar(coeffsDict_.lookup("f"))),
61  alpha_(readScalar(coeffsDict_.lookup("alpha"))),
62  beta_(readScalar(coeffsDict_.lookup("beta")))
63 {}
64 
65 
66 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
67 
69 {}
70 
71 
72 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
73 
74 inline Foam::scalar Foam::laminarFlameSpeedModels::GuldersEGR::SuRef
75 (
76  scalar phi
77 ) const
78 {
79  if (phi > small)
80  {
81  return W_*pow(phi, eta_)*exp(-xi_*sqr(phi - 1.075));
82  }
83  else
84  {
85  return 0.0;
86  }
87 }
88 
89 
90 inline Foam::scalar Foam::laminarFlameSpeedModels::GuldersEGR::Su0pTphi
91 (
92  scalar p,
93  scalar Tu,
94  scalar phi,
95  scalar Yres
96 ) const
97 {
98  static const scalar Tref = 300.0;
99  static const scalar pRef = 1.013e5;
100 
101  return SuRef(phi)*pow((Tu/Tref), alpha_)*pow((p/pRef), beta_)*(1 - f_*Yres);
102 }
103 
104 
106 Foam::laminarFlameSpeedModels::GuldersEGR::Su0pTphi
107 (
108  const volScalarField& p,
109  const volScalarField& Tu,
110  scalar phi
111 ) const
112 {
114  (
115  new volScalarField
116  (
117  IOobject
118  (
119  "Su0",
120  p.time().timeName(),
121  p.db(),
124  false
125  ),
126  p.mesh(),
127  dimensionedScalar("Su0", dimVelocity, 0.0)
128  )
129  );
130 
131  volScalarField& Su0 = tSu0.ref();
132 
133  forAll(Su0, celli)
134  {
135  Su0[celli] = Su0pTphi(p[celli], Tu[celli], phi, 0.0);
136  }
137 
138  volScalarField::Boundary& Su0Bf = Su0.boundaryFieldRef();
139 
140  forAll(Su0Bf, patchi)
141  {
142  forAll(Su0Bf[patchi], facei)
143  {
144  Su0Bf[patchi][facei] =
145  Su0pTphi
146  (
147  p.boundaryField()[patchi][facei],
148  Tu.boundaryField()[patchi][facei],
149  phi,
150  0.0
151  );
152  }
153  }
154 
155  return tSu0;
156 }
157 
158 
160 Foam::laminarFlameSpeedModels::GuldersEGR::Su0pTphi
161 (
162  const volScalarField& p,
163  const volScalarField& Tu,
164  const volScalarField& phi,
165  const volScalarField& egr
166 ) const
167 {
169  (
170  new volScalarField
171  (
172  IOobject
173  (
174  "Su0",
175  p.time().timeName(),
176  p.db(),
179  false
180  ),
181  p.mesh(),
182  dimensionedScalar("Su0", dimVelocity, 0.0)
183  )
184  );
185 
186  volScalarField& Su0 = tSu0.ref();
187 
188  forAll(Su0, celli)
189  {
190  Su0[celli] = Su0pTphi(p[celli], Tu[celli], phi[celli], egr[celli]);
191  }
192 
193  volScalarField::Boundary& Su0Bf = Su0.boundaryFieldRef();
194 
195  forAll(Su0Bf, patchi)
196  {
197  forAll(Su0Bf[patchi], facei)
198  {
199  Su0Bf[patchi][facei] =
200  Su0pTphi
201  (
202  p.boundaryField()[patchi][facei],
203  Tu.boundaryField()[patchi][facei],
204  phi.boundaryField()[patchi][facei],
205  egr.boundaryField()[patchi][facei]
206  );
207  }
208  }
209 
210  return tSu0;
211 }
212 
213 
216 {
217  if
218  (
219  psiuReactionThermo_.composition().contains("ft")
220  && psiuReactionThermo_.composition().contains("egr")
221  )
222  {
223  return Su0pTphi
224  (
225  psiuReactionThermo_.p(),
226  psiuReactionThermo_.Tu(),
228  (
229  psiuReactionThermo_.lookup("stoichiometricAirFuelMassRatio")
230  )/
231  (
232  scalar(1)/psiuReactionThermo_.composition().Y("ft")
233  - scalar(1)
234  ),
235  psiuReactionThermo_.composition().Y("egr")
236  );
237  }
238  else
239  {
240  return Su0pTphi
241  (
242  psiuReactionThermo_.p(),
243  psiuReactionThermo_.Tu(),
244  equivalenceRatio_
245  );
246  }
247 }
248 
249 
250 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
tmp< volScalarField > operator()() const
Return the laminar flame speed [m/s].
Definition: GuldersEGR.C:215
surfaceScalarField & phi
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
const Boundary & boundaryField() const
Return const-reference to the boundary field.
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:174
dimensionedSymmTensor sqr(const dimensionedVector &dv)
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:626
Macros for easy insertion into run-time selection tables.
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
Definition: dictionary.C:692
dimensionedScalar exp(const dimensionedScalar &ds)
const dictionary & optionalSubDict(const word &) const
Find and return a sub-dictionary if found.
Definition: dictionary.C:759
addToRunTimeSelectionTable(laminarFlameSpeed, constant, dictionary)
bool readScalar(const char *buf, doubleScalar &s)
Read whole of buf as a scalar. Return true if successful.
Definition: doubleScalar.H:68
Foam::psiuReactionThermo.
const Mesh & mesh() const
Return mesh.
scalar pRef
Definition: createFields.H:19
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
label patchi
Abstract class for laminar flame speed.
Boundary & boundaryFieldRef()
Return a reference to the boundary field.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
const Time & time() const
Return time.
Definition: IOobject.C:367
A class for managing temporary objects.
Definition: PtrList.H:53
const objectRegistry & db() const
Return the local objectRegistry.
Definition: IOobject.C:361
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
Namespace for OpenFOAM.
const dimensionSet dimVelocity