EDC.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) 2017 OpenFOAM Foundation
9 Copyright (C) 2019-2020 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27\*---------------------------------------------------------------------------*/
28
29#include "EDC.H"
30
31// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
32
33template<class ReactionThermo>
35(
36 const word& modelType,
37 ReactionThermo& thermo,
39 const word& combustionProperties
40)
41:
42 laminar<ReactionThermo>(modelType, thermo, turb, combustionProperties),
43 version_
44 (
45 EDCversionNames.getOrDefault
46 (
47 "version",
48 this->coeffs(),
50 )
51 ),
52 C1_(this->coeffs().getOrDefault("C1", 0.05774)),
53 C2_(this->coeffs().getOrDefault("C2", 0.5)),
54 Cgamma_(this->coeffs().getOrDefault("Cgamma", 2.1377)),
55 Ctau_(this->coeffs().getOrDefault("Ctau", 0.4083)),
56 exp1_(this->coeffs().getOrDefault("exp1", EDCexp1[int(version_)])),
57 exp2_(this->coeffs().getOrDefault("exp2", EDCexp2[int(version_)])),
58 kappa_
59 (
61 (
62 this->thermo().phasePropertyName(typeName + ":kappa"),
63 this->mesh().time().timeName(),
64 this->mesh(),
65 IOobject::NO_READ,
66 IOobject::AUTO_WRITE
67 ),
68 this->mesh(),
70 )
71{}
72
73
74// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
75
76template<class ReactionThermo>
78{}
79
80
81// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
82
83template<class ReactionThermo>
85{
86 if (this->active())
87 {
88 tmp<volScalarField> tepsilon(this->turbulence().epsilon());
89 const volScalarField& epsilon = tepsilon();
90
91 tmp<volScalarField> tmu(this->turbulence().mu());
92 const volScalarField& mu = tmu();
93
94 tmp<volScalarField> tk(this->turbulence().k());
95 const volScalarField& k = tk();
96
98 const volScalarField& rho = trho();
99
100 scalarField tauStar(epsilon.size(), Zero);
101
102 if (version_ == EDCversions::v2016)
103 {
104 tmp<volScalarField> ttc(this->chemistryPtr_->tc());
105 const volScalarField& tc = ttc();
106
107 forAll(tauStar, i)
108 {
109 const scalar nu = mu[i]/(rho[i] + SMALL);
110
111 const scalar Da =
112 max(min(sqrt(nu/(epsilon[i] + SMALL))/tc[i], 10), 1e-10);
113
114 const scalar ReT = sqr(k[i])/(nu*epsilon[i] + SMALL);
115 const scalar CtauI = min(C1_/(Da*sqrt(ReT + 1)), 2.1377);
116
117 const scalar CgammaI =
118 max(min(C2_*sqrt(Da*(ReT + 1)), 5), 0.4082);
119
120 const scalar gammaL =
121 CgammaI*pow025(nu*epsilon[i]/(sqr(k[i]) + SMALL));
122
123 tauStar[i] = CtauI*sqrt(nu/(epsilon[i] + SMALL));
124
125 if (gammaL >= 1)
126 {
127 kappa_[i] = 1;
128 }
129 else
130 {
131 kappa_[i] =
132 max
133 (
134 min
135 (
136 pow(gammaL, exp1_)/(1 - pow(gammaL, exp2_)),
137 1
138 ),
139 0
140 );
141 }
142 }
143 }
144 else
145 {
146 forAll(tauStar, i)
147 {
148 const scalar nu = mu[i]/(rho[i] + SMALL);
149 const scalar gammaL =
150 Cgamma_*pow025(nu*epsilon[i]/(sqr(k[i]) + SMALL));
151
152 tauStar[i] = Ctau_*sqrt(nu/(epsilon[i] + SMALL));
153 if (gammaL >= 1)
154 {
155 kappa_[i] = 1;
156 }
157 else
158 {
159 kappa_[i] =
160 max
161 (
162 min
163 (
164 pow(gammaL, exp1_)/(1 - pow(gammaL, exp2_)),
165 1
166 ),
167 0
168 );
169 }
170 }
171 }
172
173 Info<< "Chemistry time solved max/min : "
174 << gMax(tauStar) << " / " << gMin(tauStar) << endl;
175
176 this->chemistryPtr_->solve(tauStar);
177 }
178}
179
180
181template<class ReactionThermo>
184{
185 return kappa_*laminar<ReactionThermo>::R(Y);
186}
187
188
189template<class ReactionThermo>
192{
194 (
196 (
198 (
199 this->thermo().phasePropertyName(typeName + ":Qdot"),
200 this->mesh().time().timeName(),
201 this->mesh(),
204 false
205 ),
206 this->mesh(),
208 )
209 );
210
211 if (this->active())
212 {
213 tQdot.ref() = kappa_*this->chemistryPtr_->Qdot();
214 }
215
216 return tQdot;
217}
218
219
220template<class ReactionThermo>
222{
224 {
225 version_ =
226 (
227 EDCversionNames.getOrDefault
228 (
229 "version",
230 this->coeffs(),
232 )
233 );
234 C1_ = this->coeffs().getOrDefault("C1", 0.05774);
235 C2_ = this->coeffs().getOrDefault("C2", 0.5);
236 Cgamma_ = this->coeffs().getOrDefault("Cgamma", 2.1377);
237 Ctau_ = this->coeffs().getOrDefault("Ctau", 0.4083);
238 exp1_ = this->coeffs().getOrDefault("exp1", EDCexp1[int(version_)]);
239 exp2_ = this->coeffs().getOrDefault("exp2", EDCexp2[int(version_)]);
240
241 return true;
242 }
243
244 return false;
245}
246
247
248// ************************************************************************* //
label k
#define R(A, B, C, D, E, F, K, M)
compressible::turbulenceModel & turb
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:170
Eddy Dissipation Concept (EDC) turbulent combustion model.
Definition: EDC.H:136
virtual void correct()
Correct combustion rate.
Definition: EDC.C:84
virtual ~EDC()
Destructor.
Definition: EDC.C:77
virtual tmp< volScalarField > Qdot() const
Heat release rate [kg/m/s3].
Definition: EDC.C:191
virtual bool read()
Update properties from given dictionary.
Definition: EDC.C:221
Laminar combustion model.
Definition: laminar.H:59
Abstract base class for turbulence models (RAS, LES and laminar).
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
A class for managing temporary objects.
Definition: tmp.H:65
T & ref() const
Definition: tmpI.H:227
A class for handling words, derived from Foam::string.
Definition: word.H:68
PtrList< volScalarField > & Y
const volScalarField & mu
dynamicFvMesh & mesh
scalar epsilon
compressible::turbulenceModel & turbulence
word timeName
Definition: getTimeIndex.H:3
const Enum< EDCversions > EDCversionNames
const scalar EDCexp2[]
Definition: EDC.H:126
const scalar EDCexp1[]
Definition: EDC.H:125
const EDCversions EDCdefaultVersion
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
const dimensionSet dimless
Dimensionless.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const dimensionSet dimEnergy
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:53
messageStream Info
Information stream (stdout output on master, null elsewhere)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
dimensionedScalar sqrt(const dimensionedScalar &ds)
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
Type gMin(const FieldField< Field, Type > &f)
const dimensionSet dimVolume(pow3(dimLength))
Definition: dimensionSets.H:60
Type gMax(const FieldField< Field, Type > &f)
dimensionedScalar pow025(const dimensionedScalar &ds)
volScalarField & nu
tmp< volScalarField > trho
volScalarField & e
Definition: createFields.H:11
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333