InsertionBoundary Class Referenceabstract

Boundary structure for boundaries used for insertion of particles. More...

#include <InsertionBoundary.h>

+ Inheritance diagram for InsertionBoundary:

Public Member Functions

 InsertionBoundary ()
 
 InsertionBoundary (const InsertionBoundary &other)
 Copy constructor (with deep copy) More...
 
 ~InsertionBoundary () override
 Destructor: delete the particle that has to be copied at every insertion. More...
 
virtual BaseParticlegenerateParticle (RNG &random)
 Sets the properties of the InsertionBoundary for a single particle type ‍/ virtual void set(BaseParticle* particleToCopy, unsigned int maxFailed, Vec3D velMin, Vec3D velMax, double radMin, double radMax)=0;. More...
 
virtual void placeParticle (BaseParticle *p, RNG &random)=0
 Purely virtual function that generates the extrinsic properties (position, velocity) of a particle. More...
 
void checkBoundaryBeforeTimeStep (DPMBase *md) override
 Fills the boundary with particles. More...
 
void insertParticles (DPMBase *md)
 Fill a certain domain with particles. More...
 
unsigned int getNumberOfParticlesInserted () const
 Gets the number of particles inserted by the boundary. More...
 
Mdouble getMassOfParticlesInserted () const
 Gets the mass of particles inserted by the boundary. More...
 
Mdouble getVolumeOfParticlesInserted () const
 Gets the volume of particles inserted by the boundary. More...
 
void reset ()
 resets particle property counter variables. More...
 
void activate ()
 Turns on the InsertionBoundary. More...
 
void deactivate ()
 Turns off the InsertionBoundary. More...
 
bool isActivated ()
 Returns whether the InsertionBoundary is activated. More...
 
unsigned int getMaxFailed () const
 Gets the number of times that the boundary may fail to insert a particle. More...
 
void setParticleToCopy (std::vector< BaseParticle * > particleToCopy)
 Sets multiple different particles that will be inserted through the insertion boundary. More...
 
void setParticleToCopy (BaseParticle *particleToCopy)
 Sets the particle that will be inserted through the insertion boundary. More...
 
std::vector< BaseParticle * > getParticleToCopy ()
 Gets the particles that will be inserted through the insertion boundary. More...
 
void read (std::istream &is) override
 Reads the boundary's id_ and maxFailed_. More...
 
void write (std::ostream &os) const override
 Writes the boundary's id_ and maxFailed_. More...
 
Mdouble getVolumeFlowRate () const
 Gets the volume flow rate of the insertion routine. More...
 
void setVolumeFlowRate (Mdouble volumeFlowRate_)
 Sets the volume flow rate of the insertion routine. More...
 
Mdouble getInitialVolume () const
 Gets the initialVolume() . More...
 
void setInitialVolume (Mdouble initialVolume)
 Gets the Volume which should be inserted by the insertion routine. More...
 
void setPSD (const PSD psd)
 Sets the range of particle radii that may be generated from a user defined PSD. More...
 
void setPSD (std::vector< PSD > psd, std::vector< Mdouble > probability)
 Sets the ranges of particle radii that may be generated from user defined PSDs. More...
 
std::vector< PSDgetPSD ()
 Gets the particle size distributions set by the user. More...
 
void setVariableVolumeFlowRate (const std::vector< Mdouble > &variableCumulativeVolumeFlowRate, Mdouble samplingInterval)
 Sets a variable volume flow rate. More...
 
bool insertParticle (Mdouble time)
 Checks the inserted total volume and returns if a particle is still allowed to be inserted. More...
 
bool getCheckParticleForInteraction () const
 Gets the variable that checks if a particle has an interaction. More...
 
void setCheckParticleForInteraction (bool checkParticleForInteraction)
 Sets the variable that checks if a particle has an interaction. More...
 
void setManualInsertion (bool manualInsertion)
 Set the flag for a manual PSD insertion routine. More...
 
- Public Member Functions inherited from BaseBoundary
 BaseBoundary ()
 default constructor. More...
 
 BaseBoundary (const BaseBoundary &b)
 copy constructor More...
 
 ~BaseBoundary () override
 destructor More...
 
virtual BaseBoundarycopy () const =0
 Used to create a copy of the object NB: purely virtual function. More...
 
virtual void createPeriodicParticle (BaseParticle *p UNUSED, ParticleHandler &pH UNUSED)
 Creates a periodic particle in case of periodic boundaries in serial build. More...
 
virtual void createPeriodicParticles (ParticleHandler &pH UNUSED)
 Creates periodic copies of given particle in case of periodic boundaries. More...
 
virtual void checkBoundaryAfterParticlesMove (ParticleHandler &pH)
 Virtual function that does things to particles, each time step after particles have moved. More...
 
virtual void actionsBeforeTimeLoop ()
 Virtual function that does something after DPMBase::setupInitialConditions but before the first time step. More...
 
virtual void modifyGhostAfterCreation (BaseParticle *particle, int i)
 
virtual void writeVTK (std::fstream &file)
 
void setHandler (BoundaryHandler *handler)
 Sets the boundary's BoundaryHandler. More...
 
BoundaryHandlergetHandler () const
 Returns the boundary's BoundaryHandler. More...
 
- Public Member Functions inherited from BaseObject
 BaseObject ()=default
 Default constructor. More...
 
 BaseObject (const BaseObject &p)=default
 Copy constructor, copies all the objects BaseObject contains. More...
 
virtual ~BaseObject ()=default
 virtual destructor More...
 
virtual std::string getName () const =0
 A purely virtual function. More...
 
virtual void moveInHandler (unsigned int index)
 Except that it is virtual, it does the same thing as setIndex() does. More...
 
void setIndex (unsigned int index)
 Allows one to assign an index to an object in the handler/container. More...
 
void setId (unsigned long id)
 Assigns a unique identifier to each object in the handler (container) which remains constant even after the object is deleted from the container/handler. More...
 
unsigned int getIndex () const
 Returns the index of the object in the handler. More...
 
unsigned int getId () const
 Returns the unique identifier of any particular object. More...
 
void setGroupId (unsigned groupId)
 
unsigned getGroupId () const
 

Public Attributes

std::vector< BaseParticle * > particleToCopy_
 read Distribution class from file. ‍/ friend std::istream& operator>>(std::istream& is, InsertionBoundary::Distribution& type); More...
 
unsigned int maxFailed_
 Number of times that the wall may fail to insert a particle. More...
 
unsigned int numberOfParticlesInserted_
 Number of particles that are already inserted. More...
 
Mdouble massInserted_
 Total mass of particles inserted. More...
 
