SprayParcelI.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 
26 #include "SprayParcel.H"
27 
28 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
29 
30 template<class ParcelType>
32 :
33  ParcelType::constantProperties(),
34  sigma0_(this->dict_, 0.0),
35  mu0_(this->dict_, 0.0)
36 {}
37 
38 
39 template<class ParcelType>
41 (
42  const constantProperties& cp
43 )
44 :
45  ParcelType::constantProperties(cp),
46  sigma0_(cp.sigma0_),
47  mu0_(cp.mu0_)
48 {}
49 
50 
51 template<class ParcelType>
53 (
54  const dictionary& parentDict
55 )
56 :
57  ParcelType::constantProperties(parentDict),
58  sigma0_(this->dict_, "sigma0"),
59  mu0_(this->dict_, "mu0")
60 {}
61 
62 
63 template<class ParcelType>
65 (
66  const label parcelTypeId,
67  const scalar rhoMin,
68  const scalar rho0,
69  const scalar minParcelMass,
70  const scalar youngsModulus,
71  const scalar poissonsRatio,
72  const scalar T0,
73  const scalar TMin,
74  const scalar TMax,
75  const scalar Cp0,
76  const scalar epsilon0,
77  const scalar f0,
78  const scalar Pr,
79  const scalar pMin,
80  const Switch& constantVolume,
81  const scalar sigma0,
82  const scalar mu0
83 )
84 :
85  ParcelType::constantProperties
86  (
87  parcelTypeId,
88  rhoMin,
89  rho0,
90  minParcelMass,
91  youngsModulus,
92  poissonsRatio,
93  T0,
94  TMin,
95  TMax,
96  Cp0,
97  epsilon0,
98  f0,
99  Pr,
100  pMin,
101  constantVolume
102  ),
103  sigma0_(this->dict_, sigma0),
104  mu0_(this->dict_, mu0)
105 {}
106 
107 
108 template<class ParcelType>
110 (
111  const polyMesh& mesh,
112  const barycentric& coordinates,
113  const label celli,
114  const label tetFacei,
115  const label tetPti,
116  const label facei
117 )
118 :
119  ParcelType(mesh, coordinates, celli, tetFacei, tetPti, facei),
120  d0_(this->d()),
121  mass0_(this->mass()),
122  position0_(this->position(mesh)),
123  sigma_(0.0),
124  mu_(0.0),
125  liquidCore_(0.0),
126  KHindex_(0.0),
127  y_(0.0),
128  yDot_(0.0),
129  tc_(0.0),
130  ms_(0.0),
131  injector_(-1),
132  tMom_(great)
133 {}
134 
135 
136 template<class ParcelType>
138 (
139  const polyMesh& mesh,
140  const vector& position,
141  const label celli
142 )
143 :
144  ParcelType(mesh, position, celli),
145  d0_(this->d()),
146  mass0_(this->mass()),
147  position0_(this->position(mesh)),
148  sigma_(0.0),
149  mu_(0.0),
150  liquidCore_(0.0),
151  KHindex_(0.0),
152  y_(0.0),
153  yDot_(0.0),
154  tc_(0.0),
155  ms_(0.0),
156  injector_(-1),
157  tMom_(great)
158 {}
159 
160 
161 // * * * * * * * * * constantProperties Member Functions * * * * * * * * * * //
162 
163 template<class ParcelType>
164 inline Foam::scalar
166 {
167  return sigma0_.value();
168 }
169 
170 
171 template<class ParcelType>
172 inline Foam::scalar
174 {
175  return mu0_.value();
176 }
177 
178 
179 // * * * * * * * * * * SprayParcel Member Functions * * * * * * * * * * * * //
180 
181 template<class ParcelType>
182 inline Foam::scalar Foam::SprayParcel<ParcelType>::d0() const
183 {
184  return d0_;
185 }
186 
187 
188 template<class ParcelType>
189 inline Foam::scalar Foam::SprayParcel<ParcelType>::mass0() const
190 {
191  return mass0_;
192 }
193 
194 
195 template<class ParcelType>
197 {
198  return position0_;
199 }
200 
201 
202 template<class ParcelType>
203 inline Foam::scalar Foam::SprayParcel<ParcelType>::sigma() const
204 {
205  return sigma_;
206 }
207 
208 
209 template<class ParcelType>
210 inline Foam::scalar Foam::SprayParcel<ParcelType>::mu() const
211 {
212  return mu_;
213 }
214 
215 
216 template<class ParcelType>
217 inline Foam::scalar Foam::SprayParcel<ParcelType>::liquidCore() const
218 {
219  return liquidCore_;
220 }
221 
222 
223 template<class ParcelType>
224 inline Foam::scalar Foam::SprayParcel<ParcelType>::KHindex() const
225 {
226  return KHindex_;
227 }
228 
229 
230 template<class ParcelType>
231 inline Foam::scalar Foam::SprayParcel<ParcelType>::y() const
232 {
233  return y_;
234 }
235 
236 
237 template<class ParcelType>
238 inline Foam::scalar Foam::SprayParcel<ParcelType>::yDot() const
239 {
240  return yDot_;
241 }
242 
243 
244 template<class ParcelType>
245 inline Foam::scalar Foam::SprayParcel<ParcelType>::tc() const
246 {
247  return tc_;
248 }
249 
250 
251 template<class ParcelType>
252 inline Foam::scalar Foam::SprayParcel<ParcelType>::ms() const
253 {
254  return ms_;
255 }
256 
257 
258 template<class ParcelType>
260 {
261  return injector_;
262 }
263 
264 
265 template<class ParcelType>
266 inline Foam::scalar Foam::SprayParcel<ParcelType>::tMom() const
267 {
268  return tMom_;
269 }
270 
271 
272 template<class ParcelType>
273 inline Foam::scalar& Foam::SprayParcel<ParcelType>::d0()
274 {
275  return d0_;
276 }
277 
278 
279 template<class ParcelType>
281 {
282  return mass0_;
283 }
284 
285 
286 template<class ParcelType>
288 {
289  return position0_;
290 }
291 
292 
293 template<class ParcelType>
295 {
296  return sigma_;
297 }
298 
299 
300 template<class ParcelType>
301 inline Foam::scalar& Foam::SprayParcel<ParcelType>::mu()
302 {
303  return mu_;
304 }
305 
306 
307 template<class ParcelType>
309 {
310  return liquidCore_;
311 }
312 
313 
314 template<class ParcelType>
316 {
317  return KHindex_;
318 }
319 
320 
321 template<class ParcelType>
322 inline Foam::scalar& Foam::SprayParcel<ParcelType>::y()
323 {
324  return y_;
325 }
326 
327 
328 template<class ParcelType>
330 {
331  return yDot_;
332 }
333 
334 
335 template<class ParcelType>
336 inline Foam::scalar& Foam::SprayParcel<ParcelType>::tc()
337 {
338  return tc_;
339 }
340 
341 
342 template<class ParcelType>
343 inline Foam::scalar& Foam::SprayParcel<ParcelType>::ms()
344 {
345  return ms_;
346 }
347 
348 
349 template<class ParcelType>
351 {
352  return injector_;
353 }
354 
355 
356 template<class ParcelType>
358 {
359  return tMom_;
360 }
361 
362 
363 // ************************************************************************* //
Class to hold reacting particle constant properties.
Definition: SprayParcel.H:82
scalar mu0() const
Return const access to the initial dynamic viscosity.
Definition: SprayParcelI.H:173
scalar sigma0() const
Return const access to the initial surface tension.
Definition: SprayParcelI.H:165
SprayParcel(const polyMesh &mesh, const barycentric &coordinates, const label celli, const label tetFacei, const label tetPti, const label facei)
Construct from mesh, coordinates and topology.
Definition: SprayParcelI.H:110
scalar liquidCore() const
Return const access to liquid core.
Definition: SprayParcelI.H:217
scalar d0_
Initial droplet diameter [m].
Definition: SprayParcel.H:149
label injector() const
Return const access to injector id.
Definition: SprayParcelI.H:259
scalar tMom() const
Return const access to momentum relaxation time.
Definition: SprayParcelI.H:266
vector position0_
Injection position.
Definition: SprayParcel.H:155
label injector_
Injector id.
Definition: SprayParcel.H:182
scalar mu() const
Return const access to the liquid dynamic viscosity.
Definition: SprayParcelI.H:210
scalar mu_
Liquid dynamic viscosity [Pa.s].
Definition: SprayParcel.H:161
scalar sigma_
Liquid surface tension [N/m].
Definition: SprayParcel.H:158
scalar yDot() const
Return const access to rate of change of spherical deviation.
Definition: SprayParcelI.H:238
scalar d0() const
Return const access to initial droplet diameter.
Definition: SprayParcelI.H:182
scalar tMom_
Momentum relaxation time (needed for calculating parcel acc.)
Definition: SprayParcel.H:185
scalar KHindex() const
Return const access to Kelvin-Helmholtz breakup index.
Definition: SprayParcelI.H:224
scalar tc() const
Return const access to atomisation characteristic time.
Definition: SprayParcelI.H:245
scalar y_
Spherical deviation.
Definition: SprayParcel.H:170
scalar mass0() const
Return const access to initial mass [kg].
Definition: SprayParcelI.H:189
scalar tc_
Characteristic time (used in atomisation and/or breakup model)
Definition: SprayParcel.H:176
scalar KHindex_
Index for KH Breakup.
Definition: SprayParcel.H:167
scalar ms() const
Return const access to stripped parcel mass.
Definition: SprayParcelI.H:252
scalar y() const
Return const access to spherical deviation.
Definition: SprayParcelI.H:231
scalar mass0_
Initial mass [kg].
Definition: SprayParcel.H:152
scalar sigma() const
Return const access to the liquid surface tension.
Definition: SprayParcelI.H:203
scalar ms_
Stripped parcel mass due to breakup.
Definition: SprayParcel.H:179
scalar yDot_
Rate of change of spherical deviation.
Definition: SprayParcel.H:173
scalar liquidCore_
Part of liquid core ( >0.5=liquid, <0.5=droplet )
Definition: SprayParcel.H:164
const vector & position0() const
Return const access to initial droplet position.
Definition: SprayParcelI.H:196
A simple wrapper around bool so that it can be read as a word: true/false, on/off,...
Definition: Switch.H:61
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
const dimensionedScalar mu0
Magnetic constant/permeability of free space: default SI units: [H/m].
const dimensionedScalar epsilon0
Electric constant: default SI units: [F/m].
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
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
scalar rho0
scalar T0
Definition: createFields.H:22