CollidingParcel.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::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/m2]
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 protected:
123 
124  // Protected data
125 
126  //- Force on particle due to collisions [N]
127  vector f_;
128 
129  //- Angular momentum of Parcel in global reference frame [kg m2/s]
131 
132  //- Torque on particle due to collisions in global
133  // reference frame [Nm]
134  vector torque_;
135 
136  //- Particle collision records
137  collisionRecordList collisionRecords_;
138 
139 
140 public:
141 
142  // Static data members
143 
144  //- Runtime type information
145  TypeName("CollidingParcel");
146 
147  //- String representation of properties
149  (
150  ParcelType,
151  " (fx fy fz)"
152  + " (angularMomentumx angularMomentumy angularMomentumz)"
153  + " (torquex torquey torquez)"
154  + " collisionRecordsPairAccessed"
155  + " collisionRecordsPairOrigProcOfOther"
156  + " collisionRecordsPairOrigIdOfOther"
157  + " (collisionRecordsPairData)"
158  + " collisionRecordsWallAccessed"
159  + " collisionRecordsWallPRel"
160  + " (collisionRecordsWallData)"
161  );
162 
163 
164  // Constructors
165 
166  //- Construct from mesh, coordinates and topology
167  // Other properties initialised as null
168  inline CollidingParcel
169  (
170  const polyMesh& mesh,
171  const barycentric& coordinates,
172  const label celli,
173  const label tetFacei,
174  const label tetPti
175  );
176 
177  //- Construct from a position and a cell, searching for the rest of the
178  // required topology. Other properties are initialised as null.
179  inline CollidingParcel
180  (
181  const polyMesh& mesh,
182  const vector& position,
183  const label celli
184  );
185 
186  //- Construct from components
187  inline CollidingParcel
188  (
189  const polyMesh& mesh,
190  const barycentric& coordinates,
191  const label celli,
192  const label tetFacei,
193  const label tetPti,
194  const label typeId,
195  const scalar nParticle0,
196  const scalar d0,
197  const scalar dTarget0,
198  const vector& U0,
199  const vector& f0,
200  const vector& angularMomentum0,
201  const vector& torque0,
202  const typename ParcelType::constantProperties& constProps
203  );
204 
205  //- Construct from Istream
207  (
208  const polyMesh& mesh,
209  Istream& is,
210  bool readFields = true
211  );
212 
213  //- Construct as a copy
215 
216  //- Construct as a copy
218 
219  //- Construct and return a (basic particle) clone
220  virtual autoPtr<particle> clone() const
221  {
222  return autoPtr<particle>(new CollidingParcel(*this));
223  }
224 
225  //- Construct and return a (basic particle) clone
226  virtual autoPtr<particle> clone(const polyMesh& mesh) const
227  {
228  return autoPtr<particle>(new CollidingParcel(*this, mesh));
229  }
230 
231  //- Factory class to read-construct particles used for
232  // 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<CollidingParcel<ParcelType>> operator()(Istream& is) const
245  {
247  (
248  new CollidingParcel<ParcelType>(mesh_, is, true)
249  );
250  }
251  };
252 
253 
254  // Member Functions
255 
256  // Access
257 
258  //- Return const access to force
259  inline const vector& f() const;
260 
261  //- Return const access to angular momentum
262  inline const vector& angularMomentum() const;
263 
264  //- Return const access to torque
265  inline const vector& torque() const;
266 
267  //- Return const access to the collision records
268  inline const collisionRecordList& collisionRecords() const;
269 
270  //- Return access to force
271  inline vector& f();
272 
273  //- Return access to angular momentum
274  inline vector& angularMomentum();
275 
276  //- Return access to torque
277  inline vector& torque();
278 
279  //- Return access to collision records
280  inline collisionRecordList& collisionRecords();
281 
282  //- Particle angular velocity
283  inline vector omega() const;
284 
285 
286  // Tracking
287 
288  //- Move the parcel
289  template<class TrackData>
290  bool move(TrackData& td, const scalar trackTime);
291 
292  //- Transform the physical properties of the particle
293  // according to the given transformation tensor
294  virtual void transformProperties(const tensor& T);
295 
296  //- Transform the physical properties of the particle
297  // according to the given separation vector
298  virtual void transformProperties(const vector& separation);
299 
300 
301  // I-O
302 
303  //- Read
304  template<class CloudType>
305  static void readFields(CloudType& c);
306 
307  //- Write
308  template<class CloudType>
309  static void writeFields(const CloudType& c);
310 
311 
312  // Ostream Operator
313 
314  friend Ostream& operator<< <ParcelType>
315  (
316  Ostream&,
318  );
319 };
320 
321 
322 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
323 
324 } // End namespace Foam
325 
326 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
327 
328 #include "CollidingParcelI.H"
329 
330 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
331 
332 #ifdef NoRepository
333  #include "CollidingParcel.C"
334 #endif
335 
336 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
337 
338 #endif
339 
340 // ************************************************************************* //
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
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:137
Wrapper around kinematic parcel types to add collision modelling.
Class to hold thermo particle constant properties.
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:739
static void writeFields(const CloudType &c)
Write.
vectorFieldCompactIOField wallDataFieldCompactIOField
scalar poissonsRatio() const
Return const access to Poisson&#39;s ratio.
scalar youngsModulus() const
Return const access to Young&#39;s Modulus.
collisionRecordList collisionRecords_
Particle collision records.
dynamicFvMesh & mesh
AddToPropertyList(ParcelType, " (fx fy fz)"+" (angularMomentumx angularMomentumy angularMomentumz)"+" (torquex torquey torquez)"+" collisionRecordsPairAccessed"+" collisionRecordsPairOrigProcOfOther"+" collisionRecordsPairOrigIdOfOther"+" (collisionRecordsPairData)"+" collisionRecordsWallAccessed"+" collisionRecordsWallPRel"+" (collisionRecordsWallData)")
String representation of properties.
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:53
const vector & angularMomentum() const
Return const access to angular momentum.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
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.
const dimensionedScalar c
Speed of light in a vacuum.
CollidingParcel(const polyMesh &mesh, const barycentric &coordinates, const label celli, const label tetFacei, const label tetPti)
Construct from mesh, coordinates and topology.
bool move(TrackData &td, const scalar trackTime)
Move the parcel.
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 void transformProperties(const tensor &T)
Transform the physical properties of the particle.
virtual autoPtr< particle > clone() const
Construct and return a (basic particle) clone.
Namespace for OpenFOAM.
static void readFields(CloudType &c)
Read.