KinematicParcel.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-2018 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::KinematicParcel
26 
27 Description
28  Kinematic parcel class with rotational motion (as spherical
29  particles only) and one/two-way coupling with the continuous
30  phase.
31 
32  Sub-models include:
33  - drag
34  - turbulent dispersion
35  - wall interactions
36 
37 SourceFiles
38  KinematicParcelI.H
39  KinematicParcel.C
40  KinematicParcelIO.C
41 
42 \*---------------------------------------------------------------------------*/
43 
44 #ifndef KinematicParcel_H
45 #define KinematicParcel_H
46 
47 #include "particle.H"
48 #include "IOstream.H"
49 #include "autoPtr.H"
50 #include "interpolation.H"
51 #include "demandDrivenEntry.H"
52 
53 // #include "ParticleForceList.H" // TODO
54 
55 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
56 
57 namespace Foam
58 {
59 
60 template<class ParcelType>
61 class KinematicParcel;
62 
63 // Forward declaration of friend functions
64 
65 template<class ParcelType>
66 Ostream& operator<<
67 (
68  Ostream&,
70 );
71 
72 /*---------------------------------------------------------------------------*\
73  Class KinematicParcel Declaration
74 \*---------------------------------------------------------------------------*/
75 
76 template<class ParcelType>
77 class KinematicParcel
78 :
79  public ParcelType
80 {
81  // Private data
82 
83  //- Size in bytes of the fields
84  static const std::size_t sizeofFields_;
85 
86  //- Number of particle tracking attempts before we assume that it stalls
87  static label maxTrackAttempts;
88 
89 
90 public:
91 
92  //- Class to hold kinematic particle constant properties
93  class constantProperties
94  {
95  protected:
96 
97  // Protected data
98 
99  //- Constant properties dictionary
100  const dictionary dict_;
101 
102 
103  private:
104 
105  // Private data
106 
107  //- Parcel type id - used for post-processing to flag the type
108  // of parcels issued by this cloud
109  demandDrivenEntry<label> parcelTypeId_;
110 
111  //- Minimum density [kg/m3]
113 
114  //- Particle density [kg/m3] (constant)
116 
117  //- Minimum parcel mass [kg]
118  demandDrivenEntry<scalar> minParcelMass_;
119 
120 
121  public:
122 
123  // Constructors
124 
125  //- Null constructor
127 
128  //- Copy constructor
130 
131  //- Construct from dictionary
132  constantProperties(const dictionary& parentDict);
133 
134 
135  // Member functions
136 
137  //- Return const access to the constant properties dictionary
138  inline const dictionary& dict() const;
139 
140  //- Return const access to the parcel type id
141  inline label parcelTypeId() const;
142 
143  //- Return const access to the minimum density
144  inline scalar rhoMin() const;
145 
146  //- Return const access to the particle density
147  inline scalar rho0() const;
148 
149  //- Return const access to the minimum parcel mass
150  inline scalar minParcelMass() const;
151  };
152 
154  class trackingData
155  :
156  public ParcelType::trackingData
157  {
158  public:
160  enum trackPart
161  {
162  tpVelocityHalfStep,
163  tpLinearTrack,
164  tpRotationalTrack
165  };
166 
167 
168  private:
169 
170  // Private data
171 
172  // Interpolators for continuous phase fields
173 
174  //- Density interpolator
175  autoPtr<interpolation<scalar>> rhoInterp_;
176 
177  //- Velocity interpolator
179 
180  //- Dynamic viscosity interpolator
182 
183 
184  // Cached continuous phase properties
185 
186  //- Density [kg/m3]
187  scalar rhoc_;
188 
189  //- Velocity [m/s]
190  vector Uc_;
191 
192  //- Viscosity [Pa.s]
193  scalar muc_;
194 
195 
196  //- Local gravitational or other body-force acceleration
197  const vector& g_;
198 
199  // label specifying which part of the integration
200  // algorithm is taking place
201  trackPart part_;
202 
203 
204  public:
205 
206  // Constructors
207 
208  //- Construct from components
209  template <class TrackCloudType>
210  inline trackingData
211  (
212  const TrackCloudType& cloud,
213  trackPart part = tpLinearTrack
214  );
215 
216 
217  // Member functions
218 
219  //- Return conat access to the interpolator for continuous
220  // phase density field
221  inline const interpolation<scalar>& rhoInterp() const;
222 
223  //- Return conat access to the interpolator for continuous
224  // phase velocity field
225  inline const interpolation<vector>& UInterp() const;
226 
227  //- Return conat access to the interpolator for continuous
228  // phase dynamic viscosity field
229  inline const interpolation<scalar>& muInterp() const;
230 
231  //- Return the continuous phase density
232  inline scalar rhoc() const;
233 
234  //- Access the continuous phase density
235  inline scalar& rhoc();
236 
237  //- Return the continuous phase velocity
238  inline const vector& Uc() const;
239 
240  //- Access the continuous phase velocity
241  inline vector& Uc();
242 
243  //- Return the continuous phase viscosity
244  inline scalar muc() const;
245 
246  //- Access the continuous phase viscosity
247  inline scalar& muc();
248 
249  // Return const access to the gravitational acceleration vector
250  inline const vector& g() const;
251 
252  //- Return the part of the tracking operation taking place
253  inline trackPart part() const;
254 
255  //- Return access to the part of the tracking operation taking place
256  inline trackPart& part();
257  };
258 
259 
260 protected:
261 
262  // Protected data
263 
264  // Parcel properties
265 
266  //- Active flag - tracking inactive when active = false
267  bool active_;
268 
269  //- Parcel type id
270  label typeId_;
271 
272  //- Number of particles in Parcel
273  scalar nParticle_;
274 
275  //- Diameter [m]
276  scalar d_;
277 
278  //- Target diameter [m]
279  scalar dTarget_;
280 
281  //- Velocity of Parcel [m/s]
282  vector U_;
283 
284  //- Density [kg/m3]
285  scalar rho_;
286 
287  //- Age [s]
288  scalar age_;
289 
290  //- Time spent in turbulent eddy [s]
291  scalar tTurb_;
292 
293  //- Turbulent velocity fluctuation [m/s]
294  vector UTurb_;
295 
296 
297  // Protected Member Functions
298 
299  //- Calculate new particle velocity
300  template<class TrackCloudType>
301  const vector calcVelocity
302  (
303  TrackCloudType& cloud,
304  trackingData& td,
305  const scalar dt, // timestep
306  const scalar Re, // Reynolds number
307  const scalar mu, // local carrier viscosity
308  const scalar mass, // mass
309  const vector& Su, // explicit particle momentum source
310  vector& dUTrans, // momentum transfer to carrier
311  scalar& Spu // linearised drag coefficient
312  ) const;
313 
314 
315 public:
316 
317  // Static data members
318 
319  //- Runtime type information
320  TypeName("KinematicParcel");
321 
322  //- String representation of properties
324  (
325  ParcelType,
326  " active"
327  + " typeId"
328  + " nParticle"
329  + " d"
330  + " dTarget "
331  + " (Ux Uy Uz)"
332  + " rho"
333  + " age"
334  + " tTurb"
335  + " (UTurbx UTurby UTurbz)"
336  );
337 
338 
339  // Constructors
340 
341  //- Construct from mesh, coordinates and topology
342  // Other properties initialised as null
343  inline KinematicParcel
344  (
345  const polyMesh& mesh,
346  const barycentric& coordinates,
347  const label celli,
348  const label tetFacei,
349  const label tetPti
350  );
351 
352  //- Construct from a position and a cell, searching for the rest of the
353  // required topology. Other properties are initialised as null.
354  inline KinematicParcel
355  (
356  const polyMesh& mesh,
357  const vector& position,
358  const label celli
359  );
360 
361  //- Construct from components
362  inline KinematicParcel
363  (
364  const polyMesh& mesh,
365  const barycentric& coordinates,
366  const label celli,
367  const label tetFacei,
368  const label tetPti,
369  const label typeId,
370  const scalar nParticle0,
371  const scalar d0,
372  const scalar dTarget0,
373  const vector& U0,
374  const constantProperties& constProps
375  );
376 
377  //- Construct from Istream
379  (
380  const polyMesh& mesh,
381  Istream& is,
382  bool readFields = true
383  );
384 
385  //- Construct as a copy
387 
388  //- Construct as a copy
390 
391  //- Construct and return a (basic particle) clone
392  virtual autoPtr<particle> clone() const
393  {
394  return autoPtr<particle>(new KinematicParcel(*this));
395  }
396 
397  //- Construct and return a (basic particle) clone
398  virtual autoPtr<particle> clone(const polyMesh& mesh) const
399  {
400  return autoPtr<particle>(new KinematicParcel(*this, mesh));
401  }
402 
403  //- Factory class to read-construct particles used for
404  // parallel transfer
405  class iNew
406  {
407  const polyMesh& mesh_;
408 
409  public:
411  iNew(const polyMesh& mesh)
412  :
413  mesh_(mesh)
414  {}
416  autoPtr<KinematicParcel<ParcelType>> operator()(Istream& is) const
417  {
419  (
420  new KinematicParcel<ParcelType>(mesh_, is, true)
421  );
422  }
423  };
424 
425 
426  // Member Functions
427 
428  // Access
429 
430  //- Return const access to active flag
431  inline bool active() const;
432 
433  //- Return const access to type id
434  inline label typeId() const;
435 
436  //- Return const access to number of particles
437  inline scalar nParticle() const;
438 
439  //- Return const access to diameter
440  inline scalar d() const;
441 
442  //- Return const access to target diameter
443  inline scalar dTarget() const;
444 
445  //- Return const access to velocity
446  inline const vector& U() const;
447 
448  //- Return const access to density
449  inline scalar rho() const;
450 
451  //- Return const access to the age
452  inline scalar age() const;
453 
454  //- Return const access to time spent in turbulent eddy
455  inline scalar tTurb() const;
456 
457  //- Return const access to turbulent velocity fluctuation
458  inline const vector& UTurb() const;
459 
460 
461  // Edit
462 
463  //- Return const access to active flag
464  inline bool& active();
465 
466  //- Return access to type id
467  inline label& typeId();
468 
469  //- Return access to number of particles
470  inline scalar& nParticle();
471 
472  //- Return access to diameter
473  inline scalar& d();
474 
475  //- Return access to target diameter
476  inline scalar& dTarget();
477 
478  //- Return access to velocity
479  inline vector& U();
480 
481  //- Return access to density
482  inline scalar& rho();
483 
484  //- Return access to the age
485  inline scalar& age();
486 
487  //- Return access to time spent in turbulent eddy
488  inline scalar& tTurb();
489 
490  //- Return access to turbulent velocity fluctuation
491  inline vector& UTurb();
492 
493 
494  // Helper functions
495 
496  //- Cell owner mass
497  inline scalar massCell(const trackingData& td) const;
498 
499  //- Particle mass
500  inline scalar mass() const;
501 
502  //- Particle moment of inertia around diameter axis
503  inline scalar momentOfInertia() const;
504 
505  //- Particle volume
506  inline scalar volume() const;
507 
508  //- Particle volume for a given diameter
509  inline static scalar volume(const scalar d);
510 
511  //- Particle projected area
512  inline scalar areaP() const;
513 
514  //- Projected area for given diameter
515  inline static scalar areaP(const scalar d);
516 
517  //- Particle surface area
518  inline scalar areaS() const;
519 
520  //- Surface area for given diameter
521  inline static scalar areaS(const scalar d);
522 
523  //- Reynolds number
524  inline scalar Re(const trackingData& td) const;
525 
526  //- Reynolds number for given conditions
527  inline static scalar Re
528  (
529  const scalar rhoc,
530  const vector& U,
531  const vector& Uc,
532  const scalar d,
533  const scalar muc
534  );
535 
536  //- Weber number
537  inline scalar We
538  (
539  const trackingData& td,
540  const scalar sigma
541  ) const;
542 
543  //- Weber number for given conditions
544  inline static scalar We
545  (
546  const scalar rhoc,
547  const vector& U,
548  const vector& Uc,
549  const scalar d,
550  const scalar sigma
551  );
552 
553  //- Eotvos number
554  inline scalar Eo
555  (
556  const trackingData& td,
557  const scalar sigma
558  ) const;
559 
560  //- Eotvos number for given conditions
561  inline static scalar Eo
562  (
563  const vector& g,
564  const scalar rho,
565  const scalar rhoc,
566  const vector& U,
567  const scalar d,
568  const scalar sigma
569  );
570 
571 
572  // Main calculation loop
573 
574  //- Set cell values
575  template<class TrackCloudType>
576  void setCellValues(TrackCloudType& cloud, trackingData& td);
577 
578  //- Apply dispersion to the carrier phase velocity and update
579  // parcel turbulence parameters
580  template<class TrackCloudType>
581  void calcDispersion
582  (
583  TrackCloudType& cloud,
584  trackingData& td,
585  const scalar dt
586  );
587 
588  //- Correct cell values using latest transfer information
589  template<class TrackCloudType>
591  (
592  TrackCloudType& cloud,
593  trackingData& td,
594  const scalar dt
595  );
596 
597  //- Update parcel properties over the time interval
598  template<class TrackCloudType>
599  void calc
600  (
601  TrackCloudType& cloud,
602  trackingData& td,
603  const scalar dt
604  );
605 
606 
607  // Tracking
608 
609  //- Move the parcel
610  template<class TrackCloudType>
611  bool move
612  (
613  TrackCloudType& cloud,
614  trackingData& td,
615  const scalar trackTime
616  );
617 
618 
619  // Patch interactions
620 
621  //- Overridable function to handle the particle hitting a patch
622  // Executed before other patch-hitting functions
623  template<class TrackCloudType>
624  bool hitPatch(TrackCloudType& cloud, trackingData& td);
625 
626  //- Overridable function to handle the particle hitting a
627  // processorPatch
628  template<class TrackCloudType>
629  void hitProcessorPatch(TrackCloudType& cloud, trackingData& td);
630 
631  //- Overridable function to handle the particle hitting a wallPatch
632  template<class TrackCloudType>
633  void hitWallPatch(TrackCloudType& cloud, trackingData& td);
634 
635  //- Transform the physical properties of the particle
636  // according to the given transformation tensor
637  virtual void transformProperties(const tensor& T);
638 
639  //- Transform the physical properties of the particle
640  // according to the given separation vector
641  virtual void transformProperties(const vector& separation);
642 
643 
644  // I-O
645 
646  //- Read
647  template<class TrackCloudType>
648  static void readFields(TrackCloudType& c);
649 
650  //- Write
651  template<class TrackCloudType>
652  static void writeFields(const TrackCloudType& c);
653 
654 
655  // Ostream Operator
656 
657  friend Ostream& operator<< <ParcelType>
658  (
659  Ostream&,
661  );
662 };
663 
664 
665 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
666 
667 } // End namespace Foam
668 
669 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
670 
671 #include "KinematicParcelI.H"
673 
674 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
675 
676 #ifdef NoRepository
677  #include "KinematicParcel.C"
678 #endif
679 
680 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
681 
682 #endif
683 
684 // ************************************************************************* //
scalar massCell(const trackingData &td) const
Cell owner mass.
const vector & U() const
Return const access to velocity.
static void writeFields(const TrackCloudType &c)
Write.
scalar tTurb() const
Return const access to time spent in turbulent eddy.
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
scalar rhoMin() const
Return const access to the minimum density.
zeroField Su
Definition: alphaSuSp.H:1
const vector & UTurb() const
Return const access to turbulent velocity fluctuation.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
Kinematic parcel class with rotational motion (as spherical particles only) and one/two-way coupling ...
TypeName("KinematicParcel")
Runtime type information.
scalar nParticle() const
Return const access to number of particles.
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
Class to hold kinematic particle constant properties.
scalar mass() const
Particle mass.
scalar d() const
Return const access to diameter.
scalar tTurb_
Time spent in turbulent eddy [s].
scalar age_
Age [s].
scalar rho() const
Return const access to density.
bool cp(const fileName &src, const fileName &dst, const bool followLink=true)
Copy, recursively if necessary, the source to the destination.
Definition: POSIX.C:737
const dimensionedScalar sigma
Stefan-Boltzmann constant: default SI units: [W/m2/K4].
scalar rho0() const
Return const access to the particle density.
void setCellValues(TrackCloudType &cloud, trackingData &td)
Set cell values.
scalar areaP() const
Particle projected area.
Factory class to read-construct particles used for.
KinematicParcel(const polyMesh &mesh, const barycentric &coordinates, const label celli, const label tetFacei, const label tetPti)
Construct from mesh, coordinates and topology.
scalar dTarget() const
Return const access to target diameter.
label parcelTypeId() const
Return const access to the parcel type id.
scalar dTarget_
Target diameter [m].
scalar age() const
Return const access to the age.
void cellValueSourceCorrection(TrackCloudType &cloud, trackingData &td, const scalar dt)
Correct cell values using latest transfer information.
dynamicFvMesh & mesh
scalar Eo(const trackingData &td, const scalar sigma) const
Eotvos number.
scalar rho_
Density [kg/m3].
virtual void transformProperties(const tensor &T)
Transform the physical properties of the particle.
scalar momentOfInertia() const
Particle moment of inertia around diameter axis.
bool active_
Active flag - tracking inactive when active = false.
const dictionary dict_
Constant properties dictionary.
A cloud is a collection of lagrangian particles.
Definition: cloud.H:51
label typeId_
Parcel type id.
vector UTurb_
Turbulent velocity fluctuation [m/s].
static void readFields(TrackCloudType &c)
Read.
void calc(TrackCloudType &cloud, trackingData &td, const scalar dt)
Update parcel properties over the time interval.
bool active() const
Return const access to active flag.
const vector calcVelocity(TrackCloudType &cloud, trackingData &td, const scalar dt, const scalar Re, const scalar mu, const scalar mass, const vector &Su, vector &dUTrans, scalar &Spu) const
Calculate new particle velocity.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
scalar nParticle_
Number of particles in Parcel.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
scalar minParcelMass() const
Return const access to the minimum parcel mass.
const dimensionedScalar mu
Atomic mass unit.
vector U_
Velocity of Parcel [m/s].
virtual autoPtr< particle > clone() const
Construct and return a (basic particle) clone.
PtrList< coordinateSystem > coordinates(solidRegions.size())
const dimensionedScalar c
Speed of light in a vacuum.
scalar volume() const
Particle volume.
scalar Re(const trackingData &td) const
Reynolds number.
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
bool move(TrackCloudType &cloud, trackingData &td, const scalar trackTime)
Move the parcel.
volScalarField & p
const dimensionedVector & g
label typeId() const
Return const access to type id.
AddToPropertyList(ParcelType, " active"+" typeId"+" nParticle"+" d"+" dTarget "+" (Ux Uy Uz)"+" rho"+" age"+" tTurb"+" (UTurbx UTurby UTurbz)")
String representation of properties.
void hitProcessorPatch(TrackCloudType &cloud, trackingData &td)
Overridable function to handle the particle hitting a.
scalar areaS() const
Particle surface area.
bool hitPatch(TrackCloudType &cloud, trackingData &td)
Overridable function to handle the particle hitting a patch.
scalar d_
Diameter [m].
scalar We(const trackingData &td, const scalar sigma) const
Weber number.
void calcDispersion(TrackCloudType &cloud, trackingData &td, const scalar dt)
Apply dispersion to the carrier phase velocity and update.
Namespace for OpenFOAM.
const dictionary & dict() const
Return const access to the constant properties dictionary.
void hitWallPatch(TrackCloudType &cloud, trackingData &td)
Overridable function to handle the particle hitting a wallPatch.