Gulders.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-2020 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 "Gulders.H"
28 
29 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
30 
31 namespace Foam
32 {
33 namespace laminarFlameSpeedModels
34 {
35  defineTypeNameAndDebug(Gulders, 0);
36 
38  (
39  laminarFlameSpeed,
40  Gulders,
41  dictionary
42  );
43 }
44 }
45 
46 
47 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
48 
50 (
51  const dictionary& dict,
52  const psiuReactionThermo& ct
53 )
54 :
55  laminarFlameSpeed(dict, ct),
56 
57  coeffsDict_(dict.optionalSubDict(typeName + "Coeffs").subDict(fuel_)),
58  W_(coeffsDict_.lookup<scalar>("W")),
59  eta_(coeffsDict_.lookup<scalar>("eta")),
60  xi_(coeffsDict_.lookup<scalar>("xi")),
61  f_(coeffsDict_.lookup<scalar>("f")),
62  alpha_(coeffsDict_.lookup<scalar>("alpha")),
63  beta_(coeffsDict_.lookup<scalar>("beta"))
64 {}
65 
66 
67 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
68 
70 {}
71 
72 
73 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
74 
75 inline Foam::scalar Foam::laminarFlameSpeedModels::Gulders::SuRef
76 (
77  scalar phi
78 ) const
79 {
80  if (phi > small)
81  {
82  return W_*pow(phi, eta_)*exp(-xi_*sqr(phi - 1.075));
83  }
84  else
85  {
86  return 0.0;
87  }
88 }
89 
90 
91 inline Foam::scalar Foam::laminarFlameSpeedModels::Gulders::Su0pTphi
92 (
93  scalar p,
94  scalar Tu,
95  scalar phi,
96  scalar Yres
97 ) const
98 {
99  static const scalar Tref = 300.0;
100  static const scalar pRef = 1.013e5;
101 
102  return SuRef(phi)*pow((Tu/Tref), alpha_)*pow((p/pRef), beta_)*(1 - f_*Yres);
103 }
104 
105 
106 Foam::tmp<Foam::volScalarField> Foam::laminarFlameSpeedModels::Gulders::Su0pTphi
107 (
108  const volScalarField& p,
109  const volScalarField& Tu,
110  scalar phi
111 ) const
112 {
114  (
116  (
117  "Su0",
118  p.mesh(),
120  )
121  );
122 
123  volScalarField& Su0 = tSu0.ref();
124 
125  forAll(Su0, celli)
126  {
127  Su0[celli] = Su0pTphi(p[celli], Tu[celli], phi, 0.0);
128  }
129 
131 
132  forAll(Su0Bf, patchi)
133  {
134  forAll(Su0Bf[patchi], facei)
135  {
136  Su0Bf[patchi][facei] =
137  Su0pTphi
138  (
139  p.boundaryField()[patchi][facei],
140  Tu.boundaryField()[patchi][facei],
141  phi,
142  0.0
143  );
144  }
145  }
146 
147  return tSu0;
148 }
149 
150 
151 Foam::tmp<Foam::volScalarField> Foam::laminarFlameSpeedModels::Gulders::Su0pTphi
152 (
153  const volScalarField& p,
154  const volScalarField& Tu,
155  const volScalarField& phi
156 ) const
157 {
159  (
161  (
162  "Su0",
163  p.mesh(),
165  )
166  );
167 
168  volScalarField& Su0 = tSu0.ref();
169 
170  forAll(Su0, celli)
171  {
172  Su0[celli] = Su0pTphi(p[celli], Tu[celli], phi[celli], 0.0);
173  }
174 
176 
177  forAll(Su0Bf, patchi)
178  {
179  forAll(Su0Bf[patchi], facei)
180  {
181  Su0Bf[patchi][facei] =
182  Su0pTphi
183  (
184  p.boundaryField()[patchi][facei],
185  Tu.boundaryField()[patchi][facei],
186  phi.boundaryField()[patchi][facei],
187  0.0
188  );
189  }
190  }
191 
192  return tSu0;
193 }
194 
195 
198 {
199  if (psiuReactionThermo_.composition().contains("ft"))
200  {
201  const volScalarField& ft = psiuReactionThermo_.composition().Y("ft");
202 
203  return Su0pTphi
204  (
205  psiuReactionThermo_.p(),
206  psiuReactionThermo_.Tu(),
208  (
209  "stoichiometricAirFuelMassRatio",
210  dimless,
211  psiuReactionThermo_.properties()
212  )*ft/max(1 - ft, small)
213  );
214  }
215  else
216  {
217  return Su0pTphi
218  (
219  psiuReactionThermo_.p(),
220  psiuReactionThermo_.Tu(),
221  equivalenceRatio_
222  );
223  }
224 }
225 
226 
227 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
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:181
dimensionedSymmTensor sqr(const dimensionedVector &dv)
static tmp< GeometricField< scalar, fvPatchField, volMesh > > New(const word &name, const Internal &, const PtrList< fvPatchField< scalar >> &)
Return a temporary field constructed from name,.
const dimensionSet dimless
Macros for easy insertion into run-time selection tables.
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
Definition: dictionary.C:1002
dimensionedScalar exp(const dimensionedScalar &ds)
const dictionary & optionalSubDict(const word &) const
Find and return a sub-dictionary if found.
Definition: dictionary.C:1060
addToRunTimeSelectionTable(laminarFlameSpeed, constant, dictionary)
volScalarField::Internal & SuRef
Definition: alphaSuSp.H:34
phi
Definition: correctPhi.H:3
tmp< volScalarField > operator()() const
Return the laminar flame speed [m/s].
Definition: Gulders.C:197
Base-class for combustion fluid thermodynamic properties based on compressibility.
virtual ~Gulders()
Destructor.
Definition: Gulders.C:69
const dimensionSet dimVelocity
const Mesh & mesh() const
Return mesh.
scalar pRef
Definition: createFields.H:6
Gulders(const dictionary &, const psiuReactionThermo &)
Construct from dictionary and psiuReactionThermo.
Definition: Gulders.C:50
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.
A class for managing temporary objects.
Definition: PtrList.H:53
Namespace for OpenFOAM.