34 template<
class CloudType>
37 const word injectionMethod =
38 this->coeffDict().template lookupOrDefault<word>
44 if (injectionMethod ==
"point" || injectionMethod ==
word::null)
46 injectionMethod_ = imPoint;
50 else if (injectionMethod ==
"disc")
52 injectionMethod_ = imDisc;
55 this->coeffDict().template lookup<scalar>(
"dInner",
dimLength);
57 this->coeffDict().template lookup<scalar>(
"dOuter",
dimLength);
62 <<
"injectionMethod must be either 'point' or 'disc'"
68 template<
class CloudType>
72 this->coeffDict().template lookupOrDefault<word>
78 if (flowType ==
"constantVelocity" || flowType ==
word::null)
80 flowType_ = ftConstantVelocity;
87 this->owner().db().time().userUnits(),
93 else if (flowType ==
"pressureDrivenVelocity")
95 flowType_ = ftPressureDrivenVelocity;
102 this->owner().db().time().userUnits(),
108 else if (flowType ==
"flowRateAndDischarge")
110 flowType_ = ftFlowRateAndDischarge;
113 this->coeffDict().template lookup<scalar>(
"dInner",
dimLength);
115 this->coeffDict().template lookup<scalar>(
"dOuter",
dimLength);
122 this->owner().db().time().userUnits(),
131 <<
"flowType must be either 'constantVelocity', "
132 <<
"'pressureDrivenVelocity' or 'flowRateAndDischarge'"
140 template<
class CloudType>
149 injectionMethod_(imPoint),
150 flowType_(ftConstantVelocity),
156 this->owner().db().time().userUnits(),
166 this->owner().db().time().userUnits(),
173 injectorTetFace_(-1),
175 duration_(this->readDuration(
dict, owner)),
176 massFlowRate_(this->readMassFlowRate(
dict, owner, duration_)),
177 parcelsPerSecond_(this->readParcelsPerSecond(
dict, owner)),
183 this->owner().db().time().userUnits(),
193 this->owner().db().time().userUnits(),
203 this->coeffDict().subDict(
"sizeDistribution"),
205 owner.
rndGen().generator()
214 setInjectionMethod();
222 template<
class CloudType>
229 injectionMethod_(im.injectionMethod_),
230 flowType_(im.flowType_),
231 position_(im.position_, false),
232 direction_(im.direction_, false),
233 injectorCoordinates_(im.injectorCoordinates_),
234 injectorCell_(im.injectorCell_),
235 injectorTetFace_(im.injectorTetFace_),
236 injectorTetPt_(im.injectorTetPt_),
237 duration_(im.duration_),
238 massFlowRate_(im.massFlowRate_, false),
239 parcelsPerSecond_(im.parcelsPerSecond_, false),
240 thetaInner_(im.thetaInner_, false),
241 thetaOuter_(im.thetaOuter_, false),
242 sizeDistribution_(im.sizeDistribution_, false),
245 Umag_(im.Umag_, false),
247 Pinj_(im.Pinj_, false)
253 template<
class CloudType>
260 template<
class CloudType>
263 if (injectionMethod_ == imPoint && position_->constant())
266 this->findCellAtPosition
269 injectorCoordinates_,
278 template<
class CloudType>
281 return this->SOI_ + duration_;
285 template<
class CloudType>
292 if (t1 >= 0 && t0 < duration_)
294 return parcelsPerSecond_->integral(
max(t0, 0),
min(t1, duration_));
303 template<
class CloudType>
310 if (t1 >= 0 && t0 < duration_)
312 return massFlowRate_->integral(
max(t0, 0),
min(t1, duration_));
321 template<
class CloudType>
336 const scalar t = time - this->SOI_;
338 switch (injectionMethod_)
343 if (position_->constant())
346 celli = injectorCell_;
347 tetFacei = injectorTetFace_;
348 tetPti = injectorTetPt_;
352 this->findCellAtPosition
366 const scalar beta =
twoPi*this->globalScalar01(
rndGen);
367 const scalar frac = this->globalScalar01(
rndGen);
372 const scalar d =
sqrt((1 - frac)*
sqr(dInner_) + frac*
sqr(dOuter_));
374 this->findCellAtPosition
393 template<
class CloudType>
399 typename CloudType::parcelType::trackingData& td,
407 const scalar t = time - this->SOI_;
414 scalar theta = vGreat;
416 switch (injectionMethod_)
425 tanVec = t1*
cos(beta) + t2*
sin(beta);
429 (1 - frac)*
sqr(thetaInner_->value(t))
430 + frac*
sqr(thetaOuter_->value(t))
436 const scalar r =
mag(parcel.position(
mesh) - position_->value(t));
437 const scalar frac = (2*r - dInner_)/(dOuter_ - dInner_);
440 (1 - frac)*thetaInner_->value(t)
441 + frac*thetaOuter_->value(t);
461 case ftConstantVelocity:
463 parcel.U() = Umag_->value(t)*dirVec;
466 case ftPressureDrivenVelocity:
468 const scalar pAmbient = this->owner().pAmbient();
469 const scalar
rho = parcel.rho();
470 const scalar Umag =
::sqrt(2*(Pinj_->value(t) - pAmbient)/
rho);
471 parcel.U() = Umag*dirVec;
474 case ftFlowRateAndDischarge:
476 const scalar
A = 0.25*
pi*(
sqr(dOuter_) -
sqr(dInner_));
478 massFlowRate_->value(t)/(parcel.rho()*Cd_->value(t)*
A);
479 parcel.U() = Umag*dirVec;
489 parcel.d() = sizeDistribution_->sample();
493 template<
class CloudType>
This injector injects particles in a number of cones. The user specifies a position and a direction t...
virtual scalar nParcelsToInject(const scalar time0, const scalar time1)
Number of parcels to introduce relative to SOI.
virtual ~ConeInjection()
Destructor.
virtual void topoChange()
Set injector locations when mesh is updated.
ConeInjection(const dictionary &dict, CloudType &owner, const word &modelName)
Construct from dictionary.
virtual void setProperties(const label parcelI, const label nParcels, const scalar time, typename CloudType::parcelType::trackingData &td, typename CloudType::parcelType &parcel)
Set the parcel properties.
virtual void setPositionAndCell(const label parcelI, const label nParcels, const scalar time, barycentric &coordinates, label &celli, label &tetFacei, label &tetPti, label &facei)
Set the injection position and owner cell, tetFace and tetPt.
virtual bool fullyDescribed() const
Flag to identify whether model fully describes the parcel.
scalar timeEnd() const
Return the end-of-injection time.
virtual scalar massToInject(const scalar time0, const scalar time1)
Parcel mass to introduce relative to SOI.
Templated base class for dsmc cloud.
ParcelType parcelType
Type of parcel the cloud was instantiated for.
Run-time selectable general function of one variable.
static autoPtr< Function1< scalar > > New(const word &name, const Function1s::unitConversions &units, const dictionary &dict)
Select from dictionary.
Templated injection model class.
A list of keywords followed by any number of values (e.g. words and numbers) or sub-dictionaries.
const Type & value() const
Return const reference to value.
Base class for statistical distributions.
const polyMesh & mesh() const
Return reference to polyMesh.
Mesh consisting of general polyhedral cells.
scalar scalar01()
Return a scalar uniformly distributed between zero and one.
A class for handling words, derived from string.
static const word null
An empty word.
Foam::fvMesh mesh(Foam::IOobject(regionName, runTime.name(), runTime, Foam::IOobject::MUST_READ), false)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
const scalar twoPi(2 *pi)
static const coefficient A("A", dimPressure, 611.21)
barycentric coordinates(const polyMesh &mesh, const point &position, const label celli, const label facei, const label faceTrii, const scalar stepFraction)
Return the coordinates given the position and tet topology.
point position(const polyMesh &mesh, const barycentric &coordinates, const label celli, const label facei, const label faceTrii, const scalar stepFraction)
Return the position given the coordinates and tet topology.
errorManipArg< error, int > exit(error &err, const int errNo=1)
dimensionedScalar pos(const dimensionedScalar &ds)
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
const dimensionSet dimPressure
const dimensionSet dimless
dimensionedScalar sin(const dimensionedScalar &ds)
const dimensionSet dimLength
void mag(LagrangianPatchField< scalar > &f, const LagrangianPatchField< Type > &f1)
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
dimensionSet normalised(const dimensionSet &)
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
void sqr(LagrangianPatchField< typename outerProduct< Type, Type >::type > &f, const LagrangianPatchField< Type > &f1)
const dimensionSet dimVelocity
dimensionSet perpendicular(const dimensionSet &)
void sqrt(LagrangianPatchField< scalar > &f, const LagrangianPatchField< scalar > &f1)
tmp< DimensionedField< TypeR, GeoMesh, Field > > New(const tmp< DimensionedField< TypeR, GeoMesh, Field >> &tdf1, const word &name, const dimensionSet &dimensions)
dimensionedScalar cos(const dimensionedScalar &ds)
const unitConversion unitDegrees
randomGenerator rndGen(653213)