basicChemistryModel.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-2023 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::basicChemistryModel
26 
27 Description
28  Base class for chemistry models
29 
30 SourceFiles
31  basicChemistryModelI.H
32  basicChemistryModel.C
33 
34 \*---------------------------------------------------------------------------*/
35 
36 #ifndef basicChemistryModel_H
37 #define basicChemistryModel_H
38 
40 
41 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
42 
43 namespace Foam
44 {
45 
46 // Forward declaration of classes
47 class fvMesh;
48 
49 /*---------------------------------------------------------------------------*\
50  Class basicChemistryModel Declaration
51 \*---------------------------------------------------------------------------*/
52 
54 :
55  public IOdictionary
56 {
57 public:
58 
59  // Public Enumerations
60 
61  //- Enumeration for the type of Jacobian to be calculated by the
62  // chemistry model
63  enum class jacobianType
64  {
65  fast,
66  exact
67  };
68 
69  //- Jacobian type names
71 
72 
73 protected:
74 
75  // Protected data
76 
77  //- Reference to the mesh
78  const fvMesh& mesh_;
79 
80  //- Reference to the thermo
82 
83  //- Chemistry activation switch
85 
86  //- Initial chemical time step
87  const scalar deltaTChemIni_;
88 
89  //- Maximum chemical time step
90  const scalar deltaTChemMax_;
91 
92  //- Latest estimation of integration step
94 
95 
96 public:
97 
98  //- Runtime type information
99  TypeName("basicChemistryModel");
100 
101 
102  //- Declare run-time constructor selection tables
104  (
105  autoPtr,
107  thermo,
109  (thermo)
110  );
111 
112 
113  // Constructors
114 
115  //- Construct from thermo
117 
118  //- Disallow default bitwise copy construction
119  basicChemistryModel(const basicChemistryModel&) = delete;
120 
121 
122  // Selectors
123 
124  //- Select based on fluid reaction thermo
126  (
128  );
129 
130 
131  //- Destructor
132  virtual ~basicChemistryModel();
133 
134 
135  // Member Functions
136 
137  //- Return const access to the mesh
138  inline const fvMesh& mesh() const;
139 
140  //- Return const access to the thermo
141  inline const fluidMulticomponentThermo& thermo() const;
142 
143  //- Chemistry activation switch
144  inline Switch chemistry() const;
145 
146  //- The number of species
147  virtual label nSpecie() const = 0;
148 
149  //- The number of reactions
150  virtual label nReaction() const = 0;
151 
152  //- Return the latest estimation of integration step
153  inline const volScalarField::Internal& deltaTChem() const;
154 
155  //- Return reaction rates of the species [kg/m^3/s]
156  virtual const PtrList<volScalarField::Internal>& RR() const = 0;
157 
158  //- Return reaction rates of the species in reactioni [kg/m^3/s]
160  (
161  const label reactioni
162  ) const = 0;
163 
164  //- Calculates the reaction rates
165  virtual void calculate() = 0;
166 
167  //- Solve the reaction system for the given time step
168  // and return the characteristic time
169  virtual scalar solve(const scalar deltaT) = 0;
170 
171  //- Solve the reaction system for the given time step
172  // and return the characteristic time
173  virtual scalar solve(const scalarField& deltaT) = 0;
174 
175  //- Return the chemical time scale
176  virtual tmp<volScalarField> tc() const = 0;
177 
178  //- Return the heat release rate [kg/m/s^3]
179  virtual tmp<volScalarField> Qdot() const = 0;
180 
181 
182  // Member Operators
183 
184  //- Disallow default bitwise assignment
185  void operator=(const basicChemistryModel&) = delete;
186 };
187 
188 
189 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
190 
191 } // End namespace Foam
192 
193 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
194 
195 #include "basicChemistryModelI.H"
196 
197 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
198 
199 #endif
200 
201 // ************************************************************************* //
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:57
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: PtrList.H:75
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
Base class for chemistry models.
virtual tmp< volScalarField > tc() const =0
Return the chemical time scale.
const fluidMulticomponentThermo & thermo() const
Return const access to the thermo.
virtual tmp< volScalarField > Qdot() const =0
Return the heat release rate [kg/m/s^3].
virtual ~basicChemistryModel()
Destructor.
const fvMesh & mesh_
Reference to the mesh.
virtual void calculate()=0
Calculates the reaction rates.
basicChemistryModel(const fluidMulticomponentThermo &thermo)
Construct from thermo.
jacobianType
Enumeration for the type of Jacobian to be calculated by the.
virtual label nSpecie() const =0
The number of species.
virtual PtrList< volScalarField::Internal > reactionRR(const label reactioni) const =0
Return reaction rates of the species in reactioni [kg/m^3/s].
const fluidMulticomponentThermo & thermo_
Reference to the thermo.
static autoPtr< basicChemistryModel > New(const fluidMulticomponentThermo &thermo)
Select based on fluid reaction thermo.
virtual scalar solve(const scalar deltaT)=0
Solve the reaction system for the given time step.
const scalar deltaTChemIni_
Initial chemical time step.
Switch chemistry_
Chemistry activation switch.
const volScalarField::Internal & deltaTChem() const
Return the latest estimation of integration step.
TypeName("basicChemistryModel")
Runtime type information.
virtual const PtrList< volScalarField::Internal > & RR() const =0
Return reaction rates of the species [kg/m^3/s].
void operator=(const basicChemistryModel &)=delete
Disallow default bitwise assignment.
virtual label nReaction() const =0
The number of reactions.
const fvMesh & mesh() const
Return const access to the mesh.
declareRunTimeSelectionTable(autoPtr, basicChemistryModel, thermo,(const fluidMulticomponentThermo &thermo),(thermo))
Declare run-time constructor selection tables.
const scalar deltaTChemMax_
Maximum chemical time step.
Switch chemistry() const
Chemistry activation switch.
volScalarField::Internal deltaTChem_
Latest estimation of integration step.
static const NamedEnum< jacobianType, 2 > jacobianTypeNames_
Jacobian type names.
Base-class for multi-component fluid thermodynamic properties.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:99
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