streamlinesParticle.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-2022 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 Class
25  Foam::streamlinesParticle
26 
27 Description
28  Particle class that samples fields as it passes through. Used in streamlines
29  calculation.
30 
31 SourceFiles
32  streamlinesParticle.C
33 
34 \*---------------------------------------------------------------------------*/
35 
36 #ifndef streamlinesParticle_H
37 #define streamlinesParticle_H
38 
39 #include "particle.H"
40 #include "Cloud.H"
41 #include "autoPtr.H"
42 #include "interpolation.H"
43 #include "vectorList.H"
44 #include "DynamicField.H"
45 
46 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
47 
48 namespace Foam
49 {
50 
51 // Forward declaration of friend functions and operators
52 
53 class streamlinesParticle;
54 class streamlinesCloud;
55 
56 Ostream& operator<<(Ostream&, const streamlinesParticle&);
57 
58 /*---------------------------------------------------------------------------*\
59  Class streamlinesParticle Declaration
60 \*---------------------------------------------------------------------------*/
61 
63 :
64  public particle
65 {
66 public:
67 
68  class trackingData
69  :
71  {
72  public:
73 
74  // Public data
75 
76  #define DeclareTypeInterpolator(Type, nullArg) \
77  const PtrList<interpolation<Type>>& Type##Interp_;
79  #undef DeclareTypeInterpolator
80 
82 
83  bool trackForward_;
84 
85  bool trackOutside_;
86 
87  const label nSubCycle_;
88 
89  const scalar trackLength_;
90 
91  label tracki;
92 
94 
96 
98 
101  #define DeclareAllTypes(Type, nullArg) \
102  List<DynamicField<Type>>& all##Type##s_;
104  #undef DeclareAllTypes
105 
106 
107  // Constructors
108 
109  //- Construct from components
111  (
112  streamlinesCloud& cloud,
113  #define TypeInterpolatorArg(Type, nullArg) \
114  const PtrList<interpolation<Type>>& Type##Interp,
116  #undef TypeInterpolatorArg
117  const interpolation<vector>& UInterp,
118  const bool trackForward,
119  const bool trackOutside,
120  const label nSubCycle,
121  const scalar trackLength,
122  DynamicField<point>& allPositions,
123  DynamicField<label>& allTracks,
124  DynamicField<label>& allTrackParts,
125  DynamicField<scalar>& allAges
126  #define AllTypesArg(Type, nullArg) \
127  , List<DynamicField<Type>>& all##Type##s
128  FOR_ALL_FIELD_TYPES(AllTypesArg)
129  #undef AllTypesArg
130  )
131  :
132  particle::trackingData(cloud),
133  #define TypeInterpolatorInit(Type, nullArg) \
134  Type##Interp_(Type##Interp),
136  #undef TypeInterpolatorInit
137  UInterp_(UInterp),
138  trackForward_(trackForward),
139  trackOutside_(trackOutside),
140  nSubCycle_(nSubCycle),
141  trackLength_(trackLength),
142  allPositions_(allPositions),
143  allTracks_(allTracks),
144  allTrackParts_(allTrackParts),
145  allAges_(allAges)
146  #define AllTypesInit(Type, nullArg) \
147  , all##Type##s_(all##Type##s)
149  #undef AllTypesInit
150  {}
151  };
152 
153 
154 private:
155 
156  // Private Data
157 
158  //- Lifetime of particle. Particle dies when reaches 0.
159  label lifeTime_;
160 
161  //- Index of the track
162  label trackIndex_;
163 
164  //- Index of the part of the track
165  label trackPartIndex_;
166 
167  //- Age of the particle
168  scalar age_;
169 
170  //- Current compound transform
171  transformer transform_;
172 
173  //- Sampled positions
174  DynamicField<point> sampledPositions_;
175 
176  //- Sampled ages
177  DynamicField<scalar> sampledAges_;
178 
179  //- Sampled types
180  #define DeclareSampledTypes(Type, nullArg) \
181  List<DynamicField<Type>> sampled##Type##s_;
183  #undef DeclareSampledTypes
184 
185 
186  // Private Member Functions
187 
188  //- Interpolate all quantities; return interpolated velocity.
189  vector interpolateFields
190  (
191  const trackingData&,
192  const point&,
193  const label celli,
194  const label facei
195  );
196 
197  //- End the current track
198  void endTrack(trackingData&);
199 
200 
201 public:
202 
203  // Constructors
204 
205  //- Construct from components
207  (
208  const polyMesh& c,
209  const vector& position,
210  const label celli,
211  const label lifeTime,
212  const label trackIndex
213  );
214 
215  //- Construct from Istream
217  (
218  const polyMesh& c,
219  Istream& is,
220  bool readFields = true
221  );
222 
223  //- Construct copy
225 
226  //- Construct and return a clone
227  autoPtr<particle> clone() const
228  {
229  return autoPtr<particle>(new streamlinesParticle(*this));
230  }
231 
232  //- Factory class to read-construct particles used for parallel transfer
233  class iNew
234  {
235  const polyMesh& mesh_;
236 
237  public:
239  iNew(const polyMesh& mesh)
240  :
241  mesh_(mesh)
242  {}
244  autoPtr<streamlinesParticle> operator()(Istream& is) const
245  {
247  (
248  new streamlinesParticle(mesh_, is, true)
249  );
250  }
251  };
252 
253 
254  // Member Functions
255 
256  // Tracking
257 
258  //- Track all particles to their end point
259  bool move(streamlinesCloud&, trackingData&, const scalar);
260 
261  //- Overridable function to handle the particle hitting a wedge
263 
264  //- Overridable function to handle the particle hitting a
265  // symmetry plane
267 
268  //- Overridable function to handle the particle hitting a
269  // symmetry patch
271 
272  //- Overridable function to handle the particle hitting a cyclic
274 
275  //- Overridable function to handle the particle hitting a
276  // cyclicAMIPatch
277  void hitCyclicAMIPatch
278  (
279  const vector&,
280  const scalar,
282  trackingData&
283  );
284 
285  //- Overridable function to handle the particle hitting an
286  // nonConformalCyclicPolyPatch
288  (
289  const vector& displacement,
290  const scalar fraction,
291  const label patchi,
292  streamlinesCloud& cloud,
293  trackingData& td
294  );
295 
296  //- Overridable function to handle the particle hitting a
297  // processorPatch
299 
300  //- Overridable function to handle the particle hitting a wallPatch
302 
303 
304  // Transformations
305 
306  //- Transform the physical properties of the particle
307  // according to the given transformation tensor
308  virtual void transformProperties(const transformer&);
309 
310 
311  // I-O
312 
313  //- Read
315 
316  //- Write
317  static void writeFields(const Cloud<streamlinesParticle>&);
318 
319 
320  // Ostream Operator
321 
322  friend Ostream& operator<<(Ostream&, const streamlinesParticle&);
323 };
324 
325 
326 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
327 
328 } // End namespace Foam
329 
330 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
331 
332 #endif
333 
334 // ************************************************************************* //
#define TypeInterpolatorArg(Type, nullArg)
#define AllTypesInit(Type, nullArg)
void hitProcessorPatch(streamlinesCloud &, trackingData &)
Overridable function to handle the particle hitting a.
const polyMesh & mesh() const
Return the mesh database.
Definition: particleI.H:125
FvWallInfoData< WallInfo, label > label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
bool move(streamlinesCloud &, trackingData &, const scalar)
Track all particles to their end point.
Vector-tensor class used to perform translations, rotations and scaling operations in 3D space...
Definition: transformer.H:83
void hitCyclicAMIPatch(const vector &, const scalar, streamlinesCloud &, trackingData &)
Overridable function to handle the particle hitting a.
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: HashTable.H:59
void hitWallPatch(streamlinesCloud &, trackingData &)
Overridable function to handle the particle hitting a wallPatch.
FOR_ALL_FIELD_TYPES(DeclareTypeInterpolator)
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
trackingData(streamlinesCloud &cloud, #define TypeInterpolatorArg(Type, nullArg) const interpolation< vector > &UInterp, const bool trackForward, const bool trackOutside, const label nSubCycle, const scalar trackLength, DynamicField< point > &allPositions, DynamicField< label > &allTracks, DynamicField< label > &allTrackParts, DynamicField< scalar > &allAges #define AllTypesArg(Type, nullArg))
Construct from components.
Vector< scalar > vector
A scalar version of the templated Vector.
Definition: vector.H:49
const dimensionedScalar c
Speed of light in a vacuum.
Base particle class.
Definition: particle.H:81
static void writeFields(const Cloud< streamlinesParticle > &)
Write.
#define TypeInterpolatorInit(Type, nullArg)
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().x()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().y()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
#define DeclareSampledTypes(Type, nullArg)
Sampled types.
A Cloud of streamlines particles.
#define AllTypesArg(Type, nullArg)
const interpolation< vector > & UInterp_
Dynamically sized Field.
Definition: DynamicField.H:49
void hitWedgePatch(streamlinesCloud &, trackingData &)
Overridable function to handle the particle hitting a wedge.
#define DeclareAllTypes(Type, nullArg)
Particle class that samples fields as it passes through. Used in streamlines calculation.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
virtual void transformProperties(const transformer &)
Transform the physical properties of the particle.
void hitCyclicPatch(streamlinesCloud &, trackingData &)
Overridable function to handle the particle hitting a cyclic.
bool hitNonConformalCyclicPatch(const vector &displacement, const scalar fraction, const label patchi, streamlinesCloud &cloud, trackingData &td)
Overridable function to handle the particle hitting an.
label patchi
void hitSymmetryPatch(streamlinesCloud &, trackingData &)
Overridable function to handle the particle hitting a.
vector point
Point is a vector.
Definition: point.H:41
Factory class to read-construct particles used for parallel transfer.
autoPtr< particle > clone() const
Construct and return a clone.
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:70
Ostream & operator<<(Ostream &, const ensightPart &)
Abstract base class for interpolation.
static void readFields(Cloud< streamlinesParticle > &)
Read.
friend Ostream & operator<<(Ostream &, const streamlinesParticle &)
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:76
volScalarField & p
void hitSymmetryPlanePatch(streamlinesCloud &, trackingData &)
Overridable function to handle the particle hitting a.
#define DeclareTypeInterpolator(Type, nullArg)
streamlinesParticle(const polyMesh &c, const vector &position, const label celli, const label lifeTime, const label trackIndex)
Construct from components.
vector position() const
Return current particle position.
Definition: particleI.H:278
Namespace for OpenFOAM.