laminar.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 | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2013-2017 OpenFOAM Foundation
9 Copyright (C) 2019-2020 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27\*---------------------------------------------------------------------------*/
28
29#include "laminar.H"
30#include "fvmSup.H"
31#include "localEulerDdtScheme.H"
32
33// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
34
35template<class ReactionThermo>
37(
38 const word& modelType,
39 ReactionThermo& thermo,
41 const word& combustionProperties
42)
43:
44 ChemistryCombustion<ReactionThermo>
45 (
46 modelType,
47 thermo,
48 turb,
49 combustionProperties
50 ),
51 integrateReactionRate_
52 (
53 this->coeffs().getOrDefault("integrateReactionRate", true)
54 )
55{
56 if (integrateReactionRate_)
57 {
58 Info<< " using integrated reaction rate" << endl;
59 }
60 else
61 {
62 Info<< " using instantaneous reaction rate" << endl;
63 }
64}
65
66
67// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
68
69template<class ReactionThermo>
71{}
72
73
74// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
75
76template<class ReactionThermo>
79{
80 return this->chemistryPtr_->tc();
81}
82
83
84template<class ReactionThermo>
86{
87 if (this->active())
88 {
89 if (integrateReactionRate_)
90 {
92 {
93 const scalarField& rDeltaT =
94 fv::localEulerDdt::localRDeltaT(this->mesh());
95
96 scalar maxTime;
97 if (this->coeffs().readIfPresent("maxIntegrationTime", maxTime))
98 {
99 this->chemistryPtr_->solve
100 (
101 min(1.0/rDeltaT, maxTime)()
102 );
103 }
104 else
105 {
106 this->chemistryPtr_->solve((1.0/rDeltaT)());
107 }
108 }
109 else
110 {
111 this->chemistryPtr_->solve(this->mesh().time().deltaTValue());
112 }
113 }
114 else
115 {
116 this->chemistryPtr_->calculate();
117 }
118 }
119}
120
121
122template<class ReactionThermo>
125{
127
128 fvScalarMatrix& Su = tSu.ref();
129
130 if (this->active())
131 {
132 const label specieI =
133 this->thermo().composition().species().find(Y.member());
134
135 Su += this->chemistryPtr_->RR(specieI);
136 }
137
138 return tSu;
139}
140
141
142template<class ReactionThermo>
145{
147 (
149 (
151 (
152 this->thermo().phasePropertyName(typeName + ":Qdot"),
153 this->mesh().time().timeName(),
154 this->mesh(),
157 false
158 ),
159 this->mesh(),
161 )
162 );
163
164 if (this->active())
165 {
166 tQdot.ref() = this->chemistryPtr_->Qdot();
167 }
168
169 return tQdot;
170}
171
172
173template<class ReactionThermo>
175{
177 {
178 integrateReactionRate_ =
179 this->coeffs().getOrDefault("integrateReactionRate", true);
180 return true;
181 }
182
183 return false;
184}
185
186
187// ************************************************************************* //
#define R(A, B, C, D, E, F, K, M)
compressible::turbulenceModel & turb
Chemistry model wrapper for combustion models.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:170
Laminar combustion model.
Definition: laminar.H:59
tmp< volScalarField > tc() const
Return the chemical time scale.
Definition: laminar.C:78
virtual void correct()
Correct combustion rate.
Definition: laminar.C:85
virtual ~laminar()
Destructor.
Definition: laminar.C:70
virtual tmp< volScalarField > Qdot() const
Heat release rate [kg/m/s3].
Definition: laminar.C:144
virtual bool read()
Update properties from given dictionary.
Definition: laminar.C:174
Abstract base class for turbulence models (RAS, LES and laminar).
virtual bool enabled() const
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Definition: fvMatrix.H:121
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
A class for managing temporary objects.
Definition: tmp.H:65
T & ref() const
Definition: tmpI.H:227
A class for handling words, derived from Foam::string.
Definition: word.H:68
PtrList< volScalarField > & Y
dynamicFvMesh & mesh
Calculate the finiteVolume matrix for implicit and explicit sources.
zeroField Su
Definition: alphaSuSp.H:1
word timeName
Definition: getTimeIndex.H:3
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
const dimensionSet dimEnergy
fvMatrix< scalar > fvScalarMatrix
Definition: fvMatricesFwd.H:45
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:53
messageStream Info
Information stream (stdout output on master, null elsewhere)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:33
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
const dimensionSet dimVolume(pow3(dimLength))
Definition: dimensionSets.H:60
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:51
propsDict readIfPresent("fields", acceptFields)