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 "volFieldsFwd.H"
45 #include "DimensionedField.H"
46 #include "simpleMatrix.H"
47 
48 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
49 
50 namespace Foam
51 {
52 
53 // Forward declaration of classes
54 class fvMesh;
55 
56 /*---------------------------------------------------------------------------*\
57  Class solidChemistryModel Declaration
58 \*---------------------------------------------------------------------------*/
59 
60 template<class CompType, class SolidThermo>
62 :
63  public CompType,
64  public ODESystem
65 {
66  // Private Member Functions
67 
68  //- Disallow copy constructor
70 
71  //- Disallow default bitwise assignment
72  void operator=(const solidChemistryModel&);
73 
74 
75 protected:
76 
77  //- Reference to solid mass fractions
79 
80  //- Reactions
82 
83  //- Thermodynamic data of solids
85 
86  //- Number of solid components
88 
89  //- Number of solid reactions
91 
92  //- List of reaction rate per solid [kg/m3/s]
94 
95  //- List of active reacting cells
97 
98 
99  // Protected Member Functions
100 
101  //- Write access to source terms for solids
103 
104  //- Set reacting status of cell, celli
105  void setCellReacting(const label celli, const bool active);
106 
107 
108 public:
109 
110  //- Runtime type information
111  TypeName("solidChemistryModel");
112 
113 
114  // Constructors
115 
116  //- Construct from mesh and phase name
117  solidChemistryModel(const fvMesh& mesh, const word& phaseName);
118 
119 
120  //- Destructor
121  virtual ~solidChemistryModel();
122 
123 
124  // Member Functions
125 
126  //- The reactions
127  inline const PtrList<Reaction<SolidThermo>>& reactions() const;
128 
129  //- The number of reactions
130  inline label nReaction() const;
131 
132 
133  //- dc/dt = omega, rate of change in concentration, for each species
134  virtual scalarField omega
135  (
136  const scalarField& c,
137  const scalar T,
138  const scalar p,
139  const bool updateC0 = false
140  ) const = 0;
141 
142  //- Return the reaction rate for reaction r and the reference
143  // species and charateristic times
144  virtual scalar omega
145  (
146  const Reaction<SolidThermo>& r,
147  const scalarField& c,
148  const scalar T,
149  const scalar p,
150  scalar& pf,
151  scalar& cf,
152  label& lRef,
153  scalar& pr,
154  scalar& cr,
155  label& rRef
156  ) const = 0;
157 
158 
159  //- Return the reaction rate for iReaction and the reference
160  // species and charateristic times
161  virtual scalar omegaI
162  (
163  label iReaction,
164  const scalarField& c,
165  const scalar T,
166  const scalar p,
167  scalar& pf,
168  scalar& cf,
169  label& lRef,
170  scalar& pr,
171  scalar& cr,
172  label& rRef
173  ) const = 0;
174 
175  //- Calculates the reaction rates
176  virtual void calculate() = 0;
177 
178 
179  // Solid Chemistry model functions
180 
181  //- Return const access to the chemical source terms for solids
183  (
184  const label i
185  ) const;
186 
187  //- Return total solid source term
189 
190  //- Solve the reaction system for the given time step
191  // and return the characteristic time
192  virtual scalar solve(const scalar deltaT) = 0;
193 
194  //- Solve the reaction system for the given time step
195  // and return the characteristic time
196  virtual scalar solve(const scalarField& deltaT);
197 
198  //- Return the chemical time scale
199  virtual tmp<volScalarField> tc() const;
200 
201  //- Return source for enthalpy equation [kg/m/s3]
202  virtual tmp<volScalarField> Sh() const;
203 
204  //- Return the heat release, i.e. enthalpy/sec [m2/s3]
205  virtual tmp<volScalarField> dQ() const;
206 
207 
208  // ODE functions (overriding abstract functions in ODE.H)
209 
210  //- Number of ODE's to solve
211  virtual label nEqns() const = 0;
212 
213  virtual void derivatives
214  (
215  const scalar t,
216  const scalarField& c,
217  scalarField& dcdt
218  ) const = 0;
219 
220  virtual void jacobian
221  (
222  const scalar t,
223  const scalarField& c,
224  scalarField& dcdt,
225  scalarSquareMatrix& dfdc
226  ) const = 0;
227 
228  virtual void solve
229  (
230  scalarField &c,
231  scalar& T,
232  scalar& p,
233  scalar& deltaT,
234  scalar& subDeltaT
235  ) const = 0;
236 };
237 
238 
239 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
240 
241 } // End namespace Foam
242 
243 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
244 
245  #include "solidChemistryModelI.H"
246 
247 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
248 
249 #ifdef NoRepository
250  #include "solidChemistryModel.C"
251 #endif
252 
253 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
254 
255 #endif
256 
257 // ************************************************************************* //
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
virtual tmp< volScalarField > Sh() const
Return source for enthalpy equation [kg/m/s3].
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< DimensionedField< scalar, volMesh > > & RRs()
Write access to source terms for solids.
label nReaction_
Number of solid reactions.
PtrList< DimensionedField< scalar, volMesh > > RRs_
List of reaction rate per solid [kg/m3/s].
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.
const PtrList< Reaction< SolidThermo > > & reactions() const
The reactions.
TypeName("solidChemistryModel")
Runtime type information.
dynamicFvMesh & mesh
const PtrList< SolidThermo > & solidThermo_
Thermodynamic data of solids.
virtual tmp< volScalarField > tc() const
Return the chemical time scale.
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.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
label nReaction() const
The number of reactions.
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:62
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...
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:54
virtual label nEqns() const =0
Number of ODE&#39;s to solve.
Namespace for OpenFOAM.
virtual tmp< volScalarField > dQ() const
Return the heat release, i.e. enthalpy/sec [m2/s3].