ReactingMultiphaseParcelI.H
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-2020 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 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
27 
28 template<class ParcelType>
31 :
32  ParcelType::constantProperties(),
33  TDevol_(this->dict_, 0.0),
34  LDevol_(this->dict_, 0.0),
35  hRetentionCoeff_(this->dict_, 0.0)
36 {}
37 
38 
39 template<class ParcelType>
42 (
43  const constantProperties& cp
44 )
45 :
46  ParcelType::constantProperties(cp),
47  TDevol_(cp.TDevol_),
48  LDevol_(cp.LDevol_),
49  hRetentionCoeff_(cp.hRetentionCoeff_)
50 {}
51 
52 
53 template<class ParcelType>
56 (
57  const dictionary& parentDict
58 )
59 :
60  ParcelType::constantProperties(parentDict),
61  TDevol_(this->dict_, "TDevol"),
62  LDevol_(this->dict_, "LDevol"),
63  hRetentionCoeff_(this->dict_, "hRetentionCoeff")
64 {}
65 
66 
67 template<class ParcelType>
69 (
70  const polyMesh& mesh,
71  const barycentric& coordinates,
72  const label celli,
73  const label tetFacei,
74  const label tetPti
75 )
76 :
77  ParcelType(mesh, coordinates, celli, tetFacei, tetPti),
78  YGas_(0),
79  YLiquid_(0),
80  YSolid_(0),
81  canCombust_(0)
82 {}
83 
84 
85 template<class ParcelType>
87 (
88  const polyMesh& mesh,
89  const vector& position,
90  const label celli
91 )
92 :
93  ParcelType(mesh, position, celli),
94  YGas_(0),
95  YLiquid_(0),
96  YSolid_(0),
97  canCombust_(0)
98 {}
99 
100 
101 // * * * * * * * * * constantProperties Member Functions * * * * * * * * * * //
102 
103 template<class ParcelType>
104 inline Foam::scalar
106 {
107  return TDevol_.value();
108 }
109 
110 
111 template<class ParcelType>
112 inline Foam::scalar
114 {
115  return LDevol_.value();
116 }
117 
118 
119 template<class ParcelType>
120 inline Foam::scalar
123 {
124  scalar value = hRetentionCoeff_.value();
125 
126  if ((value < 0) || (value > 1))
127  {
129  << "hRetentionCoeff must be in the range 0 to 1" << nl
130  << exit(FatalError) << endl;
131  }
132 
133  return value;
134 }
135 
136 
137 // * * * * * * * * * * ThermoParcel Member Functions * * * * * * * * * * * * //
138 
139 template<class ParcelType>
141 YGas() const
142 {
143  return YGas_;
144 }
145 
146 
147 template<class ParcelType>
149 YLiquid() const
150 {
151  return YLiquid_;
152 }
153 
154 
155 template<class ParcelType>
157 YSolid() const
158 {
159  return YSolid_;
160 }
161 
162 
163 template<class ParcelType>
164 inline Foam::label
166 {
167  return canCombust_;
168 }
169 
170 
171 template<class ParcelType>
173 {
174  return YGas_;
175 }
176 
177 
178 template<class ParcelType>
180 {
181  return YLiquid_;
182 }
183 
184 
185 template<class ParcelType>
187 {
188  return YSolid_;
189 }
190 
191 
192 template<class ParcelType>
194 {
195  return canCombust_;
196 }
197 
198 
199 // ************************************************************************* //
label canCombust() const
Return const access to the canCombust flag.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:306
label canCombust_
Flag to identify if the particle can devolatilise and combust.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
scalar LDevol() const
Return const access to the latent heat of devolatilisation.
const Type & value() const
Return the value.
ReactingMultiphaseParcel(const polyMesh &mesh, const barycentric &coordinates, const label celli, const label tetFacei, const label tetPti)
Construct from mesh, position and topology.
scalarField YGas_
Mass fractions of gases [].
scalar hRetentionCoeff() const
Return const access to the fraction of enthalpy retained by.
const scalarField & YGas() const
Return const access to mass fractions of gases.
scalar TDevol() const
Return const access to the devolatilisation temperature.
scalarField YSolid_
Mass fractions of solids [].
static const char nl
Definition: Ostream.H:260
const scalarField & YLiquid() const
Return const access to mass fractions of liquids.
scalarField YLiquid_
Mass fractions of liquids [].
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:76
Class to hold reacting multiphase particle constant properties.
const scalarField & YSolid() const
Return const access to mass fractions of solids.