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-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 // * * * * * * * * * * * * * * * * 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 template<class ParcelType>
103 (
104  const polyMesh& mesh,
105  const barycentric& coordinates,
106  const label celli,
107  const label tetFacei,
108  const label tetPti,
109  const label typeId,
110  const scalar nParticle0,
111  const scalar d0,
112  const scalar dTarget0,
113  const vector& U0,
114  const vector& f0,
115  const vector& angularMomentum0,
116  const vector& torque0,
117  const scalarField& Y0,
118  const scalarField& YGas0,
119  const scalarField& YLiquid0,
120  const scalarField& YSolid0,
121  const constantProperties& constProps
122 )
123 :
124  ParcelType
125  (
126  mesh,
127  coordinates,
128  celli,
129  tetFacei,
130  tetPti,
131  typeId,
132  nParticle0,
133  d0,
134  dTarget0,
135  U0,
136  f0,
137  angularMomentum0,
138  torque0,
139  Y0,
140  constProps
141  ),
142  YGas_(YGas0),
143  YLiquid_(YLiquid0),
144  YSolid_(YSolid0),
145  canCombust_(0)
146 {}
147 
148 
149 // * * * * * * * * * constantProperties Member Functions * * * * * * * * * * //
150 
151 template<class ParcelType>
152 inline Foam::scalar
154 {
155  return TDevol_.value();
156 }
157 
158 
159 template<class ParcelType>
160 inline Foam::scalar
162 {
163  return LDevol_.value();
164 }
165 
166 
167 template<class ParcelType>
168 inline Foam::scalar
171 {
172  scalar value = hRetentionCoeff_.value();
173 
174  if ((value < 0) || (value > 1))
175  {
177  << "hRetentionCoeff must be in the range 0 to 1" << nl
178  << exit(FatalError) << endl;
179  }
180 
181  return value;
182 }
183 
184 
185 // * * * * * * * * * * ThermoParcel Member Functions * * * * * * * * * * * * //
186 
187 template<class ParcelType>
189 YGas() const
190 {
191  return YGas_;
192 }
193 
194 
195 template<class ParcelType>
197 YLiquid() const
198 {
199  return YLiquid_;
200 }
201 
202 
203 template<class ParcelType>
205 YSolid() const
206 {
207  return YSolid_;
208 }
209 
210 
211 template<class ParcelType>
212 inline Foam::label
214 {
215  return canCombust_;
216 }
217 
218 
219 template<class ParcelType>
221 {
222  return YGas_;
223 }
224 
225 
226 template<class ParcelType>
228 {
229  return YLiquid_;
230 }
231 
232 
233 template<class ParcelType>
235 {
236  return YSolid_;
237 }
238 
239 
240 template<class ParcelType>
242 {
243  return canCombust_;
244 }
245 
246 
247 // ************************************************************************* //
label canCombust() const
Return const access to the canCombust flag.
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
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:137
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
label canCombust_
Flag to identify if the particle can devolatilise and combust.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
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:265
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:74
Class to hold reacting multiphase particle constant properties.
const scalarField & YSolid() const
Return const access to mass fractions of solids.