solidChemistryModel.C
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-2018 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 \*---------------------------------------------------------------------------*/
25 
26 #include "solidChemistryModel.H"
27 #include "reactingMixture.H"
28 
29 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
30 
31 template<class CompType, class SolidThermo>
33 (
34  typename CompType::reactionThermo& thermo
35 )
36 :
37  CompType(thermo),
38  ODESystem(),
39  Ys_(this->solidThermo().composition().Y()),
40  reactions_
41  (
42  dynamic_cast<const reactingMixture<SolidThermo>& >
43  (
44  this->solidThermo()
45  )
46  ),
47  solidThermo_
48  (
49  dynamic_cast<const reactingMixture<SolidThermo>& >
50  (
51  this->solidThermo()
52  ).speciesData()
53  ),
54  nSolids_(Ys_.size()),
55  nReaction_(reactions_.size()),
56  RRs_(nSolids_),
57  reactingCells_(this->mesh().nCells(), true)
58 {
59  // create the fields for the chemistry sources
60  forAll(RRs_, fieldi)
61  {
62  RRs_.set
63  (
64  fieldi,
66  (
67  IOobject
68  (
69  "RRs." + Ys_[fieldi].name(),
70  this->mesh().time().timeName(),
71  this->mesh(),
72  IOobject::NO_READ,
73  IOobject::NO_WRITE
74  ),
75  this->mesh(),
77  )
78  );
79  }
80 }
81 
82 
83 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
84 
85 template<class CompType, class SolidThermo>
88 {}
89 
90 
91 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
92 
93 template<class CompType, class SolidThermo>
95 (
96  const scalarField& deltaT
97 )
98 {
100  return 0;
101 }
102 
103 
104 template<class CompType, class SolidThermo>
107 {
109  return volScalarField::null();
110 }
111 
112 
113 template<class CompType, class SolidThermo>
116 {
117  tmp<volScalarField> tQdot
118  (
120  (
121  "Qdot",
122  this->mesh_,
124  )
125  );
126 
127  if (this->chemistry_)
128  {
129  scalarField& Qdot = tQdot.ref();
130 
131  forAll(Ys_, i)
132  {
133  forAll(Qdot, celli)
134  {
135  scalar hf = solidThermo_[i].Hc();
136  Qdot[celli] -= hf*RRs_[i][celli];
137  }
138  }
139  }
140 
141  return tQdot;
142 }
143 
144 
145 template<class CompType, class SolidThermo>
147 (
148  const label celli,
149  const bool active
150 )
151 {
152  reactingCells_[celli] = active;
153 }
154 
155 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
Foam::reactingMixture.
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.
basicSpecieMixture & composition
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:174
solidChemistryModel(typename CompType::reactionThermo &thermo)
Construct from thermo.
virtual ~solidChemistryModel()
Destructor.
scalar Qdot
Definition: solveChemistry.H:2
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
virtual tmp< volScalarField > tc() const
Return the chemical time scale.
const dimensionSet dimVolume(pow3(dimLength))
Definition: dimensionSets.H:58
dynamicFvMesh & mesh
virtual scalar solve(const scalar deltaT)=0
Solve the reaction system for the given time step.
word timeName
Definition: getTimeIndex.H:3
virtual tmp< volScalarField > Qdot() const
Return the heat release rate [kg/m/s^3].
Fundamental solid thermodynamic properties.
Definition: solidThermo.H:48
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
const dimensionSet dimEnergy
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
PtrList< volScalarField > & Y
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:51
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:49
A class for managing temporary objects.
Definition: PtrList.H:53
#define NotImplemented
Issue a FatalErrorIn for a function not currently implemented.
Definition: error.H:366
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92