All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
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-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 "EDC.H"
27 
28 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
29 
30 template<class ReactionThermo>
32 (
33  const word& modelType,
34  const ReactionThermo& thermo,
36  const word& combustionProperties
37 )
38 :
39  laminar<ReactionThermo>(modelType, thermo, turb, combustionProperties),
40  version_
41  (
43  [
44  this->coeffs().lookupOrDefault
45  (
46  "version",
47  word(EDCversionNames[EDCdefaultVersion])
48  )
49  ]
50  ),
51  C1_(this->coeffs().lookupOrDefault("C1", 0.05774)),
52  C2_(this->coeffs().lookupOrDefault("C2", 0.5)),
53  Cgamma_(this->coeffs().lookupOrDefault("Cgamma", 2.1377)),
54  Ctau_(this->coeffs().lookupOrDefault("Ctau", 0.4083)),
55  exp1_(this->coeffs().lookupOrDefault("exp1", EDCexp1[int(version_)])),
56  exp2_(this->coeffs().lookupOrDefault("exp2", EDCexp2[int(version_)])),
57  kappa_
58  (
59  IOobject
60  (
61  this->thermo().phasePropertyName(typeName + ":kappa"),
62  this->mesh().time().timeName(),
63  this->mesh(),
64  IOobject::NO_READ,
65  IOobject::AUTO_WRITE
66  ),
67  this->mesh(),
69  )
70 {}
71 
72 
73 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
74 
75 template<class ReactionThermo>
77 {}
78 
79 
80 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
81 
82 template<class ReactionThermo>
84 {
85  tmp<volScalarField> tepsilon(this->turbulence().epsilon());
86  const volScalarField& epsilon = tepsilon();
87 
88  tmp<volScalarField> tmu(this->turbulence().mu());
89  const volScalarField& mu = tmu();
90 
91  tmp<volScalarField> tk(this->turbulence().k());
92  const volScalarField& k = tk();
93 
94  tmp<volScalarField> trho(this->rho());
95  const volScalarField& rho = trho();
96 
97  scalarField tauStar(epsilon.size(), 0);
98 
99  if (version_ == EDCversions::v2016)
100  {
101  tmp<volScalarField> ttc(this->chemistryPtr_->tc());
102  const volScalarField& tc = ttc();
103 
104  forAll(tauStar, i)
105  {
106  const scalar nu = mu[i]/(rho[i] + small);
107 
108  const scalar Da =
109  max(min(sqrt(nu/(epsilon[i] + small))/tc[i], 10), 1e-10);
110 
111  const scalar ReT = sqr(k[i])/(nu*epsilon[i] + small);
112  const scalar CtauI = min(C1_/(Da*sqrt(ReT + 1)), 2.1377);
113 
114  const scalar CgammaI =
115  max(min(C2_*sqrt(Da*(ReT + 1)), 5), 0.4082);
116 
117  const scalar gammaL =
118  CgammaI*pow025(nu*epsilon[i]/(sqr(k[i]) + small));
119 
120  tauStar[i] = CtauI*sqrt(nu/(epsilon[i] + small));
121 
122  if (gammaL >= 1)
123  {
124  kappa_[i] = 1;
125  }
126  else
127  {
128  kappa_[i] =
129  max
130  (
131  min
132  (
133  pow(gammaL, exp1_)/(1 - pow(gammaL, exp2_)),
134  1
135  ),
136  0
137  );
138  }
139  }
140  }
141  else
142  {
143  forAll(tauStar, i)
144  {
145  const scalar nu = mu[i]/(rho[i] + small);
146  const scalar gammaL =
147  Cgamma_*pow025(nu*epsilon[i]/(sqr(k[i]) + small));
148 
149  tauStar[i] = Ctau_*sqrt(nu/(epsilon[i] + small));
150  if (gammaL >= 1)
151  {
152  kappa_[i] = 1;
153  }
154  else
155  {
156  kappa_[i] =
157  max
158  (
159  min
160  (
161  pow(gammaL, exp1_)/(1 - pow(gammaL, exp2_)),
162  1
163  ),
164  0
165  );
166  }
167  }
168  }
169 
170  this->chemistryPtr_->solve(tauStar);
171 }
172 
173 
174 template<class ReactionThermo>
177 {
178  return kappa_*laminar<ReactionThermo>::R(Y);
179 }
180 
181 
182 template<class ReactionThermo>
185 {
186  return volScalarField::New
187  (
188  this->thermo().phasePropertyName(typeName + ":Qdot"),
189  kappa_*this->chemistryPtr_->Qdot()
190  );
191 }
192 
193 
194 template<class ReactionThermo>
196 {
198  {
199  version_ =
200  (
202  [
203  this->coeffs().lookupOrDefault
204  (
205  "version",
206  word(EDCversionNames[EDCdefaultVersion])
207  )
208  ]
209  );
210  C1_ = this->coeffs().lookupOrDefault("C1", 0.05774);
211  C2_ = this->coeffs().lookupOrDefault("C2", 0.5);
212  Cgamma_ = this->coeffs().lookupOrDefault("Cgamma", 2.1377);
213  Ctau_ = this->coeffs().lookupOrDefault("Ctau", 0.4083);
214  exp1_ = this->coeffs().lookupOrDefault("exp1", EDCexp1[int(version_)]);
215  exp2_ = this->coeffs().lookupOrDefault("exp2", EDCexp2[int(version_)]);
216 
217  return true;
218  }
219  else
220  {
221  return false;
222  }
223 }
224 
225 
226 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
virtual tmp< fvScalarMatrix > R(volScalarField &Y) const
Fuel consumption rate matrix.
Definition: EDC.C:176
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
virtual bool read()
Update properties from given dictionary.
Definition: EDC.C:195
dimensionedSymmTensor sqr(const dimensionedVector &dv)
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensionedScalar pow025(const dimensionedScalar &ds)
tmp< volScalarField > trho
label k
Boltzmann constant.
rhoReactionThermo & thermo
Definition: createFields.H:28
const scalar EDCexp1[]
Definition: EDC.H:123
Laminar combustion model.
Definition: laminar.H:51
dynamicFvMesh & mesh
const dimensionedScalar & e
Elementary charge.
Definition: doubleScalar.H:105
autoPtr< BasicCompressibleMomentumTransportModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const typename BasicCompressibleMomentumTransportModel::transportModel &transport)
A class for handling words, derived from string.
Definition: word.H:59
word timeName
Definition: getTimeIndex.H:3
Info<< "Reading field U\"<< endl;volVectorField U(IOobject("U", runTime.timeName(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE), mesh);volScalarField rho(IOobject("rho", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::AUTO_WRITE), thermo.rho());volVectorField rhoU(IOobject("rhoU", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::NO_WRITE), rho *U);volScalarField rhoE(IOobject("rhoE", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::NO_WRITE), rho *(e+0.5 *magSqr(U)));surfaceScalarField pos(IOobject("pos", runTime.timeName(), mesh), mesh, dimensionedScalar(dimless, 1.0));surfaceScalarField neg(IOobject("neg", runTime.timeName(), mesh), mesh, dimensionedScalar(dimless, -1.0));surfaceScalarField phi("phi", fvc::flux(rhoU));Info<< "Creating turbulence model\"<< endl;autoPtr< compressible::momentumTransportModel > turbulence(compressible::momentumTransportModel::New(rho, U, phi, thermo))
Definition: createFields.H:94
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
const NamedEnum< EDCversions, 4 > EDCversionNames
Definition: EDCs.C:49
const scalar EDCexp2[]
Definition: EDC.H:124
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
EDC(const word &modelType, const ReactionThermo &type, const compressibleMomentumTransportModel &turb, const word &combustionProperties)
Construct from components.
Definition: EDC.C:32
const dimensionedScalar & mu
Atomic mass unit.
#define R(A, B, C, D, E, F, K, M)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
PtrList< volScalarField > & Y
scalar epsilon
virtual tmp< volScalarField > Qdot() const
Heat release rate [kg/m/s^3].
Definition: EDC.C:184
virtual ~EDC()
Destructor.
Definition: EDC.C:76
virtual void correct()
Correct combustion rate.
Definition: EDC.C:83
A class for managing temporary objects.
Definition: PtrList.H:53
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
Abstract base class for turbulence models (RAS, LES and laminar).