Go to the documentation of this file.
39 template<
class CloudType>
48 generationSetName_(this->coeffDict().
lookup(
"generationCellSet")),
49 inflationSetName_(this->coeffDict().
lookup(
"inflationCellSet")),
52 duration_(this->coeffDict().getScalar(
"duration")),
72 volumeAccumulator_(0.0),
74 selfSeed_(this->coeffDict().getOrDefault(
"selfSeed",
false)),
80 this->coeffDict().subDict(
"sizeDistribution"),
86 const Time& time = owner.db().time();
88 flowRateProfile_->userTimeToTime(time);
89 growthRate_->userTimeToTime(time);
93 this->coeffDict().readEntry(
"dSeed", dSeed_);
96 cellSet generationCells(this->owner().
mesh(), generationSetName_);
98 generationCells_ = generationCells.toc();
100 cellSet inflationCells(this->owner().
mesh(), inflationSetName_);
103 inflationCells |= generationCells;
105 inflationCells_ = inflationCells.toc();
109 scalar generationVolume = 0.0;
111 forAll(generationCells_, gCI)
113 label cI = generationCells_[gCI];
118 scalar totalGenerationVolume = generationVolume;
122 fraction_ = generationVolume/totalGenerationVolume;
126 this->volumeTotal_ = fraction_*flowRateProfile_->integrate(0.0, duration_);
127 this->massTotal_ *= fraction_;
131 template<
class CloudType>
138 generationSetName_(im.generationSetName_),
139 inflationSetName_(im.inflationSetName_),
140 generationCells_(im.generationCells_),
141 inflationCells_(im.inflationCells_),
142 duration_(im.duration_),
143 flowRateProfile_(im.flowRateProfile_.clone()),
144 growthRate_(im.growthRate_.clone()),
145 newParticles_(im.newParticles_),
146 volumeAccumulator_(im.volumeAccumulator_),
147 fraction_(im.fraction_),
148 selfSeed_(im.selfSeed_),
150 sizeDistribution_(im.sizeDistribution_.clone())
156 template<
class CloudType>
163 template<
class CloudType>
168 template<
class CloudType>
171 return this->SOI_ + duration_;
175 template<
class CloudType>
187 scalar gR = growthRate_->value(time1);
189 scalar dT = time1 - time0;
193 forAll(inflationCells_, iCI)
195 label cI = inflationCells_[iCI];
203 scalar dTarget = pPtr->dTarget();
205 pPtr->d() =
min(dTarget, pPtr->d() + gR*dT);
211 newParticles_.clear();
219 if ((time0 >= 0.0) && (time0 < duration_))
221 volumeAccumulator_ +=
222 fraction_*flowRateProfile_->integrate(time0, time1);
228 label maxIterations =
max
231 (10*volumeAccumulator_)
235 label iterationNo = 0;
240 while (!generationCells_.empty() && volumeAccumulator_ > 0)
242 if (iterationNo > maxIterations)
245 <<
"Maximum particle split iterations ("
246 << maxIterations <<
") exceeded" <<
endl;
251 label cI = generationCells_
256 generationCells_.size() - 1
266 if (selfSeed_ && !cellCentresUsed.found(cI))
268 scalar dNew = sizeDistribution_().sample();
281 cellCentresUsed.
insert(cI);
294 scalar pD = pPtr->d();
297 if ((pD/pPtr->dTarget()) < rnd.
sample01<scalar>())
302 const point pP(pPtr->position());
303 const vector& pU = pPtr->U();
326 scalar
x = a/
sqrt(3.0);
327 scalar r = a/(2.0*
sqrt(6.0));
329 scalar d = a/(2.0*
sqrt(3.0));
331 scalar dNew = sizeDistribution_().sample();
342 volumeAccumulator_ -= volNew;
344 dNew = sizeDistribution_().sample();
353 volumeAccumulator_ -= volNew;
355 dNew = sizeDistribution_().sample();
364 volumeAccumulator_ -= volNew;
366 dNew = sizeDistribution_().sample();
375 volumeAccumulator_ -= volNew;
410 gatheredNewParticles,
417 newParticles_ = combinedNewParticles;
423 return newParticles_.size();
427 template<
class CloudType>
434 if ((time0 >= 0.0) && (time0 < duration_))
436 return fraction_*flowRateProfile_->integrate(time0, time1);
443 template<
class CloudType>
455 position = newParticles_[parcelI].first().first();
457 this->findCellAtPosition
468 template<
class CloudType>
477 parcel.U() = newParticles_[parcelI].first().second();
479 parcel.d() = newParticles_[parcelI].second().first();
481 parcel.dTarget() = newParticles_[parcelI].second().second();
485 template<
class CloudType>
492 template<
class CloudType>
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
A class for handling words, derived from Foam::string.
InflationInjection(const dictionary &dict, CloudType &owner, const word &modelName)
Construct from dictionary.
virtual label parcelsToInject(const scalar time0, const scalar time1)
Number of parcels to introduce relative to SOI.
static constexpr const zero Zero
Global zero (0)
Type sample01()
Return a sample whose components lie in the range [0,1].
virtual scalar userTimeToTime(const scalar theta) const
Convert the user-time (e.g. CA deg) to real-time (s).
Templated injection model class.
static bool master(const label communicator=worldComm)
Am I the master process.
static void scatter(const List< commsStruct > &comms, T &Value, const int tag, const label comm)
Scatter data. Distribute without modification. Reverse of gather.
Ostream & endl(Ostream &os)
Add newline and flush stream.
const wordList volume
Standard volume field types (scalar, vector, tensor, etc)
virtual bool fullyDescribed() const
Flag to identify whether model fully describes the parcel.
Inflation injection - creates new particles by splitting existing particles within in a set of genera...
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Top level data entry class for use in dictionaries. Provides a mechanism to specify a variable as a c...
Mesh consisting of general polyhedral cells.
Random & rndGen()
Return reference to the random object.
#define forAll(list, i)
Loop across all elements in list.
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
#define R(A, B, C, D, E, F, K, M)
const List< DynamicList< molecule * > > & cellOccupancy
const fvMesh & mesh() const
Return reference to the mesh.
AccessType combine(const UList< T > &lists, AccessOp aop=accessOp< T >())
Combines sub-lists into a single list.
virtual bool validInjection(const label parcelI)
Return flag to identify whether or not injection of parcelI is.
void deleteParticle(ParticleType &p)
Remove particle from cloud and delete.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Lookup type of boundary radiation properties.
Templated base class for dsmc cloud.
const List< DynamicList< ParcelType * > > & cellOccupancy() const
Return the cell occupancy addressing.
virtual ~InflationInjection()
Destructor.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
const scalarField & cellVolumes() const
A collection of cell labels.
Vector< scalar > vector
A scalar version of the templated Vector.
static autoPtr< distributionModel > New(const dictionary &dict, Random &rndGen)
Selector.
virtual scalar volumeToInject(const scalar time0, const scalar time1)
Volume of parcels to introduce relative to SOI.
static void gatherList(const List< commsStruct > &comms, List< T > &Values, const int tag, const label comm)
Gather data but keep individual values separate.
static int myProcNo(const label communicator=worldComm)
Number of this process (starting from masterNo() = 0)
An ordered pair of two objects of type <T> with first() and second() elements.
Type position(const Type &start, const Type &end)
Return a sample on the interval [start,end].
static bool & parRun() noexcept
Test if this a parallel run.
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
dimensionedScalar sqrt(const dimensionedScalar &ds)
virtual void setProperties(const label parcelI, const label nParcels, const scalar time, typename CloudType::parcelType &parcel)
Set the parcel properties.
bool insert(const Key &key)
Insert a new entry, not overwriting existing entries.
virtual void setPositionAndCell(const label parcelI, const label nParcels, const scalar time, vector &position, label &cellOwner, label &tetFacei, label &tetPti)
Set the injection position and owner cell, tetFace and tetPt.
scalar timeEnd() const
Return the end-of-injection time.
ParcelType parcelType
Type of parcel the cloud was instantiated for.
A 2-tuple for storing two objects of dissimilar types. The container is similar in purpose to std::pa...
#define WarningInFunction
Report a warning using Foam::Warning.
static label nProcs(const label communicator=worldComm)
Number of processes in parallel run, and 1 for serial run.
virtual void updateMesh()
Set injector locations when mesh is updated.