36 template<
class ReactionThermo,
class ThermoType>
42 this->coeffs().readIfPresent(
"YoxStream", YoxStream_);
43 this->coeffs().readIfPresent(
"YfStream", YfStream_);
44 this->coeffs().readIfPresent(
"sigma", sigma_);
45 this->coeffs().readIfPresent(
"ftCorr", ftCorr_);
46 this->coeffs().readIfPresent(
"alpha", alpha_);
47 this->coeffs().readIfPresent(
"laminarIgn", laminarIgn_);
54 const label nReactions = reactions_.
size();
56 for (label
k=0;
k < nReactions;
k++)
70 RijPtr_[
k].storePrevIter();
75 const label fuelIndex = species.
find(fuelNames_[
k]);
76 const label oxidantIndex = species.
find(oxidantNames_[
k]);
78 const scalar Wu = specieThermo_[fuelIndex].W();
79 const scalar Wox = specieThermo_[oxidantIndex].W();
83 const label specieI = lhs[i].index;
84 specieStoichCoeffs[specieI] = -lhs[i].stoichCoeff;
86 specieThermo_[specieI].hc()*lhs[i].stoichCoeff/Wu;
91 const label specieI = rhs[i].index;
92 specieStoichCoeffs[specieI] = rhs[i].stoichCoeff;
94 specieThermo_[specieI].hc()*rhs[i].stoichCoeff/Wu;
97 Info <<
"Fuel heat of combustion : " << qFuel_[
k] <<
endl;
100 (Wox*
mag(specieStoichCoeffs[oxidantIndex]))
101 / (Wu*
mag(specieStoichCoeffs[fuelIndex]));
103 Info <<
"stoichiometric oxygen-fuel ratio : " << s_[
k] <<
endl;
105 stoicRatio_[
k] = s_[
k]*YfStream_[
k]/YoxStream_[
k];
107 Info <<
"stoichiometric air-fuel ratio : " << stoicRatio_[
k] <<
endl;
109 const scalar fStoich = 1.0/(1.0 + stoicRatio_[
k]);
111 Info <<
"stoichiometric mixture fraction : " << fStoich <<
endl;
118 template<
class ReactionThermo,
class ThermoType>
122 const word& modelType,
125 const word& combustionProperties
144 RijPtr_(reactions_.size()),
145 Ci_(reactions_.size(), 1.0),
146 fuelNames_(this->coeffs().
lookup(
"fuels")),
147 oxidantNames_(this->coeffs().
lookup(
"oxidants")),
148 qFuel_(reactions_.size()),
149 stoicRatio_(reactions_.size()),
150 s_(reactions_.size()),
151 YoxStream_(reactions_.size(), 0.23),
152 YfStream_(reactions_.size(), 1.0),
153 sigma_(reactions_.size(), 0.02),
154 oxidantRes_(this->coeffs().
lookup(
"oxidantRes")),
155 ftCorr_(reactions_.size(),
Zero),
165 template<
class ReactionThermo,
class ThermoType>
169 return this->chemistryPtr_->tc();
173 template<
class ReactionThermo,
class ThermoType>
182 const label nReactions = reactions_.
size();
186 for (label
k=0;
k < nReactions;
k++)
203 const label fuelIndex = species.
find(fuelNames_[
k]);
207 Rijl.
ref() = -this->chemistryPtr_->calculateRR(
k, fuelIndex);
218 const label lIndex = lhs[l].index;
219 this->chemistryPtr_->RR(lIndex) =
225 const label rIndex = rhs[l].index;
226 this->chemistryPtr_->RR(rIndex) =
232 for (label
k=0;
k < nReactions;
k++)
234 const label fuelIndex = species.
find(fuelNames_[
k]);
235 const label oxidantIndex = species.
find(oxidantNames_[
k]);
247 s_[
k]*Yfuel - (Yox - YoxStream_[
k])
250 s_[
k]*YfStream_[
k] + YoxStream_[
k]
254 const scalar
sigma = sigma_[
k];
256 const scalar fStoich = 1.0/(1.0 + stoicRatio_[
k]) + ftCorr_[
k];
261 Yox/
max(oxidantRes_[
k], 1
e-3)
267 1.0 +
sqr(OAvailScaled)
295 min(RijkDiff, topHatFilter*RijlPtr[
k]*
pos(Yox)*
pos(Yfuel));
304 if (
debug && this->mesh_.time().writeTime())
314 scalar fuelStoic = 1.0;
317 const label lIndex = lhs[l].index;
318 if (lIndex == fuelIndex)
320 fuelStoic = lhs[l].stoichCoeff;
325 const scalar MwFuel = specieThermo_[fuelIndex].W();
330 const label lIndex = lhs[l].index;
332 const scalar stoichCoeff = lhs[l].stoichCoeff;
334 this->chemistryPtr_->RR(lIndex) +=
335 -Rijk*stoichCoeff*specieThermo_[lIndex].W()/fuelStoic/MwFuel;
342 const label rIndex = rhs[r].index;
344 const scalar stoichCoeff = rhs[r].stoichCoeff;
346 this->chemistryPtr_->RR(rIndex) +=
347 Rijk*stoichCoeff*specieThermo_[rIndex].W()/fuelStoic/MwFuel;
354 template<
class ReactionThermo,
class ThermoType>
362 auto&
Su = tSu.ref();
366 const label specieI =
369 Su += this->chemistryPtr_->RR(specieI);
376 template<
class ReactionThermo,
class ThermoType>
393 Qdot = this->chemistryPtr_->Qdot();
400 template<
class ReactionThermo,
class ThermoType>
406 this->coeffs().readIfPresent(
"Ci", Ci_);
407 this->coeffs().readIfPresent(
"sigma", sigma_);
408 this->coeffs().readIfPresent(
"oxidantRes", oxidantRes_);
409 this->coeffs().readIfPresent(
"ftCorr", ftCorr_);
410 this->coeffs().readIfPresent(
"alpha", alpha_);
411 this->coeffs().readIfPresent(
"laminarIgn", laminarIgn_);
static const word & zeroGradientType() noexcept
The type name for zeroGradient patch fields.
virtual void correct()
Correct combustion rate.
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
void size(const label n)
Older name for setAddressableSize.
const word zeroGradientType
A zeroGradient patch field type.
virtual tmp< volScalarField > Qdot() const
Heat release rate calculated from fuel consumption rate matrix.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
compressible::turbulenceModel & turb
dimensionedScalar sqrt(const dimensionedScalar &ds)
Ostream & endl(Ostream &os)
Add newline and flush stream.
const speciesTable & species() const
Return the table of species.
label k
Boltzmann constant.
Lookup type of boundary radiation properties.
PtrList< volScalarField > & Y()
Return the mass-fraction fields.
dimensionedScalar sigma("sigma", dimMass/sqr(dimTime), transportProperties)
Hold specie index and its coefficients in the reaction rate expression.
static word member(const word &name)
Return member (name without the extension)
#define forAll(list, i)
Loop across all elements in list.
psiReactionThermo & thermo
dimensionedScalar pos(const dimensionedScalar &ds)
const dimensionSet dimVolume(pow3(dimLength))
const dimensionedScalar e
Elementary charge.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for expressions::valueTypeCode::INVALID.
dimensionedScalar exp(const dimensionedScalar &ds)
Calculate the gradient of the given field.
A class for handling words, derived from Foam::string.
virtual bool write(const bool writeOnProc=true) const
Write using setting from DB.
static tmp< T > New(Args &&... args)
Construct tmp with forwarding arguments.
constexpr scalar pi(M_PI)
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Chemistry model wrapper for combustion models.
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
static tmp< GeometricField< scalar, fvPatchField, volMesh > > New(const word &name, IOobjectOption::registerOption regOpt, const Mesh &mesh, const dimensionSet &dims, const word &patchFieldType=fvPatchField< scalar >::calculatedType())
Return tmp field (NO_READ, NO_WRITE) from name, mesh, dimensions and patch type. [Takes current timeN...
Abstract base class for turbulence models (RAS, LES and laminar).
int debug
Static debugging option.
bool readIfPresent(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX) const
Find an entry if present, and assign to T val. FatalIOError if it is found and the number of tokens i...
Diffusion based turbulent combustion model for multicomponent species.
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("pos", dimless, 1.0));surfaceScalarField neg(IOobject("neg", runTime.timeName(), mesh), mesh, dimensionedScalar("neg", dimless, -1.0));surfaceScalarField phi("phi", fvc::flux(rhoU));Info<< "Creating turbulence model\"<< endl;autoPtr< compressible::turbulenceModel > turbulence(compressible::turbulenceModel::New(rho, U, phi, thermo))
const dimensionSet dimEnergy
virtual tmp< fvScalarMatrix > R(volScalarField &Y) const
Fuel consumption rate matrix.
void relax(const scalar alpha)
Relax field (for steady-state solution).
A wordList with hashed named lookup, which can be faster in some situations than using the normal lis...
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
virtual bool read()
Update properties from given dictionary.
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers...
PtrList< volScalarField > & Y
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
messageStream Info
Information stream (stdout output on master, null elsewhere)
Internal & ref(const bool updateAccessTime=true)
Same as internalFieldRef()
virtual basicSpecieMixture & composition()=0
Return the composition of the multi-component mixture.
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
A class for managing temporary objects.
label find(const word &val) const
Find index of the value (searches the hash).
Do not request registration (bool: false)
static constexpr const zero Zero
Global zero (0)