All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
particle.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::particle
26 
27 Description
28  Base particle class
29 
30 \*---------------------------------------------------------------------------*/
31 
32 #ifndef particle_H
33 #define particle_H
34 
35 #include "vector.H"
36 #include "barycentric.H"
37 #include "barycentricTensor.H"
38 #include "IDLList.H"
39 #include "pointField.H"
40 #include "faceList.H"
41 #include "OFstream.H"
42 #include "tetPointRef.H"
43 #include "FixedList.H"
45 #include "particleMacros.H"
46 #include "transformer.H"
47 
48 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
49 
50 namespace Foam
51 {
52 
53 // Forward declaration of classes
54 class particle;
55 
56 class polyPatch;
57 
58 class cyclicPolyPatch;
59 class cyclicAMIPolyPatch;
60 class processorPolyPatch;
61 class symmetryPlanePolyPatch;
62 class symmetryPolyPatch;
63 class wallPolyPatch;
64 class wedgePolyPatch;
65 
66 // Forward declaration of friend functions and operators
67 
68 Ostream& operator<<
69 (
70  Ostream&,
71  const particle&
72 );
73 
74 bool operator==(const particle&, const particle&);
75 
76 bool operator!=(const particle&, const particle&);
77 
78 /*---------------------------------------------------------------------------*\
79  Class Particle Declaration
80 \*---------------------------------------------------------------------------*/
81 
82 class particle
83 :
84  public IDLList<particle>::link
85 {
86  // Private member data
87 
88  //- Size in bytes of the position data
89  static const std::size_t sizeofPosition_;
90 
91  //- Size in bytes of the fields
92  static const std::size_t sizeofFields_;
93 
94  //- The value of nTracksBehind_ at which tracking is abandoned. See the
95  // description of nTracksBehind_.
96  static const label maxNTracksBehind_;
97 
98 
99 public:
101  class trackingData
102  {
103  public:
104 
105  // Public data
106 
107  //- Flag to indicate whether to keep particle (false = delete)
108  bool keepParticle;
109 
110  //- Processor to send the particle to. -1 indicates that this
111  // particle is not to be transferred.
113 
114  //- Patch from which to send the particle
116 
117  //- Patch to which to send the particle
119 
120  //- Patch face to which to send the particle
122 
123 
124  // Constructor
125  template <class TrackCloudType>
126  trackingData(const TrackCloudType& cloud)
127  :
128  keepParticle(false),
129  sendToProc(-1),
130  sendFromPatch(-1),
131  sendToPatch(-1),
132  sendToPatchFace(-1)
133  {}
134  };
135 
136 
137 private:
138 
139  // Private Data
140 
141  //- Reference to the polyMesh database
142  const polyMesh& mesh_;
143 
144  //- Coordinates of particle
145  barycentric coordinates_;
146 
147  //- Index of the cell it is in
148  label celli_;
149 
150  //- Index of the face that owns the decomposed tet that the
151  // particle is in
152  label tetFacei_;
153 
154  //- Index of the point on the face that defines the decomposed
155  // tet that the particle is in. Relative to the face base
156  // point.
157  label tetPti_;
158 
159  //- Face index if the particle is on a face otherwise -1
160  label facei_;
161 
162  //- Fraction of time-step completed
163  scalar stepFraction_;
164 
165  //- The step fraction less than the maximum reached so far
166  scalar stepFractionBehind_;
167 
168  //- The number of tracks carried out that ended in a step fraction less
169  // than the maximum reached so far. Once this reaches
170  // maxNTracksBehind_, tracking is abandoned for the current step.
171  //
172  // This is needed because when tetrahedra are inverted a straight
173  // trajectory can form a closed loop through regions of overlapping
174  // positive and negative space. Without this break clause, such loops
175  // can result in a valid track which never ends.
176  //
177  // Because the test is susceptible to round off error, a track of zero
178  // length will also increment the counter. As such, it is important
179  // that maxNTracksBehind_ is set large enough so that valid small
180  // tracks do not result in the track being abandoned. The largest
181  // number of valid small tracks that are likely to be performed
182  // sequentially is equal to the number of tetrahedra that can meet at
183  // a point. An estimate of this number is therefore used to set
184  // maxNTracksBehind_.
185  label nTracksBehind_;
186 
187  //- Originating processor id
188  label origProc_;
189 
190  //- Local particle id on originating processor
191  label origId_;
192 
193 
194 private:
195 
196  // Private Member Functions
197 
198  // Tetrahedra functions
199 
200  //- Get the vertices of the current tet
201  inline void stationaryTetGeometry
202  (
203  vector& centre,
204  vector& base,
205  vector& vertex1,
206  vector& vertex2
207  ) const;
208 
209  //- Get the transformation associated with the current tet. This
210  // will convert a barycentric position within the tet to a
211  // cartesian position in the global coordinate system. The
212  // conversion is x = A & y, where x is the cartesian position, y is
213  // the barycentric position and A is the transformation tensor.
214  inline barycentricTensor stationaryTetTransform() const;
215 
216  //- Get the reverse transform associated with the current tet. The
217  // conversion is detA*y = (x - centre) & T. The variables x, y and
218  // centre have the same meaning as for the forward transform. T is
219  // the transposed inverse of the forward transform tensor, A,
220  // multiplied by its determinant, detA. This separation allows
221  // the barycentric tracking algorithm to function on inverted or
222  // degenerate tetrahedra.
223  void stationaryTetReverseTransform
224  (
225  vector& centre,
226  scalar& detA,
228  ) const;
229 
230  //- Get the vertices of the current moving tet. Two values are
231  // returned for each vertex. The first is a constant, and the
232  // second is a linear coefficient of the track fraction.
233  inline void movingTetGeometry
234  (
235  const scalar endStepFraction,
236  Pair<vector>& centre,
237  Pair<vector>& base,
238  Pair<vector>& vertex1,
239  Pair<vector>& vertex2
240  ) const;
241 
242  //- Get the transformation associated with the current, moving, tet.
243  // This is of the same form as for the static case. As with the
244  // moving geometry, a linear function of the tracking fraction is
245  // returned for each component.
246  inline Pair<barycentricTensor> movingTetTransform
247  (
248  const scalar endStepFraction
249  ) const;
250 
251  //- Get the reverse transformation associated with the current,
252  // moving, tet. This is of the same form as for the static case. As
253  // with the moving geometry, a function of the tracking fraction is
254  // returned for each component. The functions are higher order than
255  // for the forward transform; the determinant is cubic, and the
256  // tensor is quadratic.
257  void movingTetReverseTransform
258  (
259  const scalar endStepFraction,
260  Pair<vector>& centre,
261  FixedList<scalar, 4>& detA,
263  ) const;
264 
265 
266  // Transformations
267 
268  //- Reflection transform. Corrects the coordinates when the particle
269  // moves between two tets which share a base vertex, but for which
270  // the other two non cell-centre vertices are reversed. All hits
271  // which retain the same face behave this way, as do face hits.
272  void reflect();
273 
274  //- Rotation transform. Corrects the coordinates when the particle
275  // moves between two tets with different base vertices, but are
276  // otherwise similarly oriented. Hits which change the face within
277  // the cell make use of both this and the reflect transform.
278  void rotate(const bool direction);
279 
280 
281  // Topology changes
282 
283  //- Change tet within a cell. Called after a triangle is hit.
284  void changeTet(const label tetTriI);
285 
286  //- Change tet face within a cell. Called by changeTet.
287  void changeFace(const label tetTriI);
288 
289  //- Change cell. Called when the particle hits an internal face.
290  void changeCell();
291 
292 
293  // Geometry changes
294 
295  //- Locate the particle at the given position
296  void locate
297  (
298  const vector& position,
299  label celli,
300  const bool boundaryFail,
301  const string boundaryMsg
302  );
303 
304 
305 protected:
306 
307  // Patch interactions
308 
309  //- Overridable function to handle the particle hitting a patch.
310  // Executed before other patch-hitting functions.
311  template<class TrackCloudType>
312  bool hitPatch(TrackCloudType&, trackingData&);
313 
314  //- Overridable function to handle the particle hitting a wedgePatch
315  template<class TrackCloudType>
316  void hitWedgePatch(TrackCloudType&, trackingData&);
317 
318  //- Overridable function to handle the particle hitting a
319  // symmetryPlanePatch
320  template<class TrackCloudType>
321  void hitSymmetryPlanePatch(TrackCloudType&, trackingData&);
322 
323  //- Overridable function to handle the particle hitting a symmetryPatch
324  template<class TrackCloudType>
325  void hitSymmetryPatch(TrackCloudType&, trackingData&);
326 
327  //- Overridable function to handle the particle hitting a cyclicPatch
328  template<class TrackCloudType>
329  void hitCyclicPatch(TrackCloudType&, trackingData&);
330 
331  //- Overridable function to handle the particle hitting a cyclicAMIPatch
332  template<class TrackCloudType>
333  void hitCyclicAMIPatch
334  (
335  const vector& displacement,
336  const scalar fraction,
337  TrackCloudType& cloud,
338  trackingData& td
339  );
340 
341  //- Overridable function to handle the particle hitting an
342  // nonConformalCyclicPolyPatch
343  template<class TrackCloudType>
345  (
346  const vector& displacement,
347  const scalar fraction,
348  const label patchi,
349  TrackCloudType& cloud,
350  trackingData& td
351  );
352 
353  //- Overridable function to handle the particle hitting a processorPatch
354  template<class TrackCloudType>
355  void hitProcessorPatch(TrackCloudType&, trackingData&);
356 
357  //- Overridable function to handle the particle hitting a wallPatch
358  template<class TrackCloudType>
359  void hitWallPatch(TrackCloudType&, trackingData&);
360 
361 
362 public:
363 
364  // Static Data Members
365 
366  //- Runtime type information
367  TypeName("particle");
368 
369  //- String representation of properties
371  (
372  "(coordinatesa coordinatesb coordinatesc coordinatesd) "
373  "celli tetFacei tetPti facei stepFraction "
374  "behind nBehind origProc origId"
375  );
376 
377  //- Cumulative particle counter - used to provide unique ID
378  static label particleCount_;
379 
380 
381  // Constructors
382 
383  //- Construct from components
384  particle
385  (
386  const polyMesh& mesh,
387  const barycentric& coordinates,
388  const label celli,
389  const label tetFacei,
390  const label tetPti
391  );
392 
393  //- Construct from a position and a cell, searching for the rest of the
394  // required topology
395  particle
396  (
397  const polyMesh& mesh,
398  const vector& position,
399  const label celli
400  );
401 
402  //- Construct from Istream
403  particle(const polyMesh& mesh, Istream&, bool readFields = true);
404 
405  //- Construct as a copy
406  particle(const particle& p);
407 
408  //- Construct as a copy with references to a new mesh
409  particle(const particle& p, const polyMesh& mesh);
410 
411  //- Construct a clone
412  virtual autoPtr<particle> clone() const
413  {
414  return autoPtr<particle>(new particle(*this));
415  }
416 
417  //- Factory class to read-construct particles used for
418  // parallel transfer
419  class iNew
420  {
421  const polyMesh& mesh_;
422 
423  public:
425  iNew(const polyMesh& mesh)
426  :
427  mesh_(mesh)
428  {}
430  autoPtr<particle> operator()(Istream& is) const
431  {
432  return autoPtr<particle>(new particle(mesh_, is, true));
433  }
434  };
435 
436 
437  //- Destructor
438  virtual ~particle()
439  {}
440 
441 
442  // Member Functions
443 
444  // Access
445 
446  //- Get unique particle creation id
447  inline label getNewParticleID() const;
448 
449  //- Return the mesh database
450  inline const polyMesh& mesh() const;
451 
452  //- Return current particle coordinates
453  inline const barycentric& coordinates() const;
454 
455  //- Return current cell particle is in
456  inline label cell() const;
457 
458  //- Return current tet face particle is in
459  inline label tetFace() const;
460 
461  //- Return current tet face particle is in
462  inline label tetPt() const;
463 
464  //- Return current face particle is on otherwise -1
465  inline label face() const;
466 
467  //- Return the fraction of time-step completed
468  inline scalar stepFraction() const;
469 
470  //- Return the fraction of time-step completed
471  inline scalar& stepFraction();
472 
473  //- Return the originating processor ID
474  inline label origProc() const;
475 
476  //- Return the originating processor ID
477  inline label& origProc();
478 
479  //- Return the particle ID on the originating processor
480  inline label origId() const;
481 
482  //- Return the particle ID on the originating processor
483  inline label& origId();
484 
485 
486  // Check
487 
488  //- Return the step fraction change within the overall time-step.
489  // Returns the start value and the change as a scalar pair. Always
490  // return Pair<scalar>(0, 1), unless sub-cycling is in effect, in
491  // which case the values will reflect the span of the sub-cycle
492  // within the time-step.
493  inline Pair<scalar> stepFractionSpan() const;
494 
495  //- Return the current fraction within the timestep. This differs
496  // from the stored step fraction due to sub-cycling.
497  inline scalar currentTimeFraction() const;
498 
499  //- Return the indices of the current tet that the
500  // particle occupies.
501  inline tetIndices currentTetIndices() const;
502 
503  //- Return the current tet transformation tensor
504  inline barycentricTensor currentTetTransform() const;
505 
506  //- Return the normal of the tri on tetFacei_ for the
507  // current tet.
508  inline vector normal() const;
509 
510  //- Is the particle on a face?
511  inline bool onFace() const;
512 
513  //- Is the particle on an internal face?
514  inline bool onInternalFace() const;
515 
516  //- Is the particle on a boundary face?
517  inline bool onBoundaryFace() const;
518 
519  //- Return the index of patch that the particle is on
520  inline label patch() const;
521 
522  //- Return current particle position
523  inline vector position() const;
524 
525 
526  // Track
527 
528  //- Set the step fraction and clear the behind data in preparation for
529  // a new track
530  inline void reset(const scalar stepFraction);
531 
532  //- Track along the displacement for a given fraction of the overall
533  // step. End when the track is complete, or when a boundary is hit.
534  // On exit, stepFraction_ will have been incremented to the current
535  // position, and facei_ will be set to the index of the boundary
536  // face that was hit, or -1 if the track completed within a cell.
537  // The proportion of the displacement still to be completed is
538  // returned.
539  scalar track
540  (
541  const vector& displacement,
542  const scalar fraction
543  );
544 
545  //- As particle::track, but stops when a new cell is reached.
546  scalar trackToCell
547  (
548  const vector& displacement,
549  const scalar fraction
550  );
551 
552  //- As particle::track, but stops when a face is hit.
553  scalar trackToFace
554  (
555  const vector& displacement,
556  const scalar fraction
557  );
558 
559  //- As particle::trackToFace, but stops when a tet triangle is hit. On
560  // exit, tetTriI is set to the index of the tet triangle that was hit,
561  // or -1 if the end position was reached.
562  scalar trackToTri
563  (
564  const vector& displacement,
565  const scalar fraction,
566  label& tetTriI
567  );
568 
569  //- As particle::trackToTri, but for stationary meshes
570  scalar trackToStationaryTri
571  (
572  const vector& displacement,
573  const scalar fraction,
574  label& tetTriI
575  );
576 
577  //- As particle::trackToTri, but for moving meshes
578  scalar trackToMovingTri
579  (
580  const vector& displacement,
581  const scalar fraction,
582  label& tetTriI
583  );
584 
585  //- Hit the current face. If the current face is internal than this
586  // crosses into the next cell. If it is a boundary face then this will
587  // interact the particle with the relevant patch.
588  template<class TrackCloudType>
589  void hitFace
590  (
591  const vector& displacement,
592  const scalar fraction,
593  TrackCloudType& cloud,
594  trackingData& td
595  );
596 
597  //- Convenience function. Combines trackToFace and hitFace
598  template<class TrackCloudType>
599  scalar trackToAndHitFace
600  (
601  const vector& displacement,
602  const scalar fraction,
603  TrackCloudType& cloud,
604  trackingData& td
605  );
606 
607  //- Get the displacement from the mesh centre. Used to correct the
608  // particle position in cases with reduced dimensionality. Returns a
609  // zero vector for three-dimensional cases.
611 
612 
613  // Patch data
614 
615  //- Get the normal and displacement of the current patch location
616  inline void patchData(vector& normal, vector& displacement) const;
617 
618 
619  // Transformations
620 
621  //- Transform the physical properties of the particle
622  // according to the given transformation tensor
623  virtual void transformProperties(const transformer&);
624 
625 
626  // Transfers
627 
628  //- Make changes prior to a parallel transfer. Runs either processor or
629  // nonConformalCyclic variant below.
631 
632  //- Make changes following a parallel transfer. Runs either processor
633  // or nonConformalCyclic variant below.
635 
636  //- Make changes prior to a transfer across a processor boundary.
637  // Stores the local patch face index (in facei_) so that the mesh face
638  // index can be determined on the other side.
640 
641  //- Make changes following a transfer across a processor boundary.
642  // Converts the stored patch index to a mesh index. Accounts for
643  // the receiving face being reversed relative to the sending face.
645 
646  //- Make changes prior to a transfer across a non conformal cyclic
647  // boundary. Stores the receiving patch face (in facei_). Breaks the
648  // topology and stores the cartesian position.
650  (
651  const label sendToPatch,
652  const label sendToPatchFace
653  );
654 
655  //- Make changes following a transfer across a non conformal cyclic
656  // boundary. Locates the particle using the stored face index and
657  // cartesian position.
658  void correctAfterNonConformalCyclicTransfer(const label sendToPatch);
659 
660 
661  // Interaction list referral
662 
663  //- Break the topology and store the cartesian position so that the
664  // particle can be referred.
666  (
667  const transformer& transform
668  );
669 
670  //- Correct the topology after referral. Locates the particle relative
671  // to a nearby cell/tet. The particle may end up outside this cell/tet
672  // and cannot therefore be tracked.
673  void correctAfterInteractionListReferral(const label celli);
674 
675 
676  // Decompose and reconstruct
677 
678  //- Return the tet point appropriate for decomposition or reconstruction
679  // to or from the given mesh.
681  (
682  const polyMesh& procMesh,
683  const label procCell,
684  const label procTetFace
685  ) const;
686 
687 
688  // Mapping
689 
690  //- Map after a topology change
691  void autoMap(const vector& position, const polyTopoChangeMap& mapper);
692 
693 
694  // I-O
695 
696  //- Read the fields associated with the owner cloud
697  template<class TrackCloudType>
698  static void readFields(TrackCloudType& c);
699 
700  //- Write the fields associated with the owner cloud
701  template<class TrackCloudType>
702  static void writeFields(const TrackCloudType& c);
703 
704  //- Write the particle position and cell
705  void writePosition(Ostream&) const;
706 
707 
708  // Friend Operators
709 
710  friend Ostream& operator<<(Ostream&, const particle&);
711 
712  friend bool operator==(const particle& pA, const particle& pB);
713 
714  friend bool operator!=(const particle& pA, const particle& pB);
715 };
716 
717 
718 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
719 
720 } // End namespace Foam
721 
722 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
723 
724 #include "particleI.H"
725 
726 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
727 
728 #ifdef NoRepository
729  #include "particleTemplates.C"
730 #endif
731 
732 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
733 
734 #endif
735 
736 // ************************************************************************* //
Template class for intrusive linked lists.
Definition: ILList.H:50
void reset(const scalar stepFraction)
Set the step fraction and clear the behind data in preparation for.
Definition: particleI.H:284
void hitFace(const vector &displacement, const scalar fraction, TrackCloudType &cloud, trackingData &td)
Hit the current face. If the current face is internal than this.
label origProc() const
Return the originating processor ID.
Definition: particleI.H:173
const polyMesh & mesh() const
Return the mesh database.
Definition: particleI.H:125
trackingData(const TrackCloudType &cloud)
Definition: particle.H:125
Vector-tensor class used to perform translations, rotations and scaling operations in 3D space...
Definition: transformer.H:83
vector deviationFromMeshCentre() const
Get the displacement from the mesh centre. Used to correct the.
Definition: particle.C:1008
void prepareForProcessorTransfer(trackingData &td)
Make changes prior to a transfer across a processor boundary.
Definition: particle.C:1068
virtual autoPtr< particle > clone() const
Construct a clone.
Definition: particle.H:411
A 1D vector of objects of type <T> with a fixed size <Size>.
Definition: FixedList.H:54
bool onBoundaryFace() const
Is the particle on a boundary face?
Definition: particleI.H:266
friend bool operator!=(const particle &pA, const particle &pB)
Macros for adding to particle property lists.
label procTetPt(const polyMesh &procMesh, const label procCell, const label procTetFace) const
Return the tet point appropriate for decomposition or reconstruction.
Definition: particle.C:1233
void correctAfterParallelTransfer(trackingData &td)
Make changes following a parallel transfer. Runs either processor.
Definition: particle.C:1047
void correctAfterNonConformalCyclicTransfer(const label sendToPatch)
Make changes following a transfer across a non conformal cyclic.
Definition: particle.C:1139
uint8_t direction
Definition: direction.H:45
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
scalar trackToFace(const vector &displacement, const scalar fraction)
As particle::track, but stops when a face is hit.
Definition: particle.C:606
scalar trackToTri(const vector &displacement, const scalar fraction, label &tetTriI)
As particle::trackToFace, but stops when a tet triangle is hit. On.
Definition: particle.C:991
void autoMap(const vector &position, const polyTopoChangeMap &mapper)
Map after a topology change.
Definition: particle.C:1259
scalar trackToAndHitFace(const vector &displacement, const scalar fraction, TrackCloudType &cloud, trackingData &td)
Convenience function. Combines trackToFace and hitFace.
void hitWallPatch(TrackCloudType &, trackingData &)
Overridable function to handle the particle hitting a wallPatch.
void hitSymmetryPlanePatch(TrackCloudType &, trackingData &)
Overridable function to handle the particle hitting a.
Pair< scalar > stepFractionSpan() const
Return the step fraction change within the overall time-step.
Definition: particleI.H:197
scalar trackToCell(const vector &displacement, const scalar fraction)
As particle::track, but stops when a new cell is reached.
Definition: particle.C:584
label sendFromPatch
Patch from which to send the particle.
Definition: particle.H:114
virtual void transformProperties(const transformer &)
Transform the physical properties of the particle.
Definition: particle.C:1023
void hitSymmetryPatch(TrackCloudType &, trackingData &)
Overridable function to handle the particle hitting a symmetryPatch.
const dimensionedScalar c
Speed of light in a vacuum.
Base particle class.
Definition: particle.H:81
label tetPt() const
Return current tet face particle is in.
Definition: particleI.H:149
friend Ostream & operator<<(Ostream &, const particle &)
vector normal() const
Return the normal of the tri on tetFacei_ for the.
Definition: particleI.H:248
void patchData(vector &normal, vector &displacement) const
Get the normal and displacement of the current patch location.
Definition: particleI.H:292
label cell() const
Return current cell particle is in.
Definition: particleI.H:137
void prepareForInteractionListReferral(const transformer &transform)
Break the topology and store the cartesian position so that the.
Definition: particle.C:1173
An ordered pair of two objects of type <T> with first() and second() elements.
Definition: contiguous.H:49
TypeName("particle")
Runtime type information.
label face() const
Return current face particle is on otherwise -1.
Definition: particleI.H:155
Intrusive doubly-linked list.
void hitProcessorPatch(TrackCloudType &, trackingData &)
Overridable function to handle the particle hitting a processorPatch.
tmp< fvMatrix< Type > > operator==(const fvMatrix< Type > &, const fvMatrix< Type > &)
void prepareForParallelTransfer(trackingData &td)
Make changes prior to a parallel transfer. Runs either processor or.
Definition: particle.C:1028
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
particle(const polyMesh &mesh, const barycentric &coordinates, const label celli, const label tetFacei, const label tetPti)
Construct from components.
Definition: particle.C:474
A cloud is a collection of lagrangian particles.
Definition: cloud.H:51
void correctAfterInteractionListReferral(const label celli)
Correct the topology after referral. Locates the particle relative.
Definition: particle.C:1196
void hitCyclicPatch(TrackCloudType &, trackingData &)
Overridable function to handle the particle hitting a cyclicPatch.
Storage and named access for the indices of a tet which is part of the decomposition of a cell...
Definition: tetIndices.H:81
barycentricTensor currentTetTransform() const
Return the current tet transformation tensor.
Definition: particleI.H:235
label origId() const
Return the particle ID on the originating processor.
Definition: particleI.H:185
bool onFace() const
Is the particle on a face?
Definition: particleI.H:254
label tetFace() const
Return current tet face particle is in.
Definition: particleI.H:143
label getNewParticleID() const
Get unique particle creation id.
Definition: particleI.H:111
Factory class to read-construct particles used for.
Definition: particle.H:418
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
void hitWedgePatch(TrackCloudType &, trackingData &)
Overridable function to handle the particle hitting a wedgePatch.
scalar stepFraction() const
Return the fraction of time-step completed.
Definition: particleI.H:161
tetIndices currentTetIndices() const
Return the indices of the current tet that the.
Definition: particleI.H:229
scalar trackToMovingTri(const vector &displacement, const scalar fraction, label &tetTriI)
As particle::trackToTri, but for moving meshes.
Definition: particle.C:794
label patchi
label sendToPatch
Patch to which to send the particle.
Definition: particle.H:117
void correctAfterProcessorTransfer(trackingData &td)
Make changes following a transfer across a processor boundary.
Definition: particle.C:1075
bool keepParticle
Flag to indicate whether to keep particle (false = delete)
Definition: particle.H:107
static void readFields(TrackCloudType &c)
Read the fields associated with the owner cloud.
void writePosition(Ostream &) const
Write the particle position and cell.
Definition: particleIO.C:87
static void writeFields(const TrackCloudType &c)
Write the fields associated with the owner cloud.
Templated 4x3 tensor derived from VectorSpace. Has 12 components. Can represent a barycentric transfo...
friend bool operator==(const particle &pA, const particle &pB)
bool hitPatch(TrackCloudType &, trackingData &)
Overridable function to handle the particle hitting a patch.
const barycentric & coordinates() const
Return current particle coordinates.
Definition: particleI.H:131
static label particleCount_
Cumulative particle counter - used to provide unique ID.
Definition: particle.H:374
void prepareForNonConformalCyclicTransfer(const label sendToPatch, const label sendToPatchFace)
Make changes prior to a transfer across a non conformal cyclic.
Definition: particle.C:1107
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
bool onInternalFace() const
Is the particle on an internal face?
Definition: particleI.H:260
void hitCyclicAMIPatch(const vector &displacement, const scalar fraction, TrackCloudType &cloud, trackingData &td)
Overridable function to handle the particle hitting a cyclicAMIPatch.
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:76
scalar trackToStationaryTri(const vector &displacement, const scalar fraction, label &tetTriI)
As particle::trackToTri, but for stationary meshes.
Definition: particle.C:659
volScalarField & p
#define DefinePropertyList(str)
virtual ~particle()
Destructor.
Definition: particle.H:437
scalar track(const vector &displacement, const scalar fraction)
Track along the displacement for a given fraction of the overall.
Definition: particle.C:560
bool operator!=(const particle &, const particle &)
Definition: particle.C:1282
bool hitNonConformalCyclicPatch(const vector &displacement, const scalar fraction, const label patchi, TrackCloudType &cloud, trackingData &td)
Overridable function to handle the particle hitting an.
label sendToProc
Processor to send the particle to. -1 indicates that this.
Definition: particle.H:111
vector position() const
Return current particle position.
Definition: particleI.H:278
label patch() const
Return the index of patch that the particle is on.
Definition: particleI.H:272
Namespace for OpenFOAM.
dimensionSet transform(const dimensionSet &)
Definition: dimensionSet.C:483
label sendToPatchFace
Patch face to which to send the particle.
Definition: particle.H:120
scalar currentTimeFraction() const
Return the current fraction within the timestep. This differs.
Definition: particleI.H:221