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
basic
Particle
Particle.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::Particle
26
27
Description
28
29
\*---------------------------------------------------------------------------*/
30
31
#ifndef Particle_H
32
#define Particle_H
33
34
#include <
OpenFOAM/vector.H
>
35
#include <
OpenFOAM/IDLList.H
>
36
#include <
OpenFOAM/labelList.H
>
37
#include <
OpenFOAM/pointField.H
>
38
#include <
OpenFOAM/faceList.H
>
39
#include <
OpenFOAM/typeInfo.H
>
40
#include <
OpenFOAM/OFstream.H
>
41
42
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
43
44
namespace
Foam
45
{
46
47
template
<
class
Particle>
48
class
Cloud;
49
50
class
wedgePolyPatch;
51
class
symmetryPolyPatch;
52
class
cyclicPolyPatch;
53
class
processorPolyPatch;
54
class
wallPolyPatch;
55
class
polyPatch;
56
57
// Forward declaration of friend functions and operators
58
59
template
<
class
ParticleType>
60
class
Particle;
61
62
template
<
class
ParticleType>
63
Ostream&
operator
<<
64
(
65
Ostream&,
66
const
Particle<ParticleType>&
67
);
68
69
70
/*---------------------------------------------------------------------------*\
71
Class Particle Declaration
72
\*---------------------------------------------------------------------------*/
73
74
template
<
class
ParticleType>
75
class
Particle
76
:
77
public
IDLList
<ParticleType>
::link
78
{
79
80
public
:
81
82
//- Class used to pass tracking data to the trackToFace function
83
class
trackData
84
{
85
86
// Private data
87
88
//- Reference to the cloud containing this particle
89
Cloud<ParticleType>
& cloud_;
90
91
92
public
:
93
94
bool
switchProcessor
;
95
bool
keepParticle
;
96
97
98
// Constructors
99
100
inline
trackData
(
Cloud<ParticleType>
&
cloud
);
101
102
103
// Member functions
104
105
//- Return a reference to the cloud
106
inline
Cloud<ParticleType>
&
cloud
();
107
};
108
109
110
protected
:
111
112
// Private data
113
114
//- Reference to the particle cloud
115
const
Cloud<ParticleType>
&
cloud_
;
116
117
//- Position of particle
118
vector
position_
;
119
120
//- Index of the cell it is in
121
label
celli_
;
122
123
//- Face index if the particle is on a face otherwise -1
124
label
facei_
;
125
126
//- Fraction of time-step completed
127
scalar
stepFraction_
;
128
129
//- Originating processor id
130
label
origProc_
;
131
132
//- Local particle id on originating processor
133
label
origId_
;
134
135
136
// Private member functions
137
138
//- Return the 'lambda' value for the position, p, on the face,
139
// where, p = from + lamda*(to - from)
140
// for non-static meshes
141
inline
scalar
lambda
142
(
143
const
vector
& from,
144
const
vector
& to,
145
const
label facei,
146
const
scalar
stepFraction
147
)
const
;
148
149
//- Return the 'lambda' value for the position, p, on the face,
150
// where, p = from + lamda*(to - from)
151
// for static meshes
152
inline
scalar
lambda
153
(
154
const
vector
& from,
155
const
vector
& to,
156
const
label facei
157
)
const
;
158
159
//- Find the faces between position and cell centre
160
void
findFaces
161
(
162
const
vector
&
position
,
163
DynamicList<label>
&
faceList
164
)
const
;
165
166
//- Find the faces between position and cell centre
167
void
findFaces
168
(
169
const
vector
&
position
,
170
const
label celli,
171
const
scalar
stepFraction
,
172
DynamicList<label>
&
faceList
173
)
const
;
174
175
176
// Patch interactions
177
178
//- Overridable function to handle the particle hitting a patch
179
// Executed before other patch-hitting functions
180
template
<
class
TrackData>
181
bool
hitPatch
182
(
183
const
polyPatch
&,
184
TrackData& td,
185
const
label patchI
186
);
187
188
//- Overridable function to handle the particle hitting a wedgePatch
189
template
<
class
TrackData>
190
void
hitWedgePatch
191
(
192
const
wedgePolyPatch
&,
193
TrackData& td
194
);
195
196
//- Overridable function to handle the particle hitting a
197
// symmetryPatch
198
template
<
class
TrackData>
199
void
hitSymmetryPatch
200
(
201
const
symmetryPolyPatch
&,
202
TrackData& td
203
);
204
205
//- Overridable function to handle the particle hitting a cyclicPatch
206
template
<
class
TrackData>
207
void
hitCyclicPatch
208
(
209
const
cyclicPolyPatch
&,
210
TrackData& td
211
);
212
213
//- Overridable function to handle the particle hitting a
214
// processorPatch
215
template
<
class
TrackData>
216
void
hitProcessorPatch
217
(
218
const
processorPolyPatch
&,
219
TrackData& td
220
);
221
222
//- Overridable function to handle the particle hitting a wallPatch
223
template
<
class
TrackData>
224
void
hitWallPatch
225
(
226
const
wallPolyPatch
&,
227
TrackData& td
228
);
229
230
//- Overridable function to handle the particle hitting a
231
// general patch
232
template
<
class
TrackData>
233
void
hitPatch
234
(
235
const
polyPatch
&,
236
TrackData& td
237
);
238
239
240
// Transformations
241
242
//- Transform the position the particle
243
// according to the given transformation tensor
244
virtual
void
transformPosition
(
const
tensor
&
T
);
245
246
//- Transform the physical properties of the particle
247
// according to the given transformation tensor
248
virtual
void
transformProperties
(
const
tensor
&
T
);
249
250
//- Transform the physical properties of the particle
251
// according to the given separation vector
252
virtual
void
transformProperties
(
const
vector
& separation);
253
254
255
// Parallel transfer
256
257
//- Convert global addressing to the processor patch
258
// local equivalents
259
template
<
class
TrackData>
260
void
prepareForParallelTransfer
(
const
label
patchi
, TrackData& td);
261
262
//- Convert processor patch addressing to the global equivalents
263
// and set the celli to the face-neighbour
264
template
<
class
TrackData>
265
void
correctAfterParallelTransfer
(
const
label
patchi
, TrackData& td);
266
267
268
public
:
269
270
friend
class
Cloud
<ParticleType>;
271
272
273
// Static data members
274
275
//- Runtime type information
276
TypeName
(
"Particle"
);
277
278
//- String representation of properties
279
static
string
propHeader
;
280
281
282
// Constructors
283
284
//- Construct from components
285
Particle
286
(
287
const
Cloud<ParticleType>
&,
288
const
vector
&
position
,
289
const
label celli
290
);
291
292
//- Construct from Istream
293
Particle
294
(
295
const
Cloud<ParticleType>
&,
296
Istream
&,
297
bool
readFields
=
true
298
);
299
300
//- Construct as a copy
301
Particle
(
const
Particle
&
p
);
302
303
//- Construct a clone
304
autoPtr<ParticleType>
clone
()
const
305
{
306
return
autoPtr<Particle>
(
new
Particle
(*
this
));
307
}
308
309
310
//- Factory class to read-construct particles used for
311
// parallel transfer
312
class
iNew
313
{
314
// Private data
315
316
//- Reference to the cloud
317
const
Cloud<ParticleType>
& cloud_;
318
319
320
public
:
321
322
iNew
(
const
Cloud<ParticleType>
&
cloud
)
323
:
324
cloud_(cloud)
325
{}
326
327
autoPtr<ParticleType>
operator()
(
Istream
& is)
const
328
{
329
return
autoPtr<ParticleType>
330
(
331
new
ParticleType
332
(
333
cloud_,
334
is,
335
true
336
)
337
);
338
}
339
};
340
341
342
//- Destructor
343
virtual
~Particle
()
344
{}
345
346
347
// Member Functions
348
349
// Access
350
351
//- Return true if particle is in cell
352
inline
bool
inCell
()
const
;
353
354
//- Return true if position is in cell i
355
inline
bool
inCell
356
(
357
const
vector
&
position
,
358
const
label celli,
359
const
scalar
stepFraction
360
)
const
;
361
362
//- Return current particle position
363
inline
const
vector
&
position
()
const
;
364
365
//- Return current particle position
366
inline
vector
&
position
();
367
368
//- Return current cell particle is in
369
inline
label&
cell
();
370
371
//- Return current cell particle is in
372
inline
label
cell
()
const
;
373
374
//- Return current face particle is on otherwise -1
375
inline
label
face
()
const
;
376
377
//- Return reference to the particle cloud
378
inline
const
Cloud<ParticleType>
&
cloud
()
const
;
379
380
//- Return the impact model to be used, soft or hard (default).
381
inline
bool
softImpact
()
const
;
382
383
//- Return the particle current time
384
inline
scalar
currentTime
()
const
;
385
386
387
// Check
388
389
//- Is the particle on the boundary/(or outside the domain)?
390
inline
bool
onBoundary
()
const
;
391
392
//- Which patch is particle on
393
inline
label
patch
(
const
label facei)
const
;
394
395
//- Which face of this patch is this particle on
396
inline
label
patchFace
397
(
398
const
label
patchi
,
399
const
label facei
400
)
const
;
401
402
//- The nearest distance to a wall that
403
// the particle can be in the n direction
404
inline
scalar
wallImpactDistance
(
const
vector
& n)
const
;
405
406
//- Return the fraction of time-step completed
407
inline
scalar&
stepFraction
();
408
409
//- Return the fraction of time-step completed
410
inline
scalar
stepFraction
()
const
;
411
412
//- Return the originating processor id
413
inline
label
origProc
()
const
;
414
415
//- Return the particle id on originating processor
416
inline
label
origId
()
const
;
417
418
419
// Track
420
421
//- Track particle to end of trajectory
422
// or until it hits the boundary.
423
// On entry 'stepFraction()' should be set to the fraction of the
424
// time-step at which the tracking starts and on exit it contains
425
// the fraction of the time-step completed.
426
// Returns the boundary face index if the track stops at the
427
// boundary, -1 otherwise.
428
template
<
class
TrackData>
429
label
track
430
(
431
const
vector
& endPosition,
432
TrackData& td
433
);
434
435
//- Calls the templated track with dummy TrackData
436
label
track
(
const
vector
& endPosition);
437
438
//- Track particle to a given position and returns 1.0 if the
439
// trajectory is completed without hitting a face otherwise
440
// stops at the face and returns the fraction of the trajectory
441
// completed.
442
// on entry 'stepFraction()' should be set to the fraction of the
443
// time-step at which the tracking starts.
444
template
<
class
TrackData>
445
scalar
trackToFace
446
(
447
const
vector
& endPosition,
448
TrackData& td
449
);
450
451
//- Calls the templated trackToFace with dummy TrackData
452
scalar
trackToFace
(
const
vector
& endPosition);
453
454
//- Return the index of the face to be used in the interpolation
455
// routine
456
inline
label
faceInterpolation
()
const
;
457
458
459
// I-O
460
461
//- Read the fields associated with the owner cloud
462
static
void
readFields
(
Cloud<ParticleType>
& c);
463
464
//- Write the fields associated with the owner cloud
465
static
void
writeFields
(
const
Cloud<ParticleType>
& c);
466
467
//- Write the particle data
468
void
write
(
Ostream
& os,
bool
writeFields
)
const
;
469
470
// Ostream Operator
471
472
friend
Ostream
& operator<< <ParticleType>
473
(
474
Ostream
&,
475
const
Particle<ParticleType>
&
476
);
477
};
478
479
480
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
481
482
}
// End namespace Foam
483
484
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
485
486
#include <
lagrangian/ParticleI.H
>
487
488
#define defineParticleTypeNameAndDebug(Type, DebugSwitch) \
489
template<> \
490
const Foam::word Particle<Type>::typeName(#Type); \
491
template<> \
492
int Particle<Type>::debug(Foam::debug::debugSwitch(#Type, DebugSwitch));
493
494
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
495
496
#ifdef NoRepository
497
# include <
lagrangian/Particle.C
>
498
#endif
499
500
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
501
502
#endif
503
504
// ************************ vim: set sw=4 sts=4 et: ************************ //