Go to the documentation of this file.
37 template<
class CloudType>
52 tag[particlei] =
p.tag();
53 position[particlei] =
p.position();
55 soi[particlei] =
p.soi();
68 ListListOps::combine<List<label>>
78 ListListOps::combine<List<point>>
88 ListListOps::combine<List<vector>>
98 ListListOps::combine<List<scalar>>
108 ListListOps::combine<List<scalar>>
117 maxTag =
max(maxTag, tag[particlei]);
120 label nInjectors = maxTag + 1;
130 const label tagi = tag[i];
131 const scalar t = soi[i];
132 injStartTime[tagi] =
min(t, injStartTime[tagi]);
133 injEndTime[tagi] =
max(t, injEndTime[tagi]);
134 injPosition[tagi].
append(position[i]);
136 injDiameter[tagi].
append(d[i]);
141 scalar sumVolume = 0;
142 startTime_.setSize(nInjectors, 0);
143 endTime_.setSize(nInjectors, 0);
144 sizeDistribution_.setSize(nInjectors);
145 position_.setSize(nInjectors);
146 U_.setSize(nInjectors);
147 volumeFlowRate_.setSize(nInjectors, 0);
149 scalar minTime = GREAT;
157 const label nParticle = diameters.size();
158 const scalar dTime = injEndTime[i] - injStartTime[i];
160 if ((nParticle > 1) && (dTime > ROOTVSMALL))
162 minTime =
min(minTime, injStartTime[i]);
164 startTime_[injectori] = injStartTime[i];
165 endTime_[injectori] = injEndTime[i];
168 position_[injectori].
setSize(resampleSize_);
169 U_[injectori].setSize(resampleSize_);
173 for (label samplei = 0; samplei < resampleSize_; ++samplei)
176 positioni[samplei] = injPosition[i][posi] + positionOffset_;
177 Ui[samplei] = injU[i][posi];
182 forAll(diameters, particlei)
184 sumPow3 +=
pow3(diameters[particlei]);
189 volumeFlowRate_[injectori] =
volume/dTime;
192 sizeDistribution_.
set
208 startTime_.setSize(injectori);
209 endTime_.setSize(injectori);
210 position_.setSize(injectori);
211 U_.setSize(injectori);
212 volumeFlowRate_.setSize(injectori);
213 sizeDistribution_.setSize(injectori);
216 forAll(startTime_, injectori)
218 startTime_[injectori] -= minTime;
219 endTime_[injectori] -= minTime;
224 this->volumeTotal_ = sumVolume;
227 Info<<
" Read " << position_.size() <<
" injectors with "
228 << tag.size() <<
" total particles" <<
endl;
234 template<
class CloudType>
240 const word& modelName
244 cloudName_(this->coeffDict().
lookup(
"cloud")),
245 startTime_(this->
template getModelProperty<scalarList>(
"startTime")),
246 endTime_(this->
template getModelProperty<scalarList>(
"endTime")),
248 positionOffset_(this->coeffDict().
lookup(
"positionOffset")),
251 this->
template getModelProperty<scalarList>(
"volumeFlowRate")
254 binWidth_(this->coeffDict().getScalar(
"binWidth")),
258 ceil(this->coeffDict().getScalar(
"parcelsPerInjector"))
262 this->coeffDict().getOrDefault(
"resampleSize", label(100))
264 applyDistributionMassTotal_
266 this->coeffDict().getBool(
"applyDistributionMassTotal")
270 this->coeffDict().getOrDefault(
"ignoreOutOfBounds",
false)
272 nParcelsInjected_(this->parcelsAddedTotal()),
273 nParcelsInjected0_(0),
274 currentInjectori_(0),
277 if (startTime_.size())
280 sizeDistribution_.setSize(startTime_.size());
281 forAll(sizeDistribution_, i)
287 sizeDistribution_.set
301 if (applyDistributionMassTotal_)
303 this->massTotal_ = this->volumeTotal_*this->owner().
constProps().rho0();
304 Info<<
" Set mass to inject from distribution: "
305 << this->massTotal_ <<
endl;
310 template<
class CloudType>
333 currentInjectori_(0),
336 forAll(sizeDistribution_, injectori)
338 if (sizeDistribution_.set(injectori))
340 sizeDistribution_.set
355 template<
class CloudType>
363 template<
class CloudType>
368 template<
class CloudType>
372 return max(endTime_);
376 template<
class CloudType>
386 nParcelsInjected0_ = 0;
388 if (startTime_.empty() || this->volumeTotal_ < ROOTVSMALL)
393 scalar targetVolume = 0;
394 forAll(startTime_, injectori)
396 if (time1 > startTime_[injectori])
398 scalar totalDuration =
399 min(time1, endTime_[injectori]) - startTime_[injectori];
401 targetVolume += volumeFlowRate_[injectori]*totalDuration;
405 const label targetParcels =
408 scalar(startTime_.size()*parcelsPerInjector_)
409 *targetVolume/this->volumeTotal_
412 const label nParcels = targetParcels - nParcelsInjected_;
418 template<
class CloudType>
427 forAll(startTime_, injectori)
429 if ((time1 > startTime_[injectori]) && (time1 <= endTime_[injectori]))
431 scalar duration =
min(time1, endTime_[injectori]) - time0;
432 volume += volumeFlowRate_[injectori]*duration;
440 template<
class CloudType>
444 const label nParcels,
453 currentInjectori_ = rnd.
globalPosition<label>(0, position_.size() - 1);
454 currentSamplei_ = rnd.
globalPosition<label>(0, resampleSize_ - 1);
456 position = position_[currentInjectori_][currentSamplei_];
459 this->findCellAtPosition
469 template<
class CloudType>
479 parcel.U() = U_[currentInjectori_][currentSamplei_];
482 parcel.d() = sizeDistribution_[currentInjectori_].sample();
486 nParcelsInjected0_++;
490 template<
class CloudType>
498 template<
class CloudType>
508 template<
class CloudType>
513 if (this->writeTime())
515 this->setModelProperty(
"startTime", startTime_);
516 this->setModelProperty(
"endTime", endTime_);
517 this->setModelProperty(
"position", position_);
518 this->setModelProperty(
"volumeFlowRate", volumeFlowRate_);
519 this->setModelProperty(
"U", U_);
520 forAll(sizeDistribution_, i)
List< vectorList > position_
List of position per injector.
virtual void info(Ostream &os)
Write injection info to stream.
virtual bool validInjection(const label parcelI)
Return flag to identify whether or not injection of parcelI is.
A class for handling words, derived from Foam::string.
scalarList endTime_
List of end time per injector.
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
const List< typename ParcelType::constantProperties > & constProps() const
Return all of the constant properties.
virtual void updateMesh()
Set injector locations when mesh is updated.
scalar parcelsPerInjector_
Target number of parcels to inject per injector.
Different types of constants.
label resampleSize_
Resample size.
static void scatterList(const List< commsStruct > &comms, List< T > &Values, const int tag, const label comm)
Scatter data. Reverse of gatherList.
void info(Ostream &os)
Write injection info to stream.
bool applyDistributionMassTotal_
Flag to apply mass calculated from distribution instead of.
static label nProcs(const label communicator=0)
Number of processes in parallel run.
static bool & parRun()
Is this a parallel run?
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.
const word dictName("blockMeshDict")
Templated injection model class.
scalarList startTime_
List of start time per injector.
InjectedParticleDistributionInjection(const dictionary &dict, CloudType &owner, const word &modelName)
Construct from dictionary.
void append(const T &val)
Append an element at the end of the list.
Ostream & endl(Ostream &os)
Add newline and flush stream.
const wordList volume
Standard volume field types (scalar, vector, tensor, etc)
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
general distribution model
Random & rndGen()
Return reference to the random object.
#define forAll(list, i)
Loop across all elements in list.
virtual void setProperties(const label parcelI, const label nParcels, const scalar time, typename CloudType::parcelType &parcel)
Set the parcel properties.
virtual label parcelsToInject(const scalar time0, const scalar time1)
Number of parcels to introduce relative to SOI.
dimensionedScalar pow3(const dimensionedScalar &ds)
messageStream Info
Information stream (uses stdout - output is on the master only)
PtrList< distributionModels::general > sizeDistribution_
List of size distribution model per injector.
word name(const complex &c)
Return string representation of complex.
const cellModel * lookup(const word &modelName)
Deprecated(2017-11) equivalent to cellModel::ptr static method.
void transfer(List< T > &list)
void initialise()
Initialise injectors.
label nParcelsInjected_
Running total of number of parcels injected.
virtual ~InjectedParticleDistributionInjection()
Destructor.
virtual dictionary writeDict(const word &dictName) const
Write data in dictionary format.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Templated base class for dsmc cloud.
Primarily stores particle properties so that it can be injected at a later time. Note that this store...
label nParcelsInjected0_
Number of parcels injected in last step (local proc only)
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
List< vectorList > U_
List of parcel velocity per injector.
static int myProcNo(const label communicator=0)
Number of this process (starting from masterNo() = 0)
scalar binWidth_
Bin width when generating particle distributions.
Interrogates an injectedParticleCloud to convert the raw particle data into a set of 'binned' injecto...
static void gatherList(const List< commsStruct > &comms, List< T > &Values, const int tag, const label comm)
Gather data but keep individual values separate.
const word cloudName_
Name of cloud used to seed the new particles.
constexpr scalar pi(M_PI)
Switch ignoreOutOfBounds_
Flag to suppress errors if particle injection site is out-of-bounds.
std::enable_if< std::is_same< bool, TypeT >::value, bool >::type set(const label i, bool val=true)
A bitSet::set() method for a list of bool.
ParcelType parcelType
Type of parcel the cloud was instantiated for.
virtual bool fullyDescribed() const
Flag to identify whether model fully describes the parcel.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
void setSize(const label newSize)
Alias for resize(const label)
scalar timeEnd() const
Return the end-of-injection time.
vector positionOffset_
Position offset to apply to input positions.
scalarList volumeFlowRate_
List of volume flow rate per injector [m3/s].
Type globalPosition(const Type &start, const Type &end)
Return a sample on the interval [start,end].
virtual scalar volumeToInject(const scalar time0, const scalar time1)
Volume of parcels to introduce relative to SOI.