Loading [MathJax]/extensions/tex2jax.js
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
ParticleHandler.h
Go to the documentation of this file.
1 //Copyright (c) 2013-2023, 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 #ifndef PARTICLE_HANDLER_H
28 #define PARTICLE_HANDLER_H
29 
30 #include "BaseHandler.h"
31 #include "Particles/BaseParticle.h"
32 
33 class SpeciesHandler;
34 
35 class BaseSpecies;
36 
38 {
39  MIN = 0, MAX = 1
40 };
41 
47 class ParticleHandler : public BaseHandler<BaseParticle>
48 {
49 public:
54 
59 
64 
68  ~ParticleHandler() override;
69 
73  void addExistingObject(BaseParticle* P) override;
74 
78  void addObject(BaseParticle* P) override;
79 
83  void addObject(int fromProcessor, BaseParticle* P);
84 
88  void addGhostObject(int fromProcessor, int toProcessor, BaseParticle* p);
89 
93  void addGhostObject(BaseParticle* P) override;
94 
98  void removeObject(unsigned int index) override;
99 
104  void removeGhostObject(unsigned int index);
105 
109  void removeLastObject();
110 
115 
119  void computeLargestParticle();
120 
125 
131 
136 
142 
147 
152 
157 
162 
167 
169 
170  Mdouble getKineticEnergy() const;
171 
173 
174  Mdouble getMass() const;
175 
176  Vec3D getMassTimesPosition() const;
177 
178  Vec3D getCentreOfMass() const;
179 
180  Vec3D getMomentum() const;
181 
182  Vec3D getAngularMomentum() const;
183 
184  Mdouble getVolume() const;
185 
186  Mdouble getMeanRadius() const;
187 
192 
197 
202 
207 
212 
217 
222 
227 
231  //Mdouble getParticleAttributeLocal(std::function<Mdouble (BaseParticle*)> attribute, AttributeType type) const;
232 
242  template<typename T>
243  typename std::enable_if<std::is_scalar<T>::value, T>::type
244  getParticleAttributeLocal(std::function<T(BaseParticle*)> attribute, AttributeType type) const
245  {
246  T attributeParticle;
247  T attributeFinal;
248 
249  if ((*this).getSize() == 0)
250  {
251  logger(WARN, "ParticleHandler is empty: returning 0.0");
252  attributeFinal = 0;
253  return 0.0;
254  }
255  else
256  {
257  //Initialise with the first particle found
258  attributeParticle = attribute(objects_[0]);
259  attributeFinal = attributeParticle;
260 
261  //Find the final attribute
262  for (BaseParticle* particle : (*this))
263  {
264  //Obtain the attribute
265  if (!(particle->isMPIParticle() || particle->isPeriodicGhostParticle() || particle->isFixed()))
266  {
267  attributeParticle = attribute(particle);
268  }
269 
270  //Decide what to do with the magnitude of the attribute
271  switch (type)
272  {
273  case AttributeType::MIN :
274  if (attributeParticle < attributeFinal)
275  {
276  attributeFinal = attributeParticle;
277  }
278  break;
279  case AttributeType::MAX :
280  if (attributeParticle > attributeFinal)
281  {
282  attributeFinal = attributeParticle;
283  }
284  break;
285  default :
286  logger(ERROR, "Attribute type is not recognised");
287  break;
288  }
289  }
290  return attributeFinal;
291  }
292  }
293 
303  template<typename T>
304  typename std::enable_if<std::is_scalar<T>::value, T>::type
305  getParticleAttribute(std::function<T(BaseParticle*)> attribute, AttributeType type) const
306  {
307 #ifdef MERCURYDPM_USE_MPI
308  T particleAttributeLocal = getParticleAttributeLocal(attribute, type);
309  T particleAttributeGlobal;
310 
311  //Communicate local to global using the approriate rule
312  MPIContainer& communicator = MPIContainer::Instance();
313  switch(type)
314  {
315  case AttributeType::MIN :
316  communicator.allReduce(particleAttributeLocal,particleAttributeGlobal,MPI_MIN);
317  break;
318  case AttributeType::MAX :
319  communicator.allReduce(particleAttributeLocal,particleAttributeGlobal,MPI_MAX);
320  break;
321  default :
322  logger(ERROR,"Attribute type is not recognised");
323  break;
324  }
325  return particleAttributeGlobal;
326 #else
327  return getParticleAttributeLocal(attribute, type);
328 #endif
329  }
330 
335 
339  void clear() override;
340 
345  unsigned int getNumberOfFixedParticles() const;
346 
350  unsigned int getNumberOfUnfixedParticles() const;
351 
355  static BaseParticle* createObject(const std::string& type);
356 
360  BaseParticle* readAndCreateObject(std::istream& is);
361 
362 // /*!
363 // * \brief Create a new particle, based on the information from old-style restart data.
364 // */
365 // void readAndCreateOldObject(std::istream &is, const std::string& type);
366 
370  void readAndAddObject(std::istream& is) override;
371 
372  void write(std::ostream& os) const;
373 
377  void checkExtrema(BaseParticle* P);
378 
383 
387  void computeAllMasses(unsigned int indSpecies);
388 
392  void computeAllMasses();
393 
398  void addedFixedParticle();
399 
404  void removedFixedParticle();
405 
409  std::string getName() const override;
410 
414  unsigned int getNumberOfRealObjects() const;
415 
420  unsigned int getNumberOfObjects() const override;
421 
425  unsigned int getNumberOfFixedObjectsLocal() const;
426 
430  unsigned int getNumberOfFixedObjects() const;
431 
435  unsigned int getNumberOfRealObjectsLocal() const;
436 
437  void actionsAfterTimeStep();
438 
439  double getLiquidFilmVolume() const;
440 
441  void saveNumberPSDtoCSV(std::string csvFileName, std::vector<double> diameterBins = {});
442 
443 private:
444 
446 
448 
449  Mdouble getMassLocal() const;
450 
451  Mdouble getVolumeLocal() const;
452 
454 
455  Mdouble getSumRadiusLocal() const;
456 
467 
472 
476  unsigned int NFixedParticles_;
477 };
478 
479 #endif
480 
double Mdouble
Definition: GeneralDefine.h:34
Logger< MERCURYDPM_LOGLEVEL > logger("MercuryKernel")
Definition of different loggers with certain modules. A user can define its own custom logger here.
@ WARN
@ ERROR
AttributeType
Definition: ParticleHandler.h:38
@ MIN
Definition: ParticleHandler.h:39
@ MAX
Definition: ParticleHandler.h:39
Container to store the pointers to all objects that one creates in a simulation.
Definition: BaseHandler.h:51
std::vector< BaseParticle * > objects_
The actual list of Object pointers.
Definition: BaseHandler.h:302
Definition: BaseParticle.h:54
BaseSpecies is the class from which all other species are derived.
Definition: BaseSpecies.h:50
This class contains all information and functions required for communication between processors.
Definition: MpiContainer.h:130
static MPIContainer & Instance()
fetch the instance to be used for communication
Definition: MpiContainer.h:134
Container to store all BaseParticle.
Definition: ParticleHandler.h:48
BaseParticle * getSmallestParticleLocal() const
Gets a pointer to the smallest BaseParticle (by interactionRadius) in this ParticleHandler of the loc...
Definition: ParticleHandler.cc:502
BaseParticle * getHighestPositionComponentParticle(int i) const
Gets a pointer to the particle with the highest coordinates in direction i in this ParticleHandler.
Definition: ParticleHandler.cc:892
~ParticleHandler() override
Destructor, it destructs the ParticleHandler and all BaseParticle it contains.
Definition: ParticleHandler.cc:107
std::enable_if< std::is_scalar< T >::value, T >::type getParticleAttribute(std::function< T(BaseParticle *)> attribute, AttributeType type) const
Computes an attribute type (min/max/..) of a particle attribute(position/velocity) in the global doma...
Definition: ParticleHandler.h:305
Mdouble getLargestInteractionRadiusLocal() const
Returns the largest interaction radius of the current domain.
Definition: ParticleHandler.cc:757
std::enable_if< std::is_scalar< T >::value, T >::type getParticleAttributeLocal(std::function< T(BaseParticle *)> attribute, AttributeType type) const
Computes an attribute type (min/max/..) of a particle attribute (position/velocity) in a local domain...
Definition: ParticleHandler.h:244
BaseParticle * getFastestParticle() const
Definition: ParticleHandler.cc:710
Mdouble getRotationalEnergy() const
Definition: ParticleHandler.cc:586
void clear() override
Empties the whole ParticleHandler by removing all BaseParticle.
Definition: ParticleHandler.cc:977
void removeGhostObject(unsigned int index)
Removes a BaseParticle from the ParticleHandler without a global check, this is only to be done for m...
Definition: ParticleHandler.cc:422
BaseParticle * getFastestParticleLocal() const
Gets a pointer to the fastest BaseParticle in this ParticleHandler.
Definition: ParticleHandler.cc:687
Vec3D getMomentum() const
Definition: ParticleHandler.cc:666
double getLiquidFilmVolume() const
Definition: ParticleHandler.cc:1380
BaseParticle * getHighestVelocityComponentParticle(int i) const
Gets a pointer to the particle with the highest velocity in direction i in this ParticleHandler.
Definition: ParticleHandler.cc:964
void readAndAddObject(std::istream &is) override
Definition: ParticleHandler.cc:1094
Mdouble getHighestPositionX() const
Function returns the highest position in the x-direction.
Definition: ParticleHandler.cc:995
static BaseParticle * createObject(const std::string &type)
Reads BaseParticle into the ParticleHandler from restart data.
Definition: ParticleHandler.cc:1025
Mdouble getKineticEnergyLocal() const
Definition: ParticleHandler.cc:544
void computeSmallestParticle()
Computes the smallest particle (by interaction radius) and sets it in smallestParticle_.
Definition: ParticleHandler.cc:451
Mdouble getSmallestInteractionRadiusLocal() const
Returns the smallest interaction radius of the current domain.
Definition: ParticleHandler.cc:721
std::string getName() const override
Returns the name of the handler, namely the string "ParticleHandler".
Definition: ParticleHandler.cc:1247
Vec3D getMassTimesPositionLocal() const
Definition: ParticleHandler.cc:627
Vec3D getMassTimesPosition() const
Definition: ParticleHandler.cc:636
Mdouble getMassLocal() const
Definition: ParticleHandler.cc:602
Vec3D getCentreOfMass() const
Definition: ParticleHandler.cc:654
BaseParticle * getLowestVelocityComponentParticleLocal(int i) const
Gets a pointer to the particle with the lowest velocity in direction i in this ParticleHandler.
Definition: ParticleHandler.cc:905
Vec3D getAngularMomentum() const
Definition: ParticleHandler.cc:675
unsigned int getNumberOfRealObjectsLocal() const
Returns the number of real objects on a local domain. MPI particles and periodic particles are neglec...
Definition: ParticleHandler.cc:1283
BaseParticle * getLowestPositionComponentParticleLocal(int i) const
Gets a pointer to the particle with the lowest coordinates in direction i in this ParticleHandler.
Definition: ParticleHandler.cc:829
void actionsAfterTimeStep()
Definition: ParticleHandler.cc:1372
BaseParticle * getLowestVelocityComponentParticle(int i) const
Gets a pointer to the particle with the lowest velocity in direction i in this ParticleHandler.
Definition: ParticleHandler.cc:928
void addExistingObject(BaseParticle *P) override
Adds a BaseParticle to the ParticleHandler.
Definition: ParticleHandler.cc:131
void addObject(BaseParticle *P) override
Adds a BaseParticle to the ParticleHandler.
Definition: ParticleHandler.cc:171
ParticleHandler & operator=(const ParticleHandler &rhs)
Assignment operator.
Definition: ParticleHandler.cc:85
Mdouble getLargestInteractionRadius() const
Returns the largest interaction radius.
Definition: ParticleHandler.cc:772
Mdouble getSmallestInteractionRadius() const
Returns the smallest interaction radius.
Definition: ParticleHandler.cc:736
unsigned int getNumberOfFixedObjects() const
Computes the number of fixed particles in the whole simulation.
Definition: ParticleHandler.cc:1358
void write(std::ostream &os) const
BaseParticle * readAndCreateObject(std::istream &is)
Create a new particle, based on the information provided in a restart file.
Definition: ParticleHandler.cc:1065
unsigned int getNumberOfRealObjects() const
Returns the number of real objects (on all processors)
Definition: ParticleHandler.cc:1304
Mdouble getVolume() const
Definition: ParticleHandler.cc:1263
BaseParticle * smallestParticle_
A pointer to the smallest BaseParticle (by interactionRadius) in this ParticleHandler.
Definition: ParticleHandler.h:471
void checkExtrema(BaseParticle *P)
Checks if the extrema of this ParticleHandler needs updating.
Definition: ParticleHandler.cc:1173
void removeObject(unsigned int index) override
Removes a BaseParticle from the ParticleHandler.
Definition: ParticleHandler.cc:394
unsigned int NFixedParticles_
Number of fixed particles.
Definition: ParticleHandler.h:476
BaseParticle * getLargestParticleLocal() const
Gets a pointer to the largest BaseParticle (by interactionRadius) in the ParticleHandler of the local...
Definition: ParticleHandler.cc:526
unsigned int getNumberOfObjects() const override
Returns the number of objects in the container. In parallel code this practice is forbidden to avoid ...
Definition: ParticleHandler.cc:1325
BaseParticle * largestParticle_
A pointer to the largest BaseParticle (by interactionRadius) in this ParticleHandler.
Definition: ParticleHandler.h:466
unsigned int getNumberOfFixedObjectsLocal() const
Computes the number of Fixed particles on a local domain.
Definition: ParticleHandler.cc:1342
BaseParticle * getHighestPositionComponentParticleLocal(int i) const
Gets a pointer to the particle with the highest coordinates in direction i in this ParticleHandler.
Definition: ParticleHandler.cc:868
Mdouble getVolumeLocal() const
Definition: ParticleHandler.cc:1252
void saveNumberPSDtoCSV(std::string csvFileName, std::vector< double > diameterBins={})
Definition: ParticleHandler.cc:1390
void checkExtremaOnDelete(BaseParticle *P)
Checks if the extrema of this ParticleHandler needs updating when a particle is deleted.
Definition: ParticleHandler.cc:1199
Mdouble getKineticEnergy() const
Definition: ParticleHandler.cc:557
void addGhostObject(int fromProcessor, int toProcessor, BaseParticle *p)
Adds a ghost particle located at fromProcessor to toProcessor.
Definition: ParticleHandler.cc:287
Mdouble getSumRadiusLocal() const
Definition: ParticleHandler.cc:791
BaseParticle * getLowestPositionComponentParticle(int i) const
Gets a pointer to the particle with the lowest coordinates in direction i in this ParticleHandler.
Definition: ParticleHandler.cc:855
void removeLastObject()
Removes the last BaseParticle from the ParticleHandler.
Definition: ParticleHandler.cc:441
BaseParticle * getHighestVelocityComponentParticleLocal(int i) const
Gets a pointer to the particle with the highest velocity in direction i in this ParticleHandler.
Definition: ParticleHandler.cc:941
Mdouble getRotationalEnergyLocal() const
Definition: ParticleHandler.cc:573
unsigned int getNumberOfFixedParticles() const
Gets the number of particles that are fixed.
Definition: ParticleHandler.cc:1009
void computeAllMasses()
Computes the mass for all BaseParticle in this ParticleHandler.
Definition: ParticleHandler.cc:1226
BaseParticle * getLargestParticle() const
Returns the pointer of the largest particle in the particle handler. When mercury is running in paral...
Definition: ParticleHandler.cc:534
BaseParticle * getSmallestParticle() const
Gets a pointer to the smallest BaseParticle (by interactionRadius) in this ParticleHandler of the loc...
Definition: ParticleHandler.cc:511
ParticleHandler()
Default constructor, it creates an empty ParticleHandler.
Definition: ParticleHandler.cc:50
unsigned int getNumberOfUnfixedParticles() const
Gets the number of particles that are not fixed.
Definition: ParticleHandler.cc:1017
void addedFixedParticle()
Increment of the number of fixed particles.
Definition: ParticleHandler.cc:1234
void computeLargestParticle()
Computes the largest particle (by interaction radius) and sets it in largestParticle_.
Definition: ParticleHandler.cc:475
void removedFixedParticle()
Decrement of the number of fixed particles.
Definition: ParticleHandler.cc:1239
Mdouble getMeanRadius() const
Definition: ParticleHandler.cc:804
Mdouble getMass() const
Definition: ParticleHandler.cc:611
Container to store all ParticleSpecies.
Definition: SpeciesHandler.h:37
Definition: Vector.h:51
double P
Uniform pressure.
Definition: TwenteMeshGluing.cpp:73
const std::complex< Mdouble > i
Definition: ExtendedMath.h:51