pyrolysisChemistryModel.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::pyrolysisChemistryModel
26 
27 Description
28  Pyrolysis chemistry model. It includes gas phase in the solid
29  reaction.
30 
31 SourceFiles
32  pyrolysisChemistryModelI.H
33  pyrolysisChemistryModel.C
34 
35 \*---------------------------------------------------------------------------*/
36 
37 #ifndef pyrolysisChemistryModel_H
38 #define pyrolysisChemistryModel_H
39 
40 #include "volFieldsFwd.H"
41 #include "DimensionedField.H"
42 #include "solidChemistryModel.H"
43 
44 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
45 
46 namespace Foam
47 {
48 
49 // Forward declaration of classes
50 class fvMesh;
51 
52 /*---------------------------------------------------------------------------*\
53  Class pyrolysisChemistryModel Declaration
54 \*---------------------------------------------------------------------------*/
55 
56 template<class CompType, class SolidThermo, class GasThermo>
58 :
59  public solidChemistryModel<CompType, SolidThermo>
60 {
61  // Private Member Functions
62 
63  //- Disallow default bitwise assignment
64  void operator=(const pyrolysisChemistryModel&);
65 
66 
67 protected:
68 
69  //- List of gas species present in reaction system
71 
72  //- Thermodynamic data of gases
74 
75  //- Number of gas species
76  label nGases_;
77 
78  //- Number of components being solved by ODE
80 
81  //- List of reaction rate per gas [kg/m3/s]
83 
84 
85  // Protected Member Functions
86 
87  //- Write access to source terms for gases
89 
90 
91 private:
92 
93  //- List of accumulative solid concentrations
94  mutable PtrList<volScalarField> Ys0_;
95 
96  //- Cell counter
97  label cellCounter_;
98 
99 
100 public:
101 
102  //- Runtime type information
103  TypeName("pyrolysis");
104 
105 
106  // Constructors
107 
108  //- Construct from mesh and phase name
109  pyrolysisChemistryModel(const fvMesh& mesh, const word& phaseName);
110 
111 
112  //- Destructor
113  virtual ~pyrolysisChemistryModel();
114 
115 
116  // Member Functions
117 
118  //- Thermodynamic data of gases
119  inline const PtrList<GasThermo>& gasThermo() const;
120 
121  //- Gases table
122  inline const speciesTable& gasTable() const;
123 
124  //- The number of solids
125  inline label nSpecie() const;
126 
127  //- The number of solids
128  inline label nGases() const;
129 
130 
131  //- dc/dt = omega, rate of change in concentration, for each species
132  virtual scalarField omega
133  (
134  const scalarField& c,
135  const scalar T,
136  const scalar p,
137  const bool updateC0 = false
138  ) const;
139 
140  //- Return the reaction rate for reaction r
141  // NOTE: Currently does not calculate reference specie and
142  // characteristic times (pf, cf,l Ref, etc.)
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;
156 
157 
158  //- Return the reaction rate for iReaction
159  // NOTE: Currently does not calculate reference specie and
160  // characteristic times (pf, cf,l Ref, etc.)
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;
174 
175 
176  //- Calculates the reaction rates
177  virtual void calculate();
178 
179 
180  // Chemistry model functions
181 
182  //- Return const access to the chemical source terms for gases
184  (
185  const label i
186  ) const;
187 
188  //- Return total gas source term
190 
191  //- Return sensible enthalpy for gas i [J/Kg]
192  virtual tmp<volScalarField> gasHs
193  (
194  const volScalarField& p,
195  const volScalarField& T,
196  const label i
197  ) const;
198 
199  //- Solve the reaction system for the given time step
200  // and return the characteristic time
201  virtual scalar solve(const scalar deltaT);
202 
203 
204  // ODE functions (overriding abstract functions in ODE.H)
205 
206  //- Number of ODE's to solve
207  virtual label nEqns() const;
208 
209  virtual void derivatives
210  (
211  const scalar t,
212  const scalarField& c,
213  scalarField& dcdt
214  ) const;
215 
216  virtual void jacobian
217  (
218  const scalar t,
219  const scalarField& c,
220  scalarField& dcdt,
221  scalarSquareMatrix& dfdc
222  ) const;
223 
224  virtual void solve
225  (
226  scalarField &c,
227  scalar& T,
228  scalar& p,
229  scalar& deltaT,
230  scalar& subDeltaT
231  ) const;
232 };
233 
234 
235 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
236 
237 } // End namespace Foam
238 
239 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
240 
241  #include "pyrolysisChemistryModelI.H"
242 
243 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
244 
245 #ifdef NoRepository
246  #include "pyrolysisChemistryModel.C"
247 #endif
248 
249 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
250 
251 #endif
252 
253 // ************************************************************************* //
virtual ~pyrolysisChemistryModel()
Destructor.
label nSpecie() const
The number of solids.
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 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
Return the reaction rate for iReaction.
PtrList< GasThermo > gasThermo_
Thermodynamic data of gases.
TypeName("pyrolysis")
Runtime type information.
virtual label nEqns() const
Number of ODE&#39;s to solve.
label nSpecie_
Number of components being solved by ODE.
PtrList< DimensionedField< scalar, volMesh > > & RRg()
Write access to source terms for gases.
speciesTable pyrolisisGases_
List of gas species present in reaction system.
pyrolysisChemistryModel(const fvMesh &mesh, const word &phaseName)
Construct from mesh and phase name.
Pyrolysis chemistry model. It includes gas phase in the solid reaction.
virtual void jacobian(const scalar t, const scalarField &c, scalarField &dcdt, scalarSquareMatrix &dfdc) const
Calculate the Jacobian of the system.
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
dynamicFvMesh & mesh
A class for handling words, derived from string.
Definition: word.H:59
label nGases_
Number of gas species.
virtual void derivatives(const scalar t, const scalarField &c, scalarField &dcdt) const
Calculate the derivatives in dydx.
virtual scalar solve(const scalar deltaT)
Solve the reaction system for the given time step.
const speciesTable & gasTable() const
Gases table.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
virtual tmp< volScalarField > gasHs(const volScalarField &p, const volScalarField &T, const label i) const
Return sensible enthalpy for gas i [J/Kg].
A wordList with hashed indices for faster lookup by name.
const PtrList< GasThermo > & gasThermo() const
Thermodynamic data of gases.
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 scalarField omega(const scalarField &c, const scalar T, const scalar p, const bool updateC0=false) const
dc/dt = omega, rate of change in concentration, for each species
label nGases() const
The number of solids.
volScalarField & p
A class for managing temporary objects.
Definition: PtrList.H:54
virtual void calculate()
Calculates the reaction rates.
Namespace for OpenFOAM.
PtrList< DimensionedField< scalar, volMesh > > RRg_
List of reaction rate per gas [kg/m3/s].