Go to the documentation of this file.
34 template<
class CloudType>
38 { vtFixedValue,
"fixedValue" },
39 { vtPatchValue,
"patchValue" },
40 { vtZeroGradient,
"zeroGradient" },
45 template<
class CloudType>
55 duration_(this->coeffDict().getScalar(
"duration")),
58 this->coeffDict().getScalar(
"parcelsPerSecond")
62 velocityTypeNames_.getOrDefault
71 velocityType_ == vtFixedValue
72 ? this->coeffDict().
template get<vector>(
"U0")
88 this->coeffDict().subDict(
"sizeDistribution"),
96 const Time& time = owner.db().time();
98 flowRateProfile_->userTimeToTime(time);
100 patchInjectionBase::updateMesh(owner.
mesh());
103 this->volumeTotal_ = flowRateProfile_->integrate(0.0, duration_);
107 template<
class CloudType>
115 duration_(im.duration_),
116 parcelsPerSecond_(im.parcelsPerSecond_),
117 velocityType_(im.velocityType_),
119 flowRateProfile_(im.flowRateProfile_.clone()),
120 sizeDistribution_(im.sizeDistribution_.clone()),
121 currentParceli_(im.currentParceli_),
122 currentFacei_(im.currentFacei_)
128 template<
class CloudType>
131 patchInjectionBase::updateMesh(this->owner().
mesh());
135 template<
class CloudType>
138 return this->SOI_ + duration_;
142 template<
class CloudType>
149 if ((time0 >= 0.0) && (time0 < duration_))
151 scalar nParcels = (time1 - time0)*parcelsPerSecond_;
154 label nParcelsToInject = floor(nParcels);
161 && (nParcels - scalar(nParcelsToInject) > rndPos)
167 return nParcelsToInject;
174 template<
class CloudType>
181 if ((time0 >= 0.0) && (time0 < duration_))
183 return flowRateProfile_->integrate(time0, time1);
190 template<
class CloudType>
194 const label nParcels,
202 currentParceli_ = parcelI;
204 currentFacei_ = patchInjectionBase::setPositionAndCell
206 this->owner().
mesh(),
216 template<
class CloudType>
220 const label nParcels,
226 switch (velocityType_)
235 if (parcelI != currentParceli_)
238 <<
"Synchronisation problem: "
239 <<
"attempting to set injected parcel " << parcelI
240 <<
" properties using cached parcel " << currentParceli_
241 <<
" properties" <<
endl;
244 const label patchFacei = currentFacei_;
249 <<
"Unable to set parcel velocity using patch value "
250 <<
"due to missing face index: patchFacei=" << patchFacei
255 const label patchi = this->patchId_;
256 parcel.U() =
U.boundaryField()[patchi][patchFacei];
261 const label celli = parcel.cell();
266 <<
"Unable to set parcel velocity using zeroGradient "
267 <<
"due to missing cell index"
272 parcel.U() =
U[celli];
278 <<
"Unhandled velocityType "
279 << velocityTypeNames_[velocityType_]
280 <<
". Available options are:"
281 << velocityTypeNames_.sortedToc()
287 parcel.d() = sizeDistribution_->sample();
291 template<
class CloudType>
298 template<
class CloudType>
virtual bool fullyDescribed() const
Flag to identify whether model fully describes the parcel.
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
A class for handling words, derived from Foam::string.
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)
virtual bool validInjection(const label parcelI)
Return flag to identify whether or not injection of parcelI is.
virtual scalar userTimeToTime(const scalar theta) const
Convert the user-time (e.g. CA deg) to real-time (s).
Templated injection model class.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Top level data entry class for use in dictionaries. Provides a mechanism to specify a variable as a c...
Random & rndGen()
Return reference to the random object.
PatchInjection(const dictionary &dict, CloudType &owner, const word &modelName)
Construct from dictionary.
const fvMesh & mesh() const
Return reference to the mesh.
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.
Templated base class for dsmc cloud.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
static const Enum< velocityType > velocityTypeNames_
Velocity type names.
word getWord(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Same as get< word >(const word&, keyType::option)
errorManip< error > abort(error &err)
virtual scalar timeEnd() const
Return the end-of-injection time.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
virtual void setProperties(const label parcelI, const label nParcels, const scalar time, typename CloudType::parcelType &parcel)
Set the parcel properties.
virtual scalar volumeToInject(const scalar time0, const scalar time1)
Volume of parcels to introduce relative to SOI.
ParcelType parcelType
Type of parcel the cloud was instantiated for.
#define WarningInFunction
Report a warning using Foam::Warning.
virtual void updateMesh()
Set injector locations when mesh is updated.
Type globalPosition(const Type &start, const Type &end)
Return a sample on the interval [start,end].