Mdouble volumeInserted_
 Total volume of particles inserted. More...
 
bool isActivated_
 The InsertionBoundary is activated by default. If the InsertionBoundary is deactivated, then it introduces no particles (useful for trying to maintain a certain insertion rate). More...
 
Mdouble volumeFlowRate_
 
Mdouble initialVolume_
 
std::vector< MdoublevariableCumulativeVolumeFlowRate_
 
Mdouble samplingInterval_
 
bool checkParticleForInteraction_
 Checks if a particle has an interaction with a wall or other particles. More...
 
std::vector< PSDparticleSizeDistributionVector_
 Defines a particle size distribution as an object of the PSD class; if particleSizeDistributionVector_ is empty, distribution_ is used instead. More...
 
Vec3D velMin_
 Minimum and maximum velocity of the particles to be inserted. More...
 
Vec3D velMax_
 
bool isManuallyInserting_
 A flag to enable a top-down class-by-class manual insertion of a PSD; default is FALSE. More...
 
std::vector< Mdoubleprobability_
 vector of probabilities in range [0,1] which determine the mixing ratio of partice size distributions. More...
 
int chosenSpecies_
 stores the chosen species for each timestep. More...
 

Detailed Description

Boundary structure for boundaries used for insertion of particles.

To cite the InsertionBoundary algorithm: A. R. Thornton, D. Krijgsman, A. Te Voortwis, V. Ogarko, S. Luding, R. Fransen, S. Gonzalez, O. Bokhove, O. Imole, and T. Weinhart. A review of recent work on the discrete particle method at the University of Twente: An introduction to the open-source package MercuryDPM. DEM6 - International Conference on DEMs, 2013.

Constructor & Destructor Documentation

◆ InsertionBoundary() [1/2]

InsertionBoundary::InsertionBoundary ( )
Todo:
think about making both possible; a discrete distribution and a continuous which is more accurate
 \brief Defines a custom particle size distribution; distribution_ will always be used, unless particleSizeDistributionVector_ is non-empty
&zwj;/

enum class Distribution { Uniform, Normal_1_5 // TODO add LogNormal distribution // LogNormal };

/*!
   \brief Default constructor: set everything to 0/nullptr.

Default constructor, sets all data members to 0 or default.

38  : BaseBoundary()
39 {
41  massInserted_ = 0;
42  volumeInserted_ = 0;
43  maxFailed_ = 0;
44  isActivated_ = true;
46  initialVolume_ = 0;
49  velMin_ = Vec3D(0.0, 0.0, 0.0);
50  velMax_ = Vec3D(0.0, 0.0, 0.0);
52  isManuallyInserting_ = false;
53  chosenSpecies_ = 0;
54 }
BaseBoundary()
default constructor.
Definition: BaseBoundary.cc:32
Mdouble samplingInterval_
Definition: InsertionBoundary.h:318
int chosenSpecies_
stores the chosen species for each timestep.
Definition: InsertionBoundary.h:349
unsigned int maxFailed_
Number of times that the wall may fail to insert a particle.
Definition: InsertionBoundary.h:270
unsigned int numberOfParticlesInserted_
Number of particles that are already inserted.
Definition: InsertionBoundary.h:275
bool isManuallyInserting_
A flag to enable a top-down class-by-class manual insertion of a PSD; default is FALSE.
Definition: InsertionBoundary.h:339
Mdouble initialVolume_
Definition: InsertionBoundary.h:307
bool isActivated_
The InsertionBoundary is activated by default. If the InsertionBoundary is deactivated,...
Definition: InsertionBoundary.h:293
bool checkParticleForInteraction_
Checks if a particle has an interaction with a wall or other particles.
Definition: InsertionBoundary.h:323
Vec3D velMin_
Minimum and maximum velocity of the particles to be inserted.
Definition: InsertionBoundary.h:334
Mdouble massInserted_
Total mass of particles inserted.
Definition: InsertionBoundary.h:280
std::vector< PSD > particleSizeDistributionVector_
Defines a particle size distribution as an object of the PSD class; if particleSizeDistributionVector...
Definition: InsertionBoundary.h:329
Mdouble volumeFlowRate_
Definition: InsertionBoundary.h:304
Vec3D velMax_
Definition: InsertionBoundary.h:334
Mdouble volumeInserted_
Total volume of particles inserted.
Definition: InsertionBoundary.h:285
Definition: Vector.h:51
const Mdouble inf
Definition: GeneralDefine.h:44

References checkParticleForInteraction_, chosenSpecies_, constants::inf, initialVolume_, isActivated_, isManuallyInserting_, massInserted_, maxFailed_, numberOfParticlesInserted_, particleSizeDistributionVector_, samplingInterval_, velMax_, velMin_, volumeFlowRate_, and volumeInserted_.

◆ InsertionBoundary() [2/2]

InsertionBoundary::InsertionBoundary ( const InsertionBoundary other)

Copy constructor (with deep copy)

Copy constructor

60  : BaseBoundary(other)
61 {
65  maxFailed_ = other.maxFailed_;
66  isActivated_ = other.isActivated_;
73  probability_ = other.probability_;
74  velMin_ = other.velMin_;
75  velMax_ = other.velMax_;
78 
79  for (int i = 0; i < other.particleToCopy_.size(); i++)
80  {
81  particleToCopy_.resize(other.particleToCopy_.size());
82  particleToCopy_[i] = other.particleToCopy_[i]->copy();
83  }
84 }
std::vector< Mdouble > probability_
vector of probabilities in range [0,1] which determine the mixing ratio of partice size distributions...
Definition: InsertionBoundary.h:344
std::vector< BaseParticle * > particleToCopy_
read Distribution class from file. ‍/ friend std::istream& operator>>(std::istream& is,...
Definition: InsertionBoundary.h:265
std::vector< Mdouble > variableCumulativeVolumeFlowRate_
Definition: InsertionBoundary.h:315
const std::complex< Mdouble > i
Definition: ExtendedMath.h:51

References checkParticleForInteraction_, chosenSpecies_, constants::i, initialVolume_, isActivated_, isManuallyInserting_, massInserted_, maxFailed_, numberOfParticlesInserted_, particleSizeDistributionVector_, particleToCopy_, probability_, samplingInterval_, variableCumulativeVolumeFlowRate_, velMax_, velMin_, volumeFlowRate_, and volumeInserted_.

◆ ~InsertionBoundary()

InsertionBoundary::~InsertionBoundary ( )
override

Destructor: delete the particle that has to be copied at every insertion.

Destructor that deletes the BaseParticle that is copied and inserted at every insertion.

91 {
92  for (auto& particleToCopy: particleToCopy_)
93  {
94  delete particleToCopy;
95  }
96 
97 }

References particleToCopy_.

Member Function Documentation

◆ activate()

void InsertionBoundary::activate ( )

Turns on the InsertionBoundary.

Turns on the InsertionBoundary by setting the boolean to TRUE.

376 {
377  isActivated_ = true;
378 }

References isActivated_.

Referenced by Chutebelt::actionsAfterTimeStep().

◆ checkBoundaryBeforeTimeStep()

void InsertionBoundary::checkBoundaryBeforeTimeStep ( DPMBase md)
overridevirtual

Fills the boundary with particles.

Is used to fill the insides of the boundary with particles until it is filled up.

Parameters
[in,out]mdthe problem's DPMBase object
Todo:
rename to something like "insertUntilMaxFailed"?

HERE

Todo:
create a definition for zero-size particles. Right now the zero-size particle is only used as a stop criterion for the manual PSD insertion

HERE

Reimplemented from BaseBoundary.

Reimplemented in RandomClusterInsertionBoundary.

185 {
186  logger(VERBOSE, "In InsertionBoundary::checkBoundaryBeforeTimeStep\n");
187 
188  if (!isActivated_)
189  {
190  return;
191  }
192 
193  /* Each timestep, the InsertionBoundary attempts to fill up a region with
194  * particles.
195  *
196  * It first calls generateParticle() to get a randomised particle, subject
197  * to a specified distribution over sizes and species. (The basic class
198  * supports size dispersity only but PolydisperseInsertionBoundary will
199  * support species dispersity.)
200  * Then it repeatedly calls placeParticle(), which gives the particle a
201  * random location (and possibly velocity) in a specified,
202  * geometry-dependent bound. Each time, it checks whether the new particle
203  * would have an interaction with another particle or a wall.
204  *
205  * If it manages to do that within maxFailed_ tries, then:
206  * * the new particle is inserted,
207  * * the failure counter is reset, and
208  * the processes is repeated with a new generateParticle().
209  *
210  * Otherwise, the processes terminates for this timestep.
211  * */
212 
213  // Keep count of how many successive times we have failed to place a new
214  // particle.
215  unsigned int failed = 0;
216  while (failed <= maxFailed_ && insertParticle(md->getNextTime())) // 'generating' loop
217  {
218  /* Generate random *intrinsic* properties for the new particle. */
219  logger(VERBOSE, "about to call generateParticle\n");
220 
221 
222 
223 
225 
226 
227  auto p0 = generateParticle(md->random);
228  // Important for particle generation with a particle size distribution as it generates a particle with zero
229  // radius. If a particle is not allowed to be inserted by the PSD criteria it will generate a particle with
230  // zero diameter. This if statement prevents inserting particles with zero radius, which would else be a problem.
232  if (p0->getRadius() == 0)
233  {
234  logger(VERBOSE, "The PSD for the specified volume is fully set");
235  // free up memory space
236  delete p0;
237  // out of the 'placing' loop
238  failed = maxFailed_ + 1;
239  continue;
240  }
241  logger(VERBOSE, "generated a particle with intrinsics %", p0);
242 
243  while (true) // 'placing' loop
244  {
245  /* Generate extrinsic properties (position and velocity) for this
246  * new particle. */
247 
248 
249 
251 
252 
253 
254 
255  placeParticle(p0, md->random);
256  logger(VERBOSE, "attempting to place particle at %, vel %", p0->getPosition(), p0->getVelocity());
257 
258 #ifdef MERCURYDPM_USE_MPI
259  /* Communicate the new particle's properties by setHandler (note
260  * that this doesn't actually add the particle to the handler). */
261  if (NUMBER_OF_PROCESSORS > 1)
262  {
263  MPIParticle particle;
264  // //Every domain generates a particle (to get the species right etc)
265 
266  //Send particle data from root to other processors to sync the particle properties
267  if (PROCESSOR_ID == 0)
268  {
270  }
271 
273 
274  //Process the received data
275  if (PROCESSOR_ID != 0)
276  {
277  copyDataFromMPIParticleToParticle(&particle, p0, &(md->particleHandler));
278  }
279  }
280 #endif
281  p0->setHandler(&md->particleHandler);
282  /* Check whether the particle has any interactions. */
284  {
285  //Note: in parallel only one of the domains will actually add the particle
287  failed = 0;
288 
290  const double volume = p0->getVolume();
291  volumeInserted_ += volume;
292  massInserted_ += p0->getSpecies()->getDensity() * volume;
293  logger(VERBOSE, "successfully placed a particle %, with position: % after % fails.", p0,
294  p0->getPosition(), failed);
295  /* JMFT: The generateParticle() routine allocates memory, so we should
296  * free it here. (Don't worry, the particle will have been copied to the
297  * particleHandler by this point iff we want it.) */
298  delete p0;
299 
300  break; // out of the 'placing' loop
301  }
302  else
303  {
304  failed++;
305  logger(VERBOSE, "failed to place a particle; have failed % times", failed);
306  }
307 
308  if (failed > maxFailed_)
309  {
310  logger(VERBOSE, "failed too many times; giving up");
312  {
313  particleSizeDistributionVector_[chosenSpecies_].decrementNParticlesPerClass();
314  particleSizeDistributionVector_[chosenSpecies_].decrementVolumePerClass(p0->getVolume());
315  }
316  break; // out of the 'placing' loop (and will leave the 'generating' loop too
317  }
318  }
319  logger(VERBOSE, "failed % times, so breaking out of InsertionBoundary loop for this timestep.", failed);
320  }
321  // logger(INFO, "volumeInserted_ = %", volumeInserted_);
322 }
#define PROCESSOR_ID
Definition: GeneralDefine.h:63
#define NUMBER_OF_PROCESSORS
For the MPI communication routines this quantity is often required. defining this macro makes the cod...
Definition: GeneralDefine.h:62
Logger< MERCURYDPM_LOGLEVEL > logger("MercuryKernel")
Definition of different loggers with certain modules. A user can define its own custom logger here.
@ VERBOSE
@ PARTICLE
Definition: MpiContainer.h:67
void copyDataFromMPIParticleToParticle(MPIParticle *bP, BaseParticle *p, ParticleHandler *particleHandler)
Copies data from an MPIParticle class to a BaseParticle and sets the particleHandler and species.
Definition: MpiDataClass.cc:105
std::enable_if<!std::is_pointer< U >::value, U * >::type copyAndAddObject(const U &object)
Creates a copy of a Object and adds it to the BaseHandler.
Definition: BaseHandler.h:379
Mdouble getNextTime() const
Returns the current simulation time.
Definition: DPMBase.cc:816
ParticleHandler particleHandler
An object of the class ParticleHandler, contains the pointers to all the particles created.
Definition: DPMBase.h:1437
RNG random
This is a random generator, often used for setting up the initial conditions etc.....
Definition: DPMBase.h:1432
virtual bool checkParticleForInteraction(const BaseParticle &P)
Checks whether a particle P has any interaction with walls or other particles.
Definition: DPMBase.cc:4805
virtual void placeParticle(BaseParticle *p, RNG &random)=0
Purely virtual function that generates the extrinsic properties (position, velocity) of a particle.
bool insertParticle(Mdouble time)
Checks the inserted total volume and returns if a particle is still allowed to be inserted.
Definition: InsertionBoundary.cc:150
virtual BaseParticle * generateParticle(RNG &random)
Sets the properties of the InsertionBoundary for a single particle type ‍/ virtual void set(BaseParti...
Definition: InsertionBoundary.cc:103
std::enable_if< std::is_scalar< T >::value, void >::type broadcast(T &t, int fromProcessor=0)
Broadcasts a scalar from the root to all other processors.
Definition: MpiContainer.h:441
static MPIContainer & Instance()
fetch the instance to be used for communication
Definition: MpiContainer.h:134
Data class to send a particle over MPI.
Definition: MpiDataClass.h:81
void copyDataFromParticleToMPIParticle(BaseParticle *p)
Definition: MpiDataClass.cc:131

References MPIContainer::broadcast(), DPMBase::checkParticleForInteraction(), checkParticleForInteraction_, chosenSpecies_, BaseHandler< T >::copyAndAddObject(), copyDataFromMPIParticleToParticle(), MPISphericalParticle::copyDataFromParticleToMPIParticle(), generateParticle(), DPMBase::getNextTime(), insertParticle(), MPIContainer::Instance(), isActivated_, isManuallyInserting_, logger, massInserted_, maxFailed_, NUMBER_OF_PROCESSORS, numberOfParticlesInserted_, PARTICLE, DPMBase::particleHandler, particleSizeDistributionVector_, placeParticle(), PROCESSOR_ID, DPMBase::random, VERBOSE, and volumeInserted_.

Referenced by GranuDrum::GranuDrum(), GranuHeap::GranuHeap(), insertParticles(), main(), FreeCooling2DinWalls::setupInitialConditions(), StressStrainControl::setupInitialConditions(), and ConstantMassFlowMaserSelfTest::setupInitialConditions().

◆ deactivate()

void InsertionBoundary::deactivate ( )

Turns off the InsertionBoundary.

Turns off the InsertionBoundary by setting the boolean to FALSE.

384 {
385  isActivated_ = false;
386 }

References isActivated_.

Referenced by Chutebelt::actionsAfterTimeStep(), and Chutebelt::setupInitialConditions().

◆ generateParticle()

BaseParticle * InsertionBoundary::generateParticle ( RNG random)
virtual

Sets the properties of the InsertionBoundary for a single particle type ‍/ virtual void set(BaseParticle* particleToCopy, unsigned int maxFailed, Vec3D velMin, Vec3D velMax, double radMin, double radMax)=0;.

     \brief Sets the properties of the InsertionBoundary for mutliple different particle types
    &zwj;/

virtual void set(std::vector<BaseParticle*> particleToCopy, unsigned int maxFailed, Vec3D velMin, Vec3D velMax, double radMin, double radMax)=0;

/*!

/*!
   \brief Virtual function that generates the intrinsic properties
   (species, radius) of one particle.
   \param[in] random  Random number generator

The default behaviour will be to return particleToCopy_, but this can be overridden by the children (to get size or species dispersity).

Reimplemented in PolydisperseInsertionBoundary, BidisperseCubeInsertionBoundary, and FixedClusterInsertionBoundary.

104 {
105  double check = random.getRandomNumber(0, 1);
106  for (int i = 0; i < probability_.size(); i++)
107  {
108  if (check < probability_[i])
109  {
110  chosenSpecies_ = i;
111  break;
112  }
113  else
114  {
115  check -= probability_[i];
116  }
117  }
119  if (particleSizeDistributionVector_[chosenSpecies_].getParticleSizeDistribution().empty())
120  {
121  particleSizeDistributionVector_[chosenSpecies_].setDistributionUniform(P->getRadius(), P->getRadius(), 50);
122  logger(INFO, "Assembling a discrete uniform particle size distribution (50 bins) with the radius set for the "
123  "InsertionBoundary.");
124  }
125  // manual insertion routine to insert PSDs as accurate as possible into a given volume
127  {
128  logger.assert_debug(initialVolume_ > 0.0, "Use setInitialVolume to define the particle insertion volume");
129  Mdouble radius;
130  // getVolumeFlowRate() * time + initialVolume_ - volumeInserted_ lead to more inaccurate results, therfore
131  // -volumeInserted was removed.
132  radius = particleSizeDistributionVector_[chosenSpecies_].insertManuallyByVolume(getVolumeFlowRate() *
133  getHandler()->getDPMBase()->getTime() +
135  P->setRadius(radius);
136  return P;
137  }
138  Mdouble radius;
139  radius = particleSizeDistributionVector_[chosenSpecies_].drawSample();
140  P->setRadius(radius);
141  return P;
142 }
double Mdouble
Definition: GeneralDefine.h:34
@ INFO
BoundaryHandler * getHandler() const
Returns the boundary's BoundaryHandler.
Definition: BaseBoundary.cc:143
Definition: BaseParticle.h:54
Mdouble getVolumeFlowRate() const
Gets the volume flow rate of the insertion routine.
Definition: InsertionBoundary.cc:617
Mdouble getRandomNumber()
This is a random generating routine can be used for initial positions.
Definition: RNG.cc:143
double P
Uniform pressure.
Definition: TwenteMeshGluing.cpp:73
void check(double real, double ideal, double error, std::string errorMessage)
Definition: TestHelpers.cc:37

References helpers::check(), chosenSpecies_, BaseBoundary::getHandler(), RNG::getRandomNumber(), getVolumeFlowRate(), constants::i, INFO, initialVolume_, isManuallyInserting_, logger, Global_Physical_Variables::P, particleSizeDistributionVector_, particleToCopy_, and probability_.

Referenced by checkBoundaryBeforeTimeStep(), RandomClusterInsertionBoundary::checkBoundaryBeforeTimeStep(), and CubicCell::setupInitialConditions().

◆ getCheckParticleForInteraction()

bool InsertionBoundary::getCheckParticleForInteraction ( ) const

Gets the variable that checks if a particle has an interaction.

Gets the variable that checks if a particle has an interaction with a wall or another particle.

Returns
if TRUE the particle has an interaction and if FALSE the particle has no interaction.
719 {
721 }

References checkParticleForInteraction_.

◆ getInitialVolume()

Mdouble InsertionBoundary::getInitialVolume ( ) const

Gets the initialVolume() .

Gets the volume to be inserted by the insertion routine.

Returns
A double which corresponds to volume to be inserted by the insertion routines
635 {
636  return initialVolume_;
637 }

References initialVolume_.

Referenced by NozzleDemo::actionsAfterTimeStep(), Chutebelt::actionsAfterTimeStep(), and NozzleSelfTest::actionsAfterTimeStep().

◆ getMassOfParticlesInserted()

Mdouble InsertionBoundary::getMassOfParticlesInserted ( ) const

Gets the mass of particles inserted by the boundary.

Returns the mass of particles inserted in the boundary

Returns
the mass of particles inserted
348 {
349  return massInserted_;
350 }

References massInserted_.

Referenced by Chutebelt::actionsAfterTimeStep(), BoundariesSelfTest::actionsAfterTimeStep(), and FluxBoundarySelfTest::actionsAfterTimeStep().

◆ getMaxFailed()

unsigned int InsertionBoundary::getMaxFailed ( ) const

Gets the number of times that the boundary may fail to insert a particle.

Return maxFailed_ (see InsertionBoundary::set).

Returns
the maximum number of particle insertion trials
402 {
403  return maxFailed_;
404 }

References maxFailed_.

◆ getNumberOfParticlesInserted()

unsigned int InsertionBoundary::getNumberOfParticlesInserted ( ) const

Gets the number of particles inserted by the boundary.

Returns the number of particles inserted in the boundary

Returns
the number of particles inserted
339 {
341 }

References numberOfParticlesInserted_.

Referenced by T_protectiveWall::actionsAfterTimeStep(), BoundariesSelfTest::actionsAfterTimeStep(), FluxBoundarySelfTest::actionsAfterTimeStep(), insertParticles(), and Chute::write().

◆ getParticleToCopy()

std::vector< BaseParticle * > InsertionBoundary::getParticleToCopy ( )

Gets the particles that will be inserted through the insertion boundary.

returns pointer to the particle copies which are to be inserted

459 {
460  if (particleToCopy_.empty())
461  {
462  logger(ERROR, "particleToCopy not set");
463  }
464  return particleToCopy_;
465 }
@ ERROR

References ERROR, logger, and particleToCopy_.

◆ getPSD()

std::vector< PSD > InsertionBoundary::getPSD ( )

Gets the particle size distributions set by the user.

gets the user defined particle size distributions

Returns
a vector of PSD class objects containing containing the differnt user defined particle size distributions
701 {
703 }

References particleSizeDistributionVector_.

◆ getVolumeFlowRate()

Mdouble InsertionBoundary::getVolumeFlowRate ( ) const

Gets the volume flow rate of the insertion routine.

Gets the volumetric flow rate of the insertion routine.

Returns
A double which corresponds to the flow rate of the insertion routine
618 {
619  return volumeFlowRate_;
620 }

References volumeFlowRate_.

Referenced by generateParticle(), and insertParticle().

◆ getVolumeOfParticlesInserted()

Mdouble InsertionBoundary::getVolumeOfParticlesInserted ( ) const

Gets the volume of particles inserted by the boundary.

Returns the volume of particles inserted in the boundary

Returns
the volume of particles inserted
357 {
358  return volumeInserted_;
359 }

References volumeInserted_.

Referenced by T_protectiveWall::actionsAfterTimeStep(), NozzleDemo::actionsAfterTimeStep(), BoundariesSelfTest::actionsAfterTimeStep(), FluxBoundarySelfTest::actionsAfterTimeStep(), NozzleSelfTest::actionsAfterTimeStep(), protectiveWall::actionsAfterTimeStep(), and GranuHeap::GranuHeap().

◆ insertParticle()

bool InsertionBoundary::insertParticle ( Mdouble  time)

Checks the inserted total volume and returns if a particle is still allowed to be inserted.

Checks the inserted total volume by the flowrate and initialVolume and returns if a particle is allowed to be inserted.

Returns
TRUE if the inserted volume is lower than the total volume and FALSE if the inserted volume is higher than the total volume.
151 {
152  // check if the flow rate limit has been reached
154  {
156  }
157  else
158  {
159  const Mdouble iMax = (Mdouble) variableCumulativeVolumeFlowRate_.size() - 2;
160  const Mdouble i = std::min(time / samplingInterval_, iMax);
161  if (i == iMax)
162  {
163  static unsigned count = 0;
164  if (count == 0)
165  {
166  logger(WARN, "Reached end of volume flowrate function");
167  }
168  ++count;
169  }
170  const size_t id = i;
171  const Mdouble allowedVolume = variableCumulativeVolumeFlowRate_[id] +
173  variableCumulativeVolumeFlowRate_[id]) * (i - id);
174  return volumeInserted_ < allowedVolume;
175  }
176 }
@ WARN

References getVolumeFlowRate(), constants::i, initialVolume_, logger, samplingInterval_, variableCumulativeVolumeFlowRate_, volumeInserted_, and WARN.

Referenced by checkBoundaryBeforeTimeStep(), and RandomClusterInsertionBoundary::checkBoundaryBeforeTimeStep().

◆ insertParticles()

void InsertionBoundary::insertParticles ( DPMBase md)

Fill a certain domain with particles.

\detail calls the function checkBoundaryBeforeTimeStep() to fill a domain with particles; also reports how many particles where inserted at the end of the routine.

329 {
331  logger(INFO, "Inserted % particles", getNumberOfParticlesInserted());
332 }
unsigned int getNumberOfParticlesInserted() const
Gets the number of particles inserted by the boundary.
Definition: InsertionBoundary.cc:338
void checkBoundaryBeforeTimeStep(DPMBase *md) override
Fills the boundary with particles.
Definition: InsertionBoundary.cc:184

References checkBoundaryBeforeTimeStep(), getNumberOfParticlesInserted(), INFO, and logger.

Referenced by DPMBase::fillDomainWithParticles(), and RotatingDrumWet::setupInitialConditions().

◆ isActivated()

bool InsertionBoundary::isActivated ( )

Returns whether the InsertionBoundary is activated.

checks the activation status of the InsertionBoundary by checking the respective boolean variable.

Returns
TRUE for an activated InsertionBoundary and FALSE for a deactivated InsertionBoundary
393 {
394  return isActivated_;
395 }

References isActivated_.

◆ placeParticle()

virtual void InsertionBoundary::placeParticle ( BaseParticle p,
RNG random 
)
pure virtual

Purely virtual function that generates the extrinsic properties (position, velocity) of a particle.

Parameters
[in]pThe particle to be placed
[in]randomRandom number generator

This should be implemented by the children such as CubeInsertionBoundary, as the implementation will be geometry-dependent.

Implemented in RandomClusterInsertionBoundary, PolydisperseInsertionBoundary, HopperInsertionBoundary, CubeInsertionBoundary, ChuteInsertionBoundary, BaseClusterInsertionBoundary, and FixedClusterInsertionBoundary.

Referenced by checkBoundaryBeforeTimeStep().

◆ read()

void InsertionBoundary::read ( std::istream &  is)
overridevirtual

Reads the boundary's id_ and maxFailed_.

reads the boundary's id_ and maxFailed_ from the given istream

Parameters
[in,out]isstream the data members are read from
Todo:
make theses reads non-optional

Implements BaseBoundary.

Reimplemented in PolydisperseInsertionBoundary.

472 {
473  BaseBoundary::read(is);
474  std::string dummy, type;
475  is >> dummy >> maxFailed_ >> dummy;
476  if (dummy == "volumeFlowRate")
477  {
478  is >> volumeFlowRate_ >> dummy;
479  }
480  is >> massInserted_;
481  is >> dummy >> volumeInserted_;
482  is >> dummy >> numberOfParticlesInserted_;
483  is >> dummy >> isActivated_;
484  size_t psdVectorSize;
485  is >> dummy >> psdVectorSize;
487  particleSizeDistributionVector_.resize(psdVectorSize);
488  for (auto& particleSizeDistributionVector: particleSizeDistributionVector_)
489  {
490  size_t psdSize;
491  DistributionElements radiusAndProbability{};
492  std::vector<DistributionElements> particleSizeDistribution{};
493  is >> dummy >> psdSize;
494  particleSizeDistribution.clear();
495  particleSizeDistribution.reserve(psdSize);
496  for (size_t i = 0; i < psdSize; i++)
497  {
498  is >> radiusAndProbability.internalVariable;
499  is >> radiusAndProbability.probability;
500  particleSizeDistribution.push_back(radiusAndProbability);
501  }
502  particleSizeDistributionVector.setParticleSizeDistribution(particleSizeDistribution);
503  }
504  if (psdVectorSize > 1)
505  {
506  is >> dummy;
507  Mdouble psdRatio;
508  probability_.clear();
509  probability_.reserve(psdVectorSize);
510  for (size_t i = 0; i < psdVectorSize; i++)
511  {
512  is >> psdRatio;
513  probability_.push_back(psdRatio);
514  }
515  }
516 
518  helpers::readOptionalVariable(is, "checkParticleForInteraction", checkParticleForInteraction_);
519  helpers::readOptionalVariable(is, "initialVolume", initialVolume_);
520  if (helpers::readOptionalVariable(is, "samplingInterval", samplingInterval_))
521  {
522  size_t n;
523  Mdouble flowRate;
524  is >> dummy >> n;
525  //variableCumulativeVolumeFlowRate_.clear();
527  for (size_t i = 0; i < n; ++i)
528  {
529  is >> flowRate;
530  variableCumulativeVolumeFlowRate_.push_back(flowRate);
531  }
532  }
533  is >> dummy;
534  if (dummy != "noParticleToCopy")
535  {
536  size_t particleToCopySize;
537  for (auto& particleToCopy: particleToCopy_)
538  {
539  delete particleToCopy;
540  }
541  is >> particleToCopySize;
542  particleToCopy_.resize(particleToCopySize);
543  for (auto& particleToCopy: particleToCopy_)
544  {
545  particleToCopy = getHandler()->getDPMBase()->particleHandler.readAndCreateObject(is);
546  // The .restart file records the index of the particle's species, but
547  // doesn't record the pointer, i.e. the memory address of the species within
548  // the speciesHandler. The latter needs to be reset now.
549  particleToCopy->setSpecies(getHandler()->getDPMBase()->speciesHandler.getObject(
550  particleToCopy->getIndSpecies()));
551  }
552  }
553 }
const unsigned n
Definition: CG3DPackingUnitTest.cpp:32
void read(std::istream &is) override=0
Reads the object's id_ from given istream NB: purely virtual function, overriding the version of Base...
Definition: BaseBoundary.cc:61
DPMBase * getDPMBase()
Gets the problem that is solved using this handler.
Definition: BaseHandler.h:725
T * getObject(const unsigned int id)
Gets a pointer to the Object at the specified index in the BaseHandler.
Definition: BaseHandler.h:613
void setSpecies(const ParticleSpecies *species)
Definition: BaseParticle.cc:818
class of DistributionElements which stores internalVariables and probabilities of a distribution....
Definition: DistributionElements.h:34
Mdouble internalVariable
Definition: DistributionElements.h:36
BaseParticle * readAndCreateObject(std::istream &is)
Create a new particle, based on the information provided in a restart file.
Definition: ParticleHandler.cc:1065
bool readOptionalVariable(std::istream &is, const std::string &name, T &variable)
Reads optional variables in the restart file.
Definition: FileIOHelpers.h:82

References checkParticleForInteraction_, BaseHandler< T >::getDPMBase(), BaseBoundary::getHandler(), BaseHandler< T >::getObject(), constants::i, initialVolume_, DistributionElements::internalVariable, isActivated_, massInserted_, maxFailed_, n, numberOfParticlesInserted_, DPMBase::particleHandler, particleSizeDistributionVector_, particleToCopy_, probability_, BaseBoundary::read(), ParticleHandler::readAndCreateObject(), helpers::readOptionalVariable(), samplingInterval_, BaseParticle::setSpecies(), variableCumulativeVolumeFlowRate_, volumeFlowRate_, and volumeInserted_.

Referenced by BaseClusterInsertionBoundary::read(), ChuteInsertionBoundary::read(), CubeInsertionBoundary::read(), HopperInsertionBoundary::read(), and PolydisperseInsertionBoundary::read().

◆ reset()

void InsertionBoundary::reset ( )

resets particle property counter variables.

set all particle property counter variables to zero. reset() does not activate or deactivate the InsertionBoundary.

366 {
368  massInserted_ = 0;
369  volumeInserted_ = 0;
370 }

References massInserted_, numberOfParticlesInserted_, and volumeInserted_.

◆ setCheckParticleForInteraction()

void InsertionBoundary::setCheckParticleForInteraction ( bool  checkParticleForInteraction)

Sets the variable that checks if a particle has an interaction.

Sets the distribution type from the Distribution class.

Parameters
[in]checkParticleForInteractionboolean which determines if a particle has an interaction.
728 {
729  checkParticleForInteraction_ = checkParticleForInteraction;
730 }

References checkParticleForInteraction_.

Referenced by GranuDrum::GranuDrum().

◆ setInitialVolume()

void InsertionBoundary::setInitialVolume ( Mdouble  initialVolume)

◆ setManualInsertion()

void InsertionBoundary::setManualInsertion ( bool  isManuallyInserting)

Set the flag for a manual PSD insertion routine.

sets the isManuallyInserting_ to TRUE, resulting in a top-down class-by-class insertion routine to insert PSDs as accurate as possible.

710 {
711  isManuallyInserting_ = isManuallyInserting;
712 }

References isManuallyInserting_.

Referenced by PSDManualInsertionSelfTest::setupInitialConditions().

◆ setParticleToCopy() [1/2]

void InsertionBoundary::setParticleToCopy ( BaseParticle particleToCopy)

Sets the particle that will be inserted through the insertion boundary.

Sets the vector of pointers to particles which will be inserted by the insertion boundary.

Parameters
[in]particleToCopypointer to the particle to be inserted
436 {
437  if (particleToCopy == nullptr)
438  {
439  logger(ERROR, "Setting particleToCopy to be a null pointer?");
440  }
441  else
442  {
443  if (!particleToCopy_.empty())
444  {
445  for (auto& particleToCopy: particleToCopy_)
446  {
447  delete particleToCopy;
448  }
449  }
450  particleToCopy_.resize(1);
451  particleToCopy_[0] = particleToCopy->copy();
452  }
453 }
virtual BaseParticle * copy() const =0
Particle copy method. It calls to copy constructor of this Particle, useful for polymorphism.

References BaseParticle::copy(), ERROR, logger, and particleToCopy_.

◆ setParticleToCopy() [2/2]

void InsertionBoundary::setParticleToCopy ( std::vector< BaseParticle * >  particleToCopy)

Sets multiple different particles that will be inserted through the insertion boundary.

Sets the vector of pointers to particles which will be inserted by the insertion boundary. This is mainly used to insert particles with different intrinsic properties (such as PSD, mechanical properties, etc.)

Parameters
[in]particleToCopyvector of pointers to the particles which are to be inserted
412 {
413  if (particleToCopy.empty())
414  {
415  logger(ERROR, "Setting particleToCopy to be empty?");
416  }
417  if (!particleToCopy_.empty())
418  {
419  for (auto& ParticleToCopy: particleToCopy_)
420  {
421  delete ParticleToCopy;
422  }
423  }
424  for (int i = 0; i < particleToCopy.size(); i++)
425  {
426  particleToCopy_.resize(particleToCopy.size());
427  particleToCopy_[i] = particleToCopy[i]->copy();
428  }
429 }

References ERROR, constants::i, logger, and particleToCopy_.

Referenced by RandomClusterInsertionBoundary::set(), FixedClusterInsertionBoundary::set(), HopperInsertionBoundary::set(), ChuteInsertionBoundary::set(), and CubeInsertionBoundary::set().

◆ setPSD() [1/2]

◆ setPSD() [2/2]

void InsertionBoundary::setPSD ( std::vector< PSD psd,
std::vector< Mdouble psdRatio 
)

Sets the ranges of particle radii that may be generated from user defined PSDs.

Sets the ranges of particle radii that may be generated by different PSDs defined by the user.

Todo:
TP: Consider std::move instead of a set function. This would result in a speedup + no one has to set the PSD for insertionBoundaries again as it is set when a PSD is inserted.
686 {
687  particleSizeDistributionVector_.resize(psd.size());
689  logger.assert_always(std::accumulate(psdRatio.begin(), psdRatio.end(), 0.0) == 1.0, "Please make sure that the sum of "
690  "psdRatios adds up to unity");
691  probability_.resize(psdRatio.size());
692  probability_ = psdRatio;
693 
694 }

References logger, particleSizeDistributionVector_, and probability_.

◆ setVariableVolumeFlowRate()

void InsertionBoundary::setVariableVolumeFlowRate ( const std::vector< Mdouble > &  variableCumulativeVolumeFlowRate,
Mdouble  samplingInterval 
)

Sets a variable volume flow rate.

See also
variableCumulativeVolumeFlowRate_

Sets a variable volume flow rate taken at fixed sampling intervals; the values are cumulative; thus, we need to ensure the volume inserted before time t=n*samplingInterval is less than variableCumulativeVolumeFlowRate[n].

See also
variableCumulativeVolumeFlowRate_
660 {
661  logger.assert_debug(samplingInterval > 0, "sampling interval needs to be positive");
662  const Mdouble endTime = variableCumulativeVolumeFlowRate.size() * samplingInterval;
663  logger(INFO, "variable flowrate is defined up to %", endTime);
664  logger.assert_always(getHandler()->getDPMBase()->getTimeMax() < endTime,
665  "variable flowrate is defined up to %, but tMax is set to %", endTime,
666  getHandler()->getDPMBase()->getTimeMax());
667  variableCumulativeVolumeFlowRate_ = variableCumulativeVolumeFlowRate;
668  samplingInterval_ = samplingInterval;
669 }

References BaseBoundary::getHandler(), INFO, logger, samplingInterval_, and variableCumulativeVolumeFlowRate_.

◆ setVolumeFlowRate()

void InsertionBoundary::setVolumeFlowRate ( Mdouble  volumeFlowRate)

Sets the volume flow rate of the insertion routine.

Sets the volumetric flow rate of the insertion routine.

626 {
627  volumeFlowRate_ = volumeFlowRate;
628 }

References volumeFlowRate_.

◆ write()

void InsertionBoundary::write ( std::ostream &  os) const
overridevirtual

Writes the boundary's id_ and maxFailed_.

adds the boundary's id_ and maxFailed_ to the given ostream

Parameters
[in,out]isstream the data members are to be added to

Implements BaseBoundary.

Reimplemented in PolydisperseInsertionBoundary.

560 {
561  logger(VERBOSE, "In InsertionBoundary::write\n");
563  os << " maxFailed " << maxFailed_;
564  if (std::isfinite(volumeFlowRate_))
565  os << " volumeFlowRate " << volumeFlowRate_;
566  os << " massInserted " << massInserted_;
567  os << " volumeInserted " << volumeInserted_;
568  os << " numberOfParticlesInserted " << numberOfParticlesInserted_;
569  os << " isActivated " << isActivated_;
570  os << " psdCount " << particleSizeDistributionVector_.size();
571  for (auto& particleSizeDistribution: particleSizeDistributionVector_)
572  {
573  os << " psd " << particleSizeDistribution.getParticleSizeDistribution().size();
574  for (auto p: particleSizeDistribution.getParticleSizeDistribution())
575  {
576  os << " " << p.internalVariable
577  << " " << p.probability;
578  }
579  }
580  if (!probability_.empty())
581  {
582  os << " psdRatio";
583  for (auto& psdRatio: probability_)
584  {
585  os << " " << psdRatio;
586  }
587  }
589  {
590  os << " checkParticleForInteraction " << checkParticleForInteraction_;
591  }
592  os << " initialVolume " << initialVolume_;
593  os << " samplingInterval " << samplingInterval_;
594  os << " variableCumulativeVolumeFlowRate " << variableCumulativeVolumeFlowRate_.size();
595  for (const auto flowRate: variableCumulativeVolumeFlowRate_)
596  {
597  os << ' ' << flowRate;
598  }
599  if (!particleToCopy_.empty())
600  {
601  os << " particleToCopy " << particleToCopy_.size();
602  for (auto& particleToCopy: particleToCopy_)
603  {
604  os << " " << *particleToCopy;
605  }
606  }
607  else
608  {
609  os << " noParticleToCopy";
610  }
611 }
void write(std::ostream &os) const override=0
Adds object's id_ to given ostream NB: purely virtual function, overriding the version of BaseObject.
Definition: BaseBoundary.cc:70

References checkParticleForInteraction_, initialVolume_, isActivated_, logger, massInserted_, maxFailed_, numberOfParticlesInserted_, particleSizeDistributionVector_, particleToCopy_, probability_, samplingInterval_, variableCumulativeVolumeFlowRate_, VERBOSE, volumeFlowRate_, volumeInserted_, and BaseBoundary::write().

Referenced by BaseClusterInsertionBoundary::write(), ChuteInsertionBoundary::write(), CubeInsertionBoundary::write(), HopperInsertionBoundary::write(), and PolydisperseInsertionBoundary::write().

Member Data Documentation

◆ checkParticleForInteraction_

bool InsertionBoundary::checkParticleForInteraction_

Checks if a particle has an interaction with a wall or other particles.

Referenced by checkBoundaryBeforeTimeStep(), getCheckParticleForInteraction(), InsertionBoundary(), read(), setCheckParticleForInteraction(), and write().

◆ chosenSpecies_

int InsertionBoundary::chosenSpecies_

stores the chosen species for each timestep.

Referenced by checkBoundaryBeforeTimeStep(), generateParticle(), and InsertionBoundary().

◆ initialVolume_

Mdouble InsertionBoundary::initialVolume_

◆ isActivated_

bool InsertionBoundary::isActivated_

The InsertionBoundary is activated by default. If the InsertionBoundary is deactivated, then it introduces no particles (useful for trying to maintain a certain insertion rate).

Todo:
JMFT: This is currently not saved to .restart files.

Referenced by activate(), FixedClusterInsertionBoundary::checkBoundaryBeforeTimeStep(), checkBoundaryBeforeTimeStep(), RandomClusterInsertionBoundary::checkBoundaryBeforeTimeStep(), deactivate(), InsertionBoundary(), isActivated(), read(), and write().

◆ isManuallyInserting_

bool InsertionBoundary::isManuallyInserting_

A flag to enable a top-down class-by-class manual insertion of a PSD; default is FALSE.

Referenced by checkBoundaryBeforeTimeStep(), generateParticle(), InsertionBoundary(), and setManualInsertion().

◆ massInserted_

◆ maxFailed_

◆ numberOfParticlesInserted_

◆ particleSizeDistributionVector_

std::vector<PSD> InsertionBoundary::particleSizeDistributionVector_

Defines a particle size distribution as an object of the PSD class; if particleSizeDistributionVector_ is empty, distribution_ is used instead.

Referenced by checkBoundaryBeforeTimeStep(), generateParticle(), getPSD(), InsertionBoundary(), read(), setPSD(), and write().

◆ particleToCopy_

std::vector<BaseParticle*> InsertionBoundary::particleToCopy_

read Distribution class from file. ‍/ friend std::istream& operator>>(std::istream& is, InsertionBoundary::Distribution& type);

     \brief write Distribution class to file.
    &zwj;/

friend std::ostream& operator<<(std::ostream& os, InsertionBoundary::Distribution type);

/*!

protected:

/*!
   \brief Particle that will be inserted through the insertion boundary.

Referenced by generateParticle(), FixedClusterInsertionBoundary::generateParticle(), getParticleToCopy(), InsertionBoundary(), read(), setParticleToCopy(), write(), and ~InsertionBoundary().

◆ probability_

std::vector<Mdouble> InsertionBoundary::probability_

vector of probabilities in range [0,1] which determine the mixing ratio of partice size distributions.

Referenced by generateParticle(), InsertionBoundary(), read(), setPSD(), and write().

◆ samplingInterval_

Mdouble InsertionBoundary::samplingInterval_
See also
variableCumulativeVolumeFlowRate_

Referenced by InsertionBoundary(), insertParticle(), read(), setVariableVolumeFlowRate(), and write().

◆ variableCumulativeVolumeFlowRate_

std::vector<Mdouble> InsertionBoundary::variableCumulativeVolumeFlowRate_

Defines a variable volume flow rate, taken at fixed sampling intervals; the values are cumulative; thus, we need to ensure the volume inserted before time t=n*samplingInterval is less than variableCumulativeVolumeFlowRate[n].

By default, this vector is empty; in that case, a constant volume flow rate will be used.

See also
volumeFlowRate_.

Referenced by InsertionBoundary(), insertParticle(), read(), setVariableVolumeFlowRate(), and write().

◆ velMax_

◆ velMin_

◆ volumeFlowRate_

Mdouble InsertionBoundary::volumeFlowRate_

The inflow can be controlled by setting a volume flow rate and an initial volume thus, this ensures the volume V inserted before time t is less than V = initialVolume_+volumeFlowRate_*t

The default value is volumeFlowRate_=inf, i.e. the volume is not controlled.

See also
To set a variable flow rate instead, see variableCumulativeVolumeFlowRate_.

Referenced by getVolumeFlowRate(), InsertionBoundary(), read(), setInitialVolume(), setVolumeFlowRate(), and write().

◆ volumeInserted_


The documentation for this class was generated from the following files: