All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
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-2020 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 
39 #include "fluidReactionThermo.H"
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 protected:
58 
59  // Protected data
60 
61  //- Reference to the mesh database
62  const fvMesh& mesh_;
63 
64  //- Reference to the thermo
66 
67  //- Chemistry activation switch
69 
70  //- Initial chemical time step
71  const scalar deltaTChemIni_;
72 
73  //- Maximum chemical time step
74  const scalar deltaTChemMax_;
75 
76  //- Latest estimation of integration step
78 
79 
80  // Protected Member Functions
81 
82  //- Correct function - updates due to mesh changes
83  void correct();
84 
85 
86 public:
87 
88  //- Runtime type information
89  TypeName("basicChemistryModel");
90 
91 
92  //- Declare run-time constructor selection tables
94  (
95  autoPtr,
97  thermo,
98  (const fluidReactionThermo& thermo),
99  (thermo)
100  );
101 
102 
103  // Constructors
104 
105  //- Construct from thermo
107 
108  //- Disallow default bitwise copy construction
109  basicChemistryModel(const basicChemistryModel&) = delete;
110 
111 
112  // Selectors
113 
114  //- Select based on fluid reaction thermo
116  (
117  const fluidReactionThermo& thermo
118  );
119 
120 
121  //- Destructor
122  virtual ~basicChemistryModel();
123 
124 
125  // Member Functions
126 
127  //- Return const access to the mesh database
128  inline const fvMesh& mesh() const;
129 
130  //- Return const access to the thermo
131  inline const fluidReactionThermo& thermo() const;
132 
133  //- Chemistry activation switch
134  inline Switch chemistry() const;
135 
136  //- The number of species
137  virtual label nSpecie() const = 0;
138 
139  //- The number of reactions
140  virtual label nReaction() const = 0;
141 
142  //- Return the latest estimation of integration step
143  inline const volScalarField::Internal& deltaTChem() const;
144 
145 
146  // Functions to be derived in derived classes
147 
148  // Fields
149 
150  //- Return const access to chemical source terms [kg/m^3/s]
151  virtual const volScalarField::Internal& RR
152  (
153  const label i
154  ) const = 0;
155 
156  //- Return access to chemical source terms [kg/m^3/s]
158  (
159  const label i
160  ) = 0;
161 
162  //- Return reaction rate of the speciei in reactioni
164  (
165  const label reactioni,
166  const label speciei
167  ) const = 0;
168 
169 
170  // Chemistry solution
171 
172  //- Calculates the reaction rates
173  virtual void calculate() = 0;
174 
175  //- Solve the reaction system for the given time step
176  // and return the characteristic time
177  virtual scalar solve(const scalar deltaT) = 0;
178 
179  //- Solve the reaction system for the given time step
180  // and return the characteristic time
181  virtual scalar solve(const scalarField& deltaT) = 0;
182 
183  //- Return the chemical time scale
184  virtual tmp<volScalarField> tc() const = 0;
185 
186  //- Return the heat release rate [kg/m/s^3]
187  virtual tmp<volScalarField> Qdot() const = 0;
188 
189 
190  // Member Operators
191 
192  //- Disallow default bitwise assignment
193  void operator=(const basicChemistryModel&) = delete;
194 };
195 
196 
197 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
198 
199 } // End namespace Foam
200 
201 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
202 
203 #include "basicChemistryModelI.H"
204 
205 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
206 
207 #endif
208 
209 // ************************************************************************* //
const fvMesh & mesh() const
Return const access to the mesh database.
virtual label nSpecie() const =0
The number of species.
const scalar deltaTChemMax_
Maximum chemical time step.
volScalarField::Internal deltaTChem_
Latest estimation of integration step.
declareRunTimeSelectionTable(autoPtr, basicChemistryModel, thermo,(const fluidReactionThermo &thermo),(thermo))
Declare run-time constructor selection tables.
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
const fvMesh & mesh_
Reference to the mesh database.
A simple wrapper around bool so that it can be read as a word: true/false, on/off, yes/no, y/n, t/f, or none/any.
Definition: Switch.H:60
Base-class for multi-component fluid thermodynamic properties.
const fluidReactionThermo & thermo_
Reference to the thermo.
Switch chemistry() const
Chemistry activation switch.
void correct()
Correct function - updates due to mesh changes.
void operator=(const basicChemistryModel &)=delete
Disallow default bitwise assignment.
virtual tmp< volScalarField > tc() const =0
Return the chemical time scale.
basicChemistryModel(const fluidReactionThermo &thermo)
Construct from thermo.
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:53
virtual ~basicChemistryModel()
Destructor.
TypeName("basicChemistryModel")
Runtime type information.
static autoPtr< basicChemistryModel > New(const fluidReactionThermo &thermo)
Select based on fluid reaction thermo.
Switch chemistry_
Chemistry activation switch.
virtual void calculate()=0
Calculates the reaction rates.
const fluidReactionThermo & thermo() const
Return const access to the thermo.
const scalar deltaTChemIni_
Initial chemical time step.
virtual label nReaction() const =0
The number of reactions.
virtual const volScalarField::Internal & RR(const label i) const =0
Return const access to chemical source terms [kg/m^3/s].
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
virtual tmp< volScalarField::Internal > calculateRR(const label reactioni, const label speciei) const =0
Return reaction rate of the speciei in reactioni.
const volScalarField::Internal & deltaTChem() const
Return the latest estimation of integration step.
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
A class for managing temporary objects.
Definition: PtrList.H:53
virtual tmp< volScalarField > Qdot() const =0
Return the heat release rate [kg/m/s^3].
virtual scalar solve(const scalar deltaT)=0
Solve the reaction system for the given time step.
Namespace for OpenFOAM.
Base class for chemistry models.