singleStepCombustion.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 | Copyright (C) 2011-2016 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 namespace Foam
30 {
31 namespace combustionModels
32 {
33 
34 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
35 
36 template<class CombThermoType, class ThermoType>
37 singleStepCombustion<CombThermoType, ThermoType>::singleStepCombustion
38 (
39  const word& modelType,
40  const fvMesh& mesh,
41  const word& phaseName
42 )
43 :
44  CombThermoType(modelType, mesh, phaseName),
45  singleMixturePtr_(NULL),
46  wFuel_
47  (
48  IOobject
49  (
50  IOobject::groupName("wFuel", phaseName),
51  this->mesh().time().timeName(),
52  this->mesh(),
55  ),
56  this->mesh(),
58  ),
59  semiImplicit_(readBool(this->coeffs_.lookup("semiImplicit")))
60 {
62  {
63  singleMixturePtr_ =
65  (
66  this->thermo()
67  );
68  }
69  else
70  {
72  << "Inconsistent thermo package for " << this->type() << " model:\n"
73  << " " << this->thermo().type() << nl << nl
74  << "Please select a thermo package based on "
75  << "singleStepReactingMixture" << exit(FatalError);
76  }
77 
78  if (semiImplicit_)
79  {
80  Info<< "Combustion mode: semi-implicit" << endl;
81  }
82  else
83  {
84  Info<< "Combustion mode: explicit" << endl;
85  }
86 }
87 
88 
89 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
90 
91 template<class CombThermoType, class ThermoType>
93 {}
94 
95 
96 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
97 
98 template<class CombThermoType, class ThermoType>
100 (
102 ) const
103 {
104  const label specieI =
105  this->thermoPtr_->composition().species()[Y.member()];
106 
107  volScalarField wSpecie
108  (
109  wFuel_*singleMixturePtr_->specieStoichCoeffs()[specieI]
110  );
111 
112  if (semiImplicit_)
113  {
114  const label fNorm = singleMixturePtr_->specieProd()[specieI];
115  const volScalarField fres(singleMixturePtr_->fres(specieI));
116  wSpecie /= max(fNorm*(Y - fres), scalar(1e-2));
117 
118  return -fNorm*wSpecie*fres + fNorm*fvm::Sp(wSpecie, Y);
119  }
120  else
121  {
122  return wSpecie + fvm::Sp(0.0*wSpecie, Y);
123  }
124 }
125 
126 
127 template<class CombThermoType, class ThermoType>
130 {
131  const label fuelI = singleMixturePtr_->fuelIndex();
132  volScalarField& YFuel =
133  const_cast<volScalarField&>(this->thermoPtr_->composition().Y(fuelI));
134 
135  return -singleMixturePtr_->qFuel()*(R(YFuel) & YFuel);
136 }
137 
138 
139 template<class CombThermoType, class ThermoType>
142 {
144  (
145  new volScalarField
146  (
147  IOobject
148  (
149  IOobject::groupName("dQ", this->phaseName_),
150  this->mesh_.time().timeName(),
151  this->mesh_,
154  false
155  ),
156  this->mesh_,
158  )
159  );
160 
161  if (this->active())
162  {
163  volScalarField& dQ = tdQ.ref();
164  dQ.ref() = this->mesh().V()*Sh()();
165  }
166  return tdQ;
167 }
168 
169 
170 template<class CombThermoType, class ThermoType>
172 {
173  if (CombThermoType::read())
174  {
175  return true;
176  }
177  else
178  {
179  return false;
180  }
181 }
182 
183 
184 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
185 
186 } // End namespace combustionModels
187 } // End namespace Foam
188 
189 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
bool isA(const Type &t)
Check if a dynamic_cast to typeid is possible.
Definition: typeInfo.H:134
virtual tmp< volScalarField > dQ() const
Heat release rate calculated from fuel consumption rate matrix.
virtual tmp< volScalarField > Sh() const
Sensible enthalpy source term.
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
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
const double e
Elementary charge.
Definition: doubleFloat.H:78
error FatalError
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
Single step reacting mixture.
tmp< fvMatrix< Type > > Sp(const DimensionedField< scalar, volMesh > &, const GeometricField< Type, fvPatchField, volMesh > &)
scalar Sh
Definition: solveChemistry.H:2
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
virtual tmp< fvScalarMatrix > R(volScalarField &Y) const
Fuel consumption rate matrix.
word member() const
Return member (name without the extension)
Definition: IOobject.C:254
bool readBool(Istream &)
Definition: boolIO.C:60
psiReactionThermo & thermo
Definition: createFields.H:31
const dimensionSet dimVolume(pow3(dimLength))
Definition: dimensionSets.H:58
bool read(const char *, int32_t &)
Definition: int32IO.C:85
dynamicFvMesh & mesh
dQ
Definition: YEEqn.H:14
A class for handling words, derived from string.
Definition: word.H:59
static word groupName(Name name, const word &group)
word timeName
Definition: getTimeIndex.H:3
static const char nl
Definition: Ostream.H:262
const dimensionSet dimEnergy
Internal & ref()
Return a reference to the dimensioned internal field.
virtual bool read()
Update properties from given dictionary.
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:461
#define R(A, B, C, D, E, F, K, M)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
PtrList< volScalarField > & Y
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:51
messageStream Info
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:49
A class for managing temporary objects.
Definition: PtrList.H:54
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:174
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
Calculate the matrix for implicit and explicit sources.
Namespace for OpenFOAM.