zoneCombustion.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) 2016-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 "zoneCombustion.H"
30
31// * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
32
33template<class ReactionThermo>
36(
37 const tmp<fvScalarMatrix>& tR
38) const
39{
40 fvScalarMatrix& R = tR.ref();
41 scalarField& Su = R.source();
42 scalarField filteredField(Su.size(), Zero);
43
44 forAll(zoneNames_, zonei)
45 {
46 const labelList& cells = this->mesh().cellZones()[zoneNames_[zonei]];
47
48 forAll(cells, i)
49 {
50 filteredField[cells[i]] = Su[cells[i]];
51 }
52 }
53
54 Su = filteredField;
55
56 if (R.hasDiag())
57 {
58 scalarField& Sp = R.diag();
59
60 forAll(zoneNames_, zonei)
61 {
62 const labelList& cells =
63 this->mesh().cellZones()[zoneNames_[zonei]];
64
65 forAll(cells, i)
66 {
67 filteredField[cells[i]] = Sp[cells[i]];
68 }
69 }
70
71 Sp = filteredField;
72 }
73
74 return tR;
75}
76
77
78template<class ReactionThermo>
81(
82 const tmp<volScalarField>& tS
83) const
84{
85 scalarField& S = tS.ref();
86 scalarField filteredField(S.size(), Zero);
87
88 forAll(zoneNames_, zonei)
89 {
90 const labelList& cells = this->mesh().cellZones()[zoneNames_[zonei]];
91
92 forAll(cells, i)
93 {
94 filteredField[cells[i]] = S[cells[i]];
95 }
96 }
97
98 S = filteredField;
99
100 return tS;
101}
102
103
104// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
105
106template<class ReactionThermo>
108(
109 const word& modelType,
110 ReactionThermo& thermo,
112 const word& combustionProperties
113)
114:
115 CombustionModel<ReactionThermo>
116 (
117 modelType,
118 thermo,
119 turb,
120 combustionProperties
121 ),
122 combustionModelPtr_
123 (
124 CombustionModel<ReactionThermo>::New
125 (
126 thermo,
127 turb,
128 "zoneCombustionProperties"
129 )
130 ),
131 zoneNames_()
132{
133 this->coeffs().readEntry("zones", zoneNames_);
134}
135
136
137// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
138
139template<class ReactionThermo>
141{}
142
143
144// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
145
146template<class ReactionThermo>
148{
149 return combustionModelPtr_->thermo();
150}
151
152
153template<class ReactionThermo>
154const ReactionThermo&
156{
157 return combustionModelPtr_->thermo();
158}
159
160
161template<class ReactionThermo>
163{
164 combustionModelPtr_->correct();
165}
166
167
168template<class ReactionThermo>
171(
173) const
174{
175 return filter(combustionModelPtr_->R(Y));
176}
177
178
179template<class ReactionThermo>
182{
183 return filter(combustionModelPtr_->Qdot());
184}
185
186
187template<class ReactionThermo>
189{
191 {
192 combustionModelPtr_->read();
193 return true;
194 }
195
196 return false;
197}
198
199
200// ************************************************************************* //
#define R(A, B, C, D, E, F, K, M)
compressible::turbulenceModel & turb
Combustion models for templated thermodynamics.
const dictionary & coeffs() const
Return const dictionary of the model.
Zone-filtered combustion model.
virtual void correct()
Correct combustion rate.
virtual ReactionThermo & thermo()
Return access to the thermo package.
virtual tmp< volScalarField > Qdot() const
Heat release rate [kg/m/s3].
virtual bool read()
Update properties from given dictionary.
Abstract base class for turbulence models (RAS, LES and laminar).
bool readEntry(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX, bool mandatory=true) const
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
A class for handling words, derived from Foam::string.
Definition: word.H:68
PtrList< volScalarField > & Y
dynamicFvMesh & mesh
const cellShapeList & cells
zeroField Su
Definition: alphaSuSp.H:1
zeroField Sp
Definition: alphaSuSp.H:2
List< label > labelList
A List of labels.
Definition: List.H:66
fvMatrix< scalar > fvScalarMatrix
Definition: fvMatricesFwd.H:45
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh > > &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333
Functor wrapper of allow/deny lists of wordRe for filtering.
Definition: wordRes.H:174