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"),
74 concentration_->userTimeToTime(time);
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>
118 template<
class CloudType>
121 return this->SOI_ + duration_;
125 template<
class CloudType>
134 scalar flowRateIn = 0.0;
137 flowRateIn =
max(0.0, -
sum(phip));
144 flowRateIn =
max(0.0, -
sum(phip/rhop));
147 reduce(flowRateIn, sumOp<scalar>());
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;
168 Random& rnd = this->owner().rndGen();
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;
210 this->massTotal_ =
volume*this->owner().constProps().rho0();
216 template<
class CloudType>
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>
Top level data entry class for use in dictionaries. Provides a mechanism to specify a variable as a c...
const Type & value() const noexcept
Return const reference to value.
scalar massTotal_
Total mass to inject [kg].
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &f1)
const Type & lookupObject(const word &name, const bool recursive=false) const
Lookup and return const reference to the object of the given Type. Fatal if not found or the wrong ty...
A list of keyword definitions, which are a keyword followed by a number of values (eg...
Templated injection model class.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
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.
virtual void updateMesh(const polyMesh &mesh)
Update patch geometry and derived info for injection locations.
virtual void updateMesh()
Set injector locations when mesh is updated.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
const CloudType & owner() const
Return const access to the owner cloud.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
const dimensionSet dimVolume(pow3(dimLength))
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.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
virtual bool validInjection(const label parcelI)
Return flag to identify whether or not injection of parcelI is.
const Time & time() const noexcept
Return time registry.
const objectRegistry & db() const noexcept
Return the local objectRegistry.
virtual bool fullyDescribed() const
Flag to identify whether model fully describes the parcel.
ParcelType parcelType
Type of parcel the cloud was instantiated for.
const wordList volume
Standard volume field types (scalar, vector, tensor, etc)
Type globalPosition(const Type &start, const Type &end)
Return a sample on the interval [start,end].
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 Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
Patch injection, by using patch flow rate to determine concentration and velocity.
const fvMesh & mesh() const
Return reference to the mesh.
virtual label parcelsToInject(const scalar time0, const scalar time1)
Number of parcels to introduce relative to SOI.
scalar timeEnd() const
Return the end-of-injection time.
virtual void setProperties(const label parcelI, const label nParcels, const scalar time, typename CloudType::parcelType &parcel)
Set the parcel properties.
label setPositionAndCell(const fvMesh &mesh, const scalar fraction01, Random &rnd, vector &position, label &cellOwner, label &tetFacei, label &tetPti)
Set the injection position and owner cell, tetFace and tetPt.
PatchFlowRateInjection(const dictionary &dict, CloudType &owner, const word &modelName)
Construct from dictionary.
const dimensionedScalar c
Speed of light in a vacuum.
A library of runtime-selectable doubly-truncated probability distribution models. Returns random samp...
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
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...
virtual scalar volumeToInject(const scalar time0, const scalar time1)
Volume of parcels to introduce relative to SOI.
Templated base class for dsmc cloud.
virtual scalar flowRate() const
Return the total volumetric flow rate across the patch [m3/s].
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 ~PatchFlowRateInjection()
Destructor.