Go to the documentation of this file.
36 template<
class CloudType>
46 phiName_(this->coeffDict().template getOrDefault<word>(
"phi",
"phi")),
47 rhoName_(this->coeffDict().template getOrDefault<word>(
"rho",
"rho")),
48 duration_(this->coeffDict().getScalar(
"duration")),
60 this->coeffDict().getScalar(
"parcelConcentration")
66 this->coeffDict().subDict(
"sizeDistribution"),
72 const Time& time = owner.db().time();
74 concentration_->userTimeToTime(time);
76 patchInjectionBase::updateMesh(owner.
mesh());
80 this->volumeTotal_ = 0.0;
81 this->massTotal_ = 0.0;
85 template<
class CloudType>
93 phiName_(im.phiName_),
94 rhoName_(im.rhoName_),
95 duration_(im.duration_),
96 concentration_(im.concentration_.clone()),
97 parcelConcentration_(im.parcelConcentration_),
98 sizeDistribution_(im.sizeDistribution_.clone())
104 template<
class CloudType>
111 template<
class CloudType>
114 patchInjectionBase::updateMesh(this->owner().
mesh());
118 template<
class CloudType>
121 return this->SOI_ + duration_;
125 template<
class CloudType>
135 scalar flowRateIn = 0.0;
138 flowRateIn =
max(0.0, -
sum(phip));
146 flowRateIn =
max(0.0, -
sum(phip/rhop));
155 template<
class CloudType>
162 if ((time0 >= 0.0) && (time0 < duration_))
164 scalar dt = time1 - time0;
166 scalar
c = concentration_->
value(0.5*(time0 + time1));
168 scalar nParcels = parcelConcentration_*
c*flowRate()*dt;
172 label nParcelsToInject = floor(nParcels);
180 nParcels - scalar(nParcelsToInject)
188 return nParcelsToInject;
195 template<
class CloudType>
204 if ((time0 >= 0.0) && (time0 < duration_))
206 scalar
c = concentration_->
value(0.5*(time0 + time1));
208 volume =
c*(time1 - time0)*flowRate();
211 this->volumeTotal_ =
volume;
218 template<
class CloudType>
230 patchInjectionBase::setPositionAndCell
232 this->owner().
mesh(),
242 template<
class CloudType>
252 parcel.U() = this->owner().U()[parcel.cell()];
255 parcel.d() = sizeDistribution_->sample();
259 template<
class CloudType>
266 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].
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
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 scalar userTimeToTime(const scalar theta) const
Convert the user-time (e.g. CA deg) to real-time (s).
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.
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.
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].