Go to the documentation of this file.
37 template<
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;
93 template<
class CloudType>
141 position += SMALL*(cellCentres[celli] - position);
172 <<
"Cannot find parcel injection cell. "
173 <<
"Parcel position = " <<
p0 <<
nl
186 template<
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
229 template<
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;
253 time0_ = this->owner().db().time().value();
262 template<
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)
287 template<
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_);
357 SOI_ = owner.db().time().userTimeToTime(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'"
387 template<
class CloudType>
414 template<
class CloudType>
419 template<
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;
436 template<
class CloudType>
437 template<
class TrackCloudType>
440 TrackCloudType&
cloud,
441 typename CloudType::parcelType::trackingData& td
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);
566 template<
class CloudType>
567 template<
class TrackCloudType>
570 TrackCloudType&
cloud,
571 typename CloudType::parcelType::trackingData& td,
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);
668 template<
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_);
virtual void info(Ostream &os)
Write injection info to stream.
virtual bool findCellAtPosition(label &celli, label &tetFacei, label &tetPti, vector &position, bool errorOnNotFound=true)
Find the cell that contains the supplied position.
Selector class for relaxation factors, solver type and solution.
A class for handling words, derived from Foam::string.
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
static constexpr const zero Zero
Global zero (0)
label nGeometricD() const
Return the number of valid geometric dimensions in the mesh.
InjectionModel(CloudType &owner)
Construct null from owner.
virtual void updateMesh()
Update mesh.
Templated injection model class.
Base class for cloud sub-models.
virtual bool prepareForNextTimeStep(const scalar time, label &newParcels, scalar &newVolumeFraction)
Determine properties for next time step/injection interval.
Ostream & endl(Ostream &os)
Add newline and flush stream.
const Type & value() const
Return const reference to value.
const Vector< label > & solutionD() const
Return the vector of solved-for directions in mesh.
scalar massInjected_
Total mass injected to date [kg].
virtual scalar setNumberOfParticles(const label parcels, const scalar volumeFraction, const scalar diameter, const scalar rho)
Set number of particles to inject given parcel properties.
autoPtr< Function1< scalar > > massFlowRate_
Mass flow rate profile for steady calculations.
scalar timeStep0_
Time at start of injection time step [s].
label injectorID_
Optional injector ID.
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Top level data entry class for use in dictionaries. Provides a mechanism to specify a variable as a c...
Mesh consisting of general polyhedral cells.
label findNearestCell(const point &location) const
Find the cell with the nearest cell centre to location.
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
dimensionedScalar pow3(const dimensionedScalar &ds)
const fvMesh & mesh() const
Return reference to the mesh.
messageStream Info
Information stream (stdout output on master, null elsewhere)
void injectSteadyState(TrackCloudType &cloud, typename CloudType::parcelType::trackingData &td, const scalar trackTime)
Main injection loop - steady-state.
scalar time0_
Continuous phase time at start of injection time step [s].
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Templated base class for dsmc cloud.
scalar minParticlesPerParcel_
virtual void postInjectCheck(const label parcelsAdded, const scalar massAdded)
Post injection checks.
const volVectorField & C() const
Return cell centres as volVectorField.
virtual scalar averageParcelMass()
Return the average parcel mass over the injection period.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
void findCellFacePt(const point &p, label &celli, label &tetFacei, label &tetPti) const
Find the cell, tetFacei and tetPti for point p.
scalar massTotal_
Total mass to inject [kg].
OBJstream os(runTime.globalPath()/outputName)
label parcelsAddedTotal_
Running counter of total number of parcels added.
errorManipArg< error, int > exit(error &err, const int errNo=1)
A cloud is a registry collection of lagrangian particles.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
static int myProcNo(const label communicator=worldComm)
Number of this process (starting from masterNo() = 0)
constexpr scalar pi(M_PI)
void inject(TrackCloudType &cloud, typename CloudType::parcelType::trackingData &td)
Main injection loop.
CloudType::parcelType parcelType
Convenience typedef for parcelType.
scalar SOI_
Start of injection [s].
const Time & time() const
Return the top-level database.
label nInjections_
Number of injections counter.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
const volScalarField & p0
parcelBasis parcelBasis_
Parcel basis enumeration.
dimensionedScalar pos(const dimensionedScalar &ds)