33template<
class ReactionThermo>
36 const word& modelType,
39 const word& combustionProperties
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_)])),
62 this->
thermo().phasePropertyName(typeName +
":kappa"),
76template<
class ReactionThermo>
83template<
class ReactionThermo>
109 const scalar
nu =
mu[i]/(
rho[i] + SMALL);
115 const scalar CtauI =
min(C1_/(Da*
sqrt(ReT + 1)), 2.1377);
117 const scalar CgammaI =
120 const scalar gammaL =
136 pow(gammaL, exp1_)/(1 -
pow(gammaL, exp2_)),
148 const scalar
nu =
mu[i]/(
rho[i] + SMALL);
149 const scalar gammaL =
164 pow(gammaL, exp1_)/(1 -
pow(gammaL, exp2_)),
173 Info<<
"Chemistry time solved max/min : "
176 this->chemistryPtr_->solve(tauStar);
181template<
class ReactionThermo>
189template<
class ReactionThermo>
199 this->
thermo().phasePropertyName(typeName +
":Qdot"),
213 tQdot.
ref() = kappa_*this->chemistryPtr_->Qdot();
220template<
class ReactionThermo>
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_)]);
#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,...
Eddy Dissipation Concept (EDC) turbulent combustion model.
virtual void correct()
Correct combustion rate.
virtual ~EDC()
Destructor.
virtual tmp< volScalarField > Qdot() const
Heat release rate [kg/m/s3].
virtual bool read()
Update properties from given dictionary.
Laminar combustion model.
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.
A class for handling words, derived from Foam::string.
PtrList< volScalarField > & Y
const volScalarField & mu
compressible::turbulenceModel & turbulence
const Enum< EDCversions > EDCversionNames
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.
const dimensionSet dimless
Dimensionless.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const dimensionSet dimEnergy
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
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.
dimensionedScalar sqrt(const dimensionedScalar &ds)
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
static constexpr const zero Zero
Global zero (0)
Type gMin(const FieldField< Field, Type > &f)
const dimensionSet dimVolume(pow3(dimLength))
Type gMax(const FieldField< Field, Type > &f)
dimensionedScalar pow025(const dimensionedScalar &ds)
tmp< volScalarField > trho
#define forAll(list, i)
Loop across all elements in list.