Go to the documentation of this file.
37 template<
class CloudType>
47 phiName_(this->coeffDict().template getOrDefault<word>(
"phi",
"phi")),
48 rhoName_(this->coeffDict().template getOrDefault<word>(
"rho",
"rho")),
49 duration_(this->coeffDict().getScalar(
"duration")),
61 this->coeffDict().getScalar(
"parcelConcentration")
67 this->coeffDict().subDict(
"sizeDistribution"),
72 duration_ = owner.db().time().userTimeToTime(duration_);
74 patchInjectionBase::updateMesh(owner.
mesh());
78 this->volumeTotal_ = 0.0;
79 this->massTotal_ = 0.0;
83 template<
class CloudType>
91 phiName_(im.phiName_),
92 rhoName_(im.rhoName_),
93 duration_(im.duration_),
94 concentration_(im.concentration_),
95 parcelConcentration_(im.parcelConcentration_),
96 sizeDistribution_(im.sizeDistribution_.clone())
102 template<
class CloudType>
109 template<
class CloudType>
112 patchInjectionBase::updateMesh(this->owner().
mesh());
116 template<
class CloudType>
119 return this->SOI_ + duration_;
123 template<
class CloudType>
133 scalar flowRateIn = 0.0;
136 flowRateIn =
max(0.0, -
sum(phip));
144 flowRateIn =
max(0.0, -
sum(phip/rhop));
153 template<
class CloudType>
160 if ((time0 >= 0.0) && (time0 < duration_))
162 scalar dt = time1 - time0;
164 scalar
c = concentration_.
value(0.5*(time0 + time1));
166 scalar nParcels = parcelConcentration_*
c*flowRate()*dt;
170 label nParcelsToInject = floor(nParcels);
178 nParcels - scalar(nParcelsToInject)
186 return nParcelsToInject;
193 template<
class CloudType>
202 if ((time0 >= 0.0) && (time0 < duration_))
204 scalar
c = concentration_.
value(0.5*(time0 + time1));
206 volume =
c*(time1 - time0)*flowRate();
209 this->volumeTotal_ =
volume;
216 template<
class CloudType>
228 patchInjectionBase::setPositionAndCell
230 this->owner().
mesh(),
240 template<
class CloudType>
250 parcel.U() = this->owner().U()[parcel.cell()];
253 parcel.d() = sizeDistribution_->sample();
257 template<
class CloudType>
264 template<
class CloudType>
Patch injection, by using patch flow rate to determine concentration and velocity.
virtual scalar flowRate() const
Return the total volumetric flow rate across the patch [m3/s].
A class for handling words, derived from Foam::string.
const List< typename ParcelType::constantProperties > & constProps() const
Return all of the constant properties.
const dimensionSet dimVelocity
virtual void setProperties(const label parcelI, const label nParcels, const scalar time, typename CloudType::parcelType &parcel)
Set the parcel properties.
Templated injection model class.
virtual ~PatchFlowRateInjection()
Destructor.
const wordList volume
Standard volume field types (scalar, vector, tensor, etc)
const Type & value() const
Return const reference to value.
virtual void updateMesh()
Set injector locations when mesh is updated.
Mesh consisting of general polyhedral cells.
Random & rndGen()
Return reference to the random object.
scalar timeEnd() const
Return the end-of-injection time.
const dimensionSet dimArea(sqr(dimLength))
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 fvMesh & mesh() const
Return reference to the mesh.
virtual bool validInjection(const label parcelI)
Return flag to identify whether or not injection of parcelI is.
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.
virtual scalar volumeToInject(const scalar time0, const scalar time1)
Volume of parcels to introduce relative to SOI.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
reduce(hasMovingMesh, orOp< bool >())
word getWord(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Same as get< word >(const word&, keyType::option)
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
PatchFlowRateInjection(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.
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
const dimensionedScalar c
Speed of light in a vacuum.
ParcelType parcelType
Type of parcel the cloud was instantiated for.
virtual bool fullyDescribed() const
Flag to identify whether model fully describes the parcel.
Type globalPosition(const Type &start, const Type &end)
Return a sample on the interval [start,end].