zoneCombustion.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) 2016-2018 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 "zoneCombustion.H"
27 
28 // * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
29 
30 template<class ReactionThermo>
33 (
34  const tmp<fvScalarMatrix>& tR
35 ) const
36 {
37  fvScalarMatrix& R = tR.ref();
38  scalarField& Su = R.source();
39  scalarField filteredField(Su.size(), 0);
40 
41  forAll(zoneNames_, zonei)
42  {
43  const labelList& cells = this->mesh().cellZones()[zoneNames_[zonei]];
44 
45  forAll(cells, i)
46  {
47  filteredField[cells[i]] = Su[cells[i]];
48  }
49  }
50 
51  Su = filteredField;
52 
53  if (R.hasDiag())
54  {
55  scalarField& Sp = R.diag();
56 
57  forAll(zoneNames_, zonei)
58  {
59  const labelList& cells =
60  this->mesh().cellZones()[zoneNames_[zonei]];
61 
62  forAll(cells, i)
63  {
64  filteredField[cells[i]] = Sp[cells[i]];
65  }
66  }
67 
68  Sp = filteredField;
69  }
70 
71  return tR;
72 }
73 
74 
75 template<class ReactionThermo>
78 (
79  const tmp<volScalarField>& tS
80 ) const
81 {
82  scalarField& S = tS.ref();
83  scalarField filteredField(S.size(), 0);
84 
85  forAll(zoneNames_, zonei)
86  {
87  const labelList& cells = this->mesh().cellZones()[zoneNames_[zonei]];
88 
89  forAll(cells, i)
90  {
91  filteredField[cells[i]] = S[cells[i]];
92  }
93  }
94 
95  S = filteredField;
96 
97  return tS;
98 }
99 
100 
101 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
102 
103 template<class ReactionThermo>
105 (
106  const word& modelType,
107  ReactionThermo& thermo,
108  const compressibleTurbulenceModel& turb,
109  const word& combustionProperties
110 )
111 :
113  (
114  modelType,
115  thermo,
116  turb,
117  combustionProperties
118  ),
119  combustionModelPtr_
120  (
122  (
123  thermo,
124  turb,
125  "zoneCombustionProperties"
126  )
127  ),
128  zoneNames_(this->coeffs().lookup("zones"))
129 {}
130 
131 
132 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
133 
134 template<class ReactionThermo>
136 {}
137 
138 
139 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
140 
141 template<class ReactionThermo>
143 {
144  return combustionModelPtr_->thermo();
145 }
146 
147 
148 template<class ReactionThermo>
149 const ReactionThermo&
151 {
152  return combustionModelPtr_->thermo();
153 }
154 
155 
156 template<class ReactionThermo>
158 {
159  combustionModelPtr_->correct();
160 }
161 
162 
163 template<class ReactionThermo>
166 (
167  volScalarField& Y
168 ) const
169 {
170  return filter(combustionModelPtr_->R(Y));
171 }
172 
173 
174 template<class ReactionThermo>
177 {
178  return filter(combustionModelPtr_->Qdot());
179 }
180 
181 
182 template<class ReactionThermo>
184 {
186  {
187  combustionModelPtr_->read();
188  return true;
189  }
190  else
191  {
192  return false;
193  }
194 }
195 
196 
197 // ************************************************************************* //
virtual tmp< fvScalarMatrix > R(volScalarField &Y) const
Fuel consumption rate matrix.
fvMatrix< scalar > fvScalarMatrix
Definition: fvMatricesFwd.H:42
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
virtual ReactionThermo & thermo()
Return access to the thermo package.
virtual bool read()
Update properties from given dictionary.
rhoReactionThermo & thermo
Definition: createFields.H:28
virtual tmp< volScalarField > Qdot() const
Heat release rate [kg/m/s3].
dynamicFvMesh & mesh
A class for handling words, derived from string.
Definition: word.H:59
List< label > labelList
A List of labels.
Definition: labelList.H:56
volScalarField scalarField(fieldObject, mesh)
Zone-filtered combustion model.
Abstract base class for turbulence models (RAS, LES and laminar).
Combustion models for templated thermodynamics.
A class for managing temporary objects.
Definition: PtrList.H:53
virtual void correct()
Correct combustion rate.