KinematicCloud.C
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 \*---------------------------------------------------------------------------*/
25 
26 #include "KinematicCloud.H"
27 #include "integrationScheme.H"
28 #include "interpolation.H"
29 #include "subCycleTime.H"
30 
31 #include "InjectionModelList.H"
32 #include "DispersionModel.H"
33 #include "PatchInteractionModel.H"
35 #include "SurfaceFilmModel.H"
36 
37 // * * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * //
38 
39 template<class CloudType>
41 {
42  dispersionModel_.reset
43  (
45  (
46  subModelProperties_,
47  *this
48  ).ptr()
49  );
50 
51  patchInteractionModel_.reset
52  (
54  (
55  subModelProperties_,
56  *this
57  ).ptr()
58  );
59 
60  stochasticCollisionModel_.reset
61  (
63  (
64  subModelProperties_,
65  *this
66  ).ptr()
67  );
68 
69  surfaceFilmModel_.reset
70  (
72  (
73  subModelProperties_,
74  *this
75  ).ptr()
76  );
77 
78  UIntegrator_.reset
79  (
81  (
82  "U",
83  solution_.integrationSchemes()
84  ).ptr()
85  );
86 }
87 
88 
89 template<class CloudType>
90 template<class TrackCloudType>
92 (
93  TrackCloudType& cloud,
94  typename parcelType::trackingData& td
95 )
96 {
97  if (solution_.steadyState())
98  {
99  cloud.storeState();
100 
101  cloud.preEvolve();
102 
103  evolveCloud(cloud, td);
104 
105  if (solution_.coupled())
106  {
107  cloud.relaxSources(cloud.cloudCopy());
108  }
109  }
110  else
111  {
112  cloud.preEvolve();
113 
114  evolveCloud(cloud, td);
115 
116  if (solution_.coupled())
117  {
118  cloud.scaleSources();
119  }
120  }
121 
122  cloud.info();
123 
124  cloud.postEvolve();
125 
126  if (solution_.steadyState())
127  {
128  cloud.restoreState();
129  }
130 }
131 
132 
133 template<class CloudType>
135 {
136  if (cellOccupancyPtr_.empty())
137  {
138  cellOccupancyPtr_.reset
139  (
140  new List<DynamicList<parcelType*>>(mesh_.nCells())
141  );
142  }
143  else if (cellOccupancyPtr_().size() != mesh_.nCells())
144  {
145  // If the size of the mesh has changed, reset the
146  // cellOccupancy size
147 
148  cellOccupancyPtr_().setSize(mesh_.nCells());
149  }
150 
151  List<DynamicList<parcelType*>>& cellOccupancy = cellOccupancyPtr_();
152 
153  forAll(cellOccupancy, cO)
154  {
155  cellOccupancy[cO].clear();
156  }
157 
158  forAllIter(typename KinematicCloud<CloudType>, *this, iter)
159  {
160  cellOccupancy[iter().cell()].append(&iter());
161  }
162 }
163 
164 
165 template<class CloudType>
167 {
168  // Only build the cellOccupancy if the pointer is set, i.e. it has
169  // been requested before.
170 
171  if (cellOccupancyPtr_.valid())
172  {
173  buildCellOccupancy();
174  }
175 }
176 
177 
178 template<class CloudType>
179 template<class TrackCloudType>
181 (
182  TrackCloudType& cloud,
183  typename parcelType::trackingData& td
184 )
185 {
186  if (solution_.coupled())
187  {
188  cloud.resetSourceTerms();
189  }
190 
191  if (solution_.transient())
192  {
193  label preInjectionSize = this->size();
194 
195  this->surfaceFilm().inject(cloud);
196 
197  // Update the cellOccupancy if the size of the cloud has changed
198  // during the injection.
199  if (preInjectionSize != this->size())
200  {
201  updateCellOccupancy();
202  preInjectionSize = this->size();
203  }
204 
205  injectors_.inject(cloud, td);
206 
207 
208  // Assume that motion will update the cellOccupancy as necessary
209  // before it is required.
210  cloud.motion(cloud, td);
211 
212  stochasticCollision().update(td, solution_.trackTime());
213  }
214  else
215  {
216 // this->surfaceFilm().injectSteadyState(cloud);
217 
218  injectors_.injectSteadyState(cloud, td, solution_.trackTime());
219 
220  td.part() = parcelType::trackingData::tpLinearTrack;
221  CloudType::move(cloud, td, solution_.trackTime());
222  }
223 }
224 
225 
226 template<class CloudType>
228 {
229  Info<< endl;
230 
231  if (debug)
232  {
233  this->writePositions();
234  }
235 
236  this->dispersion().cacheFields(false);
237 
238  forces_.cacheFields(false);
239 
240  functions_.postEvolve();
241 
242  solution_.nextIter();
243 
244  if (this->db().time().writeTime())
245  {
246  outputProperties_.writeObject
247  (
248  IOstream::ASCII,
249  IOstream::currentVersion,
250  this->db().time().writeCompression(),
251  true
252  );
253  }
254 }
255 
256 
257 template<class CloudType>
259 {
260  CloudType::cloudReset(c);
261 
262  rndGen_ = c.rndGen_;
263 
264  forces_.transfer(c.forces_);
265 
266  functions_.transfer(c.functions_);
267 
268  injectors_.transfer(c.injectors_);
269 
270  dispersionModel_.reset(c.dispersionModel_.ptr());
271  patchInteractionModel_.reset(c.patchInteractionModel_.ptr());
272  stochasticCollisionModel_.reset(c.stochasticCollisionModel_.ptr());
273  surfaceFilmModel_.reset(c.surfaceFilmModel_.ptr());
274 
275  UIntegrator_.reset(c.UIntegrator_.ptr());
276 }
277 
278 
279 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
280 
281 template<class CloudType>
283 (
284  const word& cloudName,
285  const volScalarField& rho,
286  const volVectorField& U,
287  const volScalarField& mu,
288  const dimensionedVector& g,
289  bool readFields
290 )
291 :
292  CloudType(rho.mesh(), cloudName, false),
293  kinematicCloud(),
294  cloudCopyPtr_(nullptr),
295  mesh_(rho.mesh()),
296  particleProperties_
297  (
298  IOobject
299  (
300  cloudName + "Properties",
301  rho.mesh().time().constant(),
302  rho.mesh(),
303  IOobject::MUST_READ_IF_MODIFIED,
304  IOobject::NO_WRITE
305  )
306  ),
307  outputProperties_
308  (
309  IOobject
310  (
311  cloudName + "OutputProperties",
312  mesh_.time().timeName(),
313  "uniform"/cloud::prefix/cloudName,
314  mesh_,
315  IOobject::READ_IF_PRESENT,
316  IOobject::NO_WRITE
317  )
318  ),
319  solution_(mesh_, particleProperties_.subDict("solution")),
320  constProps_(particleProperties_),
321  subModelProperties_
322  (
323  particleProperties_.subOrEmptyDict("subModels", solution_.active())
324  ),
325  rndGen_(0),
326  cellOccupancyPtr_(),
327  cellLengthScale_(mag(cbrt(mesh_.V()))),
328  rho_(rho),
329  U_(U),
330  mu_(mu),
331  g_(g),
332  pAmbient_(0.0),
333  forces_
334  (
335  *this,
336  mesh_,
337  subModelProperties_.subOrEmptyDict
338  (
339  "particleForces",
340  solution_.active()
341  ),
342  solution_.active()
343  ),
344  functions_
345  (
346  *this,
347  particleProperties_.subOrEmptyDict("cloudFunctions"),
348  solution_.active()
349  ),
350  injectors_
351  (
352  subModelProperties_.subOrEmptyDict("injectionModels"),
353  *this
354  ),
355  dispersionModel_(nullptr),
356  patchInteractionModel_(nullptr),
357  stochasticCollisionModel_(nullptr),
358  surfaceFilmModel_(nullptr),
359  UIntegrator_(nullptr),
360  UTrans_
361  (
363  (
364  IOobject
365  (
366  this->name() + ":UTrans",
367  this->db().time().timeName(),
368  this->db(),
369  IOobject::READ_IF_PRESENT,
370  IOobject::AUTO_WRITE
371  ),
372  mesh_,
374  )
375  ),
376  UCoeff_
377  (
379  (
380  IOobject
381  (
382  this->name() + ":UCoeff",
383  this->db().time().timeName(),
384  this->db(),
385  IOobject::READ_IF_PRESENT,
386  IOobject::AUTO_WRITE
387  ),
388  mesh_,
389  dimensionedScalar("zero", dimMass, 0.0)
390  )
391  )
392 {
393  if (solution_.active())
394  {
395  setModels();
396 
397  if (readFields)
398  {
399  parcelType::readFields(*this);
400  this->deleteLostParticles();
401  }
402  }
403 
404  if (solution_.resetSourcesOnStartup())
405  {
406  resetSourceTerms();
407  }
408 }
409 
410 
411 template<class CloudType>
413 (
415  const word& name
416 )
417 :
418  CloudType(c.mesh_, name, c),
419  kinematicCloud(),
420  cloudCopyPtr_(nullptr),
421  mesh_(c.mesh_),
422  particleProperties_(c.particleProperties_),
423  outputProperties_(c.outputProperties_),
424  solution_(c.solution_),
425  constProps_(c.constProps_),
426  subModelProperties_(c.subModelProperties_),
427  rndGen_(c.rndGen_),
428  cellOccupancyPtr_(nullptr),
429  cellLengthScale_(c.cellLengthScale_),
430  rho_(c.rho_),
431  U_(c.U_),
432  mu_(c.mu_),
433  g_(c.g_),
434  pAmbient_(c.pAmbient_),
435  forces_(c.forces_),
436  functions_(c.functions_),
437  injectors_(c.injectors_),
438  dispersionModel_(c.dispersionModel_->clone()),
439  patchInteractionModel_(c.patchInteractionModel_->clone()),
440  stochasticCollisionModel_(c.stochasticCollisionModel_->clone()),
441  surfaceFilmModel_(c.surfaceFilmModel_->clone()),
442  UIntegrator_(c.UIntegrator_->clone()),
443  UTrans_
444  (
446  (
447  IOobject
448  (
449  this->name() + ":UTrans",
450  this->db().time().timeName(),
451  this->db(),
452  IOobject::NO_READ,
453  IOobject::NO_WRITE,
454  false
455  ),
456  c.UTrans_()
457  )
458  ),
459  UCoeff_
460  (
462  (
463  IOobject
464  (
465  name + ":UCoeff",
466  this->db().time().timeName(),
467  this->db(),
468  IOobject::NO_READ,
469  IOobject::NO_WRITE,
470  false
471  ),
472  c.UCoeff_()
473  )
474  )
475 {}
476 
477 
478 template<class CloudType>
480 (
481  const fvMesh& mesh,
482  const word& name,
484 )
485 :
486  CloudType(mesh, name, IDLList<parcelType>()),
487  kinematicCloud(),
488  cloudCopyPtr_(nullptr),
489  mesh_(mesh),
490  particleProperties_
491  (
492  IOobject
493  (
494  name + "Properties",
495  mesh.time().constant(),
496  mesh,
497  IOobject::NO_READ,
498  IOobject::NO_WRITE,
499  false
500  )
501  ),
502  outputProperties_
503  (
504  IOobject
505  (
506  name + "OutputProperties",
507  mesh_.time().timeName(),
508  "uniform"/cloud::prefix/name,
509  mesh_,
510  IOobject::NO_READ,
511  IOobject::NO_WRITE,
512  false
513  )
514  ),
515  solution_(mesh),
516  constProps_(),
517  subModelProperties_(dictionary::null),
518  rndGen_(0),
519  cellOccupancyPtr_(nullptr),
520  cellLengthScale_(c.cellLengthScale_),
521  rho_(c.rho_),
522  U_(c.U_),
523  mu_(c.mu_),
524  g_(c.g_),
525  pAmbient_(c.pAmbient_),
526  forces_(*this, mesh),
527  functions_(*this),
528  injectors_(*this),
529  dispersionModel_(nullptr),
530  patchInteractionModel_(nullptr),
531  stochasticCollisionModel_(nullptr),
532  surfaceFilmModel_(nullptr),
533  UIntegrator_(nullptr),
534  UTrans_(nullptr),
535  UCoeff_(nullptr)
536 {}
537 
538 
539 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
540 
541 template<class CloudType>
543 {}
544 
545 
546 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
547 
548 template<class CloudType>
550 (
551  parcelType& parcel,
552  const scalar lagrangianDt
553 )
554 {
555  parcel.rho() = constProps_.rho0();
556 }
557 
558 
559 template<class CloudType>
561 (
562  parcelType& parcel,
563  const scalar lagrangianDt,
564  const bool fullyDescribed
565 )
566 {
567  const scalar carrierDt = mesh_.time().deltaTValue();
568  parcel.stepFraction() = (carrierDt - lagrangianDt)/carrierDt;
569 
570  if (parcel.typeId() == -1)
571  {
572  parcel.typeId() = constProps_.parcelTypeId();
573  }
574 }
575 
576 
577 template<class CloudType>
579 {
580  cloudCopyPtr_.reset
581  (
582  static_cast<KinematicCloud<CloudType>*>
583  (
584  clone(this->name() + "Copy").ptr()
585  )
586  );
587 }
588 
589 
590 template<class CloudType>
592 {
593  cloudReset(cloudCopyPtr_());
594  cloudCopyPtr_.clear();
595 }
596 
597 
598 template<class CloudType>
600 {
601  UTrans().field() = Zero;
602  UCoeff().field() = 0.0;
603 }
604 
605 
606 template<class CloudType>
607 template<class Type>
609 (
611  const DimensionedField<Type, volMesh>& field0,
612  const word& name
613 ) const
614 {
615  const scalar coeff = solution_.relaxCoeff(name);
616  field = field0 + coeff*(field - field0);
617 }
618 
619 
620 template<class CloudType>
621 template<class Type>
623 (
625  const word& name
626 ) const
627 {
628  const scalar coeff = solution_.relaxCoeff(name);
629  field *= coeff;
630 }
631 
632 
633 template<class CloudType>
635 (
636  const KinematicCloud<CloudType>& cloudOldTime
637 )
638 {
639  this->relax(UTrans_(), cloudOldTime.UTrans(), "U");
640  this->relax(UCoeff_(), cloudOldTime.UCoeff(), "U");
641 }
642 
643 
644 template<class CloudType>
646 {
647  this->scale(UTrans_(), "U");
648  this->scale(UCoeff_(), "U");
649 }
650 
651 
652 template<class CloudType>
654 {
655  // force calculaion of mesh dimensions - needed for parallel runs
656  // with topology change due to lazy evaluation of valid mesh dimensions
657  label nGeometricD = mesh_.nGeometricD();
658 
659  Info<< "\nSolving " << nGeometricD << "-D cloud " << this->name() << endl;
660 
661  this->dispersion().cacheFields(true);
662  forces_.cacheFields(true);
663  updateCellOccupancy();
664 
665  pAmbient_ = constProps_.dict().template
666  lookupOrDefault<scalar>("pAmbient", pAmbient_);
667 
668  functions_.preEvolve();
669 }
670 
671 
672 template<class CloudType>
674 {
675  if (solution_.canEvolve())
676  {
677  typename parcelType::trackingData td(*this);
678 
679  solve(*this, td);
680  }
681 }
682 
683 
684 template<class CloudType>
685 template<class TrackCloudType>
687 (
688  TrackCloudType& cloud,
689  typename parcelType::trackingData& td
690 )
691 {
692  td.part() = parcelType::trackingData::tpLinearTrack;
693  CloudType::move(cloud, td, solution_.trackTime());
694 
695  updateCellOccupancy();
696 }
697 
698 
699 template<class CloudType>
701 (
702  const parcelType& p,
703  const polyPatch& pp,
704  vector& nw,
705  vector& Up
706 ) const
707 {
708  p.patchData(nw, Up);
709 
710  // If this is a wall patch, then there may be a non-zero tangential velocity
711  // component; the lid velocity in a lid-driven cavity case, for example. We
712  // want the particle to interact with this velocity, so we look it up in the
713  // velocity field and use it to set the wall-tangential component.
714  if (isA<wallPolyPatch>(pp))
715  {
716  const label patchi = pp.index();
717  const label patchFacei = pp.whichFace(p.face());
718 
719  // We only want to use the boundary condition value onlyif it is set
720  // by the boundary condition. If the boundary values are extrapolated
721  // (e.g., slip conditions) then they represent the motion of the fluid
722  // just inside the domain rather than that of the wall itself.
723  if (U_.boundaryField()[patchi].fixesValue())
724  {
725  const vector Uw1 = U_.boundaryField()[patchi][patchFacei];
726  const vector& Uw0 =
727  U_.oldTime().boundaryField()[patchi][patchFacei];
728 
729  const scalar f = p.currentTimeFraction();
730 
731  const vector Uw = Uw0 + f*(Uw1 - Uw0);
732 
733  const tensor nnw = nw*nw;
734 
735  Up = (nnw & Up) + Uw - (nnw & Uw);
736  }
737  }
738 }
739 
740 
741 template<class CloudType>
743 {
744  updateCellOccupancy();
745  injectors_.updateMesh();
746  cellLengthScale_ = mag(cbrt(mesh_.V()));
747 }
748 
749 
750 template<class CloudType>
752 {
754 
755  updateMesh();
756 }
757 
758 
759 template<class CloudType>
761 {
762  vector linearMomentum = linearMomentumOfSystem();
763  reduce(linearMomentum, sumOp<vector>());
764 
765  scalar linearKineticEnergy = linearKineticEnergyOfSystem();
766  reduce(linearKineticEnergy, sumOp<scalar>());
767 
768  Info<< "Cloud: " << this->name() << nl
769  << " Current number of parcels = "
770  << returnReduce(this->size(), sumOp<label>()) << nl
771  << " Current mass in system = "
772  << returnReduce(massInSystem(), sumOp<scalar>()) << nl
773  << " Linear momentum = "
774  << linearMomentum << nl
775  << " |Linear momentum| = "
776  << mag(linearMomentum) << nl
777  << " Linear kinetic energy = "
778  << linearKineticEnergy << nl;
779 
780  injectors_.info(Info);
781  this->surfaceFilm().info(Info);
782  this->patchInteraction().info(Info);
783 }
784 
785 
786 // ************************************************************************* //
Template class for intrusive linked lists.
Definition: ILList.H:50
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
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
DSMCCloud< dsmcParcel > CloudType
void scaleSources()
Apply scaling to (transient) cloud sources.
void setModels()
Set cloud sub-models.
const volVectorField & U_
Velocity [m/s].
void relaxSources(const KinematicCloud< CloudType > &cloudOldTime)
Apply relaxation to (steady state) cloud sources.
void storeState()
Store the current cloud state.
UEqn relax()
autoPtr< PatchInteractionModel< KinematicCloud< CloudType > > > patchInteractionModel_
Patch interaction model.
void readFields(const typename GeoFieldType::Mesh &mesh, const IOobjectList &objects, const HashSet< word > &selectedFields, LIFOStack< regIOobject *> &storedObjects)
Read the selected GeometricFields of the specified type.
Definition: ReadFields.C:244
InjectionModelList< KinematicCloud< CloudType > > injectors_
Injector models.
Cloud< basicKinematicCollidingParcel > ::particleType parcelType
Type of parcel the cloud was instantiated for.
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: HashTable.H:60
virtual ~KinematicCloud()
Destructor.
#define forAllIter(Container, container, iter)
Iterate across all elements in the container object of type.
Definition: UList.H:453
void evolveCloud(TrackCloudType &cloud, typename parcelType::trackingData &td)
Evolve the cloud.
autoPtr< volVectorField::Internal > UTrans_
Momentum.
void setParcelThermoProperties(parcelType &parcel, const scalar lagrangianDt)
Set parcel thermo properties.
dimensioned< vector > dimensionedVector
Dimensioned vector obtained from generic dimensioned type.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
void relax(DimensionedField< Type, volMesh > &field, const DimensionedField< Type, volMesh > &field0, const word &name) const
Relax field.
volVectorField::Internal & UTrans()
Return reference to momentum source.
void postEvolve()
Post-evolve.
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:243
autoPtr< SurfaceFilmModel< KinematicCloud< CloudType > > > surfaceFilmModel_
Surface film model.
Random rndGen_
Random number generator - used by some injection routines.
void updateCellOccupancy()
Update (i.e. build) the cellOccupancy if it has.
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Definition: mapPolyMesh.H:158
void buildCellOccupancy()
Build the cellOccupancy.
const volScalarField & rho_
Density [kg/m3].
Templated patch interaction model class.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
void updateMesh()
Update mesh.
const volScalarField & mu_
Dynamic viscosity [Pa.s].
dynamicFvMesh & mesh
volScalarField::Internal & UCoeff()
Return coefficient for carrier phase U equation.
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects...
Definition: DynamicList.H:56
void resetSourceTerms()
Reset the cloud source terms.
void clear()
Clear the list, i.e. set size to zero.
Definition: ListI.H:124
A class for handling words, derived from string.
Definition: word.H:59
void motion(TrackCloudType &cloud, typename parcelType::trackingData &td)
Particle motion.
void append(const T &)
Append an element at the end of the list.
Definition: ListI.H:184
rhoEqn solve()
dimensionedScalar cbrt(const dimensionedScalar &ds)
scalar pAmbient_
Averaged ambient domain pressure.
const word & constant() const
Return constant name.
Definition: TimePaths.H:124
parcelType::constantProperties constProps_
Parcel constant properties.
const fvMesh & mesh_
References to the mesh and time databases.
word timeName
Definition: getTimeIndex.H:3
static const zero Zero
Definition: zero.H:97
void checkParcelProperties(parcelType &parcel, const scalar lagrangianDt, const bool fullyDescribed)
Check parcel properties.
const dictionary subModelProperties_
Sub-models dictionary.
Base cloud calls templated on particle type.
Definition: Cloud.H:52
regionModels::surfaceFilmModel & surfaceFilm
static const char nl
Definition: Ostream.H:265
const Mesh & mesh() const
Return mesh.
labelList f(nPoints)
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
Virtual abstract base class for templated KinematicCloud.
cloudSolution solution_
Solution properties.
Templated base class for kinematic cloud.
void patchData(const parcelType &p, const polyPatch &pp, vector &normal, vector &Up) const
Calculate the patch normal and velocity to interact with,.
label patchi
const List< DynamicList< molecule * > > & cellOccupancy
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
void evolve()
Evolve the cloud.
IOdictionary particleProperties_
Dictionary of particle properties.
autoPtr< volScalarField::Internal > UCoeff_
Coefficient for carrier phase U equation.
Templated stochastic collision model class.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
label index() const
Return the index of this patch in the boundaryMesh.
void info()
Print cloud information.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
forceType forces_
Optional particle forces.
void restoreState()
Reset the current cloud to the previously stored state.
messageStream Info
dimensioned< scalar > mag(const dimensioned< Type > &)
functionType functions_
Optional cloud function objects.
autoPtr< integrationScheme > UIntegrator_
Velocity integration.
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
void scale(DimensionedField< Type, volMesh > &field, const word &name) const
Scale field.
const word cloudName(propsDict.lookup("cloudName"))
void cloudReset(KinematicCloud< CloudType > &c)
Reset state of cloud.
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:49
virtual void autoMap(const mapPolyMesh &)
Remap the cells of particles corresponding to the.
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
autoPtr< DispersionModel< KinematicCloud< CloudType > > > dispersionModel_
Dispersion model.
IOdictionary outputProperties_
Dictionary of output properties.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
const dimensionedVector & g_
Gravity.
void preEvolve()
Pre-evolve.
void solve(TrackCloudType &cloud, typename parcelType::trackingData &td)
Solve the cloud - calls all evolution functions.
label nw
Definition: createFields.H:25
autoPtr< StochasticCollisionModel< KinematicCloud< CloudType > > > stochasticCollisionModel_
Stochastic collision model.
label whichFace(const label l) const
Return label of face in patch from global face label.
Definition: polyPatch.H:377
const dimensionSet dimVelocity
scalarField cellLengthScale_
Cell length scale.