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"),
88 flowRateProfile_->userTimeToTime(time);
89 growthRate_->userTimeToTime(time);
98 generationCells_ = generationCells.
toc();
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_);
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();
213 Random& rnd = this->owner().
rndGen();
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();
275 Pair<scalar>(dSeed_, dNew)
281 cellCentresUsed.insert(cI);
286 label cPI = rnd.position(label(0),
cellOccupancy[cI].size() - 1);
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();
338 Pair<vector>(
vector(
x, 0, -r) + pP, pU),
339 Pair<scalar>(a, dNew)
342 volumeAccumulator_ -= volNew;
344 dNew = sizeDistribution_().sample();
349 Pair<vector>(
vector(-d, a/2, -r) + pP, pU),
350 Pair<scalar>(a, dNew)
353 volumeAccumulator_ -= volNew;
355 dNew = sizeDistribution_().sample();
360 Pair<vector>(
vector(-d, -a/2, -r) + pP, pU),
361 Pair<scalar>(a, dNew)
364 volumeAccumulator_ -= volNew;
366 dNew = sizeDistribution_().sample();
371 Pair<vector>(
vector(0, 0,
R) + pP, pU),
372 Pair<scalar>(a, dNew)
375 volumeAccumulator_ -= volNew;
395 List<List<vectorPairScalarPair>> gatheredNewParticles
406 List<vectorPairScalarPair> combinedNewParticles
410 gatheredNewParticles,
411 accessOp<List<vectorPairScalarPair>>()
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>
virtual scalar volumeToInject(const scalar time0, const scalar time1)
Volume of parcels to introduce relative to SOI.
Top level data entry class for use in dictionaries. Provides a mechanism to specify a variable as a c...
virtual ~InflationInjection()
Destructor.
scalar massTotal_
Total mass to inject [kg].
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.
A list of keyword definitions, which are a keyword followed by a number of values (eg...
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
virtual bool fullyDescribed() const
Flag to identify whether model fully describes the parcel.
Templated injection model class.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
virtual void updateMesh()
Set injector locations when mesh is updated.
dimensionedScalar sqrt(const dimensionedScalar &ds)
Ostream & endl(Ostream &os)
Add newline and flush stream.
static void gatherList(const UList< commsStruct > &comms, UList< T > &values, const int tag, const label comm)
Gather data, but keep individual values separate. Uses the specified communication schedule...
static bool & parRun() noexcept
Test if this a parallel run.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tf1, const word &name, const dimensionSet &dimensions, const bool initCopy=false)
Global function forwards to reuseTmpDimensionedField::New.
static int myProcNo(const label communicator=worldComm)
Rank of this process in the communicator (starting from masterNo()). Can be negative if the process i...
Lookup type of boundary radiation properties.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Tuple2< Pair< vector >, Pair< scalar > > vectorPairScalarPair
static void broadcast(Type &value, const label comm=UPstream::worldComm)
Broadcast content (contiguous or non-contiguous) to all communicator ranks. Does nothing in non-paral...
const CloudType & owner() const
Return const access to the owner cloud.
#define forAll(list, i)
Loop across all elements in list.
InflationInjection(const dictionary &dict, CloudType &owner, const word &modelName)
Construct from dictionary.
HashSet< label, Hash< label > > labelHashSet
A HashSet of labels, uses label hasher.
bool readEntry(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX, IOobjectOption::readOption readOpt=IOobjectOption::MUST_READ) const
Find entry and assign to T val. FatalIOError if it is found and the number of tokens is incorrect...
AccessType combine(const UList< T > &lists, AccessOp aop=accessOp< T >())
Combines sub-lists into a single list.
static label nProcs(const label communicator=worldComm)
Number of ranks in parallel run (for given communicator). It is 1 for serial run. ...
virtual scalar userTimeToTime(const scalar theta) const
Convert the user-time (e.g. CA deg) to real-time (s).
A class for handling words, derived from Foam::string.
Random & rndGen()
Return reference to the random object.
const Time & time() const noexcept
Return time registry.
const objectRegistry & db() const noexcept
Return the local objectRegistry.
const dictionary & coeffDict() const
Return const access to the coefficients dictionary.
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
const vectorField & cellCentres() const
virtual bool validInjection(const label parcelI)
Return flag to identify whether or not injection of parcelI is.
ParcelType parcelType
Type of parcel the cloud was instantiated for.
const wordList volume
Standard volume field types (scalar, vector, tensor, etc)
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
Inflation injection - creates new particles by splitting existing particles within in a set of genera...
const fvMesh & mesh() const
Return reference to the mesh.
const List< DynamicList< molecule * > > & cellOccupancy
vector point
Point is a vector.
#define R(A, B, C, D, E, F, K, M)
#define WarningInFunction
Report a warning using Foam::Warning.
A collection of cell labels.
A library of runtime-selectable doubly-truncated probability distribution models. Returns random samp...
static bool master(const label communicator=worldComm)
True if process corresponds to the master rank in the communicator.
const List< DynamicList< ParcelType * > > & cellOccupancy() const
Return the cell occupancy addressing.
scalar timeEnd() const
Return the end-of-injection time.
Mesh consisting of general polyhedral cells.
scalar volumeTotal_
Total volume of particles introduced by this injector [m^3] Note: scaled to ensure massTotal is achie...
List< Key > toc() const
The table of contents (the keys) in unsorted order.
void deleteParticle(ParticleType &p)
Remove particle from cloud and delete.
virtual label parcelsToInject(const scalar time0, const scalar time1)
Number of parcels to introduce relative to SOI.
Templated base class for dsmc cloud.
void reduce(T &value, const BinaryOp &bop, const int tag=UPstream::msgType(), const label comm=UPstream::worldComm)
Reduce inplace (cf. MPI Allreduce) using linear/tree communication schedule.
virtual void setProperties(const label parcelI, const label nParcels, const scalar time, typename CloudType::parcelType &parcel)
Set the parcel properties.
const scalarField & cellVolumes() const
static constexpr const zero Zero
Global zero (0)