CollidingParcel.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-2020 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::CollidingParcel
26 
27 Description
28  Wrapper around kinematic parcel types to add collision modelling
29 
30 SourceFiles
31  CollidingParcelI.H
32  CollidingParcel.C
33  CollidingParcelIO.C
34 
35 \*---------------------------------------------------------------------------*/
36 
37 #ifndef CollidingParcel_H
38 #define CollidingParcel_H
39 
40 #include "particle.H"
41 #include "CollisionRecordList.H"
42 #include "labelFieldIOField.H"
43 #include "vectorFieldIOField.H"
44 
45 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
46 
47 namespace Foam
48 {
49 
53 
54 template<class ParcelType>
55 class CollidingParcel;
56 
57 // Forward declaration of friend functions
58 
59 template<class ParcelType>
60 Ostream& operator<<
61 (
62  Ostream&,
64 );
65 
66 /*---------------------------------------------------------------------------*\
67  Class CollidingParcel Declaration
68 \*---------------------------------------------------------------------------*/
69 
70 template<class ParcelType>
71 class CollidingParcel
72 :
73  public ParcelType
74 {
75  // Private member data
76 
77  //- Size in bytes of the fields
78  static const std::size_t sizeofFields_;
79 
80 
81 public:
82 
83  //- Class to hold thermo particle constant properties
84  class constantProperties
85  :
86  public ParcelType::constantProperties
87  {
88 
89  // Private Data
90 
91  //- Young's modulus [N/m^2]
92  demandDrivenEntry<scalar> youngsModulus_;
93 
94  //- Poisson's ratio
95  demandDrivenEntry<scalar> poissonsRatio_;
96 
97 
98  public:
99 
100  // Constructors
101 
102  //- Null constructor
104 
105  //- Copy constructor
107 
108  //- Construct from dictionary
109  constantProperties(const dictionary& parentDict);
110 
111 
112  // Member Functions
113 
114  //- Return const access to Young's Modulus
115  inline scalar youngsModulus() const;
116 
117  //- Return const access to Poisson's ratio
118  inline scalar poissonsRatio() const;
119  };
120 
121 
122  //- Class to hold temporary data during tracking
123  class trackingData
124  :
125  public ParcelType::trackingData
126  {
127  public:
129  enum trackPart
130  {
131  tpVelocityHalfStep,
132  tpLinearTrack,
133  tpRotationalTrack
134  };
135 
136 
137  private:
138 
139  // Private Data
140 
141  //- Which part of the integration algorithm is taking place
142  trackPart part_;
143 
144 
145  public:
146 
147  // Constructors
148 
149  //- Construct from components
150  template <class TrackCloudType>
151  inline trackingData(const TrackCloudType& cloud);
152 
153 
154  // Member Functions
155 
156  //- Return the part of the tracking operation taking place
157  inline trackPart part() const;
158 
159  //- Access to the part of the tracking operation taking place
160  inline trackPart& part();
161  };
162 
163 
164 protected:
165 
166  // Protected data
167 
168  //- Force on particle due to collisions [N]
169  vector f_;
170 
171  //- Angular momentum of Parcel in global reference frame [kg m2/s]
173 
174  //- Torque on particle due to collisions in global
175  // reference frame [Nm]
176  vector torque_;
177 
178  //- Particle collision records
179  collisionRecordList collisionRecords_;
180 
181 
182 public:
183 
184  // Static Data Members
185 
186  //- Runtime type information
187  TypeName("CollidingParcel");
188 
189  //- String representation of properties
191  (
192  ParcelType,
193  " (fx fy fz)"
194  + " (angularMomentumx angularMomentumy angularMomentumz)"
195  + " (torquex torquey torquez)"
196  + " collisionRecordsPairAccessed"
197  + " collisionRecordsPairOrigProcOfOther"
198  + " collisionRecordsPairOrigIdOfOther"
199  + " (collisionRecordsPairData)"
200  + " collisionRecordsWallAccessed"
201  + " collisionRecordsWallPRel"
202  + " (collisionRecordsWallData)"
203  );
204 
205 
206  // Constructors
207 
208  //- Construct from mesh, coordinates and topology
209  // Other properties initialised as null
210  inline CollidingParcel
211  (
212  const polyMesh& mesh,
213  const barycentric& coordinates,
214  const label celli,
215  const label tetFacei,
216  const label tetPti
217  );
218 
219  //- Construct from a position and a cell, searching for the rest of the
220  // required topology. Other properties are initialised as null.
221  inline CollidingParcel
222  (
223  const polyMesh& mesh,
224  const vector& position,
225  const label celli
226  );
227 
228  //- Construct from components
229  inline CollidingParcel
230  (
231  const polyMesh& mesh,
232  const barycentric& coordinates,
233  const label celli,
234  const label tetFacei,
235  const label tetPti,
236  const label typeId,
237  const scalar nParticle0,
238  const scalar d0,
239  const scalar dTarget0,
240  const vector& U0,
241  const vector& f0,
242  const vector& angularMomentum0,
243  const vector& torque0,
244  const typename ParcelType::constantProperties& constProps
245  );
246 
247  //- Construct from Istream
249  (
250  const polyMesh& mesh,
251  Istream& is,
252  bool readFields = true
253  );
254 
255  //- Construct as a copy
257 
258  //- Construct as a copy
260 
261  //- Construct and return a (basic particle) clone
262  virtual autoPtr<particle> clone() const
263  {
264  return autoPtr<particle>(new CollidingParcel(*this));
265  }
266 
267  //- Construct and return a (basic particle) clone
268  virtual autoPtr<particle> clone(const polyMesh& mesh) const
269  {
270  return autoPtr<particle>(new CollidingParcel(*this, mesh));
271  }
272 
273  //- Factory class to read-construct particles used for
274  // parallel transfer
275  class iNew
276  {
277  const polyMesh& mesh_;
278 
279  public:
281  iNew(const polyMesh& mesh)
282  :
283  mesh_(mesh)
284  {}
286  autoPtr<CollidingParcel<ParcelType>> operator()(Istream& is) const
287  {
289  (
290  new CollidingParcel<ParcelType>(mesh_, is, true)
291  );
292  }
293  };
294 
295 
296  // Member Functions
297 
298  // Access
299 
300  //- Return const access to force
301  inline const vector& f() const;
302 
303  //- Return const access to angular momentum
304  inline const vector& angularMomentum() const;
305 
306  //- Return const access to torque
307  inline const vector& torque() const;
308 
309  //- Return const access to the collision records
310  inline const collisionRecordList& collisionRecords() const;
311 
312  //- Return access to force
313  inline vector& f();
314 
315  //- Return access to angular momentum
316  inline vector& angularMomentum();
317 
318  //- Return access to torque
319  inline vector& torque();
320 
321  //- Return access to collision records
322  inline collisionRecordList& collisionRecords();
323 
324  //- Particle angular velocity
325  inline vector omega() const;
326 
327 
328  // Tracking
329 
330  //- Move the parcel
331  template<class TrackCloudType>
332  bool move
333  (
334  TrackCloudType& cloud,
335  trackingData& td,
336  const scalar trackTime
337  );
338 
339  //- Transform the physical properties of the particle
340  // according to the given transformation tensor
341  virtual void transformProperties(const transformer&);
342 
343 
344  // I-O
345 
346  //- Read
347  template<class CloudType>
348  static void readFields(CloudType& c);
349 
350  //- Write
351  template<class CloudType>
352  static void writeFields(const CloudType& c);
353 
354 
355  // Ostream Operator
356 
357  friend Ostream& operator<< <ParcelType>
358  (
359  Ostream&,
361  );
362 };
363 
364 
365 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
366 
367 } // End namespace Foam
368 
369 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
370 
371 #include "CollidingParcelI.H"
373 
374 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
375 
376 #ifdef NoRepository
377  #include "CollidingParcel.C"
378 #endif
379 
380 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
381 
382 #endif
383 
384 // ************************************************************************* //
vector f_
Force on particle due to collisions [N].
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
Vector-tensor class used to perform translations and rotations in 3D space.
Definition: transformer.H:83
CollisionRecordList< vector, vector > collisionRecordList
vectorFieldCompactIOField pairDataFieldCompactIOField
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
Wrapper around kinematic parcel types to add collision modelling.
Class to hold thermo particle constant properties.
Class to hold temporary data during tracking.
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
const vector & f() const
Return const access to force.
bool cp(const fileName &src, const fileName &dst, const bool followLink=true)
Copy, recursively if necessary, the source to the destination.
Definition: POSIX.C:753
static void writeFields(const CloudType &c)
Write.
vectorFieldCompactIOField wallDataFieldCompactIOField
scalar poissonsRatio() const
Return const access to Poisson&#39;s ratio.
virtual void transformProperties(const transformer &)
Transform the physical properties of the particle.
scalar youngsModulus() const
Return const access to Young&#39;s Modulus.
collisionRecordList collisionRecords_
Particle collision records.
const dimensionedScalar & c
Speed of light in a vacuum.
dynamicFvMesh & mesh
AddToPropertyList(ParcelType, " (fx fy fz)"+" (angularMomentumx angularMomentumy angularMomentumz)"+" (torquex torquey torquez)"+" collisionRecordsPairAccessed"+" collisionRecordsPairOrigProcOfOther"+" collisionRecordsPairOrigIdOfOther"+" (collisionRecordsPairData)"+" collisionRecordsWallAccessed"+" collisionRecordsWallPRel"+" (collisionRecordsWallData)")
String representation of properties.
A cloud is a collection of lagrangian particles.
Definition: cloud.H:51
bool move(TrackCloudType &cloud, trackingData &td, const scalar trackTime)
Move the parcel.
const vector & torque() const
Return const access to torque.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
const vector & angularMomentum() const
Return const access to angular momentum.
vector angularMomentum_
Angular momentum of Parcel in global reference frame [kg m2/s].
vector torque_
Torque on particle due to collisions in global.
vector omega() const
Particle angular velocity.
PtrList< coordinateSystem > coordinates(solidRegions.size())
A Field of objects of type <T> with automated input and output using a compact storage. Behaves like IOField except when binary output in case it writes a CompactListList.
CollidingParcel(const polyMesh &mesh, const barycentric &coordinates, const label celli, const label tetFacei, const label tetPti)
Construct from mesh, coordinates and topology.
TypeName("CollidingParcel")
Runtime type information.
Factory class to read-construct particles used for.
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:74
volScalarField & p
const collisionRecordList & collisionRecords() const
Return const access to the collision records.
Templated base class for dsmc cloud.
Definition: DSMCCloud.H:69
virtual autoPtr< particle > clone() const
Construct and return a (basic particle) clone.
Namespace for OpenFOAM.
static void readFields(CloudType &c)
Read.