solidChemistryModelI.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 | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2013-2016 OpenFOAM Foundation
9 Copyright (C) 2016 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27\*---------------------------------------------------------------------------*/
28
29#include "volFields.H"
30
31// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
32
33template<class CompType, class SolidThermo>
36{
37 return RRs_;
38}
39
40
41template<class CompType, class SolidThermo>
44{
45 return reactions_;
46}
47
48
49template<class CompType, class SolidThermo>
50inline Foam::label
52nReaction() const
53{
54 return nReaction_;
55}
56
57
58template<class CompType, class SolidThermo>
61(
62 const label i
63) const
64{
65 return RRs_[i];
66}
67
68
69template<class CompType, class SolidThermo>
72{
74 (
76 (
78 (
79 "RRs",
80 this->time().timeName(),
81 this->mesh(),
82 IOobject::NO_READ,
83 IOobject::NO_WRITE
84 ),
85 this->mesh(),
87 )
88 );
89
90 if (this->chemistry_)
91 {
92 volScalarField::Internal& RRs = tRRs.ref();
93 for (label i=0; i < nSolids_; i++)
94 {
95 RRs += RRs_[i];
96 }
97 }
98 return tRRs;
99}
100
101
102template<class CompType, class SolidThermo>
105{
107 (
109 (
111 (
112 "RRsHs",
113 this->time().timeName(),
114 this->mesh(),
115 IOobject::NO_READ,
116 IOobject::NO_WRITE
117 ),
118 this->mesh(),
120 )
121 );
122
123 if (this->chemistry_)
124 {
126
127 const volScalarField& T = this->solidThermo().T();
128 const volScalarField& p = this->solidThermo().p();
129
130 for (label i=0; i < nSolids_; i++)
131 {
132 forAll(RRs, cellI)
133 {
134 RRs[cellI] +=
135 RRs_[i][cellI]*solidThermo_[i].Hs(p[cellI], T[cellI]);
136 }
137 }
138 }
139 return tRRsHs;
140}
141
142
143// ************************************************************************* //
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:170
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
Definition: PtrList.H:73
virtual const volScalarField & T() const
Temperature [K].
Definition: basicThermo.C:622
virtual volScalarField & p()
Pressure [Pa].
Definition: basicThermo.C:610
const PtrList< Reaction< SolidThermo > > & reactions() const
The reactions.
tmp< DimensionedField< scalar, volMesh > > RRsHs() const
Return net solid sensible enthalpy.
label nReaction() const
The number of reactions.
PtrList< volScalarField::Internal > & RRs()
Write access to source terms for solids.
Fundamental solid thermodynamic properties.
Definition: solidThermo.H:55
A class for managing temporary objects.
Definition: tmp.H:65
T & ref() const
Definition: tmpI.H:227
volScalarField & p
const volScalarField & T
dynamicFvMesh & mesh
word timeName
Definition: getTimeIndex.H:3
const dimensionSet dimEnergy
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:53
const dimensionSet dimVolume(pow3(dimLength))
Definition: dimensionSets.H:60
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:51
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333