33template<
class ThermoType>
37 const scalar Wu = this->speciesData()[fuelIndex_].W();
41 const label speciei =
reaction.lhs()[i].index;
42 const scalar stoichCoeff =
reaction.lhs()[i].stoichCoeff;
43 specieStoichCoeffs_[speciei] = -stoichCoeff;
44 qFuel_.value() += this->speciesData()[speciei].hc()*stoichCoeff/Wu;
49 const label speciei =
reaction.rhs()[i].index;
50 const scalar stoichCoeff =
reaction.rhs()[i].stoichCoeff;
51 specieStoichCoeffs_[speciei] = stoichCoeff;
52 qFuel_.value() -= this->speciesData()[speciei].hc()*stoichCoeff/Wu;
53 specieProd_[speciei] = -1;
56 Info<<
"Fuel heat of combustion :" << qFuel_.value() <<
endl;
60template<
class ThermoType>
63 const label O2Index = this->species().find(
"O2");
64 const scalar Wu = this->speciesData()[fuelIndex_].W();
67 (this->speciesData()[inertIndex_].W()
68 * specieStoichCoeffs_[inertIndex_]
69 + this->speciesData()[O2Index].W()
70 *
mag(specieStoichCoeffs_[O2Index]))
71 / (Wu*
mag(specieStoichCoeffs_[fuelIndex_]));
74 (this->speciesData()[O2Index].W()
75 *
mag(specieStoichCoeffs_[O2Index]))
76 / (Wu*
mag(specieStoichCoeffs_[fuelIndex_]));
78 Info<<
"stoichiometric air-fuel ratio :" << stoicRatio_.value() <<
endl;
80 Info<<
"stoichiometric oxygen-fuel ratio :" << s_.value() <<
endl;
84template<
class ThermoType>
90 scalar totalMol = 0.0;
93 label speciei =
reaction.rhs()[i].index;
94 totalMol +=
mag(specieStoichCoeffs_[speciei]);
101 const label speciei =
reaction.rhs()[i].index;
102 Xi[i] =
mag(specieStoichCoeffs_[speciei])/totalMol;
104 Wm += Xi[i]*this->speciesData()[speciei].W();
109 const label speciei =
reaction.rhs()[i].index;
110 Yprod0_[speciei] = this->speciesData()[speciei].W()/Wm*Xi[i];
113 Info<<
"Maximum products mass concentrations:" <<
nl;
118 Info<<
" " << this->species()[i] <<
": " << Yprod0_[i] <<
nl;
123 forAll(specieStoichCoeffs_, i)
125 specieStoichCoeffs_[i] =
126 specieStoichCoeffs_[i]
127 * this->speciesData()[i].W()
128 / (this->speciesData()[fuelIndex_].W()
129 *
mag(specieStoichCoeffs_[fuelIndex_]));
134template<
class ThermoType>
139 const label O2Index = this->species().find(
"O2");
146 const label speciei =
reaction.lhs()[i].index;
147 if (speciei == fuelIndex_)
149 fres_[speciei] =
max(YFuel - YO2/s_, scalar(0));
151 else if (speciei == O2Index)
153 fres_[speciei] =
max(YO2 - YFuel*s_, scalar(0));
161 const label speciei =
reaction.rhs()[i].index;
162 if (speciei != inertIndex_)
164 forAll(fres_[speciei], celli)
166 if (fres_[fuelIndex_][celli] > 0.0)
169 fres_[speciei][celli] =
171 * (1.0 + YO2[celli]/s_.value() - YFuel[celli]);
176 fres_[speciei][celli] =
180 - YO2[celli]/s_.value()*stoicRatio_.value()
181 + YFuel[celli]*stoicRatio_.value()
192template<
class ThermoType>
197 const word& phaseName
204 specieStoichCoeffs_(this->species_.size(), Zero),
205 Yprod0_(this->species_.size(), Zero),
206 fres_(Yprod0_.size()),
207 inertIndex_(this->species().find(
thermoDict.get<
word>(
"inertSpecie"))),
209 specieProd_(Yprod0_.size(), 1)
211 if (this->
size() == 1)
217 "fres_" + this->
species()[fresI],
245 <<
"Only one reaction required for single step reaction"
253template<
class ThermoType>
Defines the attributes of an object for which implicit objectRegistry management is supported,...
virtual bool read()
Re-read model coefficients if they have changed.
Simple extension of ReactionThermo to handle reaction kinetics in addition to the equilibrium thermod...
static word timeName(const scalar t, const int precision=precision_)
label size() const noexcept
The number of elements in the list.
const speciesTable & species() const
Return the table of species.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Mesh data needed to do the Finite Volume discretisation.
const Time & time() const
Return the top-level database.
Single step reacting mixture.
void fresCorrect()
Calculates the residual for all components.
void calculateMaxProducts()
Calculate maximum products at stoichiometric mixture.
PtrList< volScalarField > fres_
List of components residual.
void massAndAirStoichRatios()
Calculate air/fuel and oxygen/fuel ratio.
void calculateqFuel()
Calculate qFuel.
A class for handling words, derived from Foam::string.
PtrList< volScalarField > & Y
const dictionary & thermoDict
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
CombustionModel< rhoReactionThermo > & reaction
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
const dimensionSet dimless
Dimensionless.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const dimensionSet dimVelocity
messageStream Info
Information stream (stdout output on master, null elsewhere)
Ostream & endl(Ostream &os)
Add newline and flush stream.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
static constexpr const zero Zero
Global zero (0)
errorManipArg< error, int > exit(error &err, const int errNo=1)
constexpr char nl
The newline '\n' character (0x0a)
#define forAll(list, i)
Loop across all elements in list.