MercuryDPM  Alpha
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
ParticleHandler.cc
Go to the documentation of this file.
1 //Copyright (c) 2013-2014, The MercuryDPM Developers Team. All rights reserved.
2 //For the list of developers, see <http://www.MercuryDPM.org/Team>.
3 //
4 //Redistribution and use in source and binary forms, with or without
5 //modification, are permitted provided that the following conditions are met:
6 // * Redistributions of source code must retain the above copyright
7 // notice, this list of conditions and the following disclaimer.
8 // * Redistributions in binary form must reproduce the above copyright
9 // notice, this list of conditions and the following disclaimer in the
10 // documentation and/or other materials provided with the distribution.
11 // * Neither the name MercuryDPM nor the
12 // names of its contributors may be used to endorse or promote products
13 // derived from this software without specific prior written permission.
14 //
15 //THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND
16 //ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
17 //WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
18 //DISCLAIMED. IN NO EVENT SHALL THE MERCURYDPM DEVELOPERS TEAM BE LIABLE FOR ANY
19 //DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
20 //(INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
21 //LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
22 //ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
23 //(INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
24 //SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
25 
26 
27 #include <limits>
28 #include <Math/Helpers.h>
30 #include "ParticleHandler.h"
31 #include "DPMBase.h"
32 #include "SpeciesHandler.h"
35 #include "Particles/BaseParticle.h"
36 
37 
43 {
44  largestParticle_ = nullptr;
45  smallestParticle_ = nullptr;
46  logger(DEBUG, "ParticleHandler::ParticleHandler() finished");
47 }
48 
57 {
58  clear();
59  setDPMBase(PH.getDPMBase());
60  largestParticle_ = nullptr;
61  smallestParticle_ = nullptr;
63  if (objects_.size() != 0)
64  {
67  }
68  logger(DEBUG, "ParticleHandler::ParticleHandler(const ParticleHandler &PH) finished");
69 }
70 
78 {
79  if (this != &rhs)
80  {
81  clear();
82  largestParticle_ = nullptr;
83  smallestParticle_ = nullptr;
85  if (objects_.size() != 0)
86  {
89  }
90  }
91  logger(DEBUG, "ParticleHandler::operator = (const ParticleHandler& rhs) finished");
92  return *this;
93 }
94 
100 {
101  //First reset the pointers, such that they are not checked twice when removing particles
102  largestParticle_ = nullptr;
103  smallestParticle_ = nullptr;
104  logger(DEBUG, "ParticleHandler::~ParticleHandler() finished");
105 }
106 
121 {
122  if (P->getSpecies() == nullptr)
123  {
124  logger(WARN, "WARNING: The particle with ID % that is added in "
125  "ParticleHandler::addObject does not have a species yet."
126  "Please make sure that you have "
127  "set the species somewhere in the driver code.", P->getId());
128  }
129  //Puts the particle in the Particle list
131  if (getDPMBase() != nullptr)
132  {
133  //This places the particle in this grid
135  //This computes where the particle currently is in the grid
137  }
138  //set the particleHandler pointer
139  P->setHandler(this);
140  //compute mass of the particle
141  P->getSpecies()->computeMass(P) ;
142  //Check if this particle has new extrema
143  checkExtrema(P);
144 }
145 
153 void ParticleHandler::removeObject(unsigned const int index)
154 {
155 #ifdef CONTACT_LIST_HGRID
156  getDPMBase()->getPossibleContactList().remove_ParticlePosibleContacts(getObject(id));
157 #endif
160 }
161 
167 {
168 #ifdef CONTACT_LIST_HGRID
169  getDPMBase()->getPossibleContactList().remove_ParticlePosibleContacts(getLastObject());
170 #endif
173 }
174 
176 {
177  if (getNumberOfObjects() == 0)
178  {
179  logger(WARN, "No particles, so cannot compute the smallest particle.");
180  return;
181  }
182  Mdouble min = std::numeric_limits<Mdouble>::max();
183  smallestParticle_ = nullptr;
184  for (BaseParticle* const particle : objects_)
185  {
186  if (particle->getInteractionRadius() < min)
187  {
188  min = particle->getInteractionRadius();
189  smallestParticle_ = particle;
190  }
191  }
192 }
193 
195 {
196  if (getNumberOfObjects() == 0)
197  {
198  logger(WARN, "No particles, so cannot compute the largest particle.");
199  return;
200  }
201  Mdouble max = -std::numeric_limits<Mdouble>::max();
202  largestParticle_ = nullptr;
203  for (BaseParticle* const particle : objects_)
204  {
205  if (particle->getInteractionRadius() > max)
206  {
207  max = particle->getInteractionRadius();
208  largestParticle_ = particle;
209  }
210  }
211 }
212 
218 {
219  if (smallestParticle_ == nullptr)
220  {
221  logger(WARN, "No particles to return in getSmallestParticle()");
222  }
223  return smallestParticle_;
224 }
225 
231 {
232  if (largestParticle_ == nullptr)
233  {
234  logger(WARN, "No particles to set get_LargestParticle()");
235  }
236  return largestParticle_;
237 }
238 
243 {
244  if (getNumberOfObjects() == 0)
245  {
246  logger(WARN, "No particles to set getFastestParticle()" );
247  return nullptr;
248  }
249  BaseParticle* p = nullptr;
250  Mdouble maxSpeed = -std::numeric_limits<Mdouble>::max();
251  for (BaseParticle* const pLoop : objects_)
252  {
253  if ((pLoop->getVelocity().getLength()) > maxSpeed)
254  {
255  maxSpeed = pLoop->getVelocity().getLength();
256  p = pLoop;
257  }
258  }
259  return p;
260 }
261 
266  Mdouble sumRadius = 0;
267  for (BaseParticle* const p : objects_) {
268  sumRadius += p->getRadius();
269  }
270  return sumRadius/getNumberOfObjects();
271 }
272 
279 {
280  if (getNumberOfObjects() == 0)
281  {
282  logger(WARN, "No getLowestPositionComponentParticle(const int i) since there are no particles.");
283  return nullptr;
284  }
285  BaseParticle* p = nullptr;
286  Mdouble min = std::numeric_limits<Mdouble>::max();
287  for (BaseParticle* const pLoop : objects_)
288  {
289  if (pLoop->getPosition().getComponent(i) < min)
290  {
291  min = pLoop->getPosition().getComponent(i);
292  p = pLoop;
293  }
294  }
295  return p;
296 }
297 
304 {
305  if (getNumberOfObjects() == 0)
306  {
307  logger(WARN, "No getHighestPositionComponentParticle(const int i) since there are no particles.");
308  return nullptr;
309  }
310  BaseParticle* p = nullptr;
311  Mdouble max = -std::numeric_limits<Mdouble>::max();
312  for (BaseParticle* const pLoop : objects_)
313  {
314  if (pLoop->getPosition().getComponent(i) > max)
315  {
316  max = pLoop->getPosition().getComponent(i);
317  p = pLoop;
318  }
319  }
320 
321  return p;
322 }
323 
330 {
331  if (getNumberOfObjects() == 0)
332  {
333  logger(WARN, "No getLowestVelocityComponentParticle(const int i) since there are no particles");
334  return nullptr;
335  }
336  BaseParticle* p = nullptr;
337  Mdouble min = std::numeric_limits<Mdouble>::max();
338  for (BaseParticle* const pLoop : objects_)
339  {
340  if (pLoop->getVelocity().getComponent(i) < min)
341  {
342  min = pLoop->getVelocity().getComponent(i);
343  p = pLoop;
344  }
345  }
346  return p;
347 }
348 
355 {
356  if (!getNumberOfObjects())
357  {
358  logger(WARN, "No getHighestVelocityComponentParticle(const int i) since there are no particles");
359  return nullptr;
360  }
361  BaseParticle* p = nullptr;
362  Mdouble max = -std::numeric_limits<Mdouble>::max();
363  for (BaseParticle* const pLoop : objects_)
364  {
365  if (pLoop->getVelocity().getComponent(i) > max)
366  {
367  max = pLoop->getVelocity().getComponent(i);
368  p = pLoop;
369  }
370  }
371  return p;
372 }
373 
378 {
379  if (getNumberOfObjects() == 0)
380  {
381  logger(WARN, "No particles to set getLightestParticle()");
382  return nullptr;
383  }
384  BaseParticle* p = nullptr;
385  Mdouble minMass = std::numeric_limits<Mdouble>::max();
386  for (BaseParticle* const pLoop : objects_)
387  {
388  if (pLoop->getMass() < minMass)
389  {
390  minMass = pLoop->getMass();
391  p = pLoop;
392  }
393  }
394  return p;
395 }
396 
403 {
404  smallestParticle_ = nullptr;
405  largestParticle_ = nullptr;
407 }
408 
410 {
411  if (type == "BaseParticle"||type == "BP"||isdigit(type[0]))
412  {
413  return new BaseParticle;
414  }
415  else if (type == "LiquidFilmParticle")
416  {
417  return new LiquidFilmParticle;
418  }
419  else
420  {
421  logger(ERROR, "Particle type % not understood in restart file. Particle has not been read.", type);
422  return nullptr;
423  }
424 }
425 
429 void ParticleHandler::readObject(std::istream& is)
430 {
431  std::string type;
432  is >> type;
433  logger.log(Log::DEBUG, "ParticleHandler::readObject(is): reading type %.", type);
434  if (type == "BaseParticle")
435  {
436  BaseParticle baseParticle;
437  is >> baseParticle;
438  //set the species pointer correctly
439  baseParticle.setSpecies(getDPMBase()->speciesHandler.getObject(baseParticle.getIndSpecies()));
440  copyAndAddObject(baseParticle);
441  setId(getLastObject(),baseParticle.getId());
442  //getLastObject()->setId(baseParticle.getId()); //to ensure old id
444  }
445  else if (type == "LiquidFilmParticle")
446  {
447  LiquidFilmParticle particle;
448  is >> particle;
449  //set the species pointer correctly
450  particle.setSpecies(getDPMBase()->speciesHandler.getObject(particle.getIndSpecies()));
451  copyAndAddObject(particle);
452  setId(getLastObject(),particle.getId());
453  //getLastObject()->setId(particle.getId()); //to ensure old id
455  }
456  else if (type == "ThermalParticle")
457  {
458  ThermalParticle particle;
459  is >> particle;
460  //set the species pointer correctly
461  particle.setSpecies(getDPMBase()->speciesHandler.getObject(particle.getIndSpecies()));
462  copyAndAddObject(particle);
463  setId(getLastObject(),particle.getId());
464  //getLastObject()->setId(particle.getId()); //to ensure old id
466  }
468  else if (type == "BP")
469  {
470  BaseParticle baseParticle;
471  baseParticle.setHandler(this);//to be able to set the species
472  baseParticle.oldRead(is);
473  //set the species pointer correctly
474  baseParticle.setSpecies(getDPMBase()->speciesHandler.getObject(baseParticle.getIndSpecies()));
475  copyAndAddObject(baseParticle);
476  setId(getLastObject(),baseParticle.getId());
477  //getLastObject()->setId(baseParticle.getId()); //to ensure old id
478  }
479  else if (isdigit(type[0]))
480  {
481  readOldObject(type, is);
482  }
483  else
484  {
485  logger(ERROR, "Particle type % not understood in restart file: You need to add this particle to ParticleHandler::readObject.", type);
486  return;
487  }
488 }
489 
499 void ParticleHandler::readOldObject(std::string type, std::istream& is)
500 {
501  //read in next line
502  std::stringstream line(std::stringstream::in | std::stringstream::out);
504  logger(VERBOSE, line.str());
505  //std::cout << line.str() << std::endl;
506 
507  BaseParticle particle;
508 
509  //Declare all properties of the particle
510  unsigned int indSpecies;
511  Mdouble radius, inverseMass, inverseInertia;
512  Vec3D position, velocity, orientation, angularVelocity;
513 
514  //Read all values
515  position.X = atof(type.c_str());
516 
517  line >> position.Y >> position.Z >> velocity >> radius >> orientation >> angularVelocity >> inverseMass >> inverseInertia >> indSpecies;
518 
519  //Put the values in the particle
520  particle.setSpecies(getDPMBase()->speciesHandler.getObject(indSpecies));
521  particle.setPosition(position);
522  particle.setVelocity(velocity);
523  particle.setRadius(radius);
524  particle.setOrientation(orientation);
525  particle.setAngularVelocity(angularVelocity);
526  if (inverseMass == 0.0)
527  particle.fixParticle();
528  else
529  {
530  particle.setInertia(1./inverseInertia);
531  }
532 
533  //Put the particle in the handler
534  copyAndAddObject(particle);
535 }
536 
537 void ParticleHandler::write(std::ostream& os) const
538 {
539  os << "Particles " << getNumberOfObjects() << std::endl;
540  for (BaseParticle* it : *this)
541  os << (*it) << std::endl;
542 }
543 
549 {
550  if (P == largestParticle_)
551  {
552  //if the properties of the largest particle changes
554  }
556  {
557  largestParticle_ = P;
558  }
559 
560  if (P == smallestParticle_)
561  {
562  //if the properties of the smallest particle changes
564  }
566  {
567  smallestParticle_ = P;
568  }
569 }
570 
575 {
576  if (P == largestParticle_)
577  {
579  }
580  if (P == smallestParticle_)
581  {
583  }
584 }
585 
590 void ParticleHandler::computeAllMasses(unsigned int indSpecies)
591 {
592  for (BaseParticle* particle : objects_)
593  {
594  if (particle->getIndSpecies() == indSpecies)
595  {
596  particle->getSpecies()->computeMass(particle);
597  }
598  }
599 }
600 
602 {
603  for (BaseParticle* particle : objects_)
604  {
605  particle->getSpecies()->computeMass(particle);
606  }
607 }
608 
612 std::string ParticleHandler::getName() const
613 {
614  return "ParticleHandler";
615 }
616 
624 {
625  static unsigned fileCounter = 0;
626  const std::string fileName = getDPMBase()->getName() + "Particle_" +
627  std::to_string(fileCounter++) + ".vtu";
628  logger(INFO, "writing % (t=%, N=%)", fileName, getDPMBase()->getTime(), this->getNumberOfObjects());
629 
630  std::fstream file;
631  file.open(fileName.c_str(), std::ios::out);
632  if (file.fail())
633  {
635  logger(WARN, "Error in writeToFile: file % could not be opened", fileName);
636  }
637  //
638  file << "<?xml version=\"1.0\"?>\n\n";
639  file << "<VTKFile type=\"UnstructuredGrid\" version=\"0.1\" byte_order=\"LittleEndian\">\n";
640  file << "<UnstructuredGrid>\n";
641  file << "<Piece NumberOfPoints=\"" << this->getNumberOfObjects() << "\" NumberOfCells=\"" << 0 << "\">\n";
642  file << "<Points>\n";
643  file << " <DataArray type=\"Float32\" Name=\"Position\" NumberOfComponents=\"3\" format=\"ascii\">\n";
644  for (const auto& p: *this)
645  {
646  file << '\t' << p->getPosition() << '\n';
647  }
648  file << " </DataArray>\n";
649  file << "</Points>\n";
650  file << "<PointData Vectors=\"vector\">\n";
651  file << " <DataArray type=\"Float32\" Name=\"Velocity\" NumberOfComponents=\"3\" format=\"ascii\">\n";
652  // Add velocity
653  for (const auto& p: *this)
654  {
655  file << '\t' << p->getVelocity() << '\n';
656  }
657  file << " </DataArray>\n";
658  file << " <DataArray type=\"Float32\" Name=\"Radius\" format=\"ascii\">\n";
659  // Add radius
660  for (const auto& p: *this)
661  {
662  file << '\t' << p->getRadius() << '\n';
663  }
664  file << " </DataArray>\n";
665  file << " <DataArray type=\"Float32\" Name=\"SpeciesType\" format=\"ascii\">\n";
666  // Add species type
667  for (const auto& p: *this)
668  {
669  file << '\t' << p->getIndSpecies() << '\n';
670  }
671  file << " </DataArray>\n";
672  file << " <DataArray type=\"Float32\" Name=\"Info\" format=\"ascii\">\n";
673  // Add species type
674  for (const auto& p: *this)
675  {
676  file << '\t' << getDPMBase()->getInfo(*p) << '\n';
677  }
678  file << " </DataArray>\n";
679  //check if this type of Interaction has extra fields
680  if (getLastObject()!= nullptr)
681  {
682  for (unsigned i = 0; i<getLastObject()->getNumberOfFieldsVTK(); i++)
683  {
684  file << " <DataArray type=\"" << getLastObject()->getTypeVTK(i) << "\" Name=\"" << getLastObject()->getNameVTK(i) << "\" format=\"ascii\">\n";
685  // Add species type
686  for (const auto &p: *this)
687  {
688  for (auto f : p->getFieldVTK(i))
689  file << '\t' << f << '\n';
690  }
691  file << " </DataArray>\n";
692  }
693  }
694  file << "</PointData>\n";
695  file << "<Cells>\n";
696  file << " <DataArray type=\"Int32\" Name=\"connectivity\" format=\"ascii\">\n";
697  file << " </DataArray>\n";
698  file << " <DataArray type=\"Int32\" Name=\"offsets\" format=\"ascii\">\n";
699  file << " </DataArray>\n";
700  file << " <DataArray type=\"UInt8\" Name=\"types\" format=\"ascii\">\n";
701  file << " </DataArray>\n";
702  file << "</Cells>\n";
703  file << "</Piece>\n";
704  file << "</UnstructuredGrid>\n";
705  file << "</VTKFile>\n";
706  file.close();
707 }
BaseParticle * getLargestParticle() const
Gets a pointer to the largest BaseParticle (by interactionRadius) in this ParticleHandler.
unsigned int getId() const
Returns the unique identifier of any particular object.
Definition: BaseObject.cc:116
virtual void hGridUpdateParticle(BaseParticle *obj UNUSED)
no implementation but can be overidden in its derived classes.
Definition: DPMBase.cc:757
void checkExtremaOnDelete(BaseParticle *P)
Checks if the extrema of this ParticleHandler needs updating when a particle is deleted.
void write(std::ostream &os) const
virtual void hGridInsertParticle(BaseParticle *obj UNUSED)
no implementation but can be overidden in its derived classes.
Definition: DPMBase.cc:750
Mdouble X
the vector components
Definition: Vector.h:52
void writeVTK() const
Writes all particles into a vtk file format (unstructured grid), consisting of particle positions...
virtual void removeObject(unsigned const int index)
Removes a BaseParticle from the ParticleHandler.
Logger< MERCURY_LOGLEVEL > logger("MercuryKernel")
void checkExtrema(BaseParticle *P)
Checks if the extrema of this ParticleHandler needs updating.
void setDPMBase(DPMBase *DPMBase)
Sets the problem that is solved using this handler.
BaseParticle * getLowestVelocityComponentParticle(const int i) const
Gets a pointer to the particle with the lowest velocity in direction i in this ParticleHandler.
const ParticleSpecies * getSpecies() const
Returns a pointer to the species of this BaseInteractable.
std::string getName() const
Returns the name of the handler, namely the string "ParticleHandler".
void setHandler(ParticleHandler *handler)
Sets the pointer to the particle's ParticleHandler.
double Mdouble
BaseParticle * getLightestParticle() const
Gets a pointer to the lightest BaseParticle (by mass) in this ParticleHandler.
void computeMass(BaseParticle *p) const
Compute Particle mass function, which required a reference to the Species vector. It computes the Par...
void computeAllMasses()
Computes the mass for all BaseParticle in this ParticleHandler.
void readObject(std::istream &is)
Reads BaseParticle into the ParticleHandler from restart data.
void setSpecies(const ParticleSpecies *species)
BaseParticle * smallestParticle_
A pointer to the smallest BaseParticle (by interactionRadius) in this ParticleHandler.
void removeLastObject()
Removes the last Object from the BaseHandler.
Definition: BaseHandler.h:375
void clear()
Empties the whole ParticleHandler by removing all BaseParticle.
BaseParticle * getSmallestParticle() const
Gets a pointer to the smallest BaseParticle (by interactionRadius) in this ParticleHandler.
void getLineFromStringStream(std::istream &in, std::stringstream &out)
Reads a line from one stringstream into another, and prepares the latter for reading in...
Definition: Helpers.cc:396
virtual void removeObject(unsigned const int id)
Removes an Object from the BaseHandler.
Definition: BaseHandler.h:333
void removeLastObject()
Removes the last BaseParticle from the ParticleHandler.
BaseParticle * getObject(const unsigned int id)
Gets a pointer to the Object at the specified index in the BaseHandler.
virtual void addObject(BaseParticle *P)
Adds a BaseParticle to the ParticleHandler.
void readOldObject(std::string type, std::istream &is)
Reads BaseParticle into the ParticleHandler from old-style restart data.
ParticleHandler()
Default constructor, it creates an empty ParticleHandler.
BaseParticle * largestParticle_
A pointer to the largest BaseParticle (by interactionRadius) in this ParticleHandler.
~ParticleHandler()
Destructor, it destructs the ParticleHandler and all BaseParticle it contains.
BaseParticle * getHighestVelocityComponentParticle(const int i) const
Gets a pointer to the particle with the highest velocity in direction i in this ParticleHandler.
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.
std::vector< BaseParticle * > objects_
The actual list of Object pointers.
Definition: BaseHandler.h:216
Container to store the pointers to all objects that one creates in a simulation.
Definition: BaseHandler.h:51
unsigned int getNumberOfObjects() const
Gets the number of Object in this BaseHandler.
Container to store all BaseParticle.
Mdouble Y
Definition: Vector.h:52
std::string to_string(const T &n)
Definition: Helpers.h:203
virtual void hGridRemoveParticle(BaseParticle *obj UNUSED)
no implementation but can be overidden in its derived classes.
Definition: DPMBase.cc:764
virtual std::string getTypeVTK(unsigned i) const
virtual double getInfo(const BaseParticle &P) const
A virtual method that allows the user to overrride and set what is written into the info column in th...
Definition: DPMBase.cc:703
void computeSmallestParticle()
Computes the smallest particle (by interaction radius) and sets it in smallestParticle_.
BaseParticle * getHighestPositionComponentParticle(const int i) const
Gets a pointer to the particle with the highest coordinates in direction i in this ParticleHandler...
virtual MERCURY_DEPRECATED void oldRead(std::istream &is)
deprecated version of the read function.
void setId(BaseParticle *object, unsigned int id)
This function sets the id and ensures that nextId is a bigger value than id.
Definition: BaseHandler.h:189
void computeLargestParticle()
Computes the largest particle (by interaction radius) and sets it in largestParticle_.
Mdouble getMeanRadius() const
Gets a pointer to the fastest BaseParticle in this ParticleHandler.
virtual void addObject(T *object)
Adds a new Object to the BaseHandler.
Definition: BaseHandler.h:313
virtual std::string getNameVTK(unsigned i) const
unsigned int getIndSpecies() const
Returns the index of the Species of this BaseInteractable.
BaseParticle * getLastObject()
Gets a pointer to the last Object in this BaseHandler.
void copyContentsFromOtherHandler(const BaseHandler< BaseParticle > &BH)
Function that copies the contents (vector of pointers, maxObject_, nextId_, DPMBase_) from one handle...
Implementation of a 3D vector (by Vitaliy).
Definition: Vector.h:45
static BaseParticle * getNewObject(const std::string &type)
DPMBase * getDPMBase()
Gets the problem that is solved using this handler.
Definition: BaseHandler.h:543
virtual unsigned getNumberOfFieldsVTK() const
BaseParticle * getLowestPositionComponentParticle(const int i) const
Gets a pointer to the particle with the lowest coordinates in direction i in this ParticleHandler...
Mdouble Z
Definition: Vector.h:52
Mdouble getInteractionRadius() const
Returns the particle's interaction radius, which might be different from radius_ (e.g., when dealing with wet particles)
const std::string & getName() const
Returns the name of the file. Does not allow to change it though.
Definition: Files.cc:132
ParticleHandler operator=(const ParticleHandler &rhs)
Assignment operator.
void clear()
Empties the whole BaseHandler by removing all Objects and setting all other variables to 0...
Definition: BaseHandler.h:392
BaseParticle * getFastestParticle() const
Gets a pointer to the fastest BaseParticle in this ParticleHandler.