KinematicParcelI.H
Go to the documentation of this file.
1/*---------------------------------------------------------------------------*\
2 ========= |
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4 \\ / O peration |
5 \\ / A nd | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2011-2017 OpenFOAM Foundation
9 Copyright (C) 2016-2021 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27\*---------------------------------------------------------------------------*/
28
30
31using namespace Foam::constant::mathematical;
32
33// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
34
35template<class ParcelType>
36inline
38:
39 dict_(),
40 parcelTypeId_(dict_, -1),
41 rhoMin_(dict_, 0.0),
42 rho0_(dict_, 0.0),
43 minParcelMass_(dict_, 0.0)
44{}
45
46
47template<class ParcelType>
49(
51)
52:
53 dict_(cp.dict_),
54 parcelTypeId_(cp.parcelTypeId_),
55 rhoMin_(cp.rhoMin_),
56 rho0_(cp.rho0_),
57 minParcelMass_(cp.minParcelMass_)
58{}
59
60
61template<class ParcelType>
63(
64 const dictionary& parentDict
65)
66:
67 dict_(parentDict.subOrEmptyDict("constantProperties")),
68 parcelTypeId_(dict_, "parcelTypeId", -1),
69 rhoMin_(dict_, "rhoMin", 1e-15),
70 rho0_(dict_, "rho0", -1),
71 minParcelMass_(dict_, "minParcelMass", 1e-15)
72{}
73
74
75template<class ParcelType>
77(
78 const polyMesh& owner,
80 const label celli,
81 const label tetFacei,
82 const label tetPti
83)
84:
85 ParcelType(owner, coordinates, celli, tetFacei, tetPti),
86 active_(true),
87 typeId_(-1),
88 nParticle_(0),
89 d_(0.0),
90 dTarget_(0.0),
91 U_(Zero),
92 rho_(0.0),
93 age_(0.0),
94 tTurb_(0.0),
96{}
97
98
99template<class ParcelType>
101(
102 const polyMesh& owner,
103 const vector& position,
104 const label celli
105)
106:
107 ParcelType(owner, position, celli),
108 active_(true),
109 typeId_(-1),
110 nParticle_(0),
111 d_(0.0),
112 dTarget_(0.0),
113 U_(Zero),
114 rho_(0.0),
115 age_(0.0),
116 tTurb_(0.0),
117 UTurb_(Zero)
118{}
119
120
121template<class ParcelType>
123(
124 const polyMesh& owner,
126 const label celli,
127 const label tetFacei,
128 const label tetPti,
129 const label typeId,
130 const scalar nParticle0,
131 const scalar d0,
132 const scalar dTarget0,
133 const vector& U0,
134 const constantProperties& constProps
135)
136:
137 ParcelType(owner, coordinates, celli, tetFacei, tetPti),
138 active_(true),
139 typeId_(typeId),
140 nParticle_(nParticle0),
141 d_(d0),
142 dTarget_(dTarget0),
143 U_(U0),
144 rho_(constProps.rho0()),
145 age_(0.0),
146 tTurb_(0.0),
147 UTurb_(Zero)
148{}
149
150
151// * * * * * * * * * constantProperties Member Functions * * * * * * * * * * //
152
153template<class ParcelType>
154inline const Foam::dictionary&
156{
157 return dict_;
158}
159
160
161template<class ParcelType>
162inline Foam::label
164{
165 return parcelTypeId_.value();
166}
167
168
169template<class ParcelType>
170inline Foam::scalar
172{
173 return rhoMin_.value();
174}
175
176
177template<class ParcelType>
178inline Foam::scalar
180{
181 return rho0_.value();
182}
183
184
185template<class ParcelType>
186inline Foam::scalar
188{
189 return minParcelMass_.value();
190}
191
192
193// * * * * * * * KinematicParcel Member Functions * * * * * * * //
194
195template<class ParcelType>
197{
198 return active_;
199}
200
201
202template<class ParcelType>
204{
205 return typeId_;
206}
207
208
209template<class ParcelType>
211{
212 return nParticle_;
213}
214
215
216template<class ParcelType>
217inline Foam::scalar Foam::KinematicParcel<ParcelType>::d() const
218{
219 return d_;
220}
221
222
223template<class ParcelType>
225{
226 return dTarget_;
227}
228
229
230template<class ParcelType>
232{
233 return U_;
234}
235
236
237template<class ParcelType>
238inline Foam::scalar Foam::KinematicParcel<ParcelType>::rho() const
239{
240 return rho_;
241}
242
243
244template<class ParcelType>
245inline Foam::scalar Foam::KinematicParcel<ParcelType>::age() const
246{
247 return age_;
248}
249
250
251template<class ParcelType>
253{
254 return tTurb_;
255}
256
257
258template<class ParcelType>
260{
261 return UTurb_;
262}
263
264
265template<class ParcelType>
266inline void Foam::KinematicParcel<ParcelType>::active(const bool state)
267{
268 active_ = state;
269}
270
271
272template<class ParcelType>
274{
275 return typeId_;
276}
277
278
279template<class ParcelType>
281{
282 return nParticle_;
283}
284
285
286template<class ParcelType>
288{
289 return d_;
290}
291
292
293template<class ParcelType>
295{
296 return dTarget_;
297}
298
299
300template<class ParcelType>
302{
303 return U_;
304}
305
306
307template<class ParcelType>
309{
310 return UCorrect_;
311}
312
313
314template<class ParcelType>
316{
317 return UCorrect_;
318}
319
320
321template<class ParcelType>
323{
324 return rho_;
325}
326
327
328template<class ParcelType>
330{
331 return age_;
332}
333
334
335template<class ParcelType>
337{
338 return tTurb_;
339}
340
341
342template<class ParcelType>
344{
345 return UTurb_;
346}
347
348
349template<class ParcelType>
351(
352 const trackingData& td
353) const
354{
355 return td.rhoc()*this->mesh().cellVolumes()[this->cell()];
356}
357
358
359template<class ParcelType>
361{
362 return rho_*volume();
363}
364
365
366template<class ParcelType>
368{
369 return 0.1*mass()*sqr(d_);
370}
371
372
373template<class ParcelType>
375{
376 return volume(d_);
377}
378
379
380template<class ParcelType>
381inline Foam::scalar Foam::KinematicParcel<ParcelType>::volume(const scalar d)
382{
383 return pi/6.0*pow3(d);
384}
385
386
387template<class ParcelType>
389{
390 return areaP(d_);
391}
392
393
394template<class ParcelType>
395inline Foam::scalar Foam::KinematicParcel<ParcelType>::areaP(const scalar d)
396{
397 return 0.25*areaS(d);
398}
399
400
401template<class ParcelType>
403{
404 return areaS(d_);
405}
406
407
408template<class ParcelType>
409inline Foam::scalar Foam::KinematicParcel<ParcelType>::areaS(const scalar d)
410{
411 return pi*d*d;
412}
413
414
415template<class ParcelType>
417(
418 const trackingData& td
419) const
420{
421 return Re(td.rhoc(), U_, td.Uc(), d_, td.muc());
422}
423
424
425template<class ParcelType>
427(
428 const scalar rhoc,
429 const vector& U,
430 const vector& Uc,
431 const scalar d,
432 const scalar muc
433)
434{
435 return rhoc*mag(U - Uc)*d/max(muc, ROOTVSMALL);
436}
437
438
439template<class ParcelType>
441(
442 const trackingData& td,
443 const scalar sigma
444) const
445{
446 return We(td.rhoc(), U_, td.Uc(), d_, sigma);
447}
448
449
450template<class ParcelType>
452(
453 const scalar rhoc,
454 const vector& U,
455 const vector& Uc,
456 const scalar d,
457 const scalar sigma
458)
459{
460 return rhoc*magSqr(U - Uc)*d/max(sigma, ROOTVSMALL);
461}
462
463
464template<class ParcelType>
466(
467 const trackingData& td,
468 const scalar sigma
469) const
470{
471 return Eo(td.g(), rho_, td.rhoc(), U_, d_, sigma);
472}
473
474
475template<class ParcelType>
477(
478 const vector& g,
479 const scalar rho,
480 const scalar rhoc,
481 const vector& U,
482 const scalar d,
483 const scalar sigma
484)
485{
486 const vector dir = U/max(mag(U), ROOTVSMALL);
487 return mag(g & dir)*mag(rho - rhoc)*sqr(d)/max(sigma, ROOTVSMALL);
488}
489
490
491// ************************************************************************* //
const uniformDimensionedVectorField & g
Class to hold kinematic particle constant properties.
scalar minParcelMass() const
Return const access to the minimum parcel mass.
scalar rho0() const
Return const access to the particle density.
label parcelTypeId() const
Return const access to the parcel type id.
const dictionary & dict() const
Return const access to the constant properties dictionary.
scalar rhoMin() const
Return const access to the minimum density.
const vector & Uc() const
Return the continuous phase velocity.
scalar muc() const
Return the continuous phase viscosity.
scalar rhoc() const
Return the continuous phase density.
Kinematic parcel class with rotational motion (as spherical particles only) and one/two-way coupling ...
const vector & UCorrect() const
Return const access to correction velocity.
label typeId_
Parcel type id.
label typeId() const
Return const access to type id.
vector UTurb_
Turbulent velocity fluctuation [m/s].
scalar momentOfInertia() const
Particle moment of inertia around diameter axis.
scalar tTurb() const
Return const access to time spent in turbulent eddy.
scalar rho_
Density [kg/m3].
const vector & U() const
Return const access to velocity.
label active_
Active flag - tracking inactive when active = false.
scalar d() const
Return const access to diameter.
scalar tTurb_
Time spent in turbulent eddy [s].
scalar volume() const
Particle volume.
scalar dTarget_
Target diameter [m].
scalar nParticle() const
Return const access to number of particles.
scalar massCell(const trackingData &td) const
Cell owner mass.
scalar d_
Diameter [m].
scalar dTarget() const
Return const access to target diameter.
scalar areaP() const
Particle projected area.
scalar rho() const
Return const access to density.
scalar mass() const
Particle mass.
scalar areaS() const
Particle surface area.
bool active() const
Return const access to active flag.
scalar age() const
Return const access to the age.
vector U_
Velocity of Parcel [m/s].
const vector & UTurb() const
Return const access to turbulent velocity fluctuation.
scalar nParticle_
Number of particles in Parcel.
A cell is defined as a list of faces with extra functionality.
Definition: cell.H:57
scalar Re() const
Real part of complex number.
Definition: complexI.H:87
tmp< volScalarField > We() const
Return the bubble Webber number.
Definition: IATEsource.C:140
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
tmp< volScalarField > Eo() const
Eotvos number.
Definition: phasePair.C:142
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:81
const scalarField & cellVolumes() const
U
Definition: pEqn.H:72
dynamicFvMesh & mesh
PtrList< coordinateSystem > coordinates(solidRegions.size())
Mathematical constants.
constexpr scalar pi(M_PI)
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
scalarField Re(const UList< complex > &cf)
Extract real component.
Definition: complexField.C:159
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedScalar pow3(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
scalar rho0
const volScalarField & cp
volScalarField & e
Definition: createFields.H:11