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-2023 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 
27 
28 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
29 
30 template<class ParcelType>
33 :
34  ParcelType::constantProperties(),
35  TDevol_(this->dict_, 0.0),
36  LDevol_(this->dict_, 0.0),
37  hRetentionCoeff_(this->dict_, 0.0)
38 {}
39 
40 
41 template<class ParcelType>
44 (
45  const constantProperties& cp
46 )
47 :
48  ParcelType::constantProperties(cp),
49  TDevol_(cp.TDevol_),
50  LDevol_(cp.LDevol_),
51  hRetentionCoeff_(cp.hRetentionCoeff_)
52 {}
53 
54 
55 template<class ParcelType>
58 (
59  const dictionary& parentDict
60 )
61 :
62  ParcelType::constantProperties(parentDict),
63  TDevol_(this->dict_, "TDevol"),
64  LDevol_(this->dict_, "LDevol"),
65  hRetentionCoeff_(this->dict_, "hRetentionCoeff")
66 {}
67 
68 
69 template<class ParcelType>
71 (
72  const polyMesh& mesh,
73  const barycentric& coordinates,
74  const label celli,
75  const label tetFacei,
76  const label tetPti,
77  const label facei
78 )
79 :
80  ParcelType(mesh, coordinates, celli, tetFacei, tetPti, facei),
81  mass0_(0),
82  YGas_(0),
83  YLiquid_(0),
84  YSolid_(0),
85  canCombust_(0)
86 {}
87 
88 
89 template<class ParcelType>
91 (
92  const polyMesh& mesh,
93  const vector& position,
94  const label celli,
95  label& nLocateBoundaryHits
96 )
97 :
98  ParcelType(mesh, position, celli, nLocateBoundaryHits),
99  mass0_(0),
100  YGas_(0),
101  YLiquid_(0),
102  YSolid_(0),
103  canCombust_(0)
104 {}
105 
106 
107 // * * * * * * * * * constantProperties Member Functions * * * * * * * * * * //
108 
109 template<class ParcelType>
110 inline Foam::scalar
112 {
113  return TDevol_.value();
114 }
115 
116 
117 template<class ParcelType>
118 inline Foam::scalar
120 {
121  return LDevol_.value();
122 }
123 
124 
125 template<class ParcelType>
126 inline Foam::scalar
128 hRetentionCoeff() const
129 {
130  scalar value = hRetentionCoeff_.value();
131 
132  if ((value < 0) || (value > 1))
133  {
135  << "hRetentionCoeff must be in the range 0 to 1" << nl
136  << exit(FatalError) << endl;
137  }
138 
139  return value;
140 }
141 
142 
143 // * * * * * * * * ReactingMultiphaseParcel Member Functions * * * * * * * * //
144 
145 template<class ParcelType>
147 {
148  return mass0_;
149 }
150 
151 
152 template<class ParcelType>
154 YGas() const
155 {
156  return YGas_;
157 }
158 
159 
160 template<class ParcelType>
162 YLiquid() const
163 {
164  return YLiquid_;
165 }
166 
167 
168 template<class ParcelType>
170 YSolid() const
171 {
172  return YSolid_;
173 }
174 
175 
176 template<class ParcelType>
177 inline Foam::label
179 {
180  return canCombust_;
181 }
182 
183 
184 template<class ParcelType>
186 {
187  return mass0_;
188 }
189 
190 
191 template<class ParcelType>
193 {
194  return YGas_;
195 }
196 
197 
198 template<class ParcelType>
200 {
201  return YLiquid_;
202 }
203 
204 
205 template<class ParcelType>
207 {
208  return YSolid_;
209 }
210 
211 
212 template<class ParcelType>
214 {
215  return canCombust_;
216 }
217 
218 
219 // ************************************************************************* //
Class to hold reacting multiphase particle constant properties.
scalar hRetentionCoeff() const
Return const access to the fraction of enthalpy retained by.
scalar TDevol() const
Return const access to the devolatilisation temperature.
scalar LDevol() const
Return const access to the latent heat of devolatilisation.
scalarField YLiquid_
Mass fractions of liquids [].
scalarField YSolid_
Mass fractions of solids [].
label canCombust_
Flag to identify if the particle can devolatilise and combust.
scalarField YGas_
Mass fractions of gases [].
label canCombust() const
Return const access to the canCombust flag.
const scalarField & YGas() const
Return const access to mass fractions of gases.
scalar mass0() const
Return const access to initial mass [kg].
ReactingMultiphaseParcel(const polyMesh &mesh, const barycentric &coordinates, const label celli, const label tetFacei, const label tetPti, const label facei)
Construct from mesh, coordinates and topology.
const scalarField & YSolid() const
Return const access to mass fractions of solids.
const scalarField & YLiquid() const
Return const access to mass fractions of liquids.
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:80
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:334
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
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
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
bool cp(const fileName &src, const fileName &dst, const bool followLink=true)
Copy, recursively if necessary, the source to the destination.
Definition: POSIX.C:753
error FatalError
static const char nl
Definition: Ostream.H:266