MPPICParcel.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) 2013-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::MPPICParcel
26 
27 Description
28  Wrapper around kinematic parcel types to add MPPIC modelling
29 
30 SourceFiles
31  MPPICParcelI.H
32  MPPICParcelTrackingDataI.H
33  MPPICParcel.C
34  MPPICParcelIO.C
35 
36 \*---------------------------------------------------------------------------*/
37 
38 #ifndef MPPICParcel_H
39 #define MPPICParcel_H
40 
41 #include "particle.H"
42 #include "labelFieldIOField.H"
43 #include "vectorFieldIOField.H"
44 
45 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
46 
47 namespace Foam
48 {
49 
50 // Forward declaration of clases
51 
52 template<class ParcelType>
53 class MPPICParcel;
54 
55 template<class Type>
56 class AveragingMethod;
57 
58 // Forward declaration of friend functions
59 
60 template<class ParcelType>
61 Ostream& operator<<
62 (
63  Ostream&,
65 );
66 
67 /*---------------------------------------------------------------------------*\
68  Class MPPICParcel Declaration
69 \*---------------------------------------------------------------------------*/
70 
71 template<class ParcelType>
72 class MPPICParcel
73 :
74  public ParcelType
75 {
76  // Private data
77 
78  //- Size in bytes of the fields
79  static const std::size_t sizeofFields_;
80 
81 
82 public:
83 
84  template<class CloudType>
85  class TrackingData
86  :
87  public ParcelType::template TrackingData<CloudType>
88  {
89 
90  public:
91 
92  enum trackPart
93  {
98  };
99 
100 
101  private:
102 
103  // Private data
104 
105  // MPPIC Averages
106 
107  //- Volume average
108  autoPtr<AveragingMethod<scalar>> volumeAverage_;
109 
110  //- Radius average [ volume^(1/3) ]
111  autoPtr<AveragingMethod<scalar>> radiusAverage_;
112 
113  //- Density average
114  autoPtr<AveragingMethod<scalar>> rhoAverage_;
115 
116  //- Velocity average
118 
119  //- Magnitude velocity sqyuared average
120  autoPtr<AveragingMethod<scalar>> uSqrAverage_;
121 
122  //- Frequency average
123  autoPtr<AveragingMethod<scalar>> frequencyAverage_;
124 
125  //- Mass average
126  autoPtr<AveragingMethod<scalar>> massAverage_;
127 
128 
129  //- Label specifying the current part of the tracking process
130  trackPart part_;
131 
132 
133  public:
134 
135  //- Constructors
136 
137  //- Construct from components
138  inline TrackingData
139  (
140  CloudType& cloud,
142  );
143 
144 
145  //- Update the MPPIC averages
146  inline void updateAverages(CloudType& cloud);
147 
148 
149  //- Access
150 
151  //- Const access to the tracking part label
152  inline trackPart part() const;
153 
154  //- Non const access to the tracking part label
155  inline trackPart& part();
156  };
157 
158 
159 protected:
160 
161  // Protected data
162 
163  //- Velocity correction due to collisions [m/s]
165 
166 
167 public:
168 
169  // Static data members
170 
171  //- Runtime type information
172  TypeName("MPPICParcel");
173 
174  //- String representation of properties
176  (
177  ParcelType,
178  "(UCorrectx UCorrecty UCorrectz)"
179  );
180 
181 
182  // Constructors
183 
184  //- Construct from mesh, coordinates and topology
185  // Other properties initialised as null
186  inline MPPICParcel
187  (
188  const polyMesh& mesh,
189  const barycentric& coordinates,
190  const label celli,
191  const label tetFacei,
192  const label tetPti
193  );
194 
195  //- Construct from a position and a cell, searching for the rest of the
196  // required topology. Other properties are initialised as null.
197  inline MPPICParcel
198  (
199  const polyMesh& mesh,
200  const vector& position,
201  const label celli
202  );
203 
204  //- Construct from components
205  inline MPPICParcel
206  (
207  const polyMesh& mesh,
208  const barycentric& coordinates,
209  const label celli,
210  const label tetFacei,
211  const label tetPti,
212  const label typeId,
213  const scalar nParticle0,
214  const scalar d0,
215  const scalar dTarget0,
216  const vector& U0,
217  const vector& UCorrect0,
218  const typename ParcelType::constantProperties& constProps
219  );
220 
221  //- Construct from Istream
223  (
224  const polyMesh& mesh,
225  Istream& is,
226  bool readFields = true
227  );
228 
229  //- Construct as a copy
230  MPPICParcel(const MPPICParcel& p);
231 
232  //- Construct as a copy
233  MPPICParcel(const MPPICParcel& p, const polyMesh& mesh);
234 
235  //- Construct and return a (basic particle) clone
236  virtual autoPtr<particle> clone() const
237  {
238  return autoPtr<particle>(new MPPICParcel(*this));
239  }
240 
241  //- Construct and return a (basic particle) clone
242  virtual autoPtr<particle> clone(const polyMesh& mesh) const
243  {
244  return autoPtr<particle>(new MPPICParcel(*this, mesh));
245  }
246 
247  //- Factory class to read-construct particles used for
248  // parallel transfer
249  class iNew
250  {
251  const polyMesh& mesh_;
252 
253  public:
255  iNew(const polyMesh& mesh)
256  :
257  mesh_(mesh)
258  {}
260  autoPtr<MPPICParcel<ParcelType>> operator()(Istream& is) const
261  {
263  (
264  new MPPICParcel<ParcelType>(mesh_, is, true)
265  );
266  }
267  };
268 
269 
270  // Member Functions
271 
272  // Access
273 
274  //- Return const access to correction velocity
275  inline const vector& UCorrect() const;
276 
277  //- Return access to correction velocity
278  inline vector& UCorrect();
279 
280 
281  // Tracking
282 
283  //- Move the parcel
284  template<class TrackData>
285  bool move(TrackData& td, const scalar trackTime);
286 
287 
288  // Friend Functions
289 
290  // I-O
291 
292  //- Read
293  template<class CloudType>
294  static void readFields(CloudType& c);
295 
296  //- Write
297  template<class CloudType>
298  static void writeFields(const CloudType& c);
299 
300 
301  // Ostream operator
302 
303  friend Ostream& operator<< <ParcelType>
304  (
305  Ostream&,
307  );
308 };
309 
310 
311 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
312 
313 } // End namespace Foam
314 
315 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
316 
317 #include "MPPICParcelI.H"
319 
320 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
321 
322 #ifdef NoRepository
323  #include "MPPICParcel.C"
324 #endif
325 
326 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
327 
328 #endif
329 
330 // ************************************************************************* //
virtual autoPtr< particle > clone() const
Construct and return a (basic particle) clone.
Definition: MPPICParcel.H:235
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
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
TrackingData(CloudType &cloud, trackPart part=tpLinearTrack)
Constructors.
dynamicFvMesh & mesh
bool move(TrackData &td, const scalar trackTime)
Move the parcel.
Definition: MPPICParcel.C:58
Wrapper around kinematic parcel types to add MPPIC modelling.
Definition: MPPICParcel.H:52
A cloud is a collection of lagrangian particles.
Definition: cloud.H:51
Base class for lagrangian averaging methods.
Definition: MPPICParcel.H:55
const vector & UCorrect() const
Return const access to correction velocity.
Definition: MPPICParcelI.H:94
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
static void writeFields(const CloudType &c)
Write.
void updateAverages(CloudType &cloud)
Update the MPPIC averages.
PtrList< coordinateSystem > coordinates(solidRegions.size())
vector UCorrect_
Velocity correction due to collisions [m/s].
Definition: MPPICParcel.H:163
MPPICParcel(const polyMesh &mesh, const barycentric &coordinates, const label celli, const label tetFacei, const label tetPti)
Construct from mesh, coordinates and topology.
Definition: MPPICParcelI.H:30
const dimensionedScalar c
Speed of light in a vacuum.
static void readFields(CloudType &c)
Read.
Definition: MPPICParcelIO.C:78
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
Factory class to read-construct particles used for.
Definition: MPPICParcel.H:248
TypeName("MPPICParcel")
Runtime type information.
AddToPropertyList(ParcelType, "(UCorrectx UCorrecty UCorrectz)")
String representation of properties.
Templated base class for dsmc cloud.
Definition: DSMCCloud.H:69
Namespace for OpenFOAM.