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-2015 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  (
73  "singleStepCombustion<CombThermoType, ThermoType>::"
74  "singleStepCombustion"
75  "("
76  "const word&, "
77  "const fvMesh& "
78  "const word&"
79  ")"
80  )
81  << "Inconsistent thermo package for " << this->type() << " model:\n"
82  << " " << this->thermo().type() << nl << nl
83  << "Please select a thermo package based on "
84  << "singleStepReactingMixture" << exit(FatalError);
85  }
86 
87  if (semiImplicit_)
88  {
89  Info<< "Combustion mode: semi-implicit" << endl;
90  }
91  else
92  {
93  Info<< "Combustion mode: explicit" << endl;
94  }
95 }
96 
97 
98 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
99 
100 template<class CombThermoType, class ThermoType>
102 {}
103 
104 
105 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
106 
107 template<class CombThermoType, class ThermoType>
109 (
111 ) const
112 {
113  const label specieI =
114  this->thermoPtr_->composition().species()[Y.member()];
115 
116  volScalarField wSpecie
117  (
118  wFuel_*singleMixturePtr_->specieStoichCoeffs()[specieI]
119  );
120 
121  if (semiImplicit_)
122  {
123  const label fNorm = singleMixturePtr_->specieProd()[specieI];
124  const volScalarField fres(singleMixturePtr_->fres(specieI));
125  wSpecie /= max(fNorm*(Y - fres), scalar(1e-2));
126 
127  return -fNorm*wSpecie*fres + fNorm*fvm::Sp(wSpecie, Y);
128  }
129  else
130  {
131  return wSpecie + fvm::Sp(0.0*wSpecie, Y);
132  }
133 }
134 
135 
136 template<class CombThermoType, class ThermoType>
139 {
140  const label fuelI = singleMixturePtr_->fuelIndex();
141  volScalarField& YFuel =
142  const_cast<volScalarField&>(this->thermoPtr_->composition().Y(fuelI));
143 
144  return -singleMixturePtr_->qFuel()*(R(YFuel) & YFuel);
145 }
146 
147 
148 template<class CombThermoType, class ThermoType>
151 {
153  (
154  new volScalarField
155  (
156  IOobject
157  (
158  IOobject::groupName("dQ", this->phaseName_),
159  this->mesh_.time().timeName(),
160  this->mesh_,
163  false
164  ),
165  this->mesh_,
166  dimensionedScalar("dQ", dimEnergy/dimTime, 0.0),
167  zeroGradientFvPatchScalarField::typeName
168  )
169  );
170 
171  if (this->active())
172  {
173  volScalarField& dQ = tdQ();
174  dQ.dimensionedInternalField() = this->mesh().V()*Sh()();
175  }
176  return tdQ;
177 }
178 
179 
180 template<class CombThermoType, class ThermoType>
182 {
183  if (CombThermoType::read())
184  {
185  return true;
186  }
187  else
188  {
189  return false;
190  }
191 }
192 
193 
194 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
195 
196 } // End namespace combustionModels
197 } // End namespace Foam
198 
199 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
word member() const
Return member (name without the extension)
Definition: IOobject.C:278
virtual tmp< volScalarField > Sh() const
Sensible enthalpy source term.
bool readBool(Istream &)
Definition: boolIO.C:60
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
scalar Sh
Definition: solveChemistry.H:2
Calculate the matrix for implicit and explicit sources.
virtual tmp< fvScalarMatrix > R(volScalarField &Y) const
Fuel consumption rate matrix.
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:51
#define R(A, B, C, D, E, F, K, M)
A class for handling words, derived from string.
Definition: word.H:59
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
DimensionedInternalField & dimensionedInternalField()
Return dimensioned internal field.
const dimensionSet dimEnergy
messageStream Info
dynamicFvMesh & mesh
virtual bool read()
Update properties from given dictionary.
tmp< fvMatrix< Type > > Sp(const DimensionedField< scalar, volMesh > &, const GeometricField< Type, fvPatchField, volMesh > &)
Namespace for OpenFOAM.
virtual tmp< volScalarField > dQ() const
Heat release rate calculated from fuel consumption rate matrix.
static word groupName(Name name, const word &group)
static const char nl
Definition: Ostream.H:260
const double e
Elementary charge.
Definition: doubleFloat.H:78
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
PtrList< volScalarField > & Y
Definition: createFields.H:36
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:49
dQ
Definition: YEEqn.H:14
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:589
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
error FatalError
const dimensionSet dimVolume(pow3(dimLength))
Definition: dimensionSets.H:58
bool read(const char *, int32_t &)
Definition: int32IO.C:87
Single step reacting mixture.
psiReactionThermo & thermo
Definition: createFields.H:32
bool isA(const Type &t)
Check if a dynamic_cast to typeid is possible.
Definition: typeInfo.H:134
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
A class for managing temporary objects.
Definition: PtrList.H:118
word timeName
Definition: getTimeIndex.H:3