37 template<
class CloudType>
44 { injectionMethod::imPoint,
"point" },
45 { injectionMethod::imDisc,
"disc" },
48 template<
class CloudType>
55 { flowType::ftConstantVelocity,
"constantVelocity" },
56 { flowType::ftPressureDrivenVelocity,
"pressureDrivenVelocity" },
57 { flowType::ftFlowRateAndDischarge,
"flowRateAndDischarge" },
63 template<
class CloudType>
66 const auto&
mesh = this->owner().mesh();
71 Function1<vector>::New(
"position", this->coeffDict(), &
mesh)
74 positionVsTime_->userTimeToTime(this->owner().time());
76 if (positionVsTime_->constant())
78 position_ = positionVsTime_->value(0);
82 directionVsTime_.reset
84 Function1<vector>::New(
"direction", this->coeffDict(), &
mesh)
87 directionVsTime_->userTimeToTime(this->owner().time());
89 if (directionVsTime_->constant())
91 direction_ = directionVsTime_->value(0);
92 direction_.normalise();
94 Random&
rndGen = this->owner().rndGen();
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_;
114 template<
class CloudType>
119 case flowType::ftConstantVelocity:
121 this->coeffDict().readEntry(
"UMag", UMag_);
124 case flowType::ftPressureDrivenVelocity:
128 Function1<scalar>::New
132 &this->owner().
mesh()
135 Pinj_->userTimeToTime(this->owner().time());
138 case flowType::ftFlowRateAndDischarge:
142 Function1<scalar>::New
146 &this->owner().
mesh()
149 Cd_->userTimeToTime(this->owner().time());
155 <<
"Unhandled flow type " 156 << flowTypeNames[flowType_]
165 template<
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),
198 parcelsPerSecond_(this->coeffDict().getScalar(
"parcelsPerSecond")),
230 this->coeffDict().subDict(
"sizeDistribution"),
242 if (innerDiameter_ >= outerDiameter_)
245 <<
"Inner diameter must be less than the outer diameter:" <<
nl 246 <<
" innerDiameter: " << innerDiameter_ <<
nl 247 <<
" outerDiameter: " << outerDiameter_
253 duration_ = time.userTimeToTime(duration_);
254 flowRateProfile_->userTimeToTime(time);
255 thetaInner_->userTimeToTime(time);
256 thetaOuter_->userTimeToTime(time);
260 omegaPtr_->userTimeToTime(time);
263 setInjectionGeometry();
269 this->
volumeTotal_ = flowRateProfile_->integrate(0.0, duration_);
275 template<
class CloudType>
282 injectionMethod_(im.injectionMethod_),
283 flowType_(im.flowType_),
284 outerDiameter_(im.outerDiameter_),
285 innerDiameter_(im.innerDiameter_),
286 duration_(im.duration_),
287 positionVsTime_(im.positionVsTime_.clone()),
288 position_(im.position_),
289 injectorCell_(im.injectorCell_),
290 tetFacei_(im.tetFacei_),
292 directionVsTime_(im.directionVsTime_.clone()),
293 direction_(im.direction_),
294 omegaPtr_(im.omegaPtr_.clone()),
295 parcelsPerSecond_(im.parcelsPerSecond_),
296 flowRateProfile_(im.flowRateProfile_.clone()),
297 thetaInner_(im.thetaInner_.clone()),
298 thetaOuter_(im.thetaOuter_.clone()),
299 sizeDistribution_(im.sizeDistribution_.clone()),
300 tanVec1_(im.tanVec1_),
301 tanVec2_(im.tanVec2_),
305 Pinj_(im.Pinj_.clone())
311 template<
class CloudType>
315 if (positionVsTime_->constant())
317 position_ = positionVsTime_->value(0);
319 this->findCellAtPosition
330 template<
class CloudType>
333 return this->SOI_ + duration_;
337 template<
class CloudType>
344 if ((time0 >= 0.0) && (time0 < duration_))
346 return floor((time1 - time0)*parcelsPerSecond_);
353 template<
class CloudType>
360 if ((time0 >= 0.0) && (time0 < duration_))
362 return flowRateProfile_->integrate(time0, time1);
369 template<
class CloudType>
382 const scalar t = time - this->SOI_;
384 if (!directionVsTime_->constant())
386 direction_ = directionVsTime_->value(t);
387 direction_.normalise();
391 scalar magTangent = 0.0;
393 while(magTangent < SMALL)
397 tangent = v - (v & direction_)*direction_;
398 magTangent =
mag(tangent);
401 tanVec1_ = tangent/magTangent;
402 tanVec2_ = direction_^tanVec1_;
408 switch (injectionMethod_)
410 case injectionMethod::imPoint:
412 if (positionVsTime_->constant())
414 position = position_;
415 cellOwner = injectorCell_;
416 tetFacei = tetFacei_;
421 position = positionVsTime_->value(t);
423 this->findCellAtPosition
433 case injectionMethod::imDisc:
436 scalar dr = outerDiameter_ - innerDiameter_;
437 scalar r = 0.5*(innerDiameter_ + frac*dr);
439 position = positionVsTime_->value(t) + r*normal_;
441 this->findCellAtPosition
453 <<
"Unhandled injection method " 454 << injectionMethodNames[injectionMethod_]
461 template<
class CloudType>
473 scalar t = time - this->SOI_;
474 scalar ti = thetaInner_->value(t);
475 scalar to = thetaOuter_->value(t);
479 scalar dcorr =
cos(coneAngle);
482 vector dirVec = dcorr*direction_;
488 case flowType::ftConstantVelocity:
490 parcel.U() = UMag_*dirVec;
493 case flowType::ftPressureDrivenVelocity:
495 scalar pAmbient = this->owner().pAmbient();
496 scalar
rho = parcel.rho();
497 scalar UMag =
::sqrt(2.0*(Pinj_->value(t) - pAmbient)/
rho);
498 parcel.U() = UMag*dirVec;
501 case flowType::ftFlowRateAndDischarge:
505 scalar massFlowRate =
507 *flowRateProfile_->value(t)
508 /this->volumeTotal();
510 scalar Umag = massFlowRate/(parcel.rho()*Cd_->value(t)*(Ao - Ai));
511 parcel.U() = Umag*dirVec;
517 <<
"Unhandled injection method " 518 << flowTypeNames[flowType_]
525 const scalar omega = omegaPtr_->value(t);
527 const vector p0(parcel.position() - positionVsTime_->value(t));
528 const vector r(
p0 - (
p0 & direction_)*direction_);
529 const scalar rMag =
mag(r);
533 parcel.U() += omega*rMag*d;
537 parcel.d() = sizeDistribution_->sample();
541 template<
class CloudType>
548 template<
class CloudType>
Different types of constants.
List< ReturnType > get(const UPtrList< T > &list, const AccessOp &aop)
List of values generated by applying the access operation to each list item.
virtual void setProperties(const label parcelI, const label nParcels, const scalar time, typename CloudType::parcelType &parcel)
Set the parcel properties.
Top level data entry class for use in dictionaries. Provides a mechanism to specify a variable as a c...
errorManipArg< error, int > exit(error &err, const int errNo=1)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
error FatalError
Error stream (stdout output on all processes), with additional 'FOAM FATAL ERROR' header text and sta...
A list of keyword definitions, which are a keyword followed by a number of values (eg...
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Templated injection model class.
Unit conversion functions.
constexpr char nl
The newline '\n' character (0x0a)
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.
dimensionedScalar sqrt(const dimensionedScalar &ds)
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.
static const Enum< flowType > flowTypeNames
Vector< Cmpt > & normalise(const scalar tol=ROOTVSMALL)
Inplace normalise the vector by its magnitude.
quaternion normalised(const quaternion &q)
Return the normalised (unit) quaternion of the given quaternion.
virtual bool validInjection(const label parcelI)
Return flag to identify whether or not injection of parcelI is.
const CloudType & owner() const
Return const access to the owner cloud.
injectionMethod
Injection method enumeration.
virtual scalar volumeToInject(const scalar time0, const scalar time1)
Volume of parcels to introduce relative to SOI.
Type sample01()
Return a sample whose components lie in the range [0,1].
dimensionedScalar cos(const dimensionedScalar &ds)
virtual label parcelsToInject(const scalar time0, const scalar time1)
Number of parcels to introduce relative to SOI.
constexpr scalar twoPi(2 *M_PI)
A class for handling words, derived from Foam::string.
Random & rndGen()
Return reference to the random object.
const Time & time() const noexcept
Return time registry.
const objectRegistry & db() const noexcept
Return the local objectRegistry.
constexpr scalar pi(M_PI)
flowType
Flow type enumeration.
ParcelType parcelType
Type of parcel the cloud was instantiated for.
dimensionedScalar sin(const dimensionedScalar &ds)
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
scalar timeEnd() const
Return the end-of-injection time.
virtual void updateMesh()
Set injector locations when mesh is updated.
Type globalSample01()
Return a sample whose components lie in the range [0,1].
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
A library of runtime-selectable doubly-truncated probability distribution models. Returns random samp...
dimensionedScalar beta("beta", dimless/dimTemperature, laminarTransport)
ConeNozzleInjection(const dictionary &dict, CloudType &owner, const word &modelName)
Construct from dictionary.
scalar volumeTotal_
Total volume of particles introduced by this injector [m^3] Note: scaled to ensure massTotal is achie...
static const Enum< injectionMethod > injectionMethodNames
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
virtual bool fullyDescribed() const
Flag to identify whether model fully describes the parcel.
constexpr scalar degToRad(const scalar deg) noexcept
Conversion from degrees to radians.
Templated base class for dsmc cloud.
const volScalarField & p0
static constexpr const zero Zero
Global zero (0)