singleStepCombustion.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) 2011-2021 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 "singleStepCombustion.H"
27 #include "fvmSup.H"
28 
29 
30 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
31 
33 {
34  const scalar Wu = mixture_.Wi(fuelIndex_);
35 
36  forAll(reaction_.lhs(), i)
37  {
38  const label speciei = reaction_.lhs()[i].index;
39  const scalar stoichCoeff = reaction_.lhs()[i].stoichCoeff;
40  specieStoichCoeffs_[speciei] = -stoichCoeff;
41  qFuel_.value() +=
42  mixture_.Hf(speciei)*mixture_.Wi(speciei)*stoichCoeff/Wu;
43  }
44 
45  forAll(reaction_.rhs(), i)
46  {
47  const label speciei = reaction_.rhs()[i].index;
48  const scalar stoichCoeff = reaction_.rhs()[i].stoichCoeff;
49  specieStoichCoeffs_[speciei] = stoichCoeff;
50  qFuel_.value() -=
51  mixture_.Hf(speciei)*mixture_.Wi(speciei)*stoichCoeff/Wu;
52  specieProd_[speciei] = -1;
53  }
54 
55  Info << "Fuel heat of combustion: " << qFuel_.value() << endl;
56 }
57 
58 
60 {
61  const label O2Index = mixture_.species()["O2"];
62  const scalar Wu = mixture_.Wi(fuelIndex_);
63 
64  stoicRatio_ =
65  (
68  + mixture_.Wi(O2Index)*mag(specieStoichCoeffs_[O2Index])
70 
71  s_ = mixture_.Wi(O2Index)*mag(specieStoichCoeffs_[O2Index])
72  /(Wu*mag(specieStoichCoeffs_[fuelIndex_]));
73 
74  Info << "stoichiometric air-fuel ratio: " << stoicRatio_.value() << endl;
75  Info << "stoichiometric oxygen-fuel ratio: " << s_.value() << endl;
76 }
77 
78 
80 {
81  scalar Wm = 0.0;
82  scalar totalMol = 0.0;
83  forAll(reaction_.rhs(), i)
84  {
85  label speciei = reaction_.rhs()[i].index;
86  totalMol += mag(specieStoichCoeffs_[speciei]);
87  }
88 
89  scalarList Xi(reaction_.rhs().size());
90 
91  forAll(reaction_.rhs(), i)
92  {
93  const label speciei = reaction_.rhs()[i].index;
94  Xi[i] = mag(specieStoichCoeffs_[speciei])/totalMol;
95  Wm += Xi[i]*mixture_.Wi(speciei);
96  }
97 
98  forAll(reaction_.rhs(), i)
99  {
100  const label speciei = reaction_.rhs()[i].index;
101  Yprod0_[speciei] = mixture_.Wi(speciei)/Wm*Xi[i];
102  }
103 
104  Info << "Maximum products mass concentrations: " << nl;
105  forAll(Yprod0_, i)
106  {
107  if (Yprod0_[i] > 0)
108  {
109  Info<< " " << mixture_.species()[i] << ": " << Yprod0_[i] << nl;
110  }
111  }
112 
113  // Normalise the stoichiometric coeff to mass
115  {
119  }
120 }
121 
122 
123 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
124 
126 (
127  const word& modelType,
128  const fluidReactionThermo& thermo,
130  const word& combustionProperties
131 )
132 :
133  combustionModel(modelType, thermo, turb, combustionProperties),
134  mixture_(dynamic_cast<const basicSpecieMixture&>(this->thermo())),
135  reaction_(mixture_.species(), this->subDict("reaction")),
136  stoicRatio_(dimensionedScalar("stoicRatio", dimless, 0)),
137  s_(dimensionedScalar("s", dimless, 0)),
138  qFuel_(dimensionedScalar("qFuel", sqr(dimVelocity), 0)),
140  Yprod0_(mixture_.species().size(), 0.0),
141  fres_(Yprod0_.size()),
142  fuelIndex_(mixture_.species()[thermo.properties().lookup("fuel")]),
143  specieProd_(Yprod0_.size(), 1),
144  wFuel_
145  (
146  IOobject
147  (
148  this->thermo().phasePropertyName("wFuel"),
149  this->mesh().time().timeName(),
150  this->mesh(),
153  ),
154  this->mesh(),
156  ),
157  semiImplicit_(readBool(this->coeffs_.lookup("semiImplicit")))
158 {
159  forAll(fres_, fresI)
160  {
161  IOobject header
162  (
163  "fres_" + mixture_.species()[fresI],
164  this->mesh().time().timeName(),
165  this->mesh()
166  );
167 
168  fres_.set
169  (
170  fresI,
171  new volScalarField
172  (
173  header,
174  this->mesh(),
175  dimensionedScalar("fres" + Foam::name(fresI), dimless, 0)
176  )
177  );
178  }
179 
180  calculateqFuel();
181 
183 
185 
186  if (semiImplicit_)
187  {
188  Info<< "Combustion mode: semi-implicit" << endl;
189  }
190  else
191  {
192  Info<< "Combustion mode: explicit" << endl;
193  }
194 }
195 
196 
197 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
198 
200 {}
201 
202 
203 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
204 
207 {
208  const label specieI = mixture_.species()[Y.member()];
209 
210  volScalarField wSpecie
211  (
212  wFuel_*specieStoichCoeffs()[specieI]
213  );
214 
215  if (semiImplicit_)
216  {
217  const label fNorm = specieProd()[specieI];
218  const volScalarField fres(this->fres(specieI));
219  wSpecie /= max(fNorm*(Y - fres), scalar(1e-2));
220 
221  return -fNorm*wSpecie*fres + fNorm*fvm::Sp(wSpecie, Y);
222  }
223  else
224  {
225  return wSpecie + fvm::Sp(0.0*wSpecie, Y);
226  }
227 }
228 
229 
232 {
233  const label fuelI = fuelIndex();
234  volScalarField& YFuel =
235  const_cast<volScalarField&>(this->thermo().composition().Y(fuelI));
236 
237  return -qFuel()*(R(YFuel) & YFuel);
238 }
239 
240 
242 {
243  if (combustionModel::read())
244  {
245  return true;
246  }
247  else
248  {
249  return false;
250  }
251 }
252 
253 
255 {
256  const label O2Index = mixture_.species()["O2"];
257  const volScalarField& YFuel = mixture_.Y()[fuelIndex_];
258  const volScalarField& YO2 = mixture_.Y()[O2Index];
259 
260  // reactants
261  forAll(reaction_.lhs(), i)
262  {
263  const label speciei = reaction_.lhs()[i].index;
264  if (speciei == fuelIndex_)
265  {
266  fres_[speciei] = max(YFuel - YO2/s_, scalar(0));
267  }
268  else if (speciei == O2Index)
269  {
270  fres_[speciei] = max(YO2 - YFuel*s_, scalar(0));
271  }
272  }
273 
274  // products
275  forAll(reaction_.rhs(), i)
276  {
277  const label speciei = reaction_.rhs()[i].index;
278  if (speciei != mixture_.defaultSpecie())
279  {
280  forAll(fres_[speciei], celli)
281  {
282  if (fres_[fuelIndex_][celli] > 0.0)
283  {
284  // rich mixture
285  fres_[speciei][celli] =
286  Yprod0_[speciei]
287  * (1.0 + YO2[celli]/s_.value() - YFuel[celli]);
288  }
289  else
290  {
291  // lean mixture
292  fres_[speciei][celli] =
293  Yprod0_[speciei]
294  * (
295  1.0
296  - YO2[celli]/s_.value()*stoicRatio_.value()
297  + YFuel[celli]*stoicRatio_.value()
298  );
299  }
300  }
301  }
302  }
303 }
304 
305 
306 // ************************************************************************* //
void calculateMaxProducts()
Calculate maximum products at stoichiometric mixture.
volScalarField wFuel_
Fuel consumption rate.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
dimensionedScalar s_
Stoichiometric oxygen-fuel mass ratio.
combustionModel(const word &modelType, const fluidReactionThermo &thermo, const compressibleMomentumTransportModel &turb, const word &combustionProperties=combustionPropertiesName)
Construct from components.
dimensionedScalar stoicRatio_
Stoichiometric air-fuel mass ratio.
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
const fvMesh & mesh() const
Return const access to the mesh database.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
label fuelIndex() const
Return the fuel specie index.
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
const basicSpecieMixture & mixture_
Reference to the mixture.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
scalarList Yprod0_
Mass concentrations at stoichiometric mixture for fres.
Base-class for multi-component fluid thermodynamic properties.
virtual scalar Wi(const label speciei) const =0
Molecular weight of the given specie [kg/kmol].
singleStepCombustion(const word &modelType, const fluidReactionThermo &thermo, const compressibleMomentumTransportModel &turb, const word &combustionProperties)
Construct from components.
PtrList< volScalarField > fres_
List of components residual.
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:636
const dimensionSet dimless
bool readBool(Istream &)
Definition: boolIO.C:60
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:239
const List< specieCoeffs > & lhs() const
Return the components of the left hand side.
Definition: reactionI.H:36
virtual basicSpecieMixture & composition()=0
Return the composition of the multi-component mixture.
static word member(const word &name)
Return member (name without the extension)
Definition: IOobject.C:149
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
Definition: dictionary.C:982
virtual const IOdictionary & properties() const =0
Return the dictionary.
const dimensionSet dimTime
bool semiImplicit_
Semi-implicit (true) or explicit (false) treatment.
A class for handling words, derived from string.
Definition: word.H:59
virtual tmp< fvScalarMatrix > R(volScalarField &Y) const
Fuel consumption rate matrix.
const List< scalar > & specieStoichCoeffs() const
Return the stoichiometric coefficient for the reaction.
const Type & value() const
Return const reference to value.
word timeName
Definition: getTimeIndex.H:3
dictionary coeffs_
Dictionary of the model.
virtual bool read()
Update properties from given dictionary.
scalarList specieStoichCoeffs_
Stoichiometric coefficient for the reaction.
static const char nl
Definition: Ostream.H:260
const dimensionSet dimVelocity
const dimensionSet dimMass
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
dimensionedScalar qFuel_
Heat of combustion [J/Kg].
virtual bool read()
Update properties from given dictionary.
const List< specieCoeffs > & rhs() const
Return the components of the right hand side.
Definition: reactionI.H:42
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
const Time & time() const
Return time.
Definition: IOobject.C:318
tmp< volScalarField > fres(const label index) const
Return the list of components residual.
PtrList< volScalarField > & Y
const dimensionedScalar & qFuel() const
Return the heat of combustion [J/Kg].
const dimensionSet dimVolume
messageStream Info
dimensioned< scalar > mag(const dimensioned< Type > &)
const doubleScalar e
Elementary charge.
Definition: doubleScalar.H:105
List< int > specieProd_
List to indicate if specie is produced/consumed.
A class for managing temporary objects.
Definition: PtrList.H:53
PtrList< volScalarField > & Y()
Return the mass-fraction fields.
virtual scalar Hf(const label speciei) const =0
Enthalpy of formation [J/kg].
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
reaction reaction_
The single-step reaction.
virtual tmp< volScalarField > Qdot() const
Heat release rate [kg/m/s^3].
void fresCorrect()
Calculates the residual for all components.
void massAndAirStoichRatios()
Calculate air/fuel and oxygen/fuel ratio.
Abstract base class for turbulence models (RAS, LES and laminar).
const List< int > & specieProd() const
Return the list to indicate if specie is produced/consumed.
Calculate the matrix for implicit and explicit sources.
const speciesTable & species() const
Return the table of species.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:844
const fluidReactionThermo & thermo() const
Return const access to the thermo.
label defaultSpecie() const
Return the index of the default specie.