KinematicCloud.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-2013 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::KinematicCloud
26 
27 Description
28  Templated base class for kinematic cloud
29 
30  - cloud function objects
31 
32  - particle forces, e.g.
33  - buoyancy
34  - drag
35  - pressure gradient
36  - ...
37 
38  - sub-models:
39  - dispersion model
40  - injection model
41  - patch interaction model
42  - stochastic collision model
43  - surface film model
44 
45 SourceFiles
46  KinematicCloudI.H
47  KinematicCloud.C
48 
49 \*---------------------------------------------------------------------------*/
50 
51 #ifndef KinematicCloud_H
52 #define KinematicCloud_H
53 
54 #include "particle.H"
55 #include "Cloud.H"
56 #include "kinematicCloud.H"
57 #include "IOdictionary.H"
58 #include "autoPtr.H"
59 #include "cachedRandom.H"
60 #include "fvMesh.H"
61 #include "volFields.H"
62 #include "fvMatrices.H"
63 #include "IntegrationSchemesFwd.H"
64 #include "cloudSolution.H"
65 
66 #include "ParticleForceList.H"
68 
69 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
70 
71 namespace Foam
72 {
73 
74 // Forward declaration of classes
75 
76 template<class CloudType>
77 class InjectionModelList;
78 
79 template<class CloudType>
80 class DispersionModel;
81 
82 template<class CloudType>
84 
85 template<class CloudType>
86 class SurfaceFilmModel;
87 
88 template<class CloudType>
90 
91 
92 /*---------------------------------------------------------------------------*\
93  Class KinematicCloud Declaration
94 \*---------------------------------------------------------------------------*/
95 
96 template<class CloudType>
97 class KinematicCloud
98 :
99  public CloudType,
100  public kinematicCloud
101 {
102 public:
103 
104  // Public typedefs
105 
106  //- Type of cloud this cloud was instantiated for
107  typedef CloudType cloudType;
108 
109  //- Type of parcel the cloud was instantiated for
110  typedef typename CloudType::particleType parcelType;
111 
112  //- Convenience typedef for this cloud type
114 
115  //- Force models type
117 
118  //- Function object type
120  functionType;
121 
122 
123 private:
124 
125  // Private data
126 
127  //- Cloud copy pointer
128  autoPtr<KinematicCloud<CloudType> > cloudCopyPtr_;
129 
130 
131  // Private Member Functions
132 
133  //- Disallow default bitwise copy construct
135 
136  //- Disallow default bitwise assignment
137  void operator=(const KinematicCloud&);
138 
139 
140 protected:
141 
142  // Protected data
143 
144  //- References to the mesh and time databases
145  const fvMesh& mesh_;
146 
147  //- Dictionary of particle properties
149 
150  //- Dictionary of output properties
152 
153  //- Solution properties
155 
156  //- Parcel constant properties
157  typename parcelType::constantProperties constProps_;
158 
159  //- Sub-models dictionary
161 
162  //- Random number generator - used by some injection routines
164 
165  //- Cell occupancy information for each parcel, (demand driven)
167 
168  //- Cell length scale
170 
171 
172  // References to the carrier gas fields
173 
174  //- Density [kg/m3]
175  const volScalarField& rho_;
176 
177  //- Velocity [m/s]
178  const volVectorField& U_;
179 
180  //- Dynamic viscosity [Pa.s]
181  const volScalarField& mu_;
182 
183 
184  // Environmental properties
185 
186  //- Gravity
187  const dimensionedVector& g_;
188 
189  //- Averaged ambient domain pressure
190  scalar pAmbient_;
191 
192 
193  //- Optional particle forces
194  forceType forces_;
195 
196  //- Optional cloud function objects
197  functionType functions_;
198 
199 
200  // References to the cloud sub-models
201 
202  //- Injector models
204 
205  //- Dispersion model
208 
209  //- Patch interaction model
212 
213  //- Stochastic collision model
216 
217  //- Surface film model
220 
221 
222  // Reference to the particle integration schemes
223 
224  //- Velocity integration
226 
227 
228  // Sources
229 
230  //- Momentum
232 
233  //- Coefficient for carrier phase U equation
235 
236 
237  // Initialisation
238 
239  //- Set cloud sub-models
240  void setModels();
241 
242 
243  // Cloud evolution functions
244 
245  //- Solve the cloud - calls all evolution functions
246  template<class TrackData>
247  void solve(TrackData& td);
248 
249  //- Build the cellOccupancy
250  void buildCellOccupancy();
251 
252  //- Update (i.e. build) the cellOccupancy if it has
253  // already been used
254  void updateCellOccupancy();
255 
256  //- Evolve the cloud
257  template<class TrackData>
258  void evolveCloud(TrackData& td);
259 
260  //- Post-evolve
261  void postEvolve();
262 
263  //- Reset state of cloud
265 
266 
267 public:
268 
269  // Constructors
270 
271  //- Construct given carrier gas fields
273  (
274  const word& cloudName,
275  const volScalarField& rho,
276  const volVectorField& U,
277  const volScalarField& mu,
278  const dimensionedVector& g,
279  bool readFields = true
280  );
281 
282  //- Copy constructor with new name
284  (
286  const word& name
287  );
288 
289  //- Copy constructor with new name - creates bare cloud
291  (
292  const fvMesh& mesh,
293  const word& name,
295  );
296 
297  //- Construct and return clone based on (this) with new name
298  virtual autoPtr<Cloud<parcelType> > clone(const word& name)
299  {
301  (
302  new KinematicCloud(*this, name)
303  );
304  }
305 
306  //- Construct and return bare clone based on (this) with new name
307  virtual autoPtr<Cloud<parcelType> > cloneBare(const word& name) const
308  {
310  (
311  new KinematicCloud(this->mesh(), name, *this)
312  );
313  }
314 
315 
316  //- Destructor
317  virtual ~KinematicCloud();
318 
319 
320  // Member Functions
321 
322  // Access
323 
324  //- Return a reference to the cloud copy
325  inline const KinematicCloud& cloudCopy() const;
326 
327  //- Switch to specify if particles of the cloud can return
328  // non-zero wall distance values - true for kinematic parcels
329  virtual bool hasWallImpactDistance() const;
330 
331 
332  // References to the mesh and databases
333 
334  //- Return reference to the mesh
335  inline const fvMesh& mesh() const;
336 
337  //- Return particle properties dictionary
338  inline const IOdictionary& particleProperties() const;
339 
340  //- Return output properties dictionary
341  inline const IOdictionary& outputProperties() const;
342 
343  //- Return non-const access to the output properties dictionary
344  inline IOdictionary& outputProperties();
345 
346  //- Return const access to the solution properties
347  inline const cloudSolution& solution() const;
348 
349  //- Return access to the solution properties
350  inline cloudSolution& solution();
351 
352  //- Return the constant properties
353  inline const typename parcelType::constantProperties&
354  constProps() const;
355 
356  //- Return access to the constant properties
357  inline typename parcelType::constantProperties& constProps();
358 
359  //- Return reference to the sub-models dictionary
360  inline const dictionary& subModelProperties() const;
361 
362 
363  // Cloud data
364 
365  //- Return reference to the random object
366  inline cachedRandom& rndGen();
367 
368  //- Return the cell occupancy information for each
369  // parcel, non-const access, the caller is
370  // responsible for updating it for its own purposes
371  // if particles are removed or created.
373 
374  //- Return the cell length scale
375  inline const scalarField& cellLengthScale() const;
376 
377 
378  // References to the carrier gas fields
379 
380  //- Return carrier gas velocity
381  inline const volVectorField& U() const;
382 
383  //- Return carrier gas density
384  inline const volScalarField& rho() const;
385 
386  //- Return carrier gas dynamic viscosity
387  inline const volScalarField& mu() const;
388 
389 
390  // Environmental properties
391 
392  //- Gravity
393  inline const dimensionedVector& g() const;
394 
395  //- Return const-access to the ambient pressure
396  inline scalar pAmbient() const;
397 
398  //- Return reference to the ambient pressure
399  inline scalar& pAmbient();
400 
401 
402  //- Optional particle forces
403 // inline const typename parcelType::forceType& forces() const;
404  inline const forceType& forces() const;
405 
406  //- Return the optional particle forces
407  inline forceType& forces();
408 
409  //- Optional cloud function objects
410  inline functionType& functions();
411 
412 
413  // Sub-models
414 
415  //- Return const access to the injection model
417  injectors() const;
418 
419  //- Return reference to the injection model
421  injectors();
422 
423  //- Return const-access to the dispersion model
425  dispersion() const;
426 
427  //- Return reference to the dispersion model
429  dispersion();
430 
431  //- Return const-access to the patch interaction model
433  patchInteraction() const;
434 
435  //- Return reference to the patch interaction model
438 
439  //- Return const-access to the stochastic collision model
440  inline const
442  stochasticCollision() const;
443 
444  //- Return reference to the stochastic collision model
447 
448  //- Return const-access to the surface film model
450  surfaceFilm() const;
451 
452  //- Return reference to the surface film model
454  surfaceFilm();
455 
456 
457  // Integration schemes
458 
459  //-Return reference to velocity integration
460  inline const vectorIntegrationScheme& UIntegrator() const;
461 
462 
463  // Sources
464 
465  // Momentum
466 
467  //- Return reference to momentum source
469 
470  //- Return const reference to momentum source
472  UTrans() const;
473 
474  //- Return coefficient for carrier phase U equation
476 
477  //- Return const coefficient for carrier phase U equation
479  UCoeff() const;
480 
481  //- Return tmp momentum source term
482  inline tmp<fvVectorMatrix> SU(volVectorField& U) const;
483 
484 
485  // Check
486 
487  //- Total number of parcels
488  inline label nParcels() const;
489 
490  //- Total mass in system
491  inline scalar massInSystem() const;
492 
493  //- Total linear momentum of the system
494  inline vector linearMomentumOfSystem() const;
495 
496  //- Total linear kinetic energy in the system
497  inline scalar linearKineticEnergyOfSystem() const;
498 
499  //- Total rotational kinetic energy in the system
500  inline scalar rotationalKineticEnergyOfSystem() const;
501 
502  //- Penetration for fraction [0-1] of the current total mass
503  inline scalar penetration(const scalar fraction) const;
504 
505  //- Mean diameter Dij
506  inline scalar Dij(const label i, const label j) const;
507 
508  //- Max diameter
509  inline scalar Dmax() const;
510 
511 
512  // Fields
513 
514  //- Volume swept rate of parcels per cell
515  inline const tmp<volScalarField> vDotSweep() const;
516 
517  //- Return the particle volume fraction field
518  // Note: for particles belonging to this cloud only
519  inline const tmp<volScalarField> theta() const;
520 
521  //- Return the particle mass fraction field
522  // Note: for particles belonging to this cloud only
523  inline const tmp<volScalarField> alpha() const;
524 
525  //- Return the particle effective density field
526  // Note: for particles belonging to this cloud only
527  inline const tmp<volScalarField> rhoEff() const;
528 
529 
530  // Cloud evolution functions
531 
532  //- Set parcel thermo properties
534  (
535  parcelType& parcel,
536  const scalar lagrangianDt
537  );
538 
539  //- Check parcel properties
541  (
542  parcelType& parcel,
543  const scalar lagrangianDt,
544  const bool fullyDescribed
545  );
546 
547  //- Store the current cloud state
548  void storeState();
549 
550  //- Reset the current cloud to the previously stored state
551  void restoreState();
552 
553  //- Reset the cloud source terms
554  void resetSourceTerms();
555 
556  //- Relax field
557  template<class Type>
558  void relax
559  (
561  const DimensionedField<Type, volMesh>& field0,
562  const word& name
563  ) const;
564 
565  //- Scale field
566  template<class Type>
567  void scale
568  (
570  const word& name
571  ) const;
572 
573  //- Apply relaxation to (steady state) cloud sources
574  void relaxSources(const KinematicCloud<CloudType>& cloudOldTime);
575 
576  //- Apply scaling to (transient) cloud sources
577  void scaleSources();
578 
579  //- Pre-evolve
580  void preEvolve();
581 
582  //- Evolve the cloud
583  void evolve();
584 
585  //- Particle motion
586  template<class TrackData>
587  void motion(TrackData& td);
588 
589  //- Calculate the patch normal and velocity to interact with,
590  // accounting for patch motion if required.
591  void patchData
592  (
593  const parcelType& p,
594  const polyPatch& pp,
595  const scalar trackFraction,
596  const tetIndices& tetIs,
597  vector& normal,
598  vector& Up
599  ) const;
600 
601 
602  // Mapping
603 
604  //- Update mesh
605  void updateMesh();
606 
607  //- Remap the cells of particles corresponding to the
608  // mesh topology change with a default tracking data object
609  virtual void autoMap(const mapPolyMesh&);
610 
611 
612  // I-O
613 
614  //- Print cloud information
615  void info();
616 };
617 
618 
619 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
620 
621 } // End namespace Foam
622 
623 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
624 
625 #include "KinematicCloudI.H"
626 
627 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
628 
629 #ifdef NoRepository
630 # include "KinematicCloud.C"
631 #endif
632 
633 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
634 
635 #endif
636 
637 // ************************************************************************* //
Templated wall surface film model class.
Top level model for Integration schemes.
const vectorIntegrationScheme & UIntegrator() const
Return reference to velocity integration.
cloudSolution solution_
Solution properties.
void resetSourceTerms()
Reset the cloud source terms.
scalar Dij(const label i, const label j) const
Mean diameter Dij.
scalar massInSystem() const
Total mass in system.
parcelType::constantProperties constProps_
Parcel constant properties.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
Foam::autoPtr< IOobject > clone() const
Clone.
Definition: IOobject.H:239
void updateMesh()
Update mesh.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
functionType & functions()
Optional cloud function objects.
scalarField cellLengthScale_
Cell length scale.
void checkParcelProperties(parcelType &parcel, const scalar lagrangianDt, const bool fullyDescribed)
Check parcel properties.
const volScalarField & mu() const
Return carrier gas dynamic viscosity.
void evolveCloud(TrackData &td)
Evolve the cloud.
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Definition: mapPolyMesh.H:158
scalar linearKineticEnergyOfSystem() const
Total linear kinetic energy in the system.
Virtual abstract base class for templated KinematicCloud.
const scalarField & cellLengthScale() const
Return the cell length scale.
void solve(TrackData &td)
Solve the cloud - calls all evolution functions.
const parcelType::constantProperties & constProps() const
Return the constant properties.
const IOdictionary & particleProperties() const
Return particle properties dictionary.
void buildCellOccupancy()
Build the cellOccupancy.
A class for handling words, derived from string.
Definition: word.H:59
scalar pAmbient_
Averaged ambient domain pressure.
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
void preEvolve()
Pre-evolve.
cachedRandom rndGen_
Random number generator - used by some injection routines.
const tmp< volScalarField > vDotSweep() const
Volume swept rate of parcels per cell.
const SurfaceFilmModel< KinematicCloud< CloudType > > & surfaceFilm() const
Return const-access to the surface film model.
tmp< fvVectorMatrix > SU(volVectorField &U) const
Return tmp momentum source term.
const tmp< volScalarField > rhoEff() const
Return the particle effective density field.
const InjectionModelList< KinematicCloud< CloudType > > & injectors() const
Return const access to the injection model.
CloudType::particleType parcelType
Type of parcel the cloud was instantiated for.
const volScalarField & rho() const
Return carrier gas density.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
void relaxSources(const KinematicCloud< CloudType > &cloudOldTime)
Apply relaxation to (steady state) cloud sources.
Templated base class for kinematic cloud.
const dimensionedVector & g_
Gravity.
Namespace for OpenFOAM.
scalar Dmax() const
Max diameter.
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
Templated patch interaction model class.
List of particle forces.
autoPtr< List< DynamicList< parcelType * > > > cellOccupancyPtr_
Cell occupancy information for each parcel, (demand driven)
autoPtr< StochasticCollisionModel< KinematicCloud< CloudType > > > stochasticCollisionModel_
Stochastic collision model.
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: HashTable.H:59
void relax(DimensionedField< Type, volMesh > &field, const DimensionedField< Type, volMesh > &field0, const word &name) const
Relax field.
const tmp< volScalarField > theta() const
Return the particle volume fraction field.
CloudType cloudType
Type of cloud this cloud was instantiated for.
const fvMesh & mesh_
References to the mesh and time databases.
volScalarField & p
Definition: createFields.H:51
scalar penetration(const scalar fraction) const
Penetration for fraction [0-1] of the current total mass.
autoPtr< DimensionedField< scalar, volMesh > > UCoeff_
Coefficient for carrier phase U equation.
const dictionary subModelProperties_
Sub-models dictionary.
scalar rotationalKineticEnergyOfSystem() const
Total rotational kinetic energy in the system.
void restoreState()
Reset the current cloud to the previously stored state.
cachedRandom & rndGen()
Return reference to the random object.
scalar pAmbient() const
Return const-access to the ambient pressure.
virtual void autoMap(const mapPolyMesh &)
Remap the cells of particles corresponding to the.
const volScalarField & mu_
Dynamic viscosity [Pa.s].
void postEvolve()
Post-evolve.
List of cloud function objects.
const DispersionModel< KinematicCloud< CloudType > > & dispersion() const
Return const-access to the dispersion model.
void scaleSources()
Apply scaling to (transient) cloud sources.
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:53
KinematicCloud< CloudType > kinematicCloudType
Convenience typedef for this cloud type.
Templated stochastic collision model class.
const KinematicCloud & cloudCopy() const
Return a reference to the cloud copy.
autoPtr< DimensionedField< vector, volMesh > > UTrans_
Momentum.
const dimensionedVector & g() const
Gravity.
const PatchInteractionModel< KinematicCloud< CloudType > > & patchInteraction() const
Return const-access to the patch interaction model.
const word & name() const
Return name.
Definition: IOobject.H:260
void setParcelThermoProperties(parcelType &parcel, const scalar lagrangianDt)
Set parcel thermo properties.
forceType forces_
Optional particle forces.
const word & cloudName() const
Return the cloud type.
Definition: DSMCCloudI.H:34
void updateCellOccupancy()
Update (i.e. build) the cellOccupancy if it has.
void scale(DimensionedField< Type, volMesh > &field, const word &name) const
Scale field.
const volScalarField & rho_
Density [kg/m3].
vector linearMomentumOfSystem() const
Total linear momentum of the system.
void storeState()
Store the current cloud state.
const volVectorField & U() const
Return carrier gas velocity.
autoPtr< SurfaceFilmModel< KinematicCloud< CloudType > > > surfaceFilmModel_
Surface film model.
Random number generator.
Definition: cachedRandom.H:63
IOdictionary outputProperties_
Dictionary of output properties.
virtual ~KinematicCloud()
Destructor.
ParticleForceList< KinematicCloud< CloudType > > forceType
Force models type.
virtual bool hasWallImpactDistance() const
Switch to specify if particles of the cloud can return.
functionType functions_
Optional cloud function objects.
DimensionedField< scalar, volMesh > & UCoeff()
Return coefficient for carrier phase U equation.
const volVectorField & U_
Velocity [m/s].
A normal distribution model.
autoPtr< PatchInteractionModel< KinematicCloud< CloudType > > > patchInteractionModel_
Patch interaction model.
const dimensionedScalar c
Speed of light in a vacuum.
autoPtr< DispersionModel< KinematicCloud< CloudType > > > dispersionModel_
Dispersion model.
List< DynamicList< parcelType * > > & cellOccupancy()
Return the cell occupancy information for each.
CloudFunctionObjectList< KinematicCloud< CloudType > > functionType
Function object type.
void evolve()
Evolve the cloud.
void cloudReset(KinematicCloud< CloudType > &c)
Reset state of cloud.
IOdictionary particleProperties_
Dictionary of particle properties.
Templated base class for dsmc cloud.
Definition: DSMCCloud.H:68
label nParcels() const
Total number of parcels.
ParcelType particleType
Definition: Cloud.H:114
InjectionModelList< KinematicCloud< CloudType > > injectors_
Injector models.
This function object reads fields from the time directories and adds them to the mesh database for fu...
Definition: readFields.H:102
void setModels()
Set cloud sub-models.
const StochasticCollisionModel< KinematicCloud< CloudType > > & stochasticCollision() const
Return const-access to the stochastic collision model.
autoPtr< vectorIntegrationScheme > UIntegrator_
Velocity integration.
Storage and named access for the indices of a tet which is part of the decomposition of a cell...
Definition: tetIndices.H:73
void motion(TrackData &td)
Particle motion.
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:117
const cloudSolution & solution() const
Return const access to the solution properties.
A class for managing temporary objects.
Definition: PtrList.H:118
Stores all relevant solution info for cloud.
Definition: cloudSolution.H:51
void patchData(const parcelType &p, const polyPatch &pp, const scalar trackFraction, const tetIndices &tetIs, vector &normal, vector &Up) const
Calculate the patch normal and velocity to interact with,.
DimensionedField< vector, volMesh > & UTrans()
Return reference to momentum source.
const dictionary & subModelProperties() const
Return reference to the sub-models dictionary.
void info()
Print cloud information.
virtual autoPtr< Cloud< parcelType > > cloneBare(const word &name) const
Construct and return bare clone based on (this) with new name.
A special matrix type and solver, designed for finite volume solutions of scalar equations.
const tmp< volScalarField > alpha() const
Return the particle mass fraction field.
const forceType & forces() const
Optional particle forces.
const fvMesh & mesh() const
Return reference to the mesh.
List of injection models.
const IOdictionary & outputProperties() const
Return output properties dictionary.