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