37template<
class CloudType>
42 scalar& newVolumeFraction
47 newVolumeFraction = 0.0;
48 bool validInjection =
false;
54 return validInjection;
58 scalar t0 = timeStep0_ - SOI_;
59 scalar t1 = time - SOI_;
62 newParcels = this->parcelsToInject(t0, t1);
66 this->volumeToInject(t0, t1)
67 /(volumeTotal_ + ROOTVSMALL);
69 if (newVolumeFraction > 0)
74 validInjection =
true;
80 validInjection =
false;
86 validInjection =
false;
89 return validInjection;
93template<
class CloudType>
119 proci = Pstream::myProcNo();
126 if (proci != Pstream::myProcNo())
141 position += SMALL*(cellCentres[celli] - position);
153 proci = Pstream::myProcNo();
159 if (proci != Pstream::myProcNo())
172 <<
"Cannot find parcel injection cell. "
173 <<
"Parcel position = " <<
p0 <<
nl
186template<
class CloudType>
190 const scalar volumeFraction,
191 const scalar diameter,
196 switch (parcelBasis_)
200 scalar volumep =
pi/6.0*
pow3(diameter);
201 scalar volumeTot = massTotal_/
rho;
203 nP = (volumeFraction*volumeTot + delayedVolume_)/(parcels*volumep);
208 nP = massTotal_/(
rho*volumeTotal_);
213 nP = nParticleFixed_;
220 <<
"Unknown parcelBasis type" <<
nl
229template<
class CloudType>
232 const label parcelsAdded,
233 const scalar massAdded
238 if (allParcelsAdded > 0)
241 <<
"Cloud: " << this->owner().
name()
242 <<
" injector: " << this->modelName() <<
nl
243 <<
" Added " << allParcelsAdded <<
" new parcels" <<
nl <<
endl;
247 parcelsAddedTotal_ += allParcelsAdded;
262template<
class CloudType>
267 volumeTotal_(this->template getModelProperty<scalar>(
"volumeTotal")),
269 massFlowRate_(nullptr),
270 massInjected_(this->template getModelProperty<scalar>(
"massInjected")),
271 nInjections_(this->template getModelProperty<label>(
"nInjections")),
274 this->template getModelProperty<scalar>(
"parcelsAddedTotal")
276 parcelBasis_(pbNumber),
277 nParticleFixed_(0.0),
279 timeStep0_(this->template getModelProperty<scalar>(
"timeStep0")),
280 minParticlesPerParcel_(1),
283 ignoreOutOfBounds_(false)
287template<
class CloudType>
292 const word& modelName,
293 const word& modelType
298 volumeTotal_(this->template getModelProperty<scalar>(
"volumeTotal")),
300 massFlowRate_(nullptr),
301 massInjected_(this->template getModelProperty<scalar>(
"massInjected")),
302 nInjections_(this->template getModelProperty<scalar>(
"nInjections")),
305 this->template getModelProperty<scalar>(
"parcelsAddedTotal")
307 parcelBasis_(pbNumber),
308 nParticleFixed_(0.0),
309 time0_(owner.db().time().value()),
310 timeStep0_(this->template getModelProperty<scalar>(
"timeStep0")),
311 minParticlesPerParcel_
313 this->coeffDict().getOrDefault(
"minParticlesPerParcel", scalar(1))
316 injectorID_(this->coeffDict().getOrDefault(
"injectorID", -1)),
319 this->coeffDict().getOrDefault(
"ignoreOutOfBounds", false)
328 if (injectorID_ != -1)
330 Info<<
" injector ID: " << injectorID_ <<
endl;
333 if (owner.solution().active())
335 if (owner.solution().transient())
337 this->coeffDict().readEntry(
"massTotal", massTotal_);
338 this->coeffDict().readEntry(
"SOI", SOI_);
351 massFlowRate_->userTimeToTime(owner.
db().
time());
352 massTotal_ = massFlowRate_->value(owner.
db().
time().
value());
353 this->coeffDict().readIfPresent(
"SOI", SOI_);
359 const word parcelBasisType(this->coeffDict().getWord(
"parcelBasisType"));
361 if (parcelBasisType ==
"mass")
363 parcelBasis_ = pbMass;
365 else if (parcelBasisType ==
"number")
367 parcelBasis_ = pbNumber;
369 else if (parcelBasisType ==
"fixed")
371 parcelBasis_ = pbFixed;
372 this->coeffDict().readEntry(
"nParticle", nParticleFixed_);
374 Info<<
" Choosing nParticle to be a fixed value, massTotal "
375 <<
"variable now does not determine anything."
381 <<
"parcelBasisType must be either 'number', 'mass' or 'fixed'"
382 <<
nl <<
exit(FatalError);
387template<
class CloudType>
395 volumeTotal_(im.volumeTotal_),
396 massTotal_(im.massTotal_),
397 massFlowRate_(im.massFlowRate_.clone()),
398 massInjected_(im.massInjected_),
399 nInjections_(im.nInjections_),
400 parcelsAddedTotal_(im.parcelsAddedTotal_),
401 parcelBasis_(im.parcelBasis_),
402 nParticleFixed_(im.nParticleFixed_),
404 timeStep0_(im.timeStep0_),
405 minParticlesPerParcel_(im.minParticlesPerParcel_),
406 delayedVolume_(im.delayedVolume_),
407 injectorID_(im.injectorID_),
408 ignoreOutOfBounds_(im.ignoreOutOfBounds_)
414template<
class CloudType>
419template<
class CloudType>
423 if (this->owner().
solution().transient())
425 nTotal = parcelsToInject(0.0, timeEnd() - timeStart());
429 nTotal = parcelsToInject(0.0, 1.0);
432 return massTotal_/nTotal;
436template<
class CloudType>
437template<
class TrackCloudType>
440 TrackCloudType&
cloud,
449 const scalar time = this->owner().
db().
time().
value();
452 label parcelsAdded = 0;
453 scalar massAdded = 0.0;
454 label newParcels = 0;
455 scalar newVolumeFraction = 0.0;
456 scalar delayedVolume = 0;
458 if (prepareForNextTimeStep(time, newParcels, newVolumeFraction))
460 const scalar trackTime = this->owner().solution().trackTime();
464 const scalar deltaT =
465 max(0.0,
min(trackTime,
min(time - SOI_, timeEnd() - time0_)));
468 const scalar padTime =
max(0.0, SOI_ - time0_);
471 for (label parcelI = 0; parcelI < newParcels; parcelI++)
473 if (validInjection(parcelI))
476 scalar timeInj = time0_ + padTime + deltaT*parcelI/newParcels;
500 const scalar dt = time - timeInj;
509 cloud.setParcelThermoProperties(*pPtr, dt);
512 setProperties(parcelI, newParcels, timeInj, *pPtr);
515 cloud.checkParcelProperties(*pPtr, dt, fullyDescribed());
535 if (pPtr->nParticle() >= minParticlesPerParcel_)
538 massAdded += pPtr->nParticle()*pPtr->mass();
540 if (pPtr->move(
cloud, td, dt))
542 pPtr->typeId() = injectorID_;
543 cloud.addParticle(pPtr);
552 delayedVolume += pPtr->nParticle()*pPtr->volume();
562 postInjectCheck(parcelsAdded, massAdded);
566template<
class CloudType>
567template<
class TrackCloudType>
570 TrackCloudType&
cloud,
572 const scalar trackTime
580 const scalar time = this->owner().
db().
time().
value();
593 label parcelsAdded = 0;
594 scalar massAdded = 0.0;
597 label newParcels = parcelsToInject(0.0, 1.0);
600 for (label parcelI = 0; parcelI < newParcels; parcelI++)
603 scalar newVolumeFraction = 1.0/scalar(newParcels);
633 cloud.setParcelThermoProperties(*pPtr, 0.0);
636 setProperties(parcelI, newParcels, 0.0, *pPtr);
639 cloud.checkParcelProperties(*pPtr, 0.0, fullyDescribed());
654 pPtr->typeId() = injectorID_;
657 cloud.addParticle(pPtr);
659 massAdded += pPtr->nParticle()*pPtr->mass();
664 postInjectCheck(parcelsAdded, massAdded);
668template<
class CloudType>
671 os <<
" Injector " << this->modelName() <<
":" <<
nl
672 <<
" - parcels added = " << parcelsAddedTotal_ <<
nl
673 <<
" - mass introduced = " << massInjected_ <<
nl;
675 if (this->writeTime())
677 this->setModelProperty(
"volumeTotal", volumeTotal_);
678 this->setModelProperty(
"massInjected", massInjected_);
679 this->setModelProperty(
"nInjections", nInjections_);
680 this->setModelProperty(
"parcelsAddedTotal", parcelsAddedTotal_);
681 this->setModelProperty(
"timeStep0", timeStep0_);
reduce(hasMovingMesh, orOp< bool >())
Base class for cloud sub-models.
Templated base class for dsmc cloud.
const fvMesh & mesh() const
Return reference to the mesh.
const Mesh & mesh() const
Return mesh.
Top level data entry class for use in dictionaries. Provides a mechanism to specify a variable as a c...
const word & name() const noexcept
Return the object name.
const objectRegistry & db() const noexcept
Return the local objectRegistry.
Templated injection model class.
virtual bool prepareForNextTimeStep(const scalar time, label &newParcels, scalar &newVolumeFraction)
Determine properties for next time step/injection interval.
void inject(TrackCloudType &cloud, typename CloudType::parcelType::trackingData &td)
Main injection loop.
virtual scalar setNumberOfParticles(const label parcels, const scalar volumeFraction, const scalar diameter, const scalar rho)
Set number of particles to inject given parcel properties.
virtual scalar averageParcelMass()
Return the average parcel mass over the injection period.
void injectSteadyState(TrackCloudType &cloud, typename CloudType::parcelType::trackingData &td, const scalar trackTime)
Main injection loop - steady-state.
virtual void updateMesh()
Update mesh.
virtual bool findCellAtPosition(label &celli, label &tetFacei, label &tetPti, vector &position, bool errorOnNotFound=true)
Find the cell that contains the supplied position.
CloudType::parcelType parcelType
Convenience typedef for parcelType.
virtual void postInjectCheck(const label parcelsAdded, const scalar massAdded)
Post injection checks.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
virtual scalar userTimeToTime(const scalar theta) const
Convert the user-time (e.g. CA deg) to real-time (s).
A cloud is a registry collection of lagrangian particles.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
const Type & value() const
Return const reference to value.
InfoProxy< ensightCells > info() const
Return info proxy.
Class used to pass data into container.
const volVectorField & C() const
Return cell centres as volVectorField.
const Time & time() const
Return the top-level database.
const Time & time() const noexcept
Return time registry.
Mesh consisting of general polyhedral cells.
void findCellFacePt(const point &p, label &celli, label &tetFacei, label &tetPti) const
Find the cell, tetFacei and tetPti for point p.
label nGeometricD() const
Return the number of valid geometric dimensions in the mesh.
const Vector< label > & solutionD() const
Return the vector of solved-for directions in mesh.
label findNearestCell(const point &location) const
Find the cell with the nearest cell centre to location.
Selector class for relaxation factors, solver type and solution.
A class for handling words, derived from Foam::string.
const volScalarField & p0
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
OBJstream os(runTime.globalPath()/outputName)
constexpr scalar pi(M_PI)
dimensionedScalar pos(const dimensionedScalar &ds)
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
dimensionedScalar pow3(const dimensionedScalar &ds)
messageStream Info
Information stream (stdout output on master, null elsewhere)
Ostream & endl(Ostream &os)
Add newline and flush stream.
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
static constexpr const zero Zero
Global zero (0)
T returnReduce(const T &value, const BinaryOp &bop, const int tag=UPstream::msgType(), const label comm=UPstream::worldComm)
Reduce (copy) and return value.
errorManipArg< error, int > exit(error &err, const int errNo=1)
constexpr char nl
The newline '\n' character (0x0a)