37template<
class CloudType>
44 { injectionMethod::imPoint,
"point" },
45 { injectionMethod::imDisc,
"disc" },
48template<
class CloudType>
55 { flowType::ftConstantVelocity,
"constantVelocity" },
56 { flowType::ftPressureDrivenVelocity,
"pressureDrivenVelocity" },
57 { flowType::ftFlowRateAndDischarge,
"flowRateAndDischarge" },
63template<
class CloudType>
66 const auto&
mesh = this->owner().mesh();
74 positionVsTime_->userTimeToTime(this->owner().time());
76 if (positionVsTime_->constant())
78 position_ = positionVsTime_->value(0);
82 directionVsTime_.reset
87 directionVsTime_->userTimeToTime(this->owner().time());
89 if (directionVsTime_->constant())
91 direction_ = directionVsTime_->value(0);
92 direction_.normalise();
98 scalar magTangent = 0.0;
100 while(magTangent < SMALL)
104 tangent = v - (v & direction_)*direction_;
105 magTangent =
mag(tangent);
108 tanVec1_ = tangent/magTangent;
109 tanVec2_ = direction_^tanVec1_;
114template<
class CloudType>
119 case flowType::ftConstantVelocity:
121 this->coeffDict().readEntry(
"UMag", UMag_);
124 case flowType::ftPressureDrivenVelocity:
132 &this->owner().
mesh()
135 Pinj_->userTimeToTime(this->owner().time());
138 case flowType::ftFlowRateAndDischarge:
146 &this->owner().
mesh()
149 Cd_->userTimeToTime(this->owner().time());
155 <<
"Unhandled flow type "
156 << flowTypeNames[flowType_]
165template<
class CloudType>
170 const word& modelName
176 injectionMethodNames.get(
"injectionMethod", this->coeffDict())
178 flowType_(flowTypeNames.get(
"flowType", this->coeffDict())),
179 outerDiameter_(this->coeffDict().getScalar(
"outerDiameter")),
180 innerDiameter_(this->coeffDict().getScalar(
"innerDiameter")),
181 duration_(this->coeffDict().getScalar(
"duration")),
182 positionVsTime_(nullptr),
187 directionVsTime_(nullptr),
189 parcelsPerSecond_(this->coeffDict().getScalar(
"parcelsPerSecond")),
221 this->coeffDict().subDict(
"sizeDistribution"),
233 if (innerDiameter_ >= outerDiameter_)
236 <<
"Inner diameter must be less than the outer diameter:" <<
nl
237 <<
" innerDiameter: " << innerDiameter_ <<
nl
238 <<
" outerDiameter: " << outerDiameter_
245 flowRateProfile_->userTimeToTime(time);
246 thetaInner_->userTimeToTime(time);
247 thetaOuter_->userTimeToTime(time);
249 setInjectionGeometry();
255 this->
volumeTotal_ = flowRateProfile_->integrate(0.0, duration_);
261template<
class CloudType>
268 injectionMethod_(im.injectionMethod_),
269 flowType_(im.flowType_),
270 outerDiameter_(im.outerDiameter_),
271 innerDiameter_(im.innerDiameter_),
272 duration_(im.duration_),
273 positionVsTime_(im.positionVsTime_.clone()),
274 position_(im.position_),
275 injectorCell_(im.injectorCell_),
276 tetFacei_(im.tetFacei_),
278 directionVsTime_(im.directionVsTime_.clone()),
279 direction_(im.direction_),
280 parcelsPerSecond_(im.parcelsPerSecond_),
281 flowRateProfile_(im.flowRateProfile_.clone()),
282 thetaInner_(im.thetaInner_.clone()),
283 thetaOuter_(im.thetaOuter_.clone()),
284 sizeDistribution_(im.sizeDistribution_.clone()),
285 tanVec1_(im.tanVec1_),
286 tanVec2_(im.tanVec2_),
290 Pinj_(im.Pinj_.clone())
296template<
class CloudType>
300 if (positionVsTime_->constant())
302 position_ = positionVsTime_->value(0);
304 this->findCellAtPosition
315template<
class CloudType>
318 return this->SOI_ + duration_;
322template<
class CloudType>
329 if ((time0 >= 0.0) && (time0 < duration_))
331 return floor((time1 - time0)*parcelsPerSecond_);
338template<
class CloudType>
345 if ((time0 >= 0.0) && (time0 < duration_))
347 return flowRateProfile_->integrate(time0, time1);
354template<
class CloudType>
367 const scalar t = time - this->SOI_;
369 if (!directionVsTime_->constant())
371 direction_ = directionVsTime_->value(t);
372 direction_.normalise();
376 scalar magTangent = 0.0;
378 while(magTangent < SMALL)
382 tangent = v - (v & direction_)*direction_;
383 magTangent =
mag(tangent);
386 tanVec1_ = tangent/magTangent;
387 tanVec2_ = direction_^tanVec1_;
393 switch (injectionMethod_)
395 case injectionMethod::imPoint:
397 if (positionVsTime_->constant())
399 position = position_;
400 cellOwner = injectorCell_;
401 tetFacei = tetFacei_;
406 position = positionVsTime_->value(t);
408 this->findCellAtPosition
418 case injectionMethod::imDisc:
421 scalar dr = outerDiameter_ - innerDiameter_;
422 scalar r = 0.5*(innerDiameter_ + frac*dr);
424 position = positionVsTime_->value(t) + r*normal_;
426 this->findCellAtPosition
438 <<
"Unhandled injection method "
439 << injectionMethodNames[injectionMethod_]
446template<
class CloudType>
458 scalar t = time - this->SOI_;
459 scalar ti = thetaInner_->value(t);
460 scalar to = thetaOuter_->value(t);
464 scalar dcorr =
cos(coneAngle);
467 vector dirVec = dcorr*direction_;
473 case flowType::ftConstantVelocity:
475 parcel.U() = UMag_*dirVec;
478 case flowType::ftPressureDrivenVelocity:
480 scalar pAmbient = this->owner().pAmbient();
481 scalar
rho = parcel.rho();
482 scalar UMag =
::sqrt(2.0*(Pinj_->value(t) - pAmbient)/
rho);
483 parcel.U() = UMag*dirVec;
486 case flowType::ftFlowRateAndDischarge:
490 scalar massFlowRate =
492 *flowRateProfile_->value(t)
493 /this->volumeTotal();
495 scalar Umag = massFlowRate/(parcel.rho()*Cd_->value(t)*(Ao - Ai));
496 parcel.U() = Umag*dirVec;
502 <<
"Unhandled injection method "
503 << flowTypeNames[flowType_]
509 parcel.d() = sizeDistribution_->sample();
513template<
class CloudType>
520template<
class CloudType>
const CloudType & owner() const
Return const access to the owner cloud.
virtual scalar volumeToInject(const scalar time0, const scalar time1)
Volume of parcels to introduce relative to SOI.
virtual label parcelsToInject(const scalar time0, const scalar time1)
Number of parcels to introduce relative to SOI.
flowType
Flow type enumeration.
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.
virtual bool validInjection(const label parcelI)
Return flag to identify whether or not injection of parcelI is.
static const Enum< injectionMethod > injectionMethodNames
injectionMethod
Injection method enumeration.
virtual void setProperties(const label parcelI, const label nParcels, const scalar time, typename CloudType::parcelType &parcel)
Set the parcel properties.
static const Enum< flowType > flowTypeNames
virtual void updateMesh()
Set injector locations when mesh is updated.
virtual bool fullyDescribed() const
Flag to identify whether model fully describes the parcel.
scalar timeEnd() const
Return the end-of-injection time.
Templated base class for dsmc cloud.
ParcelType parcelType
Type of parcel the cloud was instantiated for.
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
Top level data entry class for use in dictionaries. Provides a mechanism to specify a variable as a c...
const objectRegistry & db() const noexcept
Return the local objectRegistry.
Templated injection model class.
Type sample01()
Return a sample whose components lie in the range [0,1].
Type globalSample01()
Return a sample whose components lie in the range [0,1].
virtual scalar userTimeToTime(const scalar theta) const
Convert the user-time (e.g. CA deg) to real-time (s).
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Vector< Cmpt > & normalise(const scalar tol=ROOTVSMALL)
Inplace normalise the vector by its magnitude.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
A library of runtime-selectable doubly-truncated probability distribution models. Returns random samp...
const Time & time() const noexcept
Return time registry.
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
A class for handling words, derived from Foam::string.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
constexpr scalar pi(M_PI)
constexpr scalar twoPi(2 *M_PI)
Different types of constants.
dimensionedScalar sin(const dimensionedScalar &ds)
constexpr scalar degToRad() noexcept
Multiplication factor for degrees to radians conversion.
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
static constexpr const zero Zero
Global zero (0)
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh > > &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
errorManipArg< error, int > exit(error &err, const int errNo=1)
dimensionedScalar cos(const dimensionedScalar &ds)
constexpr char nl
The newline '\n' character (0x0a)
dimensionedScalar beta("beta", dimless/dimTemperature, laminarTransport)
Unit conversion functions.