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) 2011-2017 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 "Cloud.H"
39 #include "IDLList.H"
40 #include "pointField.H"
41 #include "faceList.H"
42 #include "OFstream.H"
43 #include "tetPointRef.H"
44 #include "FixedList.H"
46 #include "particleMacros.H"
47 #include "vectorTensorTransform.H"
48 
49 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
50 
51 namespace Foam
52 {
53 
54 // Forward declaration of classes
55 class particle;
56 
57 class polyPatch;
58 
59 class cyclicPolyPatch;
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 factor by which the displacement is increased when passing
95  // through negative space. This should be slightly bigger than one.
96  // This is needed as a straight trajectory can form a closed loop
97  // through regions of overlapping positive and negative space, leading
98  // to a track which never ends. By increasing the rate of displacement
99  // through negative regions, the change in track fraction over this
100  // part of the loop no longer exactly cancels the change over the
101  // positive part, and the track therefore terminates.
102  static const scalar negativeSpaceDisplacementFactor;
103 
104 
105 public:
106 
107  template<class CloudType>
108  class TrackingData
109  {
110  // Private data
111 
112  //- Reference to the cloud containing (this) particle
113  CloudType& cloud_;
114 
115 
116  public:
117 
118  // Public data
120  typedef CloudType cloudType;
121 
122  //- Flag to switch processor
123  bool switchProcessor;
124 
125  //- Flag to indicate whether to keep particle (false = delete)
126  bool keepParticle;
127 
128 
129  // Constructor
131  :
132  cloud_(cloud)
133  {}
134 
135 
136  // Member functions
137 
138  //- Return a reference to the cloud
139  CloudType& cloud()
140  {
141  return cloud_;
142  }
143  };
144 
145 
146 private:
147 
148  // Private data
149 
150  //- Reference to the polyMesh database
151  const polyMesh& mesh_;
152 
153  //- Coordinates of particle
154  barycentric coordinates_;
155 
156  //- Index of the cell it is in
157  label celli_;
158 
159  //- Index of the face that owns the decomposed tet that the
160  // particle is in
161  label tetFacei_;
162 
163  //- Index of the point on the face that defines the decomposed
164  // tet that the particle is in. Relative to the face base
165  // point.
166  label tetPti_;
167 
168  //- Face index if the particle is on a face otherwise -1
169  label facei_;
170 
171  //- Fraction of time-step completed
172  scalar stepFraction_;
173 
174  //- Originating processor id
175  label origProc_;
176 
177  //- Local particle id on originating processor
178  label origId_;
179 
180 
181 private:
182 
183  // Private Member Functions
184 
185  // Tetrahedra functions
186 
187  //- Get the vertices of the current tet
188  void stationaryTetGeometry
189  (
190  vector& centre,
191  vector& base,
192  vector& vertex1,
193  vector& vertex2
194  ) const;
195 
196  //- Get the transformation associated with the current tet. This
197  // will convert a barycentric position within the tet to a
198  // cartesian position in the global coordinate system. The
199  // conversion is x = A & y, where x is the cartesian position, y is
200  // the barycentric position and A is the transformation tensor.
201  barycentricTensor stationaryTetTransform() const;
202 
203  //- Get the reverse transform associated with the current tet. The
204  // conversion is detA*y = (x - centre) & T. The variables x, y and
205  // centre have the same meaning as for the forward transform. T is
206  // the transposed inverse of the forward transform tensor, A,
207  // multiplied by its determinant, detA. This separation allows
208  // the barycentric tracking algorithm to function on inverted or
209  // degenerate tetrahedra.
210  void stationaryTetReverseTransform
211  (
212  vector& centre,
213  scalar& detA,
215  ) const;
216 
217  //- Get the vertices of the current moving tet. Two values are
218  // returned for each vertex. The first is a constant, and the
219  // second is a linear coefficient of the track fraction.
220  void movingTetGeometry
221  (
222  const scalar endStepFraction,
223  Pair<vector>& centre,
224  Pair<vector>& base,
225  Pair<vector>& vertex1,
226  Pair<vector>& vertex2
227  ) const;
228 
229  //- Get the transformation associated with the current, moving, tet.
230  // This is of the same form as for the static case. As with the
231  // moving geometry, a linear function of the tracking fraction is
232  // returned for each component.
233  Pair<barycentricTensor> movingTetTransform
234  (
235  const scalar endStepFraction
236  ) const;
237 
238  //- Get the reverse transformation associated with the current,
239  // moving, tet. This is of the same form as for the static case. As
240  // with the moving geometry, a function of the tracking fraction is
241  // returned for each component. The functions are higher order than
242  // for the forward transform; the determinant is cubic, and the
243  // tensor is quadratic.
244  void movingTetReverseTransform
245  (
246  const scalar endStepFraction,
247  Pair<vector>& centre,
248  FixedList<scalar, 4>& detA,
250  ) const;
251 
252 
253  // Transformations
254 
255  //- Reflection transform. Corrects the coordinates when the particle
256  // moves between two tets which share a base vertex, but for which
257  // the other two non cell-centre vertices are reversed. All hits
258  // which retain the same face behave this way, as do face hits.
259  void reflect();
260 
261  //- Rotation transform. Corrects the coordinates when the particle
262  // moves between two tets with different base vertices, but are
263  // otherwise similarly oriented. Hits which change the face within
264  // the cell make use of both this and the reflect transform.
265  void rotate(const bool direction);
266 
267 
268  // Topology changes
269 
270  //- Change tet within a cell. Called after a triangle is hit.
271  void changeTet(const label tetTriI);
272 
273  //- Change tet face within a cell. Called by changeTet.
274  void changeFace(const label tetTriI);
275 
276  //- Change cell. Called when the particle hits an internal face.
277  void changeCell();
278 
279 
280  // Geometry changes
281 
282  //- Locate the particle at the given position
283  void locate
284  (
285  const vector& position,
286  const vector* direction,
287  const label celli,
288  const bool boundaryFail,
289  const string boundaryMsg
290  );
291 
292 
293 protected:
294 
295  // Patch interactions
296 
297  //- Overridable function to handle the particle hitting a face
298  template<class TrackData>
299  void hitFace(TrackData& td);
300 
301  //- Overridable function to handle the particle hitting a
302  // patch. Executed before other patch-hitting functions.
303  // trackFraction is passed in to allow mesh motion to
304  // interpolate in time to the correct face state.
305  template<class TrackData>
306  bool hitPatch
307  (
308  const polyPatch&,
309  TrackData& td,
310  const label patchi,
311  const scalar trackFraction,
312  const tetIndices& tetIs
313  );
314 
315  //- Overridable function to handle the particle hitting a wedgePatch
316  template<class TrackData>
317  void hitWedgePatch(const wedgePolyPatch&, TrackData& td);
318 
319  //- Overridable function to handle the particle hitting a
320  // symmetryPlanePatch
321  template<class TrackData>
323  (
324  const symmetryPlanePolyPatch&,
325  TrackData& td
326  );
327 
328  //- Overridable function to handle the particle hitting a
329  // symmetryPatch
330  template<class TrackData>
331  void hitSymmetryPatch(const symmetryPolyPatch&, TrackData& td);
332 
333  //- Overridable function to handle the particle hitting a cyclicPatch
334  template<class TrackData>
335  void hitCyclicPatch(const cyclicPolyPatch&, TrackData& td);
336 
337  //- Overridable function to handle the particle hitting a cyclicAMIPatch
338  template<class TrackData>
339  void hitCyclicAMIPatch
340  (
341  const cyclicAMIPolyPatch&,
342  TrackData& td,
343  const vector& direction
344  );
345 
346  //- Overridable function to handle the particle hitting a
347  // processorPatch
348  template<class TrackData>
349  void hitProcessorPatch(const processorPolyPatch&, TrackData& td);
350 
351  //- Overridable function to handle the particle hitting a wallPatch
352  template<class TrackData>
353  void hitWallPatch
354  (
355  const wallPolyPatch&,
356  TrackData& td,
357  const tetIndices& tetIs
358  );
359 
360  //- Overridable function to handle the particle hitting a
361  // general patch
362  template<class TrackData>
363  void hitPatch(const polyPatch&, TrackData& td);
364 
365 
366 public:
367 
368  // Static data members
369 
370  //- Runtime type information
371  TypeName("particle");
372 
373  //- String representation of properties
375  (
376  "(Px Py Pz) celli tetFacei tetPti "
377  "facei stepFraction origProc origId"
378  );
379 
380  //- Cumulative particle counter - used to provode unique ID
381  static label particleCount_;
382 
383 
384  // Constructors
385 
386  //- Construct from components
387  particle
388  (
389  const polyMesh& mesh,
390  const barycentric& coordinates,
391  const label celli,
392  const label tetFacei,
393  const label tetPti
394  );
395 
396  //- Construct from a position and a cell, searching for the rest of the
397  // required topology
398  particle
399  (
400  const polyMesh& mesh,
401  const vector& position,
402  const label celli
403  );
404 
405  //- Construct from Istream
406  particle(const polyMesh& mesh, Istream&, bool readFields = true);
407 
408  //- Construct as a copy
409  particle(const particle& p);
410 
411  //- Construct as a copy with refernce to a new mesh
412  particle(const particle& p, const polyMesh& mesh);
413 
414  //- Construct a clone
415  virtual autoPtr<particle> clone() const
416  {
417  return autoPtr<particle>(new particle(*this));
418  }
419 
420  //- Factory class to read-construct particles used for
421  // parallel transfer
422  class iNew
423  {
424  const polyMesh& mesh_;
425 
426  public:
428  iNew(const polyMesh& mesh)
429  :
430  mesh_(mesh)
431  {}
433  autoPtr<particle> operator()(Istream& is) const
434  {
435  return autoPtr<particle>(new particle(mesh_, is, true));
436  }
437  };
438 
439 
440  //- Destructor
441  virtual ~particle()
442  {}
443 
444 
445  // Member Functions
446 
447  // Access
448 
449  //- Get unique particle creation id
450  inline label getNewParticleID() const;
451 
452  //- Return the mesh database
453  inline const polyMesh& mesh() const;
454 
455  //- Return current particle coordinates
456  inline const barycentric& coordinates() const;
457 
458  //- Return current cell particle is in
459  inline label cell() const;
460 
461  //- Return current tet face particle is in
462  inline label tetFace() const;
463 
464  //- Return current tet face particle is in
465  inline label tetPt() const;
466 
467  //- Return current face particle is on otherwise -1
468  inline label face() const;
469 
470  //- Return the fraction of time-step completed
471  inline scalar stepFraction() const;
472 
473  //- Return the fraction of time-step completed
474  inline scalar& stepFraction();
475 
476  //- Return the originating processor ID
477  inline label origProc() const;
478 
479  //- Return the originating processor ID
480  inline label& origProc();
481 
482  //- Return the particle ID on the originating processor
483  inline label origId() const;
484 
485  //- Return the particle ID on the originating processor
486  inline label& origId();
487 
488 
489  // Check
490 
491  //- Return the step fraction change within the overall time-step.
492  // Returns the start value and the change as a scalar pair. Always
493  // return Pair<scalar>(0, 1), unless sub-cycling is in effect, in
494  // which case the values will reflect the span of the sub-cycle
495  // within the time-step.
496  inline Pair<scalar> stepFractionSpan() const;
497 
498  //- Return the current fraction within the timestep. This differs
499  // from the stored step fraction due to sub-cycling.
500  inline scalar currentTimeFraction() const;
501 
502  //- Return the indices of the current tet that the
503  // particle occupies.
504  inline tetIndices currentTetIndices() const;
505 
506  //- Return the current tet transformation tensor
507  inline barycentricTensor currentTetTransform() const;
508 
509  //- Return the normal of the tri on tetFacei_ for the
510  // current tet.
511  inline vector normal() const;
512 
513  //- Return the normal of the tri on tetFacei_ for the
514  // current tet at the start of the timestep, i.e. based
515  // on oldPoints
516  inline vector oldNormal() const;
517 
518  //- Is the particle on a face?
519  inline bool onFace() const;
520 
521  //- Is the particle on an internal face?
522  inline bool onInternalFace() const;
523 
524  //- Is the particle on a boundary face?
525  inline bool onBoundaryFace() const;
526 
527  //- Return the index of patch that the particle is on
528  inline label patch() const;
529 
530  //- Return current particle position
531  inline vector position() const;
532 
533 
534  // Track
535 
536  //- Track along the displacement for a given fraction of the overall
537  // step. End when the track is complete, or when a boundary is hit.
538  // On exit, stepFraction_ will have been incremented to the current
539  // position, and facei_ will be set to the index of the boundary
540  // face that was hit, or -1 if the track completed within a cell.
541  // The proportion of the displacement still to be completed is
542  // returned.
543  scalar track
544  (
545  const vector& displacement,
546  const scalar fraction
547  );
548 
549  //- As particle::track, but also stops on internal faces.
550  scalar trackToFace
551  (
552  const vector& displacement,
553  const scalar fraction
554  );
555 
556  //- As particle::trackToFace, but also stops on tet triangles. On
557  // exit, tetTriI is set to the index of the tet triangle that was
558  // hit, or -1 if the end position was reached.
559  scalar trackToTri
560  (
561  const vector& displacement,
562  const scalar fraction,
563  label& tetTriI
564  );
565 
566  //- As particle::trackToTri, but for stationary meshes
567  scalar trackToStationaryTri
568  (
569  const vector& displacement,
570  const scalar fraction,
571  label& tetTriI
572  );
573 
574  //- As particle::trackToTri, but for moving meshes
575  scalar trackToMovingTri
576  (
577  const vector& displacement,
578  const scalar fraction,
579  label& tetTriI
580  );
581 
582  //- As non-templated particle::trackToFace, but with additional
583  // boundary handling.
584  template<class TrackData>
585  void trackToFace
586  (
587  const vector& displacement,
588  const scalar fraction,
589  TrackData& td
590  );
591 
592  //- Set the constrained components of the particle position to the
593  // mesh centre.
594  void constrainToMeshCentre();
595 
596 
597  // Patch data
598 
599  //- Get the normal and velocity of the current patch location
600  void patchData(vector& n, vector& U) const;
601 
602 
603  // Transformations
604 
605  //- Transform the physical properties of the particle
606  // according to the given transformation tensor
607  virtual void transformProperties(const tensor& T);
608 
609  //- Transform the physical properties of the particle
610  // according to the given separation vector
611  virtual void transformProperties(const vector& separation);
612 
613  //- The nearest distance to a wall that
614  // the particle can be in the n direction
615  virtual scalar wallImpactDistance(const vector& n) const;
616 
617 
618  // Parallel transfer
619 
620  //- Convert global addressing to the processor patch
621  // local equivalents
622  template<class TrackData>
623  void prepareForParallelTransfer(const label patchi, TrackData& td);
624 
625  //- Convert processor patch addressing to the global equivalents
626  // and set the celli to the face-neighbour
627  template<class TrackData>
628  void correctAfterParallelTransfer(const label patchi, TrackData& td);
629 
630 
631  // Interaction list referral
632 
633  //- Break the topology and store the particle position so that the
634  // particle can be referred.
636  (
638  );
639 
640  //- Correct the topology after referral. The particle may still be
641  // outside the stored tet and therefore not track-able.
642  void correctAfterInteractionListReferral(const label celli);
643 
644 
645  // Decompose and reconstruct
646 
647  //- Return the tet point approproate for decomposition or reconstruction
648  // to or from the given mesh.
650  (
651  const polyMesh& procMesh,
652  const label procCell,
653  const label procTetFace
654  ) const;
655 
656 
657  // Mapping
658 
659  //- Map after a topology change
660  void autoMap(const vector& position, const mapPolyMesh& mapper);
661 
662 
663  // I-O
664 
665  //- Read the fields associated with the owner cloud
666  template<class CloudType>
667  static void readFields(CloudType& c);
668 
669  //- Write the fields associated with the owner cloud
670  template<class CloudType>
671  static void writeFields(const CloudType& c);
672 
673  //- Write the particle position and cell
674  void writePosition(Ostream&) const;
675 
676 
677  // Friend Operators
678 
679  friend Ostream& operator<<(Ostream&, const particle&);
680 
681  friend bool operator==(const particle& pA, const particle& pB);
682 
683  friend bool operator!=(const particle& pA, const particle& pB);
684 };
685 
686 
687 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
688 
689 } // End namespace Foam
690 
691 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
692 
693 #include "particleI.H"
694 
695 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
696 
697 #ifdef NoRepository
698  #include "particleTemplates.C"
699 #endif
700 
701 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
702 
703 #endif
704 
705 // ************************************************************************* //
Symmetry patch for non-planar or multi-plane patches.
Template class for intrusive linked lists.
Definition: ILList.H:50
void hitSymmetryPlanePatch(const symmetryPlanePolyPatch &, TrackData &td)
Overridable function to handle the particle hitting a.
label origProc() const
Return the originating processor ID.
Definition: particleI.H:93
const polyMesh & mesh() const
Return the mesh database.
Definition: particleI.H:45
U
Definition: pEqn.H:83
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
void prepareForParallelTransfer(const label patchi, TrackData &td)
Convert global addressing to the processor patch.
virtual autoPtr< particle > clone() const
Construct a clone.
Definition: particle.H:414
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:192
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 approproate for decomposition or reconstruction.
Definition: particle.C:1056
Vector-tensor class used to perform translations and rotations in 3D space.
void hitWallPatch(const wallPolyPatch &, TrackData &td, const tetIndices &tetIs)
Overridable function to handle the particle hitting a wallPatch.
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 also stops on internal faces.
Definition: particle.C:622
scalar trackToTri(const vector &displacement, const scalar fraction, label &tetTriI)
As particle::trackToFace, but also stops on tet triangles. On.
Definition: particle.C:897
Pair< scalar > stepFractionSpan() const
Return the step fraction change within the overall time-step.
Definition: particleI.H:117
static void readFields(CloudType &c)
Read the fields associated with the owner cloud.
virtual void transformProperties(const tensor &T)
Transform the physical properties of the particle.
Definition: particle.C:978
bool hitPatch(const polyPatch &, TrackData &td, const label patchi, const scalar trackFraction, const tetIndices &tetIs)
Overridable function to handle the particle hitting a.
void hitCyclicAMIPatch(const cyclicAMIPolyPatch &, TrackData &td, const vector &direction)
Overridable function to handle the particle hitting a cyclicAMIPatch.
vector oldNormal() const
Return the normal of the tri on tetFacei_ for the.
Definition: particleI.H:174
void prepareForInteractionListReferral(const vectorTensorTransform &transform)
Break the topology and store the particle position so that the.
Definition: particle.C:993
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Definition: mapPolyMesh.H:158
Base particle class.
Definition: particle.H:81
label tetPt() const
Return current tet face particle is in.
Definition: particleI.H:69
Neighbour processor patch.
friend Ostream & operator<<(Ostream &, const particle &)
vector normal() const
Return the normal of the tri on tetFacei_ for the.
Definition: particleI.H:168
CloudType & cloud()
Return a reference to the cloud.
Definition: particle.H:138
void hitProcessorPatch(const processorPolyPatch &, TrackData &td)
Overridable function to handle the particle hitting a.
label cell() const
Return current cell particle is in.
Definition: particleI.H:57
An ordered pair of two objects of type <T> with first() and second() elements.
Definition: contiguous.H:49
void patchData(vector &n, vector &U) const
Get the normal and velocity of the current patch location.
Definition: particle.C:937
TypeName("particle")
Runtime type information.
label face() const
Return current face particle is on otherwise -1.
Definition: particleI.H:75
static void writeFields(const CloudType &c)
Write the fields associated with the owner cloud.
Intrusive doubly-linked list.
tmp< fvMatrix< Type > > operator==(const fvMatrix< Type > &, const fvMatrix< Type > &)
Cyclic plane patch.
particle(const polyMesh &mesh, const barycentric &coordinates, const label celli, const label tetFacei, const label tetPti)
Construct from components.
Definition: particle.C:524
A cloud is a collection of lagrangian particles.
Definition: cloud.H:51
Foam::wallPolyPatch.
Definition: wallPolyPatch.H:48
Wedge front and back plane patch.
void correctAfterInteractionListReferral(const label celli)
Correct the topology after referral. The particle may still be.
Definition: particle.C:1018
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:155
Cyclic patch for Arbitrary Mesh Interface (AMI)
bool keepParticle
Flag to indicate whether to keep particle (false = delete)
Definition: particle.H:125
virtual scalar wallImpactDistance(const vector &n) const
The nearest distance to a wall that.
Definition: particle.C:986
void hitWedgePatch(const wedgePolyPatch &, TrackData &td)
Overridable function to handle the particle hitting a wedgePatch.
void hitFace(TrackData &td)
Overridable function to handle the particle hitting a face.
label origId() const
Return the particle ID on the originating processor.
Definition: particleI.H:105
void autoMap(const vector &position, const mapPolyMesh &mapper)
Map after a topology change.
Definition: particle.C:1082
bool onFace() const
Is the particle on a face?
Definition: particleI.H:180
label tetFace() const
Return current tet face particle is in.
Definition: particleI.H:63
label getNewParticleID() const
Get unique particle creation id.
Definition: particleI.H:31
Factory class to read-construct particles used for.
Definition: particle.H:421
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
scalar stepFraction() const
Return the fraction of time-step completed.
Definition: particleI.H:81
tetIndices currentTetIndices() const
Return the indices of the current tet that the.
Definition: particleI.H:149
scalar trackToMovingTri(const vector &displacement, const scalar fraction, label &tetTriI)
As particle::trackToTri, but for moving meshes.
Definition: particle.C:769
label patchi
void writePosition(Ostream &) const
Write the particle position and cell.
Definition: particleIO.C:84
const dimensionedScalar c
Speed of light in a vacuum.
Templated 4x3 tensor derived from VectorSpace. Has 12 components. Can represent a barycentric transfo...
friend bool operator==(const particle &pA, const particle &pB)
label n
const barycentric & coordinates() const
Return current particle coordinates.
Definition: particleI.H:51
static label particleCount_
Cumulative particle counter - used to provode unique ID.
Definition: particle.H:377
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:186
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
void correctAfterParallelTransfer(const label patchi, TrackData &td)
Convert processor patch addressing to the global equivalents.
scalar trackToStationaryTri(const vector &displacement, const scalar fraction, label &tetTriI)
As particle::trackToTri, but for stationary meshes.
Definition: particle.C:658
volScalarField & p
#define DefinePropertyList(str)
bool switchProcessor
Flag to switch processor.
Definition: particle.H:122
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
virtual ~particle()
Destructor.
Definition: particle.H:440
TrackingData(CloudType &cloud)
Definition: particle.H:129
scalar track(const vector &displacement, const scalar fraction)
Track along the displacement for a given fraction of the overall.
Definition: particle.C:603
bool operator!=(const particle &, const particle &)
Definition: particle.C:1106
void constrainToMeshCentre()
Set the constrained components of the particle position to the.
Definition: particle.C:914
Templated base class for dsmc cloud.
Definition: DSMCCloud.H:69
vector position() const
Return current particle position.
Definition: particleI.H:204
void hitCyclicPatch(const cyclicPolyPatch &, TrackData &td)
Overridable function to handle the particle hitting a cyclicPatch.
label patch() const
Return the index of patch that the particle is on.
Definition: particleI.H:198
Namespace for OpenFOAM.
dimensionSet transform(const dimensionSet &)
Definition: dimensionSet.C:477
void hitSymmetryPatch(const symmetryPolyPatch &, TrackData &td)
Overridable function to handle the particle hitting a.
scalar currentTimeFraction() const
Return the current fraction within the timestep. This differs.
Definition: particleI.H:141