laminar.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) 2013-2023 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 "laminar.H"
27 #include "fvmSup.H"
28 #include "localEulerDdtScheme.H"
30 
31 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35 namespace combustionModels
36 {
39 }
40 }
41 
42 
43 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
44 
46 (
47  const word& modelType,
50  const word& combustionProperties
51 )
52 :
54  (
55  modelType,
56  thermo,
57  turb,
58  combustionProperties
59  ),
60  integrateReactionRate_
61  (
62  this->coeffs().lookupOrDefault("integrateReactionRate", true)
63  ),
64  maxIntegrationTime_
65  (
66  this->coeffs().lookupOrDefault("maxIntegrationTime", vGreat)
67  ),
68  outerCorrect_
69  (
70  this->coeffs().lookupOrDefault("outerCorrect", false)
71  ),
72  timeIndex_(-1),
73  chemistryPtr_(basicChemistryModel::New(thermo))
74 {
75  if (integrateReactionRate_)
76  {
77  Info<< " using integrated reaction rate" << endl;
78  }
79  else
80  {
81  Info<< " using instantaneous reaction rate" << endl;
82  }
83 }
84 
85 
86 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
87 
89 {}
90 
91 
92 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
93 
95 {
96  if
97  (
98  (integrateReactionRate_ || !outerCorrect_)
99  && timeIndex_ == this->mesh().time().timeIndex()
100  )
101  {
102  return;
103  }
104 
105  if (integrateReactionRate_)
106  {
107  if (fv::localEulerDdt::enabled(this->mesh()))
108  {
109  const scalarField& rDeltaT =
110  fv::localEulerDdt::localRDeltaT(this->mesh());
111 
112  chemistryPtr_->solve(min(1/rDeltaT, maxIntegrationTime_)());
113  }
114  else
115  {
116  const scalar deltaT = this->mesh().time().deltaTValue();
117 
118  chemistryPtr_->solve(min(deltaT, maxIntegrationTime_));
119  }
120  }
121  else
122  {
123  chemistryPtr_->calculate();
124  }
125 
126  timeIndex_ = this->mesh().time().timeIndex();
127 }
128 
129 
132 {
133  return chemistryPtr_->RR()[speciei];
134 }
135 
136 
139 {
141  fvScalarMatrix& Su = tSu.ref();
142 
143  const label specieI = this->thermo().species()[Y.member()];
144  Su += chemistryPtr_->RR()[specieI];
145 
146  return tSu;
147 }
148 
149 
151 {
152  return chemistryPtr_->Qdot();
153 }
154 
155 
157 {
158  if (combustionModel::read())
159  {
160  integrateReactionRate_ =
161  this->coeffs().lookupOrDefault("integrateReactionRate", true);
162  maxIntegrationTime_ =
163  this->coeffs().lookupOrDefault("maxIntegrationTime", vGreat);
164  outerCorrect_ =
165  this->coeffs().lookupOrDefault("outerCorrect", false);
166  return true;
167  }
168  else
169  {
170  return false;
171  }
172 }
173 
174 
175 // ************************************************************************* //
Macros for easy insertion into run-time selection tables.
Generic GeometricField class.
static word member(const word &name)
Return member (name without the extension)
Definition: IOobject.C:149
Base class for chemistry models.
Base class for combustion models.
virtual bool read()
Update properties from given dictionary.
Laminar combustion model.
Definition: laminar.H:55
laminar(const word &modelType, const fluidMulticomponentThermo &thermo, const compressibleMomentumTransportModel &turb, const word &combustionProperties)
Construct from components.
Definition: laminar.C:46
virtual void correct()
Correct combustion rate.
Definition: laminar.C:94
virtual ~laminar()
Destructor.
Definition: laminar.C:88
virtual tmp< volScalarField::Internal > R(const label speciei) const
Specie consumption rate field.
Definition: laminar.C:131
virtual tmp< volScalarField > Qdot() const
Heat release rate [kg/m/s^3].
Definition: laminar.C:150
virtual bool read()
Update properties from given dictionary.
Definition: laminar.C:156
Base class for single-phase compressible turbulence models.
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
Base-class for multi-component fluid thermodynamic properties.
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Definition: fvMatrix.H:118
static const volScalarField & localRDeltaT(const fvMesh &mesh)
Return the reciprocal of the local time-step.
Definition: localEulerDdt.C:46
static bool enabled(const fvMesh &mesh)
Return true if LTS is enabled.
Definition: localEulerDdt.C:37
A class for managing temporary objects.
Definition: tmp.H:55
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:181
A class for handling words, derived from string.
Definition: word.H:62
Calculate the matrix for implicit and explicit sources.
defineTypeNameAndDebug(diffusion, 0)
addToRunTimeSelectionTable(combustionModel, diffusion, dictionary)
autoPtr< CompressibleMomentumTransportModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const viscosity &viscosity)
tmp< VolField< Type > > Su(const VolField< Type > &su, const VolField< Type > &vf)
Definition: fvcSup.C:44
Namespace for OpenFOAM.
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
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
fvMatrix< scalar > fvScalarMatrix
Definition: fvMatricesFwd.H:42
messageStream Info
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
const dimensionSet dimTime
const dimensionSet dimMass
label timeIndex
Definition: getTimeIndex.H:4
PtrList< volScalarField > & Y
fluidMulticomponentThermo & thermo
Definition: createFields.H:31