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 )
96 :
97  ParcelType(mesh, position, celli),
98  mass0_(0),
99  YGas_(0),
100  YLiquid_(0),
101  YSolid_(0),
102  canCombust_(0)
103 {}
104 
105 
106 // * * * * * * * * * constantProperties Member Functions * * * * * * * * * * //
107 
108 template<class ParcelType>
109 inline Foam::scalar
111 {
112  return TDevol_.value();
113 }
114 
115 
116 template<class ParcelType>
117 inline Foam::scalar
119 {
120  return LDevol_.value();
121 }
122 
123 
124 template<class ParcelType>
125 inline Foam::scalar
127 hRetentionCoeff() const
128 {
129  scalar value = hRetentionCoeff_.value();
130 
131  if ((value < 0) || (value > 1))
132  {
134  << "hRetentionCoeff must be in the range 0 to 1" << nl
135  << exit(FatalError) << endl;
136  }
137 
138  return value;
139 }
140 
141 
142 // * * * * * * * * ReactingMultiphaseParcel Member Functions * * * * * * * * //
143 
144 template<class ParcelType>
146 {
147  return mass0_;
148 }
149 
150 
151 template<class ParcelType>
153 YGas() const
154 {
155  return YGas_;
156 }
157 
158 
159 template<class ParcelType>
161 YLiquid() const
162 {
163  return YLiquid_;
164 }
165 
166 
167 template<class ParcelType>
169 YSolid() const
170 {
171  return YSolid_;
172 }
173 
174 
175 template<class ParcelType>
176 inline Foam::label
178 {
179  return canCombust_;
180 }
181 
182 
183 template<class ParcelType>
185 {
186  return mass0_;
187 }
188 
189 
190 template<class ParcelType>
192 {
193  return YGas_;
194 }
195 
196 
197 template<class ParcelType>
199 {
200  return YLiquid_;
201 }
202 
203 
204 template<class ParcelType>
206 {
207  return YSolid_;
208 }
209 
210 
211 template<class ParcelType>
213 {
214  return canCombust_;
215 }
216 
217 
218 // ************************************************************************* //
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:160
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:80
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:306
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:251
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:260