chemistryModel.H
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-2024 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 Class
25  Foam::chemistryModel
26 
27 Description
28  Extends base chemistry model by adding a thermo package, and ODE functions.
29  Introduces chemistry equation system and evaluation of chemical source terms
30  with optional support for TDAC mechanism reduction and tabulation.
31 
32  References:
33  \verbatim
34  Contino, F., Jeanmart, H., Lucchini, T., & D’Errico, G. (2011).
35  Coupling of in situ adaptive tabulation and dynamic adaptive chemistry:
36  An effective method for solving combustion in engine simulations.
37  Proceedings of the Combustion Institute, 33(2), 3057-3064.
38 
39  Contino, F., Lucchini, T., D'Errico, G., Duynslaegher, C.,
40  Dias, V., & Jeanmart, H. (2012).
41  Simulations of advanced combustion modes using detailed chemistry
42  combined with tabulation and mechanism reduction techniques.
43  SAE International Journal of Engines,
44  5(2012-01-0145), 185-196.
45 
46  Contino, F., Foucher, F., Dagaut, P., Lucchini, T., D’Errico, G., &
47  Mounaïm-Rousselle, C. (2013).
48  Experimental and numerical analysis of nitric oxide effect on the
49  ignition of iso-octane in a single cylinder HCCI engine.
50  Combustion and Flame, 160(8), 1476-1483.
51 
52  Contino, F., Masurier, J. B., Foucher, F., Lucchini, T., D’Errico, G., &
53  Dagaut, P. (2014).
54  CFD simulations using the TDAC method to model iso-octane combustion
55  for a large range of ozone seeding and temperature conditions
56  in a single cylinder HCCI engine.
57  Fuel, 137, 179-184.
58  \endverbatim
59 
60 SourceFiles
61  chemistryModelI.H
62  chemistryModel.C
63 
64 \*---------------------------------------------------------------------------*/
65 
66 #ifndef chemistryModel_H
67 #define chemistryModel_H
68 
69 #include "odeChemistryModel.H"
70 #include "ReactionList.H"
71 #include "ODESystem.H"
72 #include "volFields.H"
73 #include "multicomponentMixture.H"
76 #include "DynamicField.H"
77 
78 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
79 
80 namespace Foam
81 {
82 
83 /*---------------------------------------------------------------------------*\
84  Class chemistryModel Declaration
85 \*---------------------------------------------------------------------------*/
86 
87 template<class ThermoType>
88 class chemistryModel
89 :
90  public odeChemistryModel
91 {
92  // Private classes
93 
94  //- Class to define scope of reaction evaluation. Runs pre-evaluate
95  // hook on all reactions on construction and post-evaluate on
96  // destruction.
97  class reactionEvaluationScope
98  {
100 
101  public:
102 
103  reactionEvaluationScope
104  (
106  )
107  :
109  {
110  forAll(chemistry_.reactions_, i)
111  {
112  chemistry_.reactions_[i].preEvaluate();
113  }
114  }
115 
116  ~reactionEvaluationScope()
117  {
118  forAll(chemistry_.reactions_, i)
119  {
120  chemistry_.reactions_[i].postEvaluate();
121  }
122  }
123  };
124 
125 
126  // Private data
127 
128  //- Switch to select performance logging
129  Switch log_;
130 
131  //- Switch to enable per-cell CPU load caching for load-balancing
132  Switch cpuLoad_;
133 
134  //- Type of the Jacobian to be calculated
135  const jacobianType jacobianType_;
136 
137  //- Reference to the multi component mixture
138  const multicomponentMixture<ThermoType>& mixture_;
139 
140  //- Thermodynamic data of the species
141  const PtrList<ThermoType>& specieThermos_;
142 
143  //- Reactions
144  const ReactionList<ThermoType> reactions_;
145 
146  //- List of reaction rate per specie [kg/m^3/s]
148 
149  //- Temporary mass fraction field
150  mutable scalarField Y_;
151 
152  //- Temporary simplified mechanism mass fraction field
154 
155  //- Temporary concentration field
156  mutable scalarField c_;
157 
158  //- Temporary simplified mechanism concentration field
160 
161  //- Specie-temperature-pressure workspace fields
162  mutable FixedList<scalarField, 5> YTpWork_;
163 
164  //- Specie-temperature-pressure workspace matrices
165  mutable FixedList<scalarSquareMatrix, 2> YTpYTpWork_;
166 
167  //- Mechanism reduction method
169 
170  //- Mechanism reduction method reference
172 
173  //- Tabulation method
174  autoPtr<chemistryTabulationMethod> tabulationPtr_;
175 
176  //- Tabulation method reference
177  chemistryTabulationMethod& tabulation_;
178 
179  //- Log file for average time spent solving the chemistry
180  autoPtr<OFstream> cpuSolveFile_;
181 
182 
183  // Private Member Functions
184 
185  //- Solve the reaction system for the given time step
186  // of given type and return the characteristic time
187  // Variable number of species added
188  template<class DeltaTType>
189  scalar solve(const DeltaTType& deltaT);
190 
191 
192 public:
193 
194  //- Runtime type information
195  TypeName("chemistryModel");
196 
197 
198  // Constructors
199 
200  //- Construct from thermo
202 
203  //- Disallow default bitwise copy construction
204  chemistryModel(const chemistryModel&) = delete;
205 
206 
207  //- Destructor
208  virtual ~chemistryModel();
209 
210 
211  // Member Functions
212 
213  // Access
214 
215  //- Return reference to the mixture
216  inline const multicomponentMixture<ThermoType>& mixture() const;
217 
218  //- The reactions
219  inline const PtrList<Reaction<ThermoType>>& reactions() const;
220 
221  //- Thermodynamic data of the species
222  inline const PtrList<ThermoType>& specieThermos() const;
223 
224 
225  // Overrides to basicChemistryModel functions
226 
227  //- The number of reactions
228  virtual inline label nReaction() const;
229 
230  //- Return reaction rates of the species [kg/m^3/s]
231  virtual inline const PtrList<volScalarField::Internal>& RR() const;
232 
233  //- Return reaction rates of the species in reactioni [kg/m^3/s]
235  (
236  const label reactioni
237  ) const;
238 
239  //- Calculates the reaction rates
240  virtual void calculate();
241 
242  //- Solve the reaction system for the given time step
243  // and return the characteristic time
244  virtual scalar solve(const scalar deltaT);
245 
246  //- Solve the reaction system for the given time step
247  // and return the characteristic time
248  virtual scalar solve(const scalarField& deltaT);
249 
250  //- Return the chemical time scale
251  virtual tmp<volScalarField> tc() const;
252 
253  //- Return the heat release rate [kg/m/s^3]
254  virtual tmp<volScalarField> Qdot() const;
255 
256 
257  // Overrides to ODESystem
258 
259  //- Calculate the ODE derivatives
260  virtual void derivatives
261  (
262  const scalar t,
263  const scalarField& YTp,
264  const label li,
265  scalarField& dYTpdt
266  ) const;
267 
268  //- Calculate the ODE jacobian
269  virtual void jacobian
270  (
271  const scalar t,
272  const scalarField& YTp,
273  const label li,
274  scalarField& dYTpdt,
276  ) const;
277 
278 
279  // ODE solution functions
280 
281  //- Solve the ODE system
282  virtual void solve
283  (
284  scalar& p,
285  scalar& T,
286  scalarField& Y,
287  const label li,
288  scalar& deltaT,
289  scalar& subDeltaT
290  ) const = 0;
291 
292 
293  // Member Operators
294 
295  //- Disallow default bitwise assignment
296  void operator=(const chemistryModel&) = delete;
297 };
298 
299 
300 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
301 
302 } // End namespace Foam
303 
304 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
305 
306 #include "chemistryModelI.H"
307 
308 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
309 
310 #ifdef NoRepository
311  #include "chemistryModel.C"
312 #endif
313 
314 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
315 
316 #endif
317 
318 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
A 1D vector of objects of type <T> with a fixed size <Size>.
Definition: FixedList.H:78
List of templated reactions.
Definition: ReactionList.H:54
A simple wrapper around bool so that it can be read as a word: true/false, on/off,...
Definition: Switch.H:61
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: autoPtr.H:51
const fluidMulticomponentThermo & thermo() const
Return const access to the thermo.
jacobianType
Enumeration for the type of Jacobian to be calculated by the.
Switch chemistry_
Chemistry activation switch.
Switch chemistry() const
Chemistry activation switch.
Extends base chemistry model by adding a thermo package, and ODE functions. Introduces chemistry equa...
TypeName("chemistryModel")
Runtime type information.
virtual void jacobian(const scalar t, const scalarField &YTp, const label li, scalarField &dYTpdt, scalarSquareMatrix &J) const
Calculate the ODE jacobian.
virtual tmp< volScalarField > tc() const
Return the chemical time scale.
const PtrList< ThermoType > & specieThermos() const
Thermodynamic data of the species.
virtual void derivatives(const scalar t, const scalarField &YTp, const label li, scalarField &dYTpdt) const
Calculate the ODE derivatives.
void operator=(const chemistryModel &)=delete
Disallow default bitwise assignment.
virtual label nReaction() const
The number of reactions.
const multicomponentMixture< ThermoType > & mixture() const
Return reference to the mixture.
virtual ~chemistryModel()
Destructor.
virtual PtrList< volScalarField::Internal > reactionRR(const label reactioni) const
Return reaction rates of the species in reactioni [kg/m^3/s].
const PtrList< Reaction< ThermoType > > & reactions() const
The reactions.
virtual tmp< volScalarField > Qdot() const
Return the heat release rate [kg/m/s^3].
virtual void calculate()
Calculates the reaction rates.
virtual const PtrList< volScalarField::Internal > & RR() const
Return reaction rates of the species [kg/m^3/s].
An abstract class for methods of chemical mechanism reduction.
An abstract class for chemistry tabulation.
Base-class for multi-component fluid thermodynamic properties.
Foam::multicomponentMixture.
Extends base chemistry model adding an ODESystem and the reduction maps needed for tabulation.
const PtrList< volScalarField > & Y() const
Return a reference to the list of mass fraction fields.
A class for managing temporary objects.
Definition: tmp.H:55
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
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
volScalarField & p
rhoEqn solve()