All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
standardChemistryModel.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) 2011-2021 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::standardChemistryModel
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
30  terms.
31 
32 SourceFiles
33  standardChemistryModelI.H
34  standardChemistryModel.C
35 
36 \*---------------------------------------------------------------------------*/
37 
38 #ifndef standardChemistryModel_H
39 #define standardChemistryModel_H
40 
41 #include "basicChemistryModel.H"
42 #include "ReactionList.H"
43 #include "ODESystem.H"
44 #include "volFields.H"
45 #include "multiComponentMixture.H"
46 
47 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
48 
49 namespace Foam
50 {
51 
52 // Forward declaration of classes
53 class fvMesh;
54 
55 /*---------------------------------------------------------------------------*\
56  Class standardChemistryModel Declaration
57 \*---------------------------------------------------------------------------*/
58 
59 template<class ThermoType>
61 :
62  public basicChemistryModel,
63  public ODESystem
64 {
65  // Private Member Functions
66 
67  //- Solve the reaction system for the given time step
68  // of given type and return the characteristic time
69  template<class DeltaTType>
70  scalar solve(const DeltaTType& deltaT);
71 
72 
73 protected:
74 
75  // Protected classes
76 
77  //- Class to define scope of reaction evaluation. Runs pre-evaluate
78  // hook on all reactions on construction and post-evaluate on
79  // destruction.
81  {
82  const standardChemistryModel<ThermoType>& chemistry_;
83 
84  public:
85 
87  (
89  )
90  :
91  chemistry_(chemistry)
92  {
93  forAll(chemistry_.reactions_, i)
94  {
95  chemistry_.reactions_[i].preEvaluate();
96  }
97  }
98 
100  {
101  forAll(chemistry_.reactions_, i)
102  {
103  chemistry_.reactions_[i].postEvaluate();
104  }
105  }
106  };
107 
108 
109  // Protected data
110 
111  //- Reference to the field of specie mass fractions
113 
114  //- Reference to the multi component mixture
116 
117  //- Thermodynamic data of the species
119 
120  //- Reactions
122 
123  //- Number of species
124  label nSpecie_;
125 
126  //- Number of reactions
128 
129  //- Temperature below which the reaction rates are assumed 0
130  scalar Treact_;
131 
132  //- List of reaction rate per specie [kg/m^3/s]
134 
135  //- Temporary concentration field
136  mutable scalarField c_;
137 
138  //- Temporary rate-of-change of concentration field
139  mutable scalarField dcdt_;
140 
141 
142  // Protected Member Functions
143 
144  //- Write access to chemical source terms
145  // (e.g. for multi-chemistry model)
147 
148 
149 public:
150 
151  //- Runtime type information
152  TypeName("standard");
153 
154 
155  // Constructors
156 
157  //- Construct from thermo
159 
160  //- Disallow default bitwise copy construction
162 
163 
164  //- Destructor
165  virtual ~standardChemistryModel();
166 
167 
168  // Member Functions
169 
170  //- Return reference to the mixture
171  inline const multiComponentMixture<ThermoType>& mixture() const;
172 
173  //- The reactions
174  inline const PtrList<Reaction<ThermoType>>& reactions() const;
175 
176  //- Thermodynamic data of the species
177  inline const PtrList<ThermoType>& specieThermos() const;
178 
179  //- The number of species
180  virtual inline label nSpecie() const;
181 
182  //- The number of reactions
183  virtual inline label nReaction() const;
184 
185  //- Temperature below which the reaction rates are assumed 0
186  inline scalar Treact() const;
187 
188  //- Temperature below which the reaction rates are assumed 0
189  inline scalar& Treact();
190 
191  //- dc/dt = omega, rate of change in concentration, for each species
192  virtual void omega
193  (
194  const scalar p,
195  const scalar T,
196  const scalarField& c,
197  const label li,
198  scalarField& dcdt
199  ) const;
200 
201  //- Calculates the reaction rates
202  virtual void calculate();
203 
204 
205  // Chemistry model functions (overriding abstract functions in
206  // basicChemistryModel.H)
207 
208  //- Return const access to the chemical source terms for specie, i
209  inline const volScalarField::Internal& RR
210  (
211  const label i
212  ) const;
213 
214  //- Return non const access to chemical source terms [kg/m^3/s]
216  (
217  const label i
218  );
219 
220  //- Return reaction rate of the speciei in reactionI
222  (
223  const label reactionI,
224  const label speciei
225  ) const;
226 
227  //- Solve the reaction system for the given time step
228  // and return the characteristic time
229  virtual scalar solve(const scalar deltaT);
230 
231  //- Solve the reaction system for the given time step
232  // and return the characteristic time
233  virtual scalar solve(const scalarField& deltaT);
234 
235  //- Return the chemical time scale
236  virtual tmp<volScalarField> tc() const;
237 
238  //- Return the heat release rate [kg/m/s^3]
239  virtual tmp<volScalarField> Qdot() const;
240 
241 
242  // ODE functions (overriding abstract functions in ODE.H)
243 
244  //- Number of ODE's to solve
245  inline virtual label nEqns() const;
246 
247  virtual void derivatives
248  (
249  const scalar t,
250  const scalarField& c,
251  const label li,
252  scalarField& dcdt
253  ) const;
254 
255  virtual void jacobian
256  (
257  const scalar t,
258  const scalarField& c,
259  const label li,
260  scalarField& dcdt,
262  ) const;
263 
264  virtual void solve
265  (
266  scalar& p,
267  scalar& T,
268  scalarField& c,
269  const label li,
270  scalar& deltaT,
271  scalar& subDeltaT
272  ) const = 0;
273 
274 
275  // Member Operators
276 
277  //- Disallow default bitwise assignment
278  void operator=(const standardChemistryModel&) = delete;
279 };
280 
281 
282 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
283 
284 } // End namespace Foam
285 
286 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
287 
288 #include "standardChemistryModelI.H"
289 
290 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
291 
292 #ifdef NoRepository
293  #include "standardChemistryModel.C"
294 #endif
295 
296 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
297 
298 #endif
299 
300 // ************************************************************************* //
virtual label nSpecie() const
The number of species.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
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
virtual tmp< volScalarField > Qdot() const
Return the heat release rate [kg/m/s^3].
Abstract base class for the systems of ordinary differential equations.
Definition: ODESystem.H:46
scalarField dcdt_
Temporary rate-of-change of concentration field.
virtual void jacobian(const scalar t, const scalarField &c, const label li, scalarField &dcdt, scalarSquareMatrix &J) const
Calculate the Jacobian of the system.
virtual void omega(const scalar p, const scalar T, const scalarField &c, const label li, scalarField &dcdt) const
dc/dt = omega, rate of change in concentration, for each species
Base-class for multi-component fluid thermodynamic properties.
Switch chemistry() const
Chemistry activation switch.
PtrList< volScalarField::Internal > RR_
List of reaction rate per specie [kg/m^3/s].
const ReactionList< ThermoType > reactions_
Reactions.
label nSpecie_
Number of species.
const dimensionedScalar c
Speed of light in a vacuum.
void operator=(const standardChemistryModel &)=delete
Disallow default bitwise assignment.
virtual void calculate()
Calculates the reaction rates.
const PtrList< Reaction< ThermoType > > & reactions() const
The reactions.
scalar Treact() const
Temperature below which the reaction rates are assumed 0.
scalar Treact_
Temperature below which the reaction rates are assumed 0.
const PtrList< ThermoType > & specieThermos_
Thermodynamic data of the species.
const PtrList< ThermoType > & specieThermos() const
Thermodynamic data of the species.
Foam::multiComponentMixture.
standardChemistryModel(const fluidReactionThermo &thermo)
Construct from thermo.
scalarField c_
Temporary concentration field.
Extends base chemistry model by adding a thermo package, and ODE functions. Introduces chemistry equa...
label nReaction_
Number of reactions.
virtual void derivatives(const scalar t, const scalarField &c, const label li, scalarField &dcdt) const
Calculate the derivatives in dydx.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
PtrList< volScalarField::Internal > & RR()
Write access to chemical source terms.
List of templated reactions.
Definition: ReactionList.H:51
virtual tmp< volScalarField::Internal > calculateRR(const label reactionI, const label speciei) const
Return reaction rate of the speciei in reactionI.
virtual tmp< volScalarField > tc() const
Return the chemical time scale.
const fluidReactionThermo & thermo() const
Return const access to the thermo.
const PtrList< volScalarField > & Y_
Reference to the field of specie mass fractions.
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:70
virtual label nEqns() const
Number of ODE&#39;s to solve.
Class to define scope of reaction evaluation. Runs pre-evaluate.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const multiComponentMixture< ThermoType > & mixture() const
Return reference to the mixture.
volScalarField & p
A class for managing temporary objects.
Definition: PtrList.H:53
virtual ~standardChemistryModel()
Destructor.
TypeName("standard")
Runtime type information.
reactionEvaluationScope(const standardChemistryModel< ThermoType > &chemistry)
Namespace for OpenFOAM.
const multiComponentMixture< ThermoType > & mixture_
Reference to the multi component mixture.
virtual label nReaction() const
The number of reactions.
Base class for chemistry models.