EDC.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) 2017-2023 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 "EDC.H"
28 
29 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
30 
31 template<>
33  {"v1981", "v1996", "v2005", "v2016"};
34 
37 
41 
42 namespace Foam
43 {
44 namespace combustionModels
45 {
48 }
49 }
50 
51 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
52 
54 (
55  const word& modelType,
58  const word& combustionProperties
59 )
60 :
61  combustionModel(modelType, thermo, turb, combustionProperties),
62  version_
63  (
65  [
66  this->coeffs().lookupOrDefault
67  (
68  "version",
70  )
71  ]
72  ),
73  C1_(this->coeffs().lookupOrDefault("C1", 0.05774)),
74  C2_(this->coeffs().lookupOrDefault("C2", 0.5)),
75  Cgamma_(this->coeffs().lookupOrDefault("Cgamma", 2.1377)),
76  Ctau_(this->coeffs().lookupOrDefault("Ctau", 0.4083)),
77  exp1_(this->coeffs().lookupOrDefault("exp1", EDCexp1[int(version_)])),
78  exp2_(this->coeffs().lookupOrDefault("exp2", EDCexp2[int(version_)])),
79  kappa_
80  (
81  IOobject
82  (
83  this->thermo().phasePropertyName(typedName("kappa")),
84  this->mesh().time().name(),
85  this->mesh(),
86  IOobject::NO_READ,
87  IOobject::AUTO_WRITE
88  ),
89  this->mesh(),
91  ),
92  outerCorrect_
93  (
94  this->coeffs().lookupOrDefault("outerCorrect", true)
95  ),
96  timeIndex_(-1),
97  chemistryPtr_(basicChemistryModel::New(thermo))
98 {}
99 
100 
101 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
102 
104 {}
105 
106 
107 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
108 
110 {
111  if (!outerCorrect_ && timeIndex_ == this->mesh().time().timeIndex())
112  {
113  return;
114  }
115 
116  tmp<volScalarField> tepsilon(this->turbulence().epsilon());
117  const volScalarField& epsilon = tepsilon();
118 
119  tmp<volScalarField> tnu(this->turbulence().nu());
120  const volScalarField& nu = tnu();
121 
122  tmp<volScalarField> tk(this->turbulence().k());
123  const volScalarField& k = tk();
124 
125  scalarField tauStar(epsilon.size(), 0);
126 
127  if (version_ == EDCversions::v2016)
128  {
129  tmp<volScalarField> ttc(chemistryPtr_->tc());
130  const volScalarField& tc = ttc();
131 
132  forAll(tauStar, i)
133  {
134  const scalar Da =
135  max(min(sqrt(nu[i]/(epsilon[i] + small))/tc[i], 10), 1e-10);
136 
137  const scalar ReT = sqr(k[i])/(nu[i]*epsilon[i] + small);
138  const scalar CtauI = min(C1_/(Da*sqrt(ReT + 1)), Ctau_);
139 
140  const scalar CgammaI =
141  max(min(C2_*sqrt(Da*(ReT + 1)), 5), Cgamma_);
142 
143  const scalar gammaL =
144  CgammaI*pow025(nu[i]*epsilon[i]/(sqr(k[i]) + small));
145 
146  tauStar[i] = CtauI*sqrt(nu[i]/(epsilon[i] + small));
147 
148  if (gammaL >= 1)
149  {
150  kappa_[i] = 1;
151  }
152  else
153  {
154  kappa_[i] =
155  max
156  (
157  min
158  (
159  pow(gammaL, exp1_)/(1 - pow(gammaL, exp2_)),
160  1
161  ),
162  0
163  );
164  }
165  }
166  }
167  else
168  {
169  forAll(tauStar, i)
170  {
171  const scalar gammaL =
172  Cgamma_*pow025(nu[i]*epsilon[i]/(sqr(k[i]) + small));
173 
174  tauStar[i] = Ctau_*sqrt(nu[i]/(epsilon[i] + small));
175  if (gammaL >= 1)
176  {
177  kappa_[i] = 1;
178  }
179  else
180  {
181  kappa_[i] =
182  max
183  (
184  min
185  (
186  pow(gammaL, exp1_)/(1 - pow(gammaL, exp2_)),
187  1
188  ),
189  0
190  );
191  }
192  }
193  }
194 
195  chemistryPtr_->solve(tauStar);
196 
197  timeIndex_ = this->mesh().time().timeIndex();
198 }
199 
200 
203 {
204  return kappa_*chemistryPtr_->RR()[speciei];
205 }
206 
207 
210 {
212  fvScalarMatrix& Su = tSu.ref();
213 
214  const label speciei = this->thermo().species()[Y.member()];
215  Su += chemistryPtr_->RR()[speciei];
216 
217  return kappa_*tSu;
218 }
219 
220 
223 {
224  return volScalarField::New
225  (
226  this->thermo().phasePropertyName(typedName("Qdot")),
227  kappa_*chemistryPtr_->Qdot()
228  );
229 }
230 
231 
233 {
234  if (combustionModel::read())
235  {
236  version_ =
237  (
239  [
240  this->coeffs().lookupOrDefault
241  (
242  "version",
244  )
245  ]
246  );
247  C1_ = this->coeffs().lookupOrDefault("C1", 0.05774);
248  C2_ = this->coeffs().lookupOrDefault("C2", 0.5);
249  Cgamma_ = this->coeffs().lookupOrDefault("Cgamma", 2.1377);
250  Ctau_ = this->coeffs().lookupOrDefault("Ctau", 0.4083);
251  exp1_ = this->coeffs().lookupOrDefault("exp1", EDCexp1[int(version_)]);
252  exp2_ = this->coeffs().lookupOrDefault("exp2", EDCexp2[int(version_)]);
253  outerCorrect_ = this->coeffs().lookupOrDefault("outerCorrect", true);
254 
255  return true;
256  }
257  else
258  {
259  return false;
260  }
261 }
262 
263 
264 // ************************************************************************* //
label k
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
Macros for easy insertion into run-time selection tables.
Generic GeometricField class.
static tmp< GeometricField< Type, PatchField, GeoMesh > > New(const word &name, const Internal &, const PtrList< PatchField< Type >> &, const HashPtrTable< Source > &=HashPtrTable< Source >())
Return a temporary field constructed from name,.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:99
static word member(const word &name)
Return member (name without the extension)
Definition: IOobject.C:149
Initialise the NamedEnum HashTable from the static list of names.
Definition: NamedEnum.H:54
Base class for chemistry models.
Base class for combustion models.
virtual bool read()
Update properties from given dictionary.
Eddy Dissipation Concept (EDC) turbulent combustion model.
Definition: EDC.H:135
virtual void correct()
Correct combustion rate.
Definition: EDC.C:109
EDC(const word &modelType, const fluidMulticomponentThermo &type, const compressibleMomentumTransportModel &turb, const word &combustionProperties)
Construct from components.
Definition: EDC.C:54
virtual ~EDC()
Destructor.
Definition: EDC.C:103
virtual tmp< volScalarField::Internal > R(const label speciei) const
Specie consumption rate field.
Definition: EDC.C:202
virtual tmp< volScalarField > Qdot() const
Heat release rate [kg/m/s^3].
Definition: EDC.C:222
virtual bool read()
Update properties from given dictionary.
Definition: EDC.C:232
Base class for single-phase compressible turbulence models.
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
Base-class for multi-component fluid thermodynamic properties.
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Definition: fvMatrix.H:118
A class for managing temporary objects.
Definition: tmp.H:55
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:181
A class for handling words, derived from string.
Definition: word.H:62
const scalar epsilon
defineTypeNameAndDebug(diffusion, 0)
const scalar EDCexp2[]
Definition: EDC.H:126
const scalar EDCexp1[]
Definition: EDC.H:125
const NamedEnum< EDCversions, 4 > EDCversionNames
Definition: EDC.C:36
const EDCversions EDCdefaultVersion
Definition: EDC.C:39
EDCversions
EDC model versions.
Definition: EDC.H:115
addToRunTimeSelectionTable(combustionModel, diffusion, dictionary)
autoPtr< CompressibleMomentumTransportModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const viscosity &viscosity)
tmp< VolField< Type > > Su(const VolField< Type > &su, const VolField< Type > &vf)
Definition: fvcSup.C:44
Namespace for OpenFOAM.
const doubleScalar e
Definition: doubleScalar.H:106
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
dimensionedSymmTensor sqr(const dimensionedVector &dv)
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
fvMatrix< scalar > fvScalarMatrix
Definition: fvMatricesFwd.H:42
const dimensionSet dimless
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
const dimensionSet dimTime
dimensionedScalar sqrt(const dimensionedScalar &ds)
word typedName(Name name)
Return the name of the object within the given type.
Definition: typeInfo.H:176
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
const dimensionSet dimMass
dimensionedScalar pow025(const dimensionedScalar &ds)
label timeIndex
Definition: getTimeIndex.H:4
autoPtr< incompressible::momentumTransportModel > turbulence(incompressible::momentumTransportModel::New(U, phi, viscosity))
PtrList< volScalarField > & Y
fluidMulticomponentThermo & thermo
Definition: createFields.H:31