MercuryDPM  Alpha
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
BaseParticle Class Reference

#include <BaseParticle.h>

+ Inheritance diagram for BaseParticle:

Public Member Functions

 BaseParticle ()
 Basic Particle constructor, creates an Particle at (0,0,0) with radius, mass and inertia equal to 1. More...
 
 BaseParticle (const BaseParticle &p)
 Particle copy constructor, which accepts as input a reference to a Particle. It creates a copy of this Particle and all it's information. Usually it is better to use the copy() function for polymorfism. More...
 
virtual ~BaseParticle ()
 Particle destructor, needs to be implemented and checked if it removes tangential spring information. More...
 
virtual BaseParticlecopy () const
 Particle copy method. It calls to copy constructor of this Particle, useful for polymorfism. More...
 
Mdouble getVolume () const
 Get Particle volume function, which required a reference to the Species vector. It returns the volume of the Particle. More...
 
void fixParticle ()
 Fix Particle function. It fixes a Particle by setting its inverse mass and inertia and velocities to zero. More...
 
bool isFixed () const
 Is fixed Particle function. It returns whether a Particle is fixed or not, by checking its inverse Mass. More...
 
void unfix ()
 Unfix Particle function, which required a reference to the Species vector. It unfixes a Particle by computing the Particles mass and inertia. More...
 
virtual void read (std::istream &is)
 Particle read function, which accepts an std::istream as input. More...
 
virtual MERCURY_DEPRECATED void oldRead (std::istream &is)
 deprecated version of the read function. More...
 
virtual void write (std::ostream &os) const
 Particle print function, which accepts an std::ostream as input. More...
 
virtual std::string getName () const
 Returns the name of the object. More...
 
void printHGrid (std::ostream &os) const
 Adds particle's HGrid level and cell coordinates to an ostream. More...
 
unsigned int getHGridLevel () const
 Returns particle's HGrid level. More...
 
BaseParticlegetHGridNextObject () const
 Returns pointer to next object in particle's HGrid level & cell. More...
 
BaseParticlegetHGridPrevObject () const
 Returns pointer to previous object in particle's HGrid level & cell. More...
 
int getHGridX () const
 Returns particle's HGrid cell X-coordinate. More...
 
int getHGridY () const
 Returns particle's HGrid cell Y-coordinate. More...
 
int getHGridZ () const
 Returns particle's HGrid cell Z-coordinate. More...
 
Mdouble getInertia () const
 Returns the particle's inertia_. More...
 
Mdouble getInvInertia () const
 Returns the particle's invInertia_. More...
 
Mdouble getInvMass () const
 Returns the particle's invMass_. More...
 
Mdouble getKineticEnergy () const
 Calculates the particle's kinetic energy. More...
 
Mdouble getMass () const
 Returns the particle's mass_. More...
 
BaseParticlegetPeriodicFromParticle () const
 Returns the 'original' particle this one's a periodic copy of. More...
 
Mdouble getRadius () const
 Returns the particle's radius_. More...
 
Mdouble getInteractionRadius () const
 Returns the particle's interaction radius, which might be different from radius_ (e.g., when dealing with wet particles) More...
 
Mdouble getWallInteractionRadius () const
 Returns the interaction radius for interaction with walls. See also BaseParticle::getInteractionRadius(). More...
 
const Vec3DgetDisplacement () const
 Returns the particle's displacement relative to the previous time step. More...
 
const Vec3DgetPreviousPosition () const
 Returns the particle's position in the previous time step. More...
 
const Vec3D getDisplacement2 (Mdouble xmin, Mdouble xmax, Mdouble ymin, Mdouble ymax, Mdouble zmin, Mdouble zmax, Mdouble t) const
 
void setInertia (const Mdouble newInertia)
 Sets the particle's inertia_ (and adjusts invInertia_ accordingly) More...
 
void setInfiniteInertia ()
 Sets the particle's inertia_ to 'infinite' (1e20) and its invInertia_ to 0. More...
 
void setPeriodicFromParticle (BaseParticle *p)
 Assigns the pointer to the 'original' particle this one's a periodic copy of. More...
 
void setHGridX (const int x)
 Sets the particle's HGrid cell X-coordinate. More...
 
void setHGridY (const int y)
 Sets the particle's HGrid cell Y-coordinate. More...
 
void setHGridZ (const int z)
 Sets the particle's HGrid cell Z-coordinate. More...
 
void setHGridLevel (const unsigned int level)
 Sets the particle's HGrid level. More...
 
void setHGridNextObject (BaseParticle *p)
 Sets the pointer to the next object in the particle's HGrid cell & level. More...
 
void setHGridPrevObject (BaseParticle *p)
 Sets the pointer to the previous object in the particle's HGrid cell & level. More...
 
void setRadius (const Mdouble radius)
 Sets the particle's radius_ (and adjusts the mass_ accordingly, based on the particle's species) More...
 
MERCURY_DEPRECATED void setMass (const Mdouble mass)
 Sets the particle's mass. More...
 
void setMassForP3Statistics (const Mdouble mass)
 Sets the particle's mass This function should not be used, but is necessary to extend the CG toolbox to non-spherical particles. More...
 
void setDisplacement (const Vec3D &disp)
 Sets the particle's displacement (= difference between current position and that of the previous time step) More...
 
void setPreviousPosition (const Vec3D &pos)
 Sets the particle's position in the previous time step. More...
 
void movePrevious (const Vec3D &posMove)
 Adds a vector to the particle's previousPosition_. More...
 
void accelerate (const Vec3D &vel)
 Increases the particle's velocity_ by the given vector. More...
 
void angularAccelerate (const Vec3D &angVel)
 Increases the particle's angularVelocity_ by the given vector. More...
 
void addDisplacement (const Vec3D &addDisp)
 Adds a vector to the particle's displacement_. More...
 
void setHandler (ParticleHandler *handler)
 Sets the pointer to the particle's ParticleHandler. More...
 
ParticleHandlergetHandler () const
 Returns pointer to the particle's ParticleHandler. More...
 
std::vector< BaseInteraction * > getInteractionWith (BaseParticle *P, Mdouble timeStamp, InteractionHandler *interactionHandler)
 Checks if particle is in interaction with given particle P, and if so, returns pointer to the associated BaseInteraction object (else returns 0). More...
 
void integrateBeforeForceComputation (double time, double timeStep)
 First step of Velocity Verlet integration. More...
 
void integrateAfterForceComputation (double time, double timeStep)
 Second step of Velocity Verlet integration. More...
 
unsigned int getParticleDimensions () const
 Returns the particle's dimensions (either 2 or 3). More...
 
MERCURY_DEPRECATED void setIndSpecies (unsigned int indSpecies)
 
void setSpecies (const ParticleSpecies *species)
 
virtual unsigned getNumberOfFieldsVTK () const
 
virtual std::string getTypeVTK (unsigned i) const
 
virtual std::string getNameVTK (unsigned i) const
 
virtual std::vector< MdoublegetFieldVTK (unsigned i) const
 
- Public Member Functions inherited from BaseInteractable
 BaseInteractable ()
 Default BaseInteractable constructor, it simply creates an empty BaseInteractable. More...
 
 BaseInteractable (const BaseInteractable &p)
 Copy constructor. It copies the BaseInteractable and all objects it contains. More...
 
virtual ~BaseInteractable ()
 Destructor, it simply destructs the BaseInteractable and all the objects it contains. More...
 
unsigned int getIndSpecies () const
 Returns the index of the Species of this BaseInteractable. More...
 
const ParticleSpeciesgetSpecies () const
 Returns a pointer to the species of this BaseInteractable. More...
 
void setSpecies (const ParticleSpecies *species)
 Sets the species of this BaseInteractable. More...
 
const Vec3DgetForce () const
 Returns the force on this BaseInteractable. More...
 
const Vec3DgetTorque () const
 Returns the torque on this BaseInteractable. More...
 
void setForce (const Vec3D &force)
 Sets the force on this BaseInteractable. More...
 
void setTorque (const Vec3D &torque)
 Sets the torque on this BaseInteractable. More...
 
void addForce (const Vec3D &addForce)
 Adds an amount to the force on this BaseInteractable. More...
 
void addTorque (const Vec3D &addTorque)
 Adds an amount to the torque on this BaseInteractable. More...
 
const Vec3DgetPosition () const
 Returns the position of this BaseInteractable. More...
 
const Vec3DgetOrientation () const
 Returns the orientation of this BaseInteractable. More...
 
void setPosition (const Vec3D &position)
 Sets the position of this BaseInteractable. More...
 
void setOrientation (const Vec3D &orientation)
 Sets the orientation of this BaseInteractable. More...
 
virtual void move (const Vec3D &move)
 Moves this BaseInteractable by adding an amount to the position. More...
 
void rotate (const Vec3D &rotate)
 Rotates this BaseInteractable. More...
 
const std::list
< BaseInteraction * > & 
getInteractions () const
 Returns a reference to the list of interactions in this BaseInteractable. More...
 
void addInteraction (BaseInteraction *I)
 Adds an interaction to this BaseInteractable. More...
 
bool removeInteraction (BaseInteraction *I)
 Removes an interaction from this BaseInteractable. More...
 
void copyInteractionsForPeriodicParticles (const BaseInteractable &p)
 Copies interactions to this BaseInteractable whenever a periodic copy made. More...
 
void setVelocity (const Vec3D &velocity)
 set the velocity of the BaseInteractable. More...
 
void setAngularVelocity (const Vec3D &angularVelocity)
 set the angular velocity of the BaseInteractble. More...
 
void addVelocity (const Vec3D &velocity)
 adds an increment to the velocity. More...
 
void addAngularVelocity (const Vec3D &angularVelocity)
 add an increment to the angular velocity. More...
 
virtual const Vec3DgetVelocity () const
 Returns the velocity of this interactable. More...
 
virtual const Vec3DgetAngularVelocity () const
 Returns the angular velocity of this interactable. More...
 
void setPrescribedPosition (const std::function< Vec3D(double)> &prescribedPosition)
 Allows the position of an infinite mass interactable to be prescribed. More...
 
void applyPrescribedPosition (double time)
 Computes the position from the user defined prescribed position function. More...
 
void setPrescribedVelocity (const std::function< Vec3D(double)> &prescribedVelocity)
 Allows the velocity of an infinite mass interactable to be prescribed. More...
 
void applyPrescribedVelocity (double time)
 Computes the velocity from the user defined prescribed velocity function. More...
 
void setPrescribedOrientation (const std::function< Vec3D(double)> &prescribedOrientation)
 Allows the orientation of the infinite mass interactbale to be prescribed. More...
 
void applyPrescribedOrientation (double time)
 Computes the orientation from the user defined prescribed orientation function. More...
 
void setPrescribedAngularVelocity (const std::function< Vec3D(double)> &prescribedAngularVelocity)
 Allows the angular velocity of the infinite mass interactable to be prescribed. More...
 
void applyPrescribedAngularVelocity (double time)
 Computes the angular velocity from the user defined prescribed angular velocity. More...
 
virtual const Vec3D getVelocityAtContact (const Vec3D &contact) const
 Returns the velocity at the contact point, use by many force laws. More...
 
void integrateBeforeForceComputation (double time, double timeStep)
 This is part of integrate routine for objects with infinite mass. More...
 
void integrateAfterForceComputation (double time, double timeStep)
 This is part of the integration routine for objects with infinite mass. More...
 
- Public Member Functions inherited from BaseObject
 BaseObject ()
 Default constructor. More...
 
 BaseObject (const BaseObject &p)
 Copy constructor, copies all the objects BaseObject contains. More...
 
virtual ~BaseObject ()
 virtual destructor More...
 
virtual void moveInHandler (const unsigned int index)
 Except that it is virtual, it does the same thing as setIndex() does. More...
 
void setIndex (const unsigned int index)
 Allows one to assign an index to an object in the handler/container. More...
 
void setId (const unsigned int 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...
 

Private Attributes

ParticleHandlerhandler_
 Pointer to the particle's ParticleHandler container. More...
 
int HGridX_
 Hgrid attributes. More...
 
int HGridY_
 
int HGridZ_
 
unsigned int HGridLevel_
 Cell position in the grid. More...
 
BaseParticleHGridNextObject_
 Grid level for the object. More...
 
BaseParticleHGridPrevObject_
 Pointer to the next Particle in the same HGrid cell. More...
 
Mdouble mass_
 Pointer to the previous Particle in the same HGrid cell. More...
 
Mdouble invMass_
 Particle mass_. More...
 
Mdouble inertia_
 Inverse Particle mass (for computation optimization) More...
 
Mdouble invInertia_
 Particle inertia_. More...
 
Mdouble radius_
 Inverse Particle inverse inertia (for computation optimization) More...
 
BaseParticleperiodicFromParticle_
 Particle radius_. More...
 
Vec3D displacement_
 Pointer to originating Particle. More...
 
Vec3D previousPosition_
 Displacement (only used in StatisticsVector, StatisticsPoint) More...
 

Friends

class ParticleSpecies
 Particle's position at previous time step. More...
 

Detailed Description

Definition at line 49 of file BaseParticle.h.

Constructor & Destructor Documentation

BaseParticle::BaseParticle ( )

Basic Particle constructor, creates an Particle at (0,0,0) with radius, mass and inertia equal to 1.

default constructor, creates an Particle at (0,0,0) with radius, mass and inertia equal to 1

Definition at line 38 of file BaseParticle.cc.

References DEBUG, displacement_, handler_, HGridLevel_, HGridNextObject_, HGridPrevObject_, HGridX_, HGridY_, HGridZ_, inertia_, invInertia_, invMass_, logger, mass_, periodicFromParticle_, radius_, and Vec3D::setZero().

Referenced by copy().

39 {
40  handler_ = nullptr;
42  radius_ = 1.0;
43  mass_ = invMass_ = 1.0;
44  inertia_ = invInertia_ = 1.0;
45  HGridNextObject_ = nullptr;
46 
47  periodicFromParticle_ = nullptr;
48 #ifdef CONTACT_LIST_HGRID
49  firstPossibleContact = nullptr;
50 #endif
51  HGridNextObject_ = nullptr;
52  HGridPrevObject_ = nullptr;
53  HGridLevel_ = 99999;
54  HGridX_ = 99999;
55  HGridY_ = 99999;
56  HGridZ_ = 99999;
57 
58  logger(DEBUG, "BaseParticle::BaseParticle() finished");
59 }
Mdouble invInertia_
Particle inertia_.
Definition: BaseParticle.h:409
Logger< MERCURY_LOGLEVEL > logger("MercuryKernel")
Mdouble mass_
Pointer to the previous Particle in the same HGrid cell.
Definition: BaseParticle.h:406
ParticleHandler * handler_
Pointer to the particle's ParticleHandler container.
Definition: BaseParticle.h:390
Mdouble invMass_
Particle mass_.
Definition: BaseParticle.h:407
void setZero()
Sets all elements to zero.
Definition: Vector.cc:52
int HGridX_
Hgrid attributes.
Definition: BaseParticle.h:400
BaseParticle * HGridNextObject_
Grid level for the object.
Definition: BaseParticle.h:402
Mdouble inertia_
Inverse Particle mass (for computation optimization)
Definition: BaseParticle.h:408
Vec3D displacement_
Pointer to originating Particle.
Definition: BaseParticle.h:413
unsigned int HGridLevel_
Cell position in the grid.
Definition: BaseParticle.h:401
BaseParticle * HGridPrevObject_
Pointer to the next Particle in the same HGrid cell.
Definition: BaseParticle.h:403
Mdouble radius_
Inverse Particle inverse inertia (for computation optimization)
Definition: BaseParticle.h:410
BaseParticle * periodicFromParticle_
Particle radius_.
Definition: BaseParticle.h:411
BaseParticle::BaseParticle ( const BaseParticle p)

Particle copy constructor, which accepts as input a reference to a Particle. It creates a copy of this Particle and all it's information. Usually it is better to use the copy() function for polymorfism.

Constructor that copies most of the properties of the given particle. Please note that not everything is copied, for example the position in the HGrid is not determined yet by the end of this constructor. It also does not copy the interactions and the pointer to the handler that handles this particle. Use with care.

Parameters
[in,out]pReference to the BaseParticle this one should become a copy of.

Definition at line 69 of file BaseParticle.cc.

References DEBUG, displacement_, getInertia(), getInvInertia(), getInvMass(), getMass(), handler_, HGridLevel_, HGridNextObject_, HGridPrevObject_, HGridX_, HGridY_, HGridZ_, inertia_, invInertia_, invMass_, logger, mass_, periodicFromParticle_, and radius_.

71 {
72  handler_ = nullptr;
74  radius_ = p.radius_;
75  mass_ = p.getMass();
76  invMass_ = p.getInvMass();
77  inertia_ = p.getInertia();
79 
80  HGridNextObject_ = nullptr;
81  HGridPrevObject_ = nullptr;
82  HGridX_ = 99999;
83  HGridY_ = 99999;
84  HGridZ_ = 99999;
86 
88 #ifdef CONTACT_LIST_HGRID
89  firstPossibleContact = nullptr;
90 #endif
91  logger(DEBUG, "BaseParticle::BaseParticle(BaseParticle &p) finished");
92 }
Mdouble invInertia_
Particle inertia_.
Definition: BaseParticle.h:409
Logger< MERCURY_LOGLEVEL > logger("MercuryKernel")
Mdouble mass_
Pointer to the previous Particle in the same HGrid cell.
Definition: BaseParticle.h:406
ParticleHandler * handler_
Pointer to the particle's ParticleHandler container.
Definition: BaseParticle.h:390
Mdouble invMass_
Particle mass_.
Definition: BaseParticle.h:407
Mdouble getInvMass() const
Returns the particle's invMass_.
int HGridX_
Hgrid attributes.
Definition: BaseParticle.h:400
BaseParticle * HGridNextObject_
Grid level for the object.
Definition: BaseParticle.h:402
Mdouble inertia_
Inverse Particle mass (for computation optimization)
Definition: BaseParticle.h:408
Vec3D displacement_
Pointer to originating Particle.
Definition: BaseParticle.h:413
Mdouble getMass() const
Returns the particle's mass_.
unsigned int HGridLevel_
Cell position in the grid.
Definition: BaseParticle.h:401
BaseParticle * HGridPrevObject_
Pointer to the next Particle in the same HGrid cell.
Definition: BaseParticle.h:403
Mdouble getInertia() const
Returns the particle's inertia_.
Mdouble getInvInertia() const
Returns the particle's invInertia_.
Mdouble radius_
Inverse Particle inverse inertia (for computation optimization)
Definition: BaseParticle.h:410
BaseParticle * periodicFromParticle_
Particle radius_.
Definition: BaseParticle.h:411
BaseInteractable()
Default BaseInteractable constructor, it simply creates an empty BaseInteractable.
BaseParticle::~BaseParticle ( )
virtual

Particle destructor, needs to be implemented and checked if it removes tangential spring information.

Destructor. It asks the ParticleHandler to check if this was the smallest or largest particle and adjust itself accordingly.

Definition at line 98 of file BaseParticle.cc.

References ParticleHandler::checkExtremaOnDelete(), DEBUG, getHandler(), BaseObject::getId(), and logger.

99 {
100  if (getHandler() != nullptr)
101  {
103  }
104  logger(DEBUG, "BaseParticle::~BaseParticle() of particle % finished.", getId());
105 }
unsigned int getId() const
Returns the unique identifier of any particular object.
Definition: BaseObject.cc:116
void checkExtremaOnDelete(BaseParticle *P)
Checks if the extrema of this ParticleHandler needs updating when a particle is deleted.
Logger< MERCURY_LOGLEVEL > logger("MercuryKernel")
ParticleHandler * getHandler() const
Returns pointer to the particle's ParticleHandler.

Member Function Documentation

void BaseParticle::accelerate ( const Vec3D vel)

Increases the particle's velocity_ by the given vector.

increases the the particle's velocity_ (BaseInteractable member) by adding the given vector.

Parameters
[in]velvector to be added to the velocity_

Definition at line 662 of file BaseParticle.cc.

References BaseInteractable::addVelocity().

Referenced by integrateAfterForceComputation(), integrateBeforeForceComputation(), and AngledPeriodicBoundary::shiftPosition().

663 {
664  addVelocity(vel);
665 }
void addVelocity(const Vec3D &velocity)
adds an increment to the velocity.
void BaseParticle::addDisplacement ( const Vec3D addDisp)

Adds a vector to the particle's displacement_.

Lets you add a vector to the particle's displacement_ vector.

Parameters
[in]addDispvector to be added.

Definition at line 681 of file BaseParticle.cc.

References displacement_.

682 {
683  displacement_ += addDisp;
684 }
Vec3D displacement_
Pointer to originating Particle.
Definition: BaseParticle.h:413
void BaseParticle::angularAccelerate ( const Vec3D angVel)

Increases the particle's angularVelocity_ by the given vector.

increases the particle's angularVelocity_ (BaseInteractable member) by adding the given vector.

Parameters
[in]angVelvector to be added to the angularVelocity_

Definition at line 672 of file BaseParticle.cc.

References BaseInteractable::addAngularVelocity().

Referenced by integrateAfterForceComputation(), integrateBeforeForceComputation(), and AngledPeriodicBoundary::shiftPosition().

673 {
674  addAngularVelocity(angVel);
675 }
void addAngularVelocity(const Vec3D &angularVelocity)
add an increment to the angular velocity.
BaseParticle * BaseParticle::copy ( ) const
virtual
void BaseParticle::fixParticle ( )

Fix Particle function. It fixes a Particle by setting its inverse mass and inertia and velocities to zero.

Fixes a BaseParticle by setting its inverse mass and inertia and velocities to zero.

Definition at line 146 of file BaseParticle.cc.

References inertia_, invInertia_, invMass_, mass_, BaseInteractable::setAngularVelocity(), and BaseInteractable::setVelocity().

Referenced by DPMBase::setFixedParticles().

147 {
148  mass_ = 1e20;
149  invMass_ = 0.0;
150  inertia_ = 1e20;
151  invInertia_ = 0.0;
152  setVelocity(Vec3D(0.0, 0.0, 0.0));
153  setAngularVelocity(Vec3D(0.0, 0.0, 0.0));
154 }
Mdouble invInertia_
Particle inertia_.
Definition: BaseParticle.h:409
void setVelocity(const Vec3D &velocity)
set the velocity of the BaseInteractable.
Mdouble mass_
Pointer to the previous Particle in the same HGrid cell.
Definition: BaseParticle.h:406
Mdouble invMass_
Particle mass_.
Definition: BaseParticle.h:407
Mdouble inertia_
Inverse Particle mass (for computation optimization)
Definition: BaseParticle.h:408
Implementation of a 3D vector (by Vitaliy).
Definition: Vector.h:45
void setAngularVelocity(const Vec3D &angularVelocity)
set the angular velocity of the BaseInteractble.
const Vec3D & BaseParticle::getDisplacement ( ) const

Returns the particle's displacement relative to the previous time step.

Returns the particle's displacement_, which is the difference between the current particle's position and its position in the previous time step.

Returns
(reference to) the particle displacement vector

Definition at line 421 of file BaseParticle.cc.

References displacement_.

Referenced by CircularPeriodicBoundary::checkBoundaryAfterParticleMoved(), and integrateBeforeForceComputation().

422 {
423  return displacement_;
424 }
Vec3D displacement_
Pointer to originating Particle.
Definition: BaseParticle.h:413
const Vec3D BaseParticle::getDisplacement2 ( Mdouble  xmin,
Mdouble  xmax,
Mdouble  ymin,
Mdouble  ymax,
Mdouble  zmin,
Mdouble  zmax,
Mdouble  t 
) const
Todo:
see .cc file.
Todo:
Rewrite, redefine (TW). Is only used in StatisticsVector.hcc, consider moving to that class.

Definition at line 439 of file BaseParticle.cc.

References BaseInteractable::getPosition(), getPreviousPosition(), Vec3D::X, Vec3D::Y, and Vec3D::Z.

440 {
442  if (xmax > xmin && fabs(disp.X) > .5 * (xmax - xmin))
443  {
444  if (disp.X > 0)
445  disp.X -= xmax - xmin;
446  else
447  disp.X += xmax - xmin;
448  }
449  if (ymax > ymin && fabs(disp.Y) > .5 * (ymax - ymin))
450  {
451  if (disp.Y > 0)
452  disp.Y -= ymax - ymin;
453  else
454  disp.Y += ymax - ymin;
455  }
456  if (zmax > zmin && fabs(disp.Z) > .5 * (zmax - zmin))
457  {
458  if (disp.Z > 0)
459  disp.Z -= zmax - zmin;
460  else
461  disp.Z += zmax - zmin;
462  }
463  disp /= t;
464  return disp;
465 }
Mdouble X
the vector components
Definition: Vector.h:52
const Vec3D & getPreviousPosition() const
Returns the particle's position in the previous time step.
const Vec3D & getPosition() const
Returns the position of this BaseInteractable.
Mdouble Y
Definition: Vector.h:52
Implementation of a 3D vector (by Vitaliy).
Definition: Vector.h:45
Mdouble Z
Definition: Vector.h:52
std::vector< Mdouble > BaseParticle::getFieldVTK ( unsigned  i) const
virtual

Reimplemented in LiquidFilmParticle.

Definition at line 864 of file BaseParticle.cc.

864  {
865  return std::vector<Mdouble>();
866 }
ParticleHandler * BaseParticle::getHandler ( ) const

Returns pointer to the particle's ParticleHandler.

Returns the particle's ParticleHandler

Returns
pointer to the particle's ParticleHandler; null if handler not yet specified

Definition at line 701 of file BaseParticle.cc.

References handler_.

Referenced by getParticleDimensions(), integrateAfterForceComputation(), integrateBeforeForceComputation(), setHandler(), setRadius(), and ~BaseParticle().

702 {
703  return handler_;
704 }
ParticleHandler * handler_
Pointer to the particle's ParticleHandler container.
Definition: BaseParticle.h:390
BaseParticle * BaseParticle::getHGridNextObject ( ) const

Returns pointer to next object in particle's HGrid level & cell.

Returns the next object in the particle's HGrid cell

Returns
pointer to the next object in the particle's HGrid cell

Definition at line 275 of file BaseParticle.cc.

References HGridNextObject_.

Referenced by Mercury2D::hGridFindContactsWithinTargetCell(), Mercury3D::hGridFindContactsWithinTargetCell(), Mercury2D::hGridFindContactsWithTargetCell(), Mercury3D::hGridFindContactsWithTargetCell(), Mercury3D::hGridHasContactsInTargetCell(), Mercury2D::hGridRemoveParticle(), and Mercury3D::hGridRemoveParticle().

276 {
277  return HGridNextObject_;
278 }
BaseParticle * HGridNextObject_
Grid level for the object.
Definition: BaseParticle.h:402
BaseParticle * BaseParticle::getHGridPrevObject ( ) const

Returns pointer to previous object in particle's HGrid level & cell.

Returns the previous object in the particle's HGrid cell

Returns
pointer to the previous object in the particle's HGrid cell

Definition at line 285 of file BaseParticle.cc.

References HGridPrevObject_.

Referenced by Mercury2D::hGridRemoveParticle(), and Mercury3D::hGridRemoveParticle().

286 {
287  return HGridPrevObject_;
288 }
BaseParticle * HGridPrevObject_
Pointer to the next Particle in the same HGrid cell.
Definition: BaseParticle.h:403
int BaseParticle::getHGridZ ( ) const

Returns particle's HGrid cell Z-coordinate.

Get the particle's HGrid cell's Z-coordinate

Returns
the particle's HGrid cell's Z-coordinate

Definition at line 323 of file BaseParticle.cc.

References HGridZ_.

Referenced by Mercury3D::hGridFindContactsWithinTargetCell(), Mercury3D::hGridFindContactsWithTargetCell(), Mercury3D::hGridFindOneSidedContacts(), Mercury3D::hGridHasContactsInTargetCell(), Mercury3D::hGridRemoveParticle(), and Mercury3D::hGridUpdateParticle().

324 {
325  return HGridZ_;
326 }
Mdouble BaseParticle::getInertia ( ) const

Returns the particle's inertia_.

Returns the particle's inertia_

Returns
the particle's inertia_

Definition at line 332 of file BaseParticle.cc.

References inertia_.

Referenced by BaseParticle().

333 {
334  return inertia_;
335 }
Mdouble inertia_
Inverse Particle mass (for computation optimization)
Definition: BaseParticle.h:408
Mdouble BaseParticle::getInteractionRadius ( ) const

Returns the particle's interaction radius, which might be different from radius_ (e.g., when dealing with wet particles)

Calculates the interaction radius of the particle (when it comes to interaction with other particles), including the effect of a possible additional 'interaction distance' besides the 'normal' radius. The interaction radius differs from the radius_, for example, when dealing with wet particles (i.e. particles with an additional liquid layer, which is dealt with in the particle's species).

Returns
the particle's interaction radius for particle-particle interaction

Definition at line 401 of file BaseParticle.cc.

References BaseSpecies::getInteractionDistance(), BaseInteractable::getSpecies(), and radius_.

Referenced by DPMBase::areInContact(), ParticleHandler::checkExtrema(), DPMBase::checkParticleForInteractionLocal(), DPMBase::checkParticleForInteractionLocalPeriodic(), ParticleHandler::computeLargestParticle(), ParticleHandler::computeSmallestParticle(), LeesEdwardsBoundary::createHorizontalPeriodicParticles(), ShearBoxBoundary::createHorizontalPeriodicParticles(), CircularPeriodicBoundary::createPeriodicParticles(), MaserBoundary::createPeriodicParticles(), AngledPeriodicBoundary::createPeriodicParticles(), PeriodicBoundary::createPeriodicParticles(), ShearBoxBoundary::createVerticalPeriodicParticles(), LeesEdwardsBoundary::createVerticalPeriodicParticles(), IntersectionOfWalls::getDistanceAndNormal(), MercuryBase::getHGridTargetMaxInteractionRadius(), MercuryBase::getHGridTargetMinInteractionRadius(), getInteractionWith(), Mercury3D::hGridFindOneSidedContacts(), Mercury2D::hGridFindOneSidedContacts(), Mercury2D::hGridHasParticleContacts(), MercuryBase::hGridNeedsRebuilding(), HGridOptimiser::initialise(), and HGrid::insertParticleToHgrid().

402 {
403  return radius_ + getSpecies()->getInteractionDistance() * 0.5;
404 }
virtual Mdouble getInteractionDistance() const =0
returns the largest separation distance at which adhesive short-range forces can occur.
const ParticleSpecies * getSpecies() const
Returns a pointer to the species of this BaseInteractable.
Mdouble radius_
Inverse Particle inverse inertia (for computation optimization)
Definition: BaseParticle.h:410
std::vector< BaseInteraction * > BaseParticle::getInteractionWith ( BaseParticle P,
Mdouble  timeStamp,
InteractionHandler interactionHandler 
)
virtual

Checks if particle is in interaction with given particle P, and if so, returns pointer to the associated BaseInteraction object (else returns 0).

Creates/updates a BaseInteraction object, treating the interaction between this particle and a given one, in case there is an overlap between the two.

Parameters
[in]Pparticle to check the interaction with
[in]timeStamptime stamp to be assigned to the interaction object (i.e., the current time)
[in,out]interactionHandlerBaseInteraction container from where the interaction is retrieved, and to which it is assigned (if it is a new interaction).
Returns
the pointer to the interaction object (if the particles overlap), or 0 (if they don't overlap).
Todo:
We should consider setting the contact point to
Author
weinhartt

Implements BaseInteractable.

Definition at line 717 of file BaseParticle.cc.

References InteractionHandler::getInteraction(), getInteractionRadius(), Vec3D::getLengthSquared(), BaseInteraction::getNormal(), BaseInteraction::getOverlap(), BaseInteractable::getPosition(), getRadius(), BaseInteraction::setContactPoint(), BaseInteraction::setDistance(), BaseInteraction::setNormal(), and BaseInteraction::setOverlap().

Referenced by DPMBase::computeInternalForces(), and FileReader::read().

718 {
719  //get the normal (from P away from the contact)
720  Vec3D branchVector = P->getPosition() - getPosition();
721  //Get the square of the distance between particle i and particle j
722  Mdouble distanceSquared = Vec3D::getLengthSquared(branchVector);
723  Mdouble sumOfInteractionRadii = P->getInteractionRadius() + getInteractionRadius();
724  std::vector<BaseInteraction*> interactions;
725  if (distanceSquared < (sumOfInteractionRadii * sumOfInteractionRadii))
726  {
727  BaseInteraction* C = interactionHandler->getInteraction(P, this, timeStamp);
728  Mdouble distance = std::sqrt(distanceSquared);
729  C->setNormal(branchVector / distance);
730  C->setOverlap(P->getRadius() + getRadius() - distance);
731  C->setDistance(distance);
732  C->setContactPoint(P->getPosition() - (P->getRadius() - 0.5 * C->getOverlap()) * C->getNormal());
734  //Mdouble ratio=P->getRadius()/(getRadius()+P->getRadius());
735  //C->setContactPoint(P->getPosition() - (P->getRadius() - ratio * C->getOverlap()) * C->getNormal());
736  interactions.push_back(C);
737  }
738  return interactions;
739 }
void setNormal(Vec3D normal)
Sets the normal vector between the two interacting objects.
void setOverlap(Mdouble overlap)
Set the overlap between the two interacting object.
BaseInteraction * getInteraction(BaseInteractable *P, BaseInteractable *I, Mdouble timeStamp)
Returns the Interaction between the BaseInteractable's P and I.
double Mdouble
void setContactPoint(Vec3D contactPoint)
Set the location of the contact point between the two interacting objects.
void setDistance(Mdouble distance)
Sets the interaction distance between the two interacting objects.
Mdouble getLengthSquared() const
Calculates the squared length of this Vec3D: .
Definition: Vector.cc:300
const Vec3D & getPosition() const
Returns the position of this BaseInteractable.
Stores information about interactions between two interactable objects; often particles but could be ...
const Vec3D & getNormal() const
Gets the normal vector between the two interacting objects.
Mdouble getRadius() const
Returns the particle's radius_.
Mdouble getOverlap() const
Returns a Mdouble with the current overlap between the two interacting objects.
Implementation of a 3D vector (by Vitaliy).
Definition: Vector.h:45
Mdouble getInteractionRadius() const
Returns the particle's interaction radius, which might be different from radius_ (e.g., when dealing with wet particles)
Mdouble BaseParticle::getInvInertia ( ) const

Returns the particle's invInertia_.

Returns the particle's inverse inertia_

Returns
the particles invInertia_

Definition at line 341 of file BaseParticle.cc.

References invInertia_.

Referenced by BaseParticle(), integrateAfterForceComputation(), and integrateBeforeForceComputation().

342 {
343  return invInertia_;
344 }
Mdouble invInertia_
Particle inertia_.
Definition: BaseParticle.h:409
Mdouble BaseParticle::getInvMass ( ) const

Returns the particle's invMass_.

Returns the particle's inverse mass_

Returns
the particle's invMass_

Definition at line 350 of file BaseParticle.cc.

References invMass_.

Referenced by BaseParticle(), ThermalInteraction< NormalForceInteraction >::computeNormalForce(), integrateAfterForceComputation(), and integrateBeforeForceComputation().

351 {
352  return invMass_;
353 }
Mdouble invMass_
Particle mass_.
Definition: BaseParticle.h:407
Mdouble BaseParticle::getKineticEnergy ( ) const

Calculates the particle's kinetic energy.

Calculates the particle's kinetic energy

Returns
the particle's kinetic energy

Definition at line 359 of file BaseParticle.cc.

References Vec3D::getLengthSquared(), BaseInteractable::getVelocity(), and mass_.

360 {
361  return 0.5 * mass_ * getVelocity().getLengthSquared();
362 }
static Mdouble getLengthSquared(const Vec3D &a)
Calculates the squared length of a Vec3D: .
Definition: Vector.cc:291
Mdouble mass_
Pointer to the previous Particle in the same HGrid cell.
Definition: BaseParticle.h:406
virtual const Vec3D & getVelocity() const
Returns the velocity of this interactable.
Mdouble BaseParticle::getMass ( ) const

Returns the particle's mass_.

Returns the mass_ of the particle

Returns
the mass_ of the particle

Definition at line 368 of file BaseParticle.cc.

References mass_.

Referenced by BaseParticle(), DPMBase::computeExternalForces(), ParticleSpecies::computeMass(), and BaseInteraction::getEffectiveMass().

369 {
370  return mass_;
371 }
Mdouble mass_
Pointer to the previous Particle in the same HGrid cell.
Definition: BaseParticle.h:406
std::string BaseParticle::getName ( ) const
virtual

Returns the name of the object.

Returns the name of the object; in this case 'BaseParticle'.

Returns
The object name.

Implements BaseObject.

Reimplemented in LiquidFilmParticle, and ThermalParticle.

Definition at line 194 of file BaseParticle.cc.

Referenced by InsertionBoundary::write().

195 {
196  return "BaseParticle";
197 }
std::string BaseParticle::getNameVTK ( unsigned  i) const
virtual

Reimplemented in LiquidFilmParticle.

Definition at line 860 of file BaseParticle.cc.

Referenced by ParticleHandler::writeVTK().

860  {
861  return "";
862 }
unsigned BaseParticle::getNumberOfFieldsVTK ( ) const
virtual

Reimplemented in LiquidFilmParticle.

Definition at line 852 of file BaseParticle.cc.

Referenced by ParticleHandler::writeVTK().

852  {
853  return 0;
854 }
unsigned int BaseParticle::getParticleDimensions ( ) const

Returns the particle's dimensions (either 2 or 3).

Returns the amount of dimensions of the particle (2 or 3, basically)

Returns
the number of dimension of the particle

Definition at line 801 of file BaseParticle.cc.

References BaseHandler< T >::getDPMBase(), getHandler(), and DPMBase::getParticleDimensions().

Referenced by ParticleSpecies::computeMass(), and getVolume().

802 {
804 }
unsigned int getParticleDimensions() const
Returns the particle dimensions.
Definition: DPMBase.cc:599
ParticleHandler * getHandler() const
Returns pointer to the particle's ParticleHandler.
DPMBase * getDPMBase()
Gets the problem that is solved using this handler.
Definition: BaseHandler.h:543
const Vec3D & BaseParticle::getPreviousPosition ( ) const

Returns the particle's position in the previous time step.

Returns the particle's position in the previous time step.

Returns
(reference to) the previous position of the particle

Definition at line 430 of file BaseParticle.cc.

References previousPosition_.

Referenced by getDisplacement2().

431 {
432  return previousPosition_;
433 }
Vec3D previousPosition_
Displacement (only used in StatisticsVector, StatisticsPoint)
Definition: BaseParticle.h:414
std::string BaseParticle::getTypeVTK ( unsigned  i) const
virtual

Reimplemented in LiquidFilmParticle.

Definition at line 856 of file BaseParticle.cc.

Referenced by ParticleHandler::writeVTK().

856  {
857  return "";
858 }
Mdouble BaseParticle::getVolume ( ) const

Get Particle volume function, which required a reference to the Species vector. It returns the volume of the Particle.

Returns the volume of the BaseParticle, which is calculated using its number of dimensions and radius.

Returns
The actual volume of this BaseParticle.

Definition at line 122 of file BaseParticle.cc.

References ERROR, getParticleDimensions(), handler_, logger, constants::pi, and radius_.

123 {
124  if (handler_ == nullptr)
125  {
126  logger(ERROR, "[BaseParticle::getVolume] no particle handler specified");
127  return 0;
128  }
129  switch (getParticleDimensions())
130  {
131  case 3:
132  return (4.0 / 3.0 * constants::pi * radius_ * radius_ * radius_);
133  case 2:
134  return (constants::pi * radius_ * radius_);
135  case 1:
136  return (2.0 * radius_);
137  default:
138  logger(ERROR, "[BaseParticle::getVolume] dimension of the particle is not set");
139  return 0;
140  }
141 }
Logger< MERCURY_LOGLEVEL > logger("MercuryKernel")
ParticleHandler * handler_
Pointer to the particle's ParticleHandler container.
Definition: BaseParticle.h:390
const Mdouble pi
Definition: ExtendedMath.h:42
Mdouble radius_
Inverse Particle inverse inertia (for computation optimization)
Definition: BaseParticle.h:410
unsigned int getParticleDimensions() const
Returns the particle's dimensions (either 2 or 3).
Mdouble BaseParticle::getWallInteractionRadius ( ) const

Returns the interaction radius for interaction with walls. See also BaseParticle::getInteractionRadius().

The interaction radius of the particle (when it comes to interaction with walls). See also BaseParticle::getInteractionRadius().

Returns
the particle's interaction radius for particle-wall interaction

Definition at line 411 of file BaseParticle.cc.

References BaseSpecies::getInteractionDistance(), BaseInteractable::getSpecies(), and radius_.

Referenced by Screw::getDistanceAndNormal(), TriangulatedWall::Face::getDistanceAndNormal(), Coil::getDistanceAndNormal(), RestrictedWall::getDistanceAndNormal(), SphericalWall::getDistanceAndNormal(), AxisymmetricIntersectionOfWalls::getDistanceAndNormal(), and InfiniteWall::getDistanceAndNormal().

412 {
414 }
virtual Mdouble getInteractionDistance() const =0
returns the largest separation distance at which adhesive short-range forces can occur.
const ParticleSpecies * getSpecies() const
Returns a pointer to the species of this BaseInteractable.
Mdouble radius_
Inverse Particle inverse inertia (for computation optimization)
Definition: BaseParticle.h:410
void BaseParticle::integrateAfterForceComputation ( double  time,
double  timeStep 
)

Second step of Velocity Verlet integration.

Second step of Velocity Verlet integration (see also http://en.wikipedia.org/wiki/Verlet_integration#Velocity_Verlet).

Parameters
[in]timecurrent time
[in]timeStepcurrent time step

Definition at line 781 of file BaseParticle.cc.

References accelerate(), angularAccelerate(), BaseInteractable::getForce(), getHandler(), getInvInertia(), getInvMass(), BaseInteractable::getTorque(), and BaseInteractable::integrateAfterForceComputation().

782 {
783  if (getInvMass() == 0.0)
784  {
786  }
787  else
788  {
789  accelerate(getForce() * getInvMass() * 0.5 * timeStep);
790  if (getHandler()->getDPMBase()->getRotation())
791  {
792  angularAccelerate(getTorque() * getInvInertia() * 0.5 * timeStep);
793  }
794  }
795 }
const Vec3D & getForce() const
Returns the force on this BaseInteractable.
ParticleHandler * getHandler() const
Returns pointer to the particle's ParticleHandler.
Mdouble getInvMass() const
Returns the particle's invMass_.
void accelerate(const Vec3D &vel)
Increases the particle's velocity_ by the given vector.
const Vec3D & getTorque() const
Returns the torque on this BaseInteractable.
Mdouble getInvInertia() const
Returns the particle's invInertia_.
void angularAccelerate(const Vec3D &angVel)
Increases the particle's angularVelocity_ by the given vector.
void integrateAfterForceComputation(double time, double timeStep)
This is part of the integration routine for objects with infinite mass.
void BaseParticle::integrateBeforeForceComputation ( double  time,
double  timeStep 
)

First step of Velocity Verlet integration.

Todo:
TW: add sum sum fi Wi (the stress gradient) to the statistics

First step of Velocity Verlet integration (see also http://en.wikipedia.org/wiki/Verlet_integration#Velocity_Verlet).

Parameters
[in]timecurrent time
[in]timeStepcurrent time step
Todo:
If the position is described by the used, one should also call BaseInteractabld::integrateBeforeForceComputation. To check if it works correctly, remove the p0.fixParticle() line from the DrivenParticleUnitTest
Author
irana

Definition at line 751 of file BaseParticle.cc.

References accelerate(), angularAccelerate(), BaseInteractable::getAngularVelocity(), getDisplacement(), BaseHandler< T >::getDPMBase(), BaseInteractable::getForce(), getHandler(), getInvInertia(), getInvMass(), BaseInteractable::getTorque(), BaseInteractable::getVelocity(), DPMBase::hGridUpdateMove(), BaseInteractable::integrateBeforeForceComputation(), BaseInteractable::move(), BaseInteractable::rotate(), and setDisplacement().

752 {
757  if (getInvMass() == 0.0)
758  {
760  }
761  else
762  {
763  accelerate(getForce() * getInvMass() * 0.5 * timeStep);
764  setDisplacement(getVelocity() * timeStep);
766  getHandler()->getDPMBase()->hGridUpdateMove(this, getDisplacement().getLength());
767  if (getHandler()->getDPMBase()->getRotation())
768  {
769  angularAccelerate(getTorque() * getInvInertia() * 0.5 * timeStep);
770  rotate(getAngularVelocity() * timeStep);
771  }
772  }
773 }
void rotate(const Vec3D &rotate)
Rotates this BaseInteractable.
const Vec3D & getForce() const
Returns the force on this BaseInteractable.
void setDisplacement(const Vec3D &disp)
Sets the particle's displacement (= difference between current position and that of the previous time...
void integrateBeforeForceComputation(double time, double timeStep)
This is part of integrate routine for objects with infinite mass.
virtual const Vec3D & getAngularVelocity() const
Returns the angular velocity of this interactable.
ParticleHandler * getHandler() const
Returns pointer to the particle's ParticleHandler.
Mdouble getInvMass() const
Returns the particle's invMass_.
void accelerate(const Vec3D &vel)
Increases the particle's velocity_ by the given vector.
const Vec3D & getDisplacement() const
Returns the particle's displacement relative to the previous time step.
const Vec3D & getTorque() const
Returns the torque on this BaseInteractable.
Mdouble getInvInertia() const
Returns the particle's invInertia_.
void angularAccelerate(const Vec3D &angVel)
Increases the particle's angularVelocity_ by the given vector.
virtual void hGridUpdateMove(BaseParticle *, Mdouble)
no implementation but can be overidden in its derived classes.
Definition: DPMBase.cc:845
virtual const Vec3D & getVelocity() const
Returns the velocity of this interactable.
virtual void move(const Vec3D &move)
Moves this BaseInteractable by adding an amount to the position.
DPMBase * getDPMBase()
Gets the problem that is solved using this handler.
Definition: BaseHandler.h:543
bool BaseParticle::isFixed ( ) const

Is fixed Particle function. It returns whether a Particle is fixed or not, by checking its inverse Mass.

Checks whether a BaseParticle is fixed or not, by checking its inverse Mass.

Returns
TRUE if particle is fixed, i.e. if the inverse mass (invMass_) is 0.

Definition at line 160 of file BaseParticle.cc.

References invMass_.

Referenced by DPMBase::computeExternalForces(), DPMBase::computeInternalForces(), ParticleSpecies::computeMass(), BaseInteraction::gatherContactStatistics(), and BaseInteraction::writeToFStat().

161 {
162  return (invMass_ == 0.0);
163 }
Mdouble invMass_
Particle mass_.
Definition: BaseParticle.h:407
void BaseParticle::movePrevious ( const Vec3D posMove)

Adds a vector to the particle's previousPosition_.

Lets you add a vector to the particle's previousPosition_ vector.

Parameters
[in]posMovethe vector to be added to the current previousPosition_ vector.

Definition at line 652 of file BaseParticle.cc.

References previousPosition_.

653 {
654  previousPosition_ += posMove;
655 }
Vec3D previousPosition_
Displacement (only used in StatisticsVector, StatisticsPoint)
Definition: BaseParticle.h:414
void BaseParticle::oldRead ( std::istream &  is)
virtual

deprecated version of the read function.

Deprecated:
Should NOT BE USED by any new users! Is expected to be obsolete by Mercury 2.0. Please use BaseParticle::read() instead.

This is the previously used version of the read function. Now just kept for legacy purposes.

Deprecated:
Should be gone in Mercury 2.0. Use BaseParticle::read() instead.
Todo:
{fix this function} IFCD: what is broken about it?

Definition at line 227 of file BaseParticle.cc.

References inertia_, invInertia_, invMass_, mass_, setIndSpecies(), BaseInteractable::setOrientation(), and BaseInteractable::setPosition().

Referenced by ParticleHandler::readObject().

228 {
230 
231  unsigned int indSpecies;
232  Vec3D orientation;
233  Vec3D position;
234  is >> invMass_ >> invInertia_ >> indSpecies;
235  setPosition(position);
236  setOrientation(orientation);
237  setIndSpecies(indSpecies);
238  if (invMass_ != 0.0)
239  mass_ = 1.0 / invMass_;
240  else
241  mass_ = 1e20;
242  if (invInertia_ != 0.0)
243  inertia_ = 1.0 / invInertia_;
244  else
245  inertia_ = 1e20;
246 }
Mdouble invInertia_
Particle inertia_.
Definition: BaseParticle.h:409
void setOrientation(const Vec3D &orientation)
Sets the orientation of this BaseInteractable.
Mdouble mass_
Pointer to the previous Particle in the same HGrid cell.
Definition: BaseParticle.h:406
Mdouble invMass_
Particle mass_.
Definition: BaseParticle.h:407
Mdouble inertia_
Inverse Particle mass (for computation optimization)
Definition: BaseParticle.h:408
MERCURY_DEPRECATED void setIndSpecies(unsigned int indSpecies)
void setPosition(const Vec3D &position)
Sets the position of this BaseInteractable.
Implementation of a 3D vector (by Vitaliy).
Definition: Vector.h:45
void BaseParticle::printHGrid ( std::ostream &  os) const

Adds particle's HGrid level and cell coordinates to an ostream.

Adds the particle's HGridLevel_ and HGRid x/y/z positions to an std::ostream.

Parameters
[in,out]osthe ostream which has the mentioned properties added.

Definition at line 253 of file BaseParticle.cc.

References HGridLevel_, HGridX_, HGridY_, and HGridZ_.

254 {
255  os << "Particle( HGRID_Level:" << HGridLevel_
256  << ", HGRID_x:" << HGridX_
257  << ", HGRID_y:" << HGridY_
258  << ", HGRID_z:" << HGridZ_
259  << ")";
260 }
int HGridX_
Hgrid attributes.
Definition: BaseParticle.h:400
unsigned int HGridLevel_
Cell position in the grid.
Definition: BaseParticle.h:401
void BaseParticle::read ( std::istream &  is)
virtual

Particle read function, which accepts an std::istream as input.

Particle read function. Has an std::istream as argument, from which it extracts the radius_, invMass_ and invInertia_, respectively. From these the mass_ and inertia_ are deduced. An additional set of properties is read through the call to the parent's method BaseInteractable::read().

Parameters
[in,out]isinput stream with particle properties.

Implements BaseInteractable.

Reimplemented in LiquidFilmParticle, and ThermalParticle.

Definition at line 207 of file BaseParticle.cc.

References inertia_, invInertia_, invMass_, mass_, radius_, and BaseInteractable::read().

Referenced by LiquidFilmParticle::read(), and ThermalParticle::read().

208 {
210  std::string dummy;
211  is >> dummy >> radius_ >> dummy >> invMass_ >> dummy >> invInertia_;
212  if (invMass_ != 0.0)
213  mass_ = 1.0 / invMass_;
214  else
215  mass_ = 1e20;
216  if (invInertia_ != 0.0)
217  inertia_ = 1.0 / invInertia_;
218  else
219  inertia_ = 1e20;
220 }
Mdouble invInertia_
Particle inertia_.
Definition: BaseParticle.h:409
Mdouble mass_
Pointer to the previous Particle in the same HGrid cell.
Definition: BaseParticle.h:406
Mdouble invMass_
Particle mass_.
Definition: BaseParticle.h:407
Mdouble inertia_
Inverse Particle mass (for computation optimization)
Definition: BaseParticle.h:408
virtual void read(std::istream &is)=0
Reads a BaseInteractable from an input stream.
Mdouble radius_
Inverse Particle inverse inertia (for computation optimization)
Definition: BaseParticle.h:410
void BaseParticle::setDisplacement ( const Vec3D disp)

Sets the particle's displacement (= difference between current position and that of the previous time step)

This is used to set the particle displacement_

Parameters
[in]dispthe displacement vector

Definition at line 633 of file BaseParticle.cc.

References displacement_.

Referenced by CircularPeriodicBoundary::checkBoundaryAfterParticleMoved(), and integrateBeforeForceComputation().

634 {
635  displacement_ = disp;
636 }
Vec3D displacement_
Pointer to originating Particle.
Definition: BaseParticle.h:413
void BaseParticle::setHandler ( ParticleHandler handler)

Sets the pointer to the particle's ParticleHandler.

Assigns the particle to a ParticleHandler, and assigns a species to it based on the particles indSpecies_ (BaseInteractable data member).

Parameters
[in]handlerpointer to the ParticleHandler

Definition at line 691 of file BaseParticle.cc.

References getHandler(), BaseInteractable::getIndSpecies(), BaseHandler< T >::getObject(), handler_, and setSpecies().

Referenced by ParticleHandler::addObject(), Chute::createBottom(), ParticleHandler::readObject(), setSpecies(), and ChuteBottom::setupInitialConditions().

692 {
693  handler_ = handler;
694  setSpecies(getHandler()->getDPMBase()->speciesHandler.getObject(getIndSpecies()));
695 }
ParticleHandler * handler_
Pointer to the particle's ParticleHandler container.
Definition: BaseParticle.h:390
void setSpecies(const ParticleSpecies *species)
ParticleHandler * getHandler() const
Returns pointer to the particle's ParticleHandler.
T * getObject(const unsigned int id)
Gets a pointer to the Object at the specified index in the BaseHandler.
Definition: BaseHandler.h:451
unsigned int getIndSpecies() const
Returns the index of the Species of this BaseInteractable.
void BaseParticle::setHGridLevel ( const unsigned int  level)

Sets the particle's HGrid level.

Sets the particle's HGridLevel_

Parameters
[in]levelthe particle's HGridLevel_

Definition at line 536 of file BaseParticle.cc.

References HGridLevel_.

Referenced by HGrid::insertParticleToHgrid().

537 {
538  HGridLevel_ = level;
539 }
unsigned int HGridLevel_
Cell position in the grid.
Definition: BaseParticle.h:401
void BaseParticle::setHGridNextObject ( BaseParticle p)

Sets the pointer to the next object in the particle's HGrid cell & level.

Sets the pointer to the next object in the particle's HGrid cell

Parameters
[in]ppointer to the next object

Definition at line 545 of file BaseParticle.cc.

References HGridNextObject_.

Referenced by Mercury2D::hGridRemoveParticle(), Mercury3D::hGridRemoveParticle(), Mercury3D::hGridUpdateParticle(), and Mercury2D::hGridUpdateParticle().

546 {
547  HGridNextObject_ = p;
548 }
BaseParticle * HGridNextObject_
Grid level for the object.
Definition: BaseParticle.h:402
void BaseParticle::setHGridPrevObject ( BaseParticle p)

Sets the pointer to the previous object in the particle's HGrid cell & level.

Sets the pointer to the previous object in the particle's HGrid cell

Parameters
[in]ppointer to the previous object

Definition at line 554 of file BaseParticle.cc.

References HGridPrevObject_.

Referenced by Mercury2D::hGridRemoveParticle(), Mercury3D::hGridRemoveParticle(), Mercury3D::hGridUpdateParticle(), and Mercury2D::hGridUpdateParticle().

555 {
556  HGridPrevObject_ = p;
557 }
BaseParticle * HGridPrevObject_
Pointer to the next Particle in the same HGrid cell.
Definition: BaseParticle.h:403
void BaseParticle::setHGridX ( const int  x)

Sets the particle's HGrid cell X-coordinate.

Set the x-index of the particle's HGrid cell position

Parameters
[in]xx-index of particle's HGrid cell

Definition at line 509 of file BaseParticle.cc.

References HGridX_.

Referenced by Mercury2D::hGridUpdateParticle(), and Mercury3D::hGridUpdateParticle().

510 {
511  HGridX_ = x;
512 }
int HGridX_
Hgrid attributes.
Definition: BaseParticle.h:400
void BaseParticle::setHGridY ( const int  y)

Sets the particle's HGrid cell Y-coordinate.

Set the y-index of the particle's HGrid cell position

Parameters
[in]yy-index of particle's HGrid cell

Definition at line 518 of file BaseParticle.cc.

References HGridY_.

Referenced by Mercury2D::hGridUpdateParticle(), and Mercury3D::hGridUpdateParticle().

519 {
520  HGridY_ = y;
521 }
void BaseParticle::setHGridZ ( const int  z)

Sets the particle's HGrid cell Z-coordinate.

Set the z-index of the particle's HGrid cell position

Parameters
[in]zz-index of particle's HGrid cell

Definition at line 527 of file BaseParticle.cc.

References HGridZ_.

Referenced by Mercury3D::hGridUpdateParticle().

528 {
529  HGridZ_ = z;
530 }
void BaseParticle::setIndSpecies ( unsigned int  indSpecies)
virtual
Deprecated:
Please use setSpecies(const ParticleSpecies*) instead.

Set the particle's species and species' index. Logs a warning if no ParticleHandler is assigned.

Parameters
[in]indSpeciesThe index of the species in the SpeciesHandler.
Todo:
TW do we have to update the species stored in the interactions here?

Reimplemented from BaseInteractable.

Definition at line 811 of file BaseParticle.cc.

References ERROR, BaseHandler< T >::getDPMBase(), BaseObject::getId(), BaseHandler< T >::getObject(), handler_, logger, BaseInteractable::setIndSpecies(), setSpecies(), and DPMBase::speciesHandler.

Referenced by oldRead().

812 {
813  if (handler_ != nullptr)
814  {
817  }
818  else
819  {
821  logger(ERROR, "setIndSpecies called on a particle with no particle handler.\n"
822  "Therefore I can't request the given species from the species handler.\n"
823  " PartID = %", getId());
824  }
825 }
unsigned int getId() const
Returns the unique identifier of any particular object.
Definition: BaseObject.cc:116
Logger< MERCURY_LOGLEVEL > logger("MercuryKernel")
ParticleHandler * handler_
Pointer to the particle's ParticleHandler container.
Definition: BaseParticle.h:390
void setSpecies(const ParticleSpecies *species)
T * getObject(const unsigned int id)
Gets a pointer to the Object at the specified index in the BaseHandler.
Definition: BaseHandler.h:451
SpeciesHandler speciesHandler
A handler to that stores the species type i.e. elastic, linear visco-elastic... et cetera...
Definition: DPMBase.h:991
virtual void setIndSpecies(unsigned int indSpecies)
Sets the index of the Species of this BaseInteractable.
DPMBase * getDPMBase()
Gets the problem that is solved using this handler.
Definition: BaseHandler.h:543
void BaseParticle::setInertia ( const Mdouble  newInertia)

Sets the particle's inertia_ (and adjusts invInertia_ accordingly)

Sets the particle's inertia_ and invInertia_.

Parameters
[in]newInertiathe new inertia_ to be set.

Definition at line 471 of file BaseParticle.cc.

References ERROR, inertia_, invInertia_, and logger.

472 {
473  if (newInertia >= 0)
474  {
475  inertia_ = newInertia;
476  invInertia_ = 1.0 / newInertia;
477  }
478  else
479  {
480  logger(ERROR, "Error in set_inertia (%)", newInertia);
481  }
482 }
Mdouble invInertia_
Particle inertia_.
Definition: BaseParticle.h:409
Logger< MERCURY_LOGLEVEL > logger("MercuryKernel")
Mdouble inertia_
Inverse Particle mass (for computation optimization)
Definition: BaseParticle.h:408
void BaseParticle::setInfiniteInertia ( )

Sets the particle's inertia_ to 'infinite' (1e20) and its invInertia_ to 0.

Sets the inertia_ to 1e20 and the invInertia_ (which is actually used in the calculations) to 0.

Definition at line 488 of file BaseParticle.cc.

References inertia_, and invInertia_.

489 {
490  inertia_ = 1e20;
491  invInertia_ = 0;
492 } //> i.e. no rotations
Mdouble invInertia_
Particle inertia_.
Definition: BaseParticle.h:409
Mdouble inertia_
Inverse Particle mass (for computation optimization)
Definition: BaseParticle.h:408
void BaseParticle::setMass ( const Mdouble  mass)

Sets the particle's mass.

Deprecated:
Please do not set the mass yourself, but use ParticleSpecies->computeMass instead. That makes sure

Sets the mass_ of the particle

Parameters
[in]massthe new particle's mass

Definition at line 590 of file BaseParticle.cc.

References ERROR, invMass_, logger, mass_, and WARN.

591 {
592  logger(WARN, "WARNING: Do not use particle->setMass, instead use "
593  "particleSpecies->computeMass, since this function can cause "
594  "inconsistencies between the mass, density and radius of this particle!");
595  if (mass >= 0.0)
596  {
597  if (invMass_ != 0.0) //InvMass=0 is a flag for a fixed particle
598  {
599  mass_ = mass;
600  invMass_ = 1.0 / mass;
601  }
602  }
603  else
604  {
605  logger(ERROR, "Error in BaseParticle::setMass, the given mass to be set must be positive.");
606  }
607 }
Logger< MERCURY_LOGLEVEL > logger("MercuryKernel")
Mdouble mass_
Pointer to the previous Particle in the same HGrid cell.
Definition: BaseParticle.h:406
Mdouble invMass_
Particle mass_.
Definition: BaseParticle.h:407
void BaseParticle::setMassForP3Statistics ( const Mdouble  mass)

Sets the particle's mass This function should not be used, but is necessary to extend the CG toolbox to non-spherical particles.

Sets the mass_ of the particle

Parameters
[in]massthe new particle's mass

Definition at line 613 of file BaseParticle.cc.

References ERROR, invMass_, logger, and mass_.

614 {
615  if(mass >= 0.0)
616  {
617  if(invMass_ != 0.0) //InvMass=0 is a flag for a fixed particle
618  {
619  mass_ = mass;
620  invMass_ = 1.0 / mass;
621  }
622  }
623  else
624  {
625  logger(ERROR, "Error in BaseParticle::setMass, the given mass to be set must be positive.");
626  }
627 }
Logger< MERCURY_LOGLEVEL > logger("MercuryKernel")
Mdouble mass_
Pointer to the previous Particle in the same HGrid cell.
Definition: BaseParticle.h:406
Mdouble invMass_
Particle mass_.
Definition: BaseParticle.h:407
void BaseParticle::setPeriodicFromParticle ( BaseParticle p)

Assigns the pointer to the 'original' particle this one's a periodic copy of.

Lets you set which particle this one is actually a periodic copy of (used in periodic boundary condition implementations).

Parameters
[in]ppointer to the 'original' particle this one is a periodic copy of.

Definition at line 500 of file BaseParticle.cc.

References periodicFromParticle_.

Referenced by LeesEdwardsBoundary::createHorizontalPeriodicParticles(), ShearBoxBoundary::createHorizontalPeriodicParticles(), CircularPeriodicBoundary::createPeriodicParticles(), MaserBoundary::createPeriodicParticles(), AngledPeriodicBoundary::createPeriodicParticles(), PeriodicBoundary::createPeriodicParticles(), LeesEdwardsBoundary::createVerticalPeriodicParticles(), and ShearBoxBoundary::createVerticalPeriodicParticles().

501 {
503 }
BaseParticle * periodicFromParticle_
Particle radius_.
Definition: BaseParticle.h:411
void BaseParticle::setPreviousPosition ( const Vec3D pos)

Sets the particle's position in the previous time step.

This is used to set the particle's previous position

Parameters
[in]posthe particle's previous position vector.

Definition at line 642 of file BaseParticle.cc.

References previousPosition_.

643 {
644  previousPosition_ = pos;
645 }
Vec3D previousPosition_
Displacement (only used in StatisticsVector, StatisticsPoint)
Definition: BaseParticle.h:414
void BaseParticle::setRadius ( const Mdouble  radius)

Sets the particle's radius_ (and adjusts the mass_ accordingly, based on the particle's species)

Sets the radius of the particle, and from that computes the new mass (using its species) and checks whether it now is either the smallest or biggest particle in its ParticleHandler.

Parameters
[in]radiusthe new radius

Definition at line 576 of file BaseParticle.cc.

References ParticleHandler::checkExtrema(), ParticleSpecies::computeMass(), getHandler(), BaseInteractable::getSpecies(), and radius_.

Referenced by Chute::createBottom(), ChuteInsertionBoundary::generateParticle(), CubeInsertionBoundary::generateParticle(), HopperInsertionBoundary::generateParticle(), helpers::loadingTest(), helpers::normalAndTangentialLoadingTest(), helpers::objectivenessTest(), FileReader::read(), and ChuteBottom::setupInitialConditions().

577 {
578  radius_ = radius;
579  if (getHandler())
580  {
581  getSpecies()->computeMass(this);
582  getHandler()->checkExtrema(this);
583  }
584 }
void checkExtrema(BaseParticle *P)
Checks if the extrema of this ParticleHandler needs updating.
const ParticleSpecies * getSpecies() const
Returns a pointer to the species of this BaseInteractable.
void computeMass(BaseParticle *p) const
Compute Particle mass function, which required a reference to the Species vector. It computes the Par...
ParticleHandler * getHandler() const
Returns pointer to the particle's ParticleHandler.
Mdouble radius_
Inverse Particle inverse inertia (for computation optimization)
Definition: BaseParticle.h:410
void BaseParticle::setSpecies ( const ParticleSpecies species)

In addition to the functionality of BaseInteractable::setSpecies, this function sets the pointer to the particleHandler, which is needed to retrieve species information.

Todo:
TW: this function should also check if the particle is the correct particle for the species type

Sets the particle's species. If this particle does not have a handler yet, this function also assigns the ParticleHandler in the same DPMBase as the SpeciesHandler of the given species as its handler.

Parameters
[in]speciespointer to the ParticleSpecies object, to be set as the particle's species.
Todo:
TW should we chaeck here if we have the right kind of species for the right kind of particle?

Definition at line 834 of file BaseParticle.cc.

References BaseHandler< T >::getDPMBase(), BaseSpecies::getHandler(), handler_, DPMBase::particleHandler, setHandler(), and BaseInteractable::setSpecies().

Referenced by MaserBoundary::addParticleToMaser(), MaserBoundary::checkBoundaryAfterParticleMoved(), Chute::createBottom(), helpers::loadingTest(), ChuteBottom::makeRoughBottom(), helpers::normalAndTangentialLoadingTest(), helpers::objectivenessTest(), FileReader::read(), DPMBase::readNextDataFile(), ParticleHandler::readObject(), setHandler(), setIndSpecies(), ChuteBottom::setupInitialConditions(), Chute::setupInitialConditions(), and ChuteWithHopper::setupInitialConditions().

835 {
838  //set pointer to the ParticleHandler handler_, which is needed to retrieve
839  //species information
840  //\todo maybe these if statements should throw warnings
841  if (handler_ == nullptr)
842  {
843  SpeciesHandler* sH = species->getHandler();
844  DPMBase* dB = sH->getDPMBase();
845  if (dB != nullptr)
846  {
848  }
849  }
850 }
Container to store all ParticleSpecies.
The DPMBase header includes quite a few header files, defining all the handlers, which are essential...
Definition: DPMBase.h:65
void setHandler(ParticleHandler *handler)
Sets the pointer to the particle's ParticleHandler.
ParticleHandler * handler_
Pointer to the particle's ParticleHandler container.
Definition: BaseParticle.h:390
SpeciesHandler * getHandler() const
Returns the pointer to the handler to which this species belongs.
Definition: BaseSpecies.cc:74
void setSpecies(const ParticleSpecies *species)
Sets the species of this BaseInteractable.
ParticleHandler particleHandler
An object of the class ParticleHandler, contains the pointers to all the particles created...
Definition: DPMBase.h:1001
DPMBase * getDPMBase()
Gets the problem that is solved using this handler.
Definition: BaseHandler.h:543
void BaseParticle::unfix ( )

Unfix Particle function, which required a reference to the Species vector. It unfixes a Particle by computing the Particles mass and inertia.

Unfixes the particle by computing the Particles mass and inertia, using the species and radius.

Definition at line 170 of file BaseParticle.cc.

References ParticleSpecies::computeMass(), BaseInteractable::getSpecies(), and invMass_.

Referenced by DPMBase::readNextDataFile().

171 {
172  invMass_ = 1.0;
173  getSpecies()->computeMass(this);
174 }
const ParticleSpecies * getSpecies() const
Returns a pointer to the species of this BaseInteractable.
Mdouble invMass_
Particle mass_.
Definition: BaseParticle.h:407
void computeMass(BaseParticle *p) const
Compute Particle mass function, which required a reference to the Species vector. It computes the Par...
void BaseParticle::write ( std::ostream &  os) const
virtual

Particle print function, which accepts an std::ostream as input.

BaseParticle print method, which accepts an os std::ostream as input. It prints human readable BaseParticle information to the std::ostream.

Parameters
[in,out]osstream to which the info is written

Implements BaseInteractable.

Reimplemented in LiquidFilmParticle, and ThermalParticle.

Definition at line 182 of file BaseParticle.cc.

References invInertia_, invMass_, radius_, and BaseInteractable::write().

Referenced by LiquidFilmParticle::write(), and ThermalParticle::write().

183 {
185  os << " radius " << radius_
186  << " invMass " << invMass_
187  << " invInertia " << invInertia_;
188 }
Mdouble invInertia_
Particle inertia_.
Definition: BaseParticle.h:409
Mdouble invMass_
Particle mass_.
Definition: BaseParticle.h:407
Mdouble radius_
Inverse Particle inverse inertia (for computation optimization)
Definition: BaseParticle.h:410
virtual void write(std::ostream &os) const =0
Write a BaseInteractable to an output stream.

Friends And Related Function Documentation

friend class ParticleSpecies
friend

Particle's position at previous time step.

Since ParticleSpecies is allowed to set the mass of a BaseParticle, it is a friend of this class.

Todo:
Is it possible to only make the function ParticleSpecies::computeMass a friend?

Definition at line 421 of file BaseParticle.h.

Member Data Documentation

Vec3D BaseParticle::displacement_
private

Pointer to originating Particle.

Definition at line 413 of file BaseParticle.h.

Referenced by addDisplacement(), BaseParticle(), getDisplacement(), and setDisplacement().

ParticleHandler* BaseParticle::handler_
private

Pointer to the particle's ParticleHandler container.

Definition at line 390 of file BaseParticle.h.

Referenced by BaseParticle(), getHandler(), getVolume(), setHandler(), setIndSpecies(), and setSpecies().

unsigned int BaseParticle::HGridLevel_
private

Cell position in the grid.

Definition at line 401 of file BaseParticle.h.

Referenced by BaseParticle(), getHGridLevel(), printHGrid(), and setHGridLevel().

BaseParticle* BaseParticle::HGridNextObject_
private

Grid level for the object.

Definition at line 402 of file BaseParticle.h.

Referenced by BaseParticle(), getHGridNextObject(), and setHGridNextObject().

BaseParticle* BaseParticle::HGridPrevObject_
private

Pointer to the next Particle in the same HGrid cell.

Definition at line 403 of file BaseParticle.h.

Referenced by BaseParticle(), getHGridPrevObject(), and setHGridPrevObject().

int BaseParticle::HGridX_
private

Hgrid attributes.

Definition at line 400 of file BaseParticle.h.

Referenced by BaseParticle(), getHGridX(), printHGrid(), and setHGridX().

int BaseParticle::HGridY_
private

Definition at line 400 of file BaseParticle.h.

Referenced by BaseParticle(), getHGridY(), printHGrid(), and setHGridY().

int BaseParticle::HGridZ_
private

Definition at line 400 of file BaseParticle.h.

Referenced by BaseParticle(), getHGridZ(), printHGrid(), and setHGridZ().

Mdouble BaseParticle::inertia_
private

Inverse Particle mass (for computation optimization)

Definition at line 408 of file BaseParticle.h.

Referenced by BaseParticle(), ParticleSpecies::computeMass(), fixParticle(), getInertia(), oldRead(), read(), setInertia(), and setInfiniteInertia().

Mdouble BaseParticle::invInertia_
private
Mdouble BaseParticle::invMass_
private

Particle mass_.

Todo:
TW: why do we need to store mass and inertia; we already have invMass and invInertia; can we take it out?

Definition at line 407 of file BaseParticle.h.

Referenced by BaseParticle(), ParticleSpecies::computeMass(), fixParticle(), getInvMass(), isFixed(), oldRead(), read(), setMass(), setMassForP3Statistics(), unfix(), and write().

Mdouble BaseParticle::mass_
private

Pointer to the previous Particle in the same HGrid cell.

Particle attributes

Definition at line 406 of file BaseParticle.h.

Referenced by BaseParticle(), ParticleSpecies::computeMass(), fixParticle(), getKineticEnergy(), getMass(), oldRead(), read(), setMass(), and setMassForP3Statistics().

BaseParticle* BaseParticle::periodicFromParticle_
private

Particle radius_.

Definition at line 411 of file BaseParticle.h.

Referenced by BaseParticle(), getPeriodicFromParticle(), and setPeriodicFromParticle().

Vec3D BaseParticle::previousPosition_
private

Displacement (only used in StatisticsVector, StatisticsPoint)

Definition at line 414 of file BaseParticle.h.

Referenced by getPreviousPosition(), movePrevious(), and setPreviousPosition().

Mdouble BaseParticle::radius_
private

Inverse Particle inverse inertia (for computation optimization)

Definition at line 410 of file BaseParticle.h.

Referenced by BaseParticle(), getInteractionRadius(), getRadius(), getVolume(), getWallInteractionRadius(), read(), setRadius(), and write().


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