solidChemistryModel.H
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2013-2016 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::solidChemistryModel
26 
27 Description
28  Extends base solid chemistry model by adding a thermo package, and ODE
29  functions.
30  Introduces chemistry equation system and evaluation of chemical source
31  terms.
32 
33 SourceFiles
34  solidChemistryModelI.H
35  solidChemistryModel.C
36 
37 \*---------------------------------------------------------------------------*/
38 
39 #ifndef solidChemistryModel_H
40 #define solidChemistryModel_H
41 
42 #include "Reaction.H"
43 #include "ODESystem.H"
44 #include "volFields.H"
45 #include "simpleMatrix.H"
46 
47 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
48 
49 namespace Foam
50 {
51 
52 // Forward declaration of classes
53 class fvMesh;
54 
55 /*---------------------------------------------------------------------------*\
56  Class solidChemistryModel Declaration
57 \*---------------------------------------------------------------------------*/
58 
59 template<class CompType, class SolidThermo>
61 :
62  public CompType,
63  public ODESystem
64 {
65  // Private Member Functions
66 
67  //- Disallow copy constructor
69 
70  //- Disallow default bitwise assignment
71  void operator=(const solidChemistryModel&);
72 
73 
74 protected:
75 
76  //- Reference to solid mass fractions
78 
79  //- Reactions
81 
82  //- Thermodynamic data of solids
84 
85  //- Number of solid components
87 
88  //- Number of solid reactions
90 
91  //- List of reaction rate per solid [kg/m3/s]
93 
94  //- List of active reacting cells
96 
97 
98  // Protected Member Functions
99 
100  //- Write access to source terms for solids
102 
103  //- Set reacting status of cell, celli
104  void setCellReacting(const label celli, const bool active);
105 
106 
107 public:
108 
109  //- Runtime type information
110  TypeName("solidChemistryModel");
111 
112 
113  // Constructors
114 
115  //- Construct from mesh and phase name
116  solidChemistryModel(const fvMesh& mesh, const word& phaseName);
117 
118 
119  //- Destructor
120  virtual ~solidChemistryModel();
121 
122 
123  // Member Functions
124 
125  //- The reactions
126  inline const PtrList<Reaction<SolidThermo>>& reactions() const;
127 
128  //- The number of reactions
129  inline label nReaction() const;
130 
131 
132  //- dc/dt = omega, rate of change in concentration, for each species
133  virtual scalarField omega
134  (
135  const scalarField& c,
136  const scalar T,
137  const scalar p,
138  const bool updateC0 = false
139  ) const = 0;
140 
141  //- Return the reaction rate for reaction r and the reference
142  // species and charateristic times
143  virtual scalar omega
144  (
145  const Reaction<SolidThermo>& r,
146  const scalarField& c,
147  const scalar T,
148  const scalar p,
149  scalar& pf,
150  scalar& cf,
151  label& lRef,
152  scalar& pr,
153  scalar& cr,
154  label& rRef
155  ) const = 0;
156 
157 
158  //- Return the reaction rate for iReaction and the reference
159  // species and charateristic times
160  virtual scalar omegaI
161  (
162  label iReaction,
163  const scalarField& c,
164  const scalar T,
165  const scalar p,
166  scalar& pf,
167  scalar& cf,
168  label& lRef,
169  scalar& pr,
170  scalar& cr,
171  label& rRef
172  ) const = 0;
173 
174  //- Calculates the reaction rates
175  virtual void calculate() = 0;
176 
177 
178  // Solid Chemistry model functions
179 
180  //- Return const access to the chemical source terms for solids
181  inline const volScalarField::Internal& RRs
182  (
183  const label i
184  ) const;
185 
186  //- Return total solid source term
187  inline tmp<volScalarField::Internal> RRs() const;
188 
189  //- Solve the reaction system for the given time step
190  // and return the characteristic time
191  virtual scalar solve(const scalar deltaT) = 0;
192 
193  //- Solve the reaction system for the given time step
194  // and return the characteristic time
195  virtual scalar solve(const scalarField& deltaT);
196 
197  //- Return the chemical time scale
198  virtual tmp<volScalarField> tc() const;
199 
200  //- Return the heat release rate [kg/m/s3]
201  virtual tmp<volScalarField> Qdot() const;
202 
203 
204  // ODE functions (overriding abstract functions in ODE.H)
205 
206  //- Number of ODE's to solve
207  virtual label nEqns() const = 0;
208 
209  virtual void derivatives
210  (
211  const scalar t,
212  const scalarField& c,
213  scalarField& dcdt
214  ) const = 0;
215 
216  virtual void jacobian
217  (
218  const scalar t,
219  const scalarField& c,
220  scalarField& dcdt,
221  scalarSquareMatrix& dfdc
222  ) const = 0;
223 
224  virtual void solve
225  (
226  scalarField &c,
227  scalar& T,
228  scalar& p,
229  scalar& deltaT,
230  scalar& subDeltaT
231  ) const = 0;
232 };
233 
234 
235 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
236 
237 } // End namespace Foam
238 
239 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
240 
241  #include "solidChemistryModelI.H"
242 
243 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
244 
245 #ifdef NoRepository
246  #include "solidChemistryModel.C"
247 #endif
248 
249 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
250 
251 #endif
252 
253 // ************************************************************************* //
List< bool > reactingCells_
List of active reacting cells.
virtual scalar omegaI(label iReaction, const scalarField &c, const scalar T, const scalar p, scalar &pf, scalar &cf, label &lRef, scalar &pr, scalar &cr, label &rRef) const =0
Return the reaction rate for iReaction and the reference.
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
Abstract base class for the systems of ordinary differential equations.
Definition: ODESystem.H:46
void setCellReacting(const label celli, const bool active)
Set reacting status of cell, celli.
PtrList< volScalarField::Internal > & RRs()
Write access to source terms for solids.
const PtrList< Reaction< SolidThermo > > & reactions() const
The reactions.
label nReaction_
Number of solid reactions.
Extends base solid chemistry model by adding a thermo package, and ODE functions. Introduces chemistr...
Simple extension of ReactionThermo to handle reaction kinetics in addition to the equilibrium thermod...
Definition: Reaction.H:53
virtual ~solidChemistryModel()
Destructor.
PtrList< volScalarField::Internal > RRs_
List of reaction rate per solid [kg/m3/s].
virtual tmp< volScalarField > tc() const
Return the chemical time scale.
TypeName("solidChemistryModel")
Runtime type information.
dynamicFvMesh & mesh
const PtrList< SolidThermo > & solidThermo_
Thermodynamic data of solids.
A class for handling words, derived from string.
Definition: word.H:59
virtual scalar solve(const scalar deltaT)=0
Solve the reaction system for the given time step.
label nSolids_
Number of solid components.
virtual scalarField omega(const scalarField &c, const scalar T, const scalar p, const bool updateC0=false) const =0
dc/dt = omega, rate of change in concentration, for each species
virtual void jacobian(const scalar t, const scalarField &c, scalarField &dcdt, scalarSquareMatrix &dfdc) const =0
Calculate the Jacobian of the system.
const PtrList< Reaction< SolidThermo > > & reactions_
Reactions.
virtual tmp< volScalarField > Qdot() const
Return the heat release rate [kg/m/s3].
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
PtrList< volScalarField > & Ys_
Reference to solid 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:63
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
const dimensionedScalar c
Speed of light in a vacuum.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
label nReaction() const
The number of reactions.
virtual void calculate()=0
Calculates the reaction rates.
virtual void derivatives(const scalar t, const scalarField &c, scalarField &dcdt) const =0
Calculate the derivatives in dydx.
volScalarField & p
A class for managing temporary objects.
Definition: PtrList.H:53
virtual label nEqns() const =0
Number of ODE&#39;s to solve.
Namespace for OpenFOAM.