Home
Downloads
Documentation
Installation
User Guide
man-pages
API Documentation
README
Release Notes
Changes
License
Support
SourceForge Project
Main Page
Related Pages
Namespaces
Classes
Files
File List
File Members
All
Classes
Namespaces
Files
Functions
Variables
Typedefs
Enumerations
Enumerator
Friends
Macros
Pages
src
lagrangian
intermediate
parcels
Templates
KinematicParcel
KinematicParcel.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 | Copyright (C) 1991-2010 OpenCFD Ltd.
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
Class
25
Foam::KinematicParcel
26
27
Description
28
Kinematic parcel class with one/two-way coupling with the continuous
29
phase.
30
31
Sub-models include:
32
- drag
33
- turbulent dispersion
34
- wall interactions
35
36
SourceFiles
37
KinematicParcelI.H
38
KinematicParcel.C
39
KinematicParcelIO.C
40
41
\*---------------------------------------------------------------------------*/
42
43
#ifndef KinematicParcel_H
44
#define KinematicParcel_H
45
46
#include <
lagrangian/Particle.H
>
47
#include <
OpenFOAM/IOstream.H
>
48
#include <
OpenFOAM/autoPtr.H
>
49
#include <
finiteVolume/interpolationCellPoint.H
>
50
#include <
OpenFOAM/contiguous.H
>
51
52
#include <
lagrangianIntermediate/KinematicCloud_.H
>
53
54
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
55
56
namespace
Foam
57
{
58
59
template
<
class
ParcelType>
60
class
KinematicParcel;
61
62
// Forward declaration of friend functions
63
64
template
<
class
ParcelType>
65
Ostream&
operator
<<
66
(
67
Ostream&,
68
const
KinematicParcel<ParcelType>&
69
);
70
71
/*---------------------------------------------------------------------------*\
72
Class KinematicParcel Declaration
73
\*---------------------------------------------------------------------------*/
74
75
template
<
class
ParcelType>
76
class
KinematicParcel
77
:
78
public
Particle
<ParcelType>
79
{
80
public
:
81
82
//- Class to hold kinematic particle constant properties
83
class
constantProperties
84
{
85
// Private data
86
87
//- Constant properties dictionary
88
const
dictionary
dict_;
89
90
//- Minimum density [kg/m3]
91
const
scalar rhoMin_;
92
93
//- Particle density [kg/m3] (constant)
94
const
scalar rho0_;
95
96
//- Minimum particle mass [kg]
97
const
scalar minParticleMass_;
98
99
100
public
:
101
102
//- Constructor
103
constantProperties
(
const
dictionary
& parentDict);
104
105
// Member functions
106
107
//- Return const access to the constant properties dictionary
108
inline
const
dictionary
&
dict
()
const
;
109
110
//- Return const access to the minimum density
111
inline
scalar
rhoMin
()
const
;
112
113
//- Return const access to the particle density
114
inline
scalar
rho0
()
const
;
115
116
//- Return const access to the minimum particle mass
117
inline
scalar
minParticleMass
()
const
;
118
};
119
120
121
//- Class used to pass kinematic tracking data to the trackToFace function
122
class
trackData
123
:
124
public
Particle
<ParcelType>
::trackData
125
{
126
// Private data
127
128
//- Reference to the cloud containing this particle
129
KinematicCloud<ParcelType>
& cloud_;
130
131
//- Particle constant properties
132
const
constantProperties
& constProps_;
133
134
135
// Interpolators for continuous phase fields
136
137
//- Density interpolator
138
const
interpolation<scalar>
& rhoInterp_;
139
140
//- Velocity interpolator
141
const
interpolation<vector>
& UInterp_;
142
143
//- Dynamic viscosity interpolator
144
const
interpolation<scalar>
& muInterp_;
145
146
//- Local gravitational or other body-force acceleration
147
const
vector
& g_;
148
149
150
public
:
151
152
// Constructors
153
154
//- Construct from components
155
inline
trackData
156
(
157
KinematicCloud<ParcelType>
&
cloud
,
158
const
constantProperties
&
constProps
,
159
const
interpolation<scalar>
&
rhoInterp
,
160
const
interpolation<vector>
&
UInterp
,
161
const
interpolation<scalar>
&
muInterp
,
162
const
vector
&
g
163
);
164
165
166
// Member functions
167
168
//- Return access to the owner cloud
169
inline
KinematicCloud<ParcelType>
&
cloud
();
170
171
//- Return const access to the constant properties
172
inline
const
constantProperties
&
constProps
()
const
;
173
174
//- Return conat access to the interpolator for continuous
175
// phase density field
176
inline
const
interpolation<scalar>
&
rhoInterp
()
const
;
177
178
//- Return conat access to the interpolator for continuous
179
// phase velocity field
180
inline
const
interpolation<vector>
&
UInterp
()
const
;
181
182
//- Return conat access to the interpolator for continuous
183
// phase dynamic viscosity field
184
inline
const
interpolation<scalar>
&
muInterp
()
const
;
185
186
// Return const access to the gravitational acceleration vector
187
inline
const
vector
&
g
()
const
;
188
};
189
190
191
protected
:
192
193
// Protected data
194
195
// Parcel properties
196
197
//- Active flag - tracking inactive when active = false
198
bool
active_
;
199
200
//- Parcel type id
201
label
typeId_
;
202
203
//- Number of particles in Parcel
204
scalar
nParticle_
;
205
206
//- Diameter [m]
207
scalar
d_
;
208
209
//- Velocity of Parcel [m/s]
210
vector
U_
;
211
212
//- Density [kg/m3]
213
scalar
rho_
;
214
215
//- Time spent in turbulent eddy [s]
216
scalar
tTurb_
;
217
218
//- Turbulent velocity fluctuation [m/s]
219
vector
UTurb_
;
220
221
222
// Cell-based quantities
223
224
//- Density [kg/m3]
225
scalar
rhoc_
;
226
227
//- Velocity [m/s]
228
vector
Uc_
;
229
230
//- Viscosity [Pa.s]
231
scalar
muc_
;
232
233
234
// Protected member functions
235
236
//- Calculate new particle velocity
237
template
<
class
TrackData>
238
const
vector
calcVelocity
239
(
240
TrackData& td,
241
const
scalar dt,
// timestep
242
const
label cellI,
// owner cell
243
const
scalar
Re
,
// Reynolds number
244
const
scalar
mu
,
// local carrier viscosity
245
const
scalar
d
,
// diameter
246
const
vector
&
U
,
// velocity
247
const
scalar
rho
,
// density
248
const
scalar
mass
,
// mass
249
const
vector
&
Su
,
// explicit particle momentum source
250
vector
& dUTrans
// momentum transfer to carrier
251
)
const
;
252
253
254
public
:
255
256
// Static data members
257
258
//- String representation of properties
259
static
string
propHeader
;
260
261
//- Runtime type information
262
TypeName
(
"KinematicParcel"
);
263
264
265
friend
class
Cloud
<ParcelType>;
266
267
268
// Constructors
269
270
//- Construct from owner, position, and cloud owner
271
// Other properties initialised as null
272
inline
KinematicParcel
273
(
274
KinematicCloud<ParcelType>
& owner,
275
const
vector
&
position
,
276
const
label cellI
277
);
278
279
//- Construct from components
280
inline
KinematicParcel
281
(
282
KinematicCloud<ParcelType>
& owner,
283
const
vector
&
position
,
284
const
label cellI,
285
const
label
typeId
,
286
const
scalar nParticle0,
287
const
scalar d0,
288
const
vector
& U0,
289
const
constantProperties& constProps
290
);
291
292
//- Construct from Istream
293
KinematicParcel
294
(
295
const
Cloud<ParcelType>
& c,
296
Istream
& is,
297
bool
readFields
=
true
298
);
299
300
//- Construct as a copy
301
KinematicParcel
(
const
KinematicParcel
&
p
);
302
303
//- Construct and return a clone
304
autoPtr<KinematicParcel>
clone
()
const
305
{
306
return
autoPtr<KinematicParcel>
(
new
KinematicParcel
(*
this
));
307
}
308
309
310
// Member Functions
311
312
// Access
313
314
//- Return const access to active flag
315
inline
bool
active
()
const
;
316
317
//- Return const access to type id
318
inline
label
typeId
()
const
;
319
320
//- Return const access to number of particles
321
inline
scalar
nParticle
()
const
;
322
323
//- Return const access to diameter
324
inline
scalar
d
()
const
;
325
326
//- Return const access to velocity
327
inline
const
vector
&
U
()
const
;
328
329
//- Return const access to density
330
inline
scalar
rho
()
const
;
331
332
//- Return const access to time spent in turbulent eddy
333
inline
scalar
tTurb
()
const
;
334
335
//- Return const access to turbulent velocity fluctuation
336
inline
const
vector
&
UTurb
()
const
;
337
338
339
// Edit
340
341
//- Return const access to active flag
342
inline
bool
&
active
();
343
344
//- Return access to type id
345
inline
label
typeId
();
346
347
//- Return access to number of particles
348
inline
scalar&
nParticle
();
349
350
//- Return access to diameter
351
inline
scalar&
d
();
352
353
//- Return access to velocity
354
inline
vector
&
U
();
355
356
//- Return access to density
357
inline
scalar&
rho
();
358
359
//- Return access to time spent in turbulent eddy
360
inline
scalar&
tTurb
();
361
362
//- Return access to turbulent velocity fluctuation
363
inline
vector
&
UTurb
();
364
365
366
// Helper functions
367
368
//- The nearest distance to a wall that
369
// the particle can be in the n direction
370
inline
scalar
wallImpactDistance
(
const
vector
& n)
const
;
371
372
//- Return the index of the face to be used in the interpolation
373
// routine
374
inline
label
faceInterpolation
()
const
;
375
376
//- Cell owner mass
377
inline
scalar
massCell
(
const
label cellI)
const
;
378
379
//- Particle mass
380
inline
scalar
mass
()
const
;
381
382
//- Particle volume
383
inline
scalar
volume
()
const
;
384
385
//- Particle volume for a given diameter
386
inline
scalar
volume
(
const
scalar
d
)
const
;
387
388
//- Particle projected area
389
inline
scalar
areaP
()
const
;
390
391
//- Projected area for given diameter
392
inline
scalar
areaP
(
const
scalar
d
)
const
;
393
394
//- Particle surface area
395
inline
scalar
areaS
()
const
;
396
397
//- Surface area for given diameter
398
inline
scalar
areaS
(
const
scalar
d
)
const
;
399
400
//- Reynolds number
401
inline
scalar
Re
402
(
403
const
vector
&
U
,
// particle velocity
404
const
scalar
d
,
// particle diameter
405
const
scalar rhoc,
// carrier density
406
const
scalar muc
// carrier dynamic viscosity
407
)
const
;
408
409
410
// Main calculation loop
411
412
//- Set cell values
413
template
<
class
TrackData>
414
void
setCellValues
415
(
416
TrackData& td,
417
const
scalar dt,
418
const
label cellI
419
);
420
421
//- Correct cell values using latest transfer information
422
template
<
class
TrackData>
423
void
cellValueSourceCorrection
424
(
425
TrackData& td,
426
const
scalar dt,
427
const
label cellI
428
);
429
430
//- Update parcel properties over the time interval
431
template
<
class
TrackData>
432
void
calc
433
(
434
TrackData& td,
435
const
scalar dt,
436
const
label cellI
437
);
438
439
440
// Tracking
441
442
//- Move the parcel
443
template
<
class
TrackData>
444
bool
move
(TrackData& td);
445
446
447
// Patch interactions
448
449
//- Overridable function to handle the particle hitting a patch
450
// Executed before other patch-hitting functions
451
template
<
class
TrackData>
452
bool
hitPatch
453
(
454
const
polyPatch
&
p
,
455
TrackData& td,
456
const
label patchI
457
);
458
459
460
//- Overridable function to handle the particle hitting a patch
461
// Executed before other patch-hitting functions without trackData
462
bool
hitPatch
463
(
464
const
polyPatch
&
p
,
465
int
& td,
466
const
label patchI
467
);
468
469
470
//- Overridable function to handle the particle hitting a
471
// processorPatch
472
template
<
class
TrackData>
473
void
hitProcessorPatch
474
(
475
const
processorPolyPatch
&,
476
TrackData& td
477
);
478
479
//- Overridable function to handle the particle hitting a
480
// processorPatch without trackData
481
void
hitProcessorPatch
482
(
483
const
processorPolyPatch
&,
484
int
&
485
);
486
487
//- Overridable function to handle the particle hitting a wallPatch
488
template
<
class
TrackData>
489
void
hitWallPatch
490
(
491
const
wallPolyPatch
&,
492
TrackData& td
493
);
494
495
//- Overridable function to handle the particle hitting a wallPatch
496
// without trackData
497
void
hitWallPatch
498
(
499
const
wallPolyPatch
&,
500
int
&
501
);
502
503
//- Overridable function to handle the particle hitting a polyPatch
504
template
<
class
TrackData>
505
void
hitPatch
506
(
507
const
polyPatch
&,
508
TrackData& td
509
);
510
511
//- Overridable function to handle the particle hitting a polyPatch
512
//- without trackData
513
void
hitPatch
514
(
515
const
polyPatch
&,
516
int
&
517
);
518
519
//- Transform the physical properties of the particle
520
// according to the given transformation tensor
521
void
transformProperties
(
const
tensor
&
T
);
522
523
//- Transform the physical properties of the particle
524
// according to the given separation vector
525
void
transformProperties
(
const
vector
& separation);
526
527
528
// I-O
529
530
//- Read
531
static
void
readFields
(
Cloud<ParcelType>
& c);
532
533
//- Write
534
static
void
writeFields
(
const
Cloud<ParcelType>
& c);
535
536
537
// Ostream Operator
538
539
friend
Ostream
& operator<< <ParcelType>
540
(
541
Ostream
&,
542
const
KinematicParcel<ParcelType>
&
543
);
544
};
545
546
547
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
548
549
}
// End namespace Foam
550
551
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
552
553
#include <
lagrangianIntermediate/KinematicParcelI.H
>
554
555
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
556
557
#define defineParcelTypeNameAndDebug(Type, DebugSwitch) \
558
template<> \
559
const Foam::word KinematicParcel<Type>::typeName(#Type); \
560
template<> \
561
int KinematicParcel<Type>::debug \
562
( \
563
Foam::debug::debugSwitch(#Type, DebugSwitch) \
564
);
565
566
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
567
568
#ifdef NoRepository
569
#include <
lagrangianIntermediate/KinematicParcel.C
>
570
#endif
571
572
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
573
574
#endif
575
576
// ************************ vim: set sw=4 sts=4 et: ************************ //
577
578