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
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 | 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
\*---------------------------------------------------------------------------*/
25
26
#include <
OpenFOAM/mathematicalConstants.H
>
27
28
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
29
30
template
<
class
ParcelType>
31
inline
Foam::KinematicParcel<ParcelType>::constantProperties::constantProperties
32
(
33
const
dictionary
& parentDict
34
)
35
:
36
dict_(parentDict.
subDict
(
"constantProperties"
)),
37
rhoMin_(
dimensionedScalar
(dict_.lookup(
"rhoMin"
)).value()),
38
rho0_(
dimensionedScalar
(dict_.lookup(
"rho0"
)).value()),
39
minParticleMass_
40
(
41
dimensionedScalar
(dict_.lookup(
"minParticleMass"
)).value()
42
)
43
{}
44
45
46
template
<
class
ParcelType>
47
inline
Foam::KinematicParcel<ParcelType>::trackData::trackData
48
(
49
KinematicCloud<ParcelType>
&
cloud
,
50
const
constantProperties
& constProps,
51
const
interpolation<scalar>
& rhoInterp,
52
const
interpolation<vector>
& UInterp,
53
const
interpolation<scalar>
& muInterp,
54
const
vector
& g
55
)
56
:
57
Particle<ParcelType>::trackData
(cloud),
58
cloud_(cloud),
59
constProps_(constProps),
60
rhoInterp_(rhoInterp),
61
UInterp_(UInterp),
62
muInterp_(muInterp),
63
g_(g)
64
{}
65
66
67
template
<
class
ParcelType>
68
inline
Foam::KinematicParcel<ParcelType>::KinematicParcel
69
(
70
KinematicCloud<ParcelType>
& owner,
71
const
vector
& position,
72
const
label cellI
73
)
74
:
75
Particle<ParcelType>
(owner, position, cellI),
76
active_(
true
),
77
typeId_(owner.parcelTypeId()),
78
nParticle_(0),
79
d_(0.0),
80
U_(vector::zero),
81
rho_(0.0),
82
tTurb_(0.0),
83
UTurb_(vector::zero),
84
rhoc_(0.0),
85
Uc_(vector::zero),
86
muc_(0.0)
87
{}
88
89
90
template
<
class
ParcelType>
91
inline
Foam::KinematicParcel<ParcelType>::KinematicParcel
92
(
93
KinematicCloud<ParcelType>
& owner,
94
const
vector
& position,
95
const
label cellI,
96
const
label typeId,
97
const
scalar nParticle0,
98
const
scalar d0,
99
const
vector
& U0,
100
const
constantProperties
& constProps
101
)
102
:
103
Particle<ParcelType>
(owner, position, cellI),
104
active_(
true
),
105
typeId_(typeId),
106
nParticle_(nParticle0),
107
d_(d0),
108
U_(U0),
109
rho_(constProps.
rho0
()),
110
tTurb_(0.0),
111
UTurb_(vector::zero),
112
rhoc_(0.0),
113
Uc_(vector::zero),
114
muc_(0.0)
115
{}
116
117
118
// * * * * * * * * * constantProperties Member Functions * * * * * * * * * * //
119
120
template
<
class
ParcelType>
121
inline
const
Foam::dictionary
&
122
Foam::KinematicParcel<ParcelType>::constantProperties::dict
()
const
123
{
124
return
dict_;
125
}
126
127
128
template
<
class
ParcelType>
129
inline
Foam::scalar
130
Foam::KinematicParcel<ParcelType>::constantProperties::rhoMin
()
const
131
{
132
return
rhoMin_;
133
}
134
135
136
template
<
class
ParcelType>
137
inline
Foam::scalar
138
Foam::KinematicParcel<ParcelType>::constantProperties::rho0
()
const
139
{
140
return
rho0_;
141
}
142
143
144
template
<
class
ParcelType>
145
inline
Foam::scalar
146
Foam::KinematicParcel<ParcelType>::constantProperties::minParticleMass
()
const
147
{
148
return
minParticleMass_;
149
}
150
151
152
// * * * * * * * * * * * trackData Member Functions * * * * * * * * * * * * //
153
154
template
<
class
ParcelType>
155
inline
Foam::KinematicCloud<ParcelType>
&
156
Foam::KinematicParcel<ParcelType>::trackData::cloud
()
157
{
158
return
cloud_
;
159
}
160
161
162
template
<
class
ParcelType>
163
inline
const
typename
Foam::KinematicParcel<ParcelType>::constantProperties
&
164
Foam::KinematicParcel<ParcelType>::trackData::constProps
()
const
165
{
166
return
constProps_;
167
}
168
169
170
template
<
class
ParcelType>
171
inline
const
Foam::interpolation<Foam::scalar>
&
172
Foam::KinematicParcel<ParcelType>::trackData::rhoInterp
()
const
173
{
174
return
rhoInterp_;
175
}
176
177
178
template
<
class
ParcelType>
179
inline
const
Foam::interpolation<Foam::vector>
&
180
Foam::KinematicParcel<ParcelType>::trackData::UInterp
()
const
181
{
182
return
UInterp_;
183
}
184
185
186
template
<
class
ParcelType>
187
inline
const
Foam::interpolation<Foam::scalar>
&
188
Foam::KinematicParcel<ParcelType>::trackData::muInterp
()
const
189
{
190
return
muInterp_;
191
}
192
193
194
template
<
class
ParcelType>
195
inline
const
Foam::vector
&
196
Foam::KinematicParcel<ParcelType>::trackData::g
()
const
197
{
198
return
g_;
199
}
200
201
202
// * * * * * * * * * * KinematicParcel Member Functions * * * * * * * * * * //
203
204
template
<
class
ParcelType>
205
inline
bool
Foam::KinematicParcel<ParcelType>::active
()
const
206
{
207
return
active_
;
208
}
209
210
211
template
<
class
ParcelType>
212
inline
Foam::label
Foam::KinematicParcel<ParcelType>::typeId
()
const
213
{
214
return
typeId_
;
215
}
216
217
218
template
<
class
ParcelType>
219
inline
Foam::scalar
Foam::KinematicParcel<ParcelType>::nParticle
()
const
220
{
221
return
nParticle_
;
222
}
223
224
225
template
<
class
ParcelType>
226
inline
Foam::scalar
Foam::KinematicParcel<ParcelType>::d
()
const
227
{
228
return
d_
;
229
}
230
231
232
template
<
class
ParcelType>
233
inline
const
Foam::vector
&
Foam::KinematicParcel<ParcelType>::U
()
const
234
{
235
return
U_
;
236
}
237
238
239
template
<
class
ParcelType>
240
inline
Foam::scalar
Foam::KinematicParcel<ParcelType>::rho
()
const
241
{
242
return
rho_
;
243
}
244
245
246
template
<
class
ParcelType>
247
inline
Foam::scalar
Foam::KinematicParcel<ParcelType>::tTurb
()
const
248
{
249
return
tTurb_
;
250
}
251
252
253
template
<
class
ParcelType>
254
inline
const
Foam::vector
&
Foam::KinematicParcel<ParcelType>::UTurb
()
const
255
{
256
return
UTurb_
;
257
}
258
259
260
template
<
class
ParcelType>
261
inline
bool
&
Foam::KinematicParcel<ParcelType>::active
()
262
{
263
return
active_
;
264
}
265
266
267
template
<
class
ParcelType>
268
inline
Foam::label
Foam::KinematicParcel<ParcelType>::typeId
()
269
{
270
return
typeId_
;
271
}
272
273
274
template
<
class
ParcelType>
275
inline
Foam::scalar&
Foam::KinematicParcel<ParcelType>::nParticle
()
276
{
277
return
nParticle_
;
278
}
279
280
281
template
<
class
ParcelType>
282
inline
Foam::scalar&
Foam::KinematicParcel<ParcelType>::d
()
283
{
284
return
d_
;
285
}
286
287
288
template
<
class
ParcelType>
289
inline
Foam::vector
&
Foam::KinematicParcel<ParcelType>::U
()
290
{
291
return
U_
;
292
}
293
294
295
template
<
class
ParcelType>
296
inline
Foam::scalar&
Foam::KinematicParcel<ParcelType>::rho
()
297
{
298
return
rho_
;
299
}
300
301
302
template
<
class
ParcelType>
303
inline
Foam::scalar&
Foam::KinematicParcel<ParcelType>::tTurb
()
304
{
305
return
tTurb_
;
306
}
307
308
309
template
<
class
ParcelType>
310
inline
Foam::vector
&
Foam::KinematicParcel<ParcelType>::UTurb
()
311
{
312
return
UTurb_
;
313
}
314
315
316
template
<
class
ParcelType>
317
inline
Foam::scalar
Foam::KinematicParcel<ParcelType>::wallImpactDistance
318
(
319
const
vector
&
320
)
const
321
{
322
return
0.5*
d_
;
323
}
324
325
326
template
<
class
ParcelType>
327
inline
Foam::label
Foam::KinematicParcel<ParcelType>::faceInterpolation
()
const
328
{
329
// Use volume-based interpolation if dealing with external faces
330
if
(this->
cloud
().internalFace(this->
face
()))
331
{
332
return
this->
face
();
333
}
334
else
335
{
336
return
-1;
337
}
338
}
339
340
341
template
<
class
ParcelType>
342
inline
Foam::scalar
Foam::KinematicParcel<ParcelType>::massCell
343
(
344
const
label cellI
345
)
const
346
{
347
return
rhoc_
*this->
cloud
().
pMesh
().
cellVolumes
()[cellI];
348
}
349
350
351
template
<
class
ParcelType>
352
inline
Foam::scalar
Foam::KinematicParcel<ParcelType>::mass
()
const
353
{
354
return
rho_
*
volume
();
355
}
356
357
358
template
<
class
ParcelType>
359
inline
Foam::scalar
Foam::KinematicParcel<ParcelType>::volume
()
const
360
{
361
return
volume
(
d_
);
362
}
363
364
365
template
<
class
ParcelType>
366
inline
Foam::scalar
367
Foam::KinematicParcel<ParcelType>::volume
(
const
scalar
d
)
const
368
{
369
return
mathematicalConstant::pi
/6.0*
pow3
(d);
370
}
371
372
373
template
<
class
ParcelType>
374
inline
Foam::scalar
Foam::KinematicParcel<ParcelType>::areaP
()
const
375
{
376
return
areaP
(
d_
);
377
}
378
379
380
template
<
class
ParcelType>
381
inline
Foam::scalar
382
Foam::KinematicParcel<ParcelType>::areaP
(
const
scalar
d
)
const
383
{
384
return
0.25*
areaS
(d);
385
}
386
387
388
template
<
class
ParcelType>
389
inline
Foam::scalar
Foam::KinematicParcel<ParcelType>::areaS
()
const
390
{
391
return
areaS
(
d_
);
392
}
393
394
395
template
<
class
ParcelType>
396
inline
Foam::scalar
397
Foam::KinematicParcel<ParcelType>::areaS
(
const
scalar
d
)
const
398
{
399
return
mathematicalConstant::pi
*d*
d
;
400
}
401
402
403
template
<
class
ParcelType>
404
inline
Foam::scalar
405
Foam::KinematicParcel<ParcelType>::Re
406
(
407
const
vector
&
U
,
408
const
scalar
d
,
409
const
scalar rhoc,
410
const
scalar muc
411
)
const
412
{
413
return
rhoc*
mag
(U -
Uc_
)*d/muc;
414
}
415
416
417
// ************************ vim: set sw=4 sts=4 et: ************************ //