MercuryDPM  Trunk
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
InteractionHandler.cc
Go to the documentation of this file.
1 //Copyright (c) 2013-2020, 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 #include <iostream>
26 #include <iomanip>
27 #include <string>
28 #include <cstdlib>
29 #include "Particles/BaseParticle.h"
30 #include "InteractionHandler.h"
31 #include "SpeciesHandler.h"
32 #include "Species/BaseSpecies.h"
34 #include "DPMBase.h"
35 #include "MpiDataClass.h"
36 #ifdef MERCURY_USE_OMP
37 #include "omp.h"
38 #endif
39 
44 {
46  logger(DEBUG, "InteractionHandler::InteractionHandler() finished");
47 }
48 
57 {
58  writeVTK_ = IH.writeVTK_;
59  //By default interactions are not copied.
60  logger(DEBUG, "InteractionHandler::InteractionHandler(const "
61  "InteractionHandler &IH) finished, please note that no interactions"
62  " have been copied.");
63 }
64 
69 {
70  if (this != &rhs)
71  {
72  clear();
73  }
74  writeVTK_ = rhs.writeVTK_;
75  logger(DEBUG, "InteractionHandler::operator =(const InteractionHandler& rhs) finished.");
76  return *this;
77 }
78 
80 {
81  logger(DEBUG, "InteractionHandler::~InteractionHandler() finished");
82 }
83 
88 {
89  //set the particleHandler pointer
90  I->setHandler(this);
91  //Puts the interaction in the Interaction list
93  I->getP()->addInteraction(I);
94  I->getI()->addInteraction(I);
95 }
96 
104 {
105  //for particle-particle collision it is assumed BaseInteractable P has a lower index then I, so we only have to check for I, not P
106  {
107  for (unsigned j = 0; j < P->getInteractions().size(); ++j) {
108  auto i = P->getInteractions()[j];
109  if (i->getI() == I) {
110  return i;
111  }
112  }
113  }
114  return nullptr;
115 }
116 
118 {
120  BaseInteraction* C = species->getNewInteraction(P, I, timeStamp);
121  addObject(C);
122  return C;
123 }
124 
126  #ifdef MERCURY_USE_OMP
127  newObjects_.clear();
128  newObjects_.resize(getDPMBase()->getNumberOfOMPThreads());
129  #endif
130 }
131 
133  for (auto objects : newObjects_) {
134  for (auto object : objects) {
135  addObject(object);
136  }
137  }
138 }
139 
149 InteractionHandler::getInteraction(BaseInteractable* const P, BaseInteractable* const I, const unsigned timeStamp)
150 {
151  const BaseSpecies* const species = getDPMBase()->speciesHandler.getMixedObject(P->getIndSpecies(),
152  I->getIndSpecies());
153 
155  if (C == nullptr) {
156  C = species->getNewInteraction(P, I, timeStamp);
157  #ifdef MERCURY_USE_OMP
158  if (omp_get_num_threads()>1) {
159  newObjects_[omp_get_thread_num()].push_back(C);
160  C->setHandler(this);
161  } else {
162  addObject(C);
163  }
164  #else
165  addObject(C);
166  #endif
167  }
168 
169  //set timeStamp
170  C->setTimeStamp(timeStamp);
171 
173  //set species of collision
174  C->setSpecies(species);
175 
176  return C;
177 }
178 
186 {
187  //NOTE: assuming the interaction type of a species is the same for all species
188  BaseSpecies* species = this->getDPMBase()->speciesHandler.getObject(0);
189  BaseInteraction* emptyInteraction = species->getEmptyInteraction();
190 
191  return emptyInteraction;
192 }
193 
201 {
202  BaseSpecies* species = this->getDPMBase()->speciesHandler.getObject(0);
203  species->deleteEmptyInteraction(interaction);
204 }
205 
214 void* InteractionHandler::createMPIInteractionDataArray(unsigned int numberOfInteractions) const
215 {
216  BaseInteraction* emptyInteraction = this->createEmptyInteraction();
217  //create a vector based on the first interaction. This is to avoid templating the whole class
218  void* historyData = emptyInteraction->createMPIInteractionDataArray(numberOfInteractions);
219  this->deleteEmptyInteraction(emptyInteraction);
220 
221  return historyData;
222 }
223 
231 {
232  BaseInteraction* emptyInteraction = this->createEmptyInteraction();
233  emptyInteraction->deleteMPIInteractionDataArray(dataArray);
234  this->deleteEmptyInteraction(emptyInteraction);
235 }
236 
248 void InteractionHandler::getInteractionDetails(void* interactionData, unsigned int index, unsigned int& identificationP,
249  unsigned int& identificationI, bool& isWallInteraction,
250  unsigned& timeStamp)
251 {
252  BaseInteraction* emptyInteraction = this->createEmptyInteraction();
253  emptyInteraction->getInteractionDetails(interactionData, index, identificationP, identificationI, isWallInteraction,
254  timeStamp);
255  this->deleteEmptyInteraction(emptyInteraction);
256 }
257 
268 {
269  BaseInteraction* c = nullptr;
270  Mdouble dNormalMax = 0.9; //minimal agreement required
271  for (const auto& i : P->getInteractions())
272  {
273  if (i->getI() == I && i->getTimeStamp() != timeStamp)
274  {
275  const Mdouble dNormal = Vec3D::dot(normal, i->getNormal());
276  if (dNormal > dNormalMax)
277  {
278  c = i;
279  dNormalMax = dNormal;
280  }
281  }
282  }
283 
285 
286  if (c)
287  {
288  c->setTimeStamp(timeStamp);
289  }
290  else
291  {
292  c = species->getNewInteraction(P, I, timeStamp);
293  c->setSpecies(species);
294  addObject(c);
295  //logger(INFO,"new interaction t=%: i=%, n=%",timeStamp,getNumberOfObjects(),normal);
296  }
297 
298  return c;
299 }
300 
320 {
321  BaseInteraction* iMain = getObject(id);
322 
323  BaseParticle* P = dynamic_cast<BaseParticle*>(iMain->getP());
324  BaseParticle* I = dynamic_cast<BaseParticle*>(iMain->getI());
325  if (P != nullptr && I != nullptr) //check that both P and I are particles (not walls)
326  {
328  if (realI != nullptr && !P->getPeriodicFromParticle()) //check that P is a real and I is a ghost particle
329  {
330  if (P->getIndex() < realI->getIndex())
331  {
332  BaseInteraction* iOther = getExistingInteraction(P, realI);
333  //You have to also check for existing interactions of the particles in reverse order since the copySwitchPointer function can revert the order of the particles
335  if (iOther == nullptr)
336  {
337  iOther = getExistingInteraction(realI, P);
338  }
339  if (iOther != nullptr) //if the interaction existed before the ghost particles were created
340  {
341  //Here we decide which of the two interactions should be kept:
342  //the interaction between the two real particles (iMain), or
343  //the interaction between the real and a ghost particle (iOther).
344  //It picks the one for which a collision has happened,
345  //i.e. the one with the newer timeStamp.
347  if (iOther->getTimeStamp() <
348  iMain->getTimeStamp()) //if the interaction has been active during the last computeForce routine, make this the new (real) interaction.
349 
350  {
351  iMain->setI(realI);
352  removeObject(iOther->getIndex());
353  return;
354  }
355  else //if the interaction has not been active during the last computeForce routine
356  {
358  return;
359  }
360  }
361  else //if the interaction has been created during the last computeForce routine, make this a new (real) interaction.
362  {
363  iMain->setI(realI);
364  return;
365  }
366  }
367  }
368  }
369  //this point is reached if either P or I are walls, or P and I are both ghost particles; in these cases, the interaction gets deleted
371 }
372 
380 void InteractionHandler::eraseOldInteractions(unsigned currentNTime)
381 {
382  //std::cout<<"void InteractionHandler::eraseOldInteractions(Mdouble lastTimeStep)"<<std::endl;
383  //std::cout<<"Current interactions="<<getNumberOfObjects()<<std::endl;
384  //Remove_if reconstructs the vector with only elements passing the check_spring_time function
385  //Erase removes the end of the vector
387  for (unsigned int id = 0; id < getNumberOfObjects(); id++)
388  {
389  if (getObject(id)->getTimeStamp() <= currentNTime)
390  {
391  getObject(id)->actionsOnErase();
392  removeObject(id);
393  --id;
394  }
395  }
396 }
397 
399 {
400  for (auto i: *this)
401  {
402  i->actionsAfterTimeStep();
403  }
404 }
405 
410 {
411  Mdouble sum = 0;
412  Mdouble n = 0;
413  for (BaseInteraction* const p : objects_)
414  {
415  if (p->getOverlap() > 0)
416  {
417  sum += p->getOverlap();
418  n++;
419  }
420  }
421  return sum / n;
422 }
423 
427 std::string InteractionHandler::getName() const
428 {
429  return "InteractionHandler";
430 }
431 
436 void InteractionHandler::write(std::ostream& os) const
437 {
438 #ifdef MERCURY_USE_MPI
439  //note: this function only prints the particles on processor 0
440  //The rest of the particles are processed in the restart file write function
441  MPIContainer& communicator = MPIContainer::Instance();
442  unsigned int totalNumberOfInteractions = getNumberOfObjects();
443  os << "Interactions " << totalNumberOfInteractions << '\n';
444  //logger(INFO,"% interactions",getNumberOfObjects());
445  for (BaseInteraction* it : *this)
446  {
447  os << (*it) << '\n';
448  }
449 #else
450  os << "Interactions " << getNumberOfObjects() << '\n';
451  for (BaseInteraction* i : objects_)
452  os << (*i) << '\n';
453 #endif
454 }
455 
457 void InteractionHandler::read(std::istream& is)
458 {
459  clear();
460  unsigned int N;
461  std::string dummy;
462  is >> dummy;
463  std::stringstream line;
465  line >> N;
466  if (N > 1e5) logger(INFO, "Reading % %", N, dummy);
467  logger(VERBOSE, "In %::read(is): reading in % objects.", getName(), N);
469  // set map
470  particleById.clear();
471  for (BaseParticle* p : getDPMBase()->particleHandler) {
472  particleById[p->getId()] = p;
473  }
474  wallById.clear();
475  for (BaseWall* w : getDPMBase()->wallHandler) {
476  wallById[w->getId()] = w;
477  }
478  for (unsigned int i = 0; i < N; i++) {
479  readAndAddObject(is);
480  }
481  particleById.clear();
482  wallById.clear();
483 }
484 
485 
490 {
491  std::string type, dummy, idType;
492  unsigned int id0, id1;
493  Mdouble doubleTimeStamp;
494  unsigned timeStamp;
495 
497  is >> type;
498  logger(DEBUG, "InteractionHandler::readAndAddObject(is): reading type %.", type);
499  Mdouble timeStampDouble;
500  is >> idType >> id0 >> id1 >> dummy >> timeStampDouble;
501  timeStamp = timeStampDouble; //in order to read old restart files
502 
504  BaseInteraction* C;
505 
506 #ifdef MERCURY_USE_MPI
507  if (idType.compare("particleIds") == 0)
508  {
509  std::vector<BaseParticle*> list0 = getDPMBase()->particleHandler.getObjectsById(id0);
510  std::vector<BaseParticle*> list1 = getDPMBase()->particleHandler.getObjectsById(id1);
511  for (int p0 = 0; p0 < list0.size(); p0++)
512  {
513  for (int p1 = 0; p1 < list1.size(); p1++)
514  {
515  C = getInteraction(list0[p0],list1[p1], timeStamp);
516  if (C != nullptr)
517  {
518  is >> *C;
519  }
520  }
521  }
522 
523  }
524  else
525  {
526  std::vector<BaseParticle*> list0 = getDPMBase()->particleHandler.getObjectsById(id0);
527  for (int p0 = 0; p0 < list0.size(); p0++)
528  {
529  C = getInteraction(list0[p0],getDPMBase()->wallHandler.getObjectById(id1), timeStamp);
530  if (C != nullptr)
531  {
532  is >> *C;
533  }
534  }
535  }
536 #else
537  //this requires particleById/wallById to be set, which is done in the read() function
538  if (idType == "particleIds")
539  C = getInteraction(particleById[id0],particleById[id1],timeStamp);
540  else
541  C = getInteraction(particleById[id0],wallById[id1], timeStamp);
542  is >> (*C);
543 #endif
544  is.ignore(std::numeric_limits<std::streamsize>::max(), '\n');
545 }
546 
548 {
549  writeVTK_ = fileType;
550 }
551 
553 {
554  return writeVTK_;
555 }
556 
558  double liquidVolume = 0;
559  for (auto i : objects_) {
560  auto j = dynamic_cast<LiquidMigrationWilletInteraction*>(i);
561  if (j and !static_cast<BaseParticle*>(j->getP())->isMPIParticle()) liquidVolume += j->getLiquidBridgeVolume();
562  }
563  return getMPISum(liquidVolume);
564 };
This class contains all information and functions required for communication between processors...
Definition: MpiContainer.h:125
InteractionHandler & operator=(const InteractionHandler &rhs)
Assignment operator.
BaseInteraction * createEmptyInteraction() const
Creates an empty interaction.
unsigned int getIndex() const
Returns the index of the object in the handler.
Definition: BaseObject.h:118
virtual void deleteEmptyInteraction(BaseInteraction *interaction) const =0
std::enable_if<!std::is_pointer< typename U::MixedSpeciesType >::value, typename U::MixedSpeciesType * >::type getMixedObject(const U *S, const U *T)
void eraseOldInteractions(unsigned)
erases interactions which have an old timestamp.
virtual void actionsOnErase()
If an interaction needs to do something before it gets erased, add it here. E.g. Liquid bridges ruptu...
static MPIContainer & Instance()
fetch the instance to be used for communication
Definition: MpiContainer.h:130
Defines the liquid bridge willet interaction between two particles or walls.
BaseSpecies is the class from which all other species are derived.
Definition: BaseSpecies.h:49
virtual void getInteractionDetails(void *interactionDataArray, unsigned int index, unsigned int &identificationP, unsigned int &identificationI, bool &isWallInteraction, unsigned &timeStamp)
void write(std::ostream &os) const
Writes the InteractionHandler to an output stream, for example a restart file.
FileType getWriteVTK() const
virtual void * createMPIInteractionDataArray(unsigned int numberOfInteractions) const
Logger< MERCURY_LOGLEVEL > logger("MercuryKernel")
double Mdouble
Definition: GeneralDefine.h:34
void addInteraction(BaseInteraction *I)
Adds an interaction to this BaseInteractable.
void * createMPIInteractionDataArray(unsigned int numberOfInteractions) const
creates an empty MPIInteractionDataArray
std::map< unsigned, BaseParticle * > particleById
void read(std::istream &is)
void setStorageCapacity(const unsigned int N)
Sets the storage capacity of this BaseHandler.
void setI(BaseInteractable *I)
Sets the second object involved in the interaction (often particle or wall).
const std::complex< Mdouble > i
Definition: ExtendedMath.h:50
std::vector< T * > getObjectsById(const unsigned int id)
Gets a vector of pointers to the objects with the specific id.
Definition: BaseHandler.h:593
void removeObjectKeepingPeriodics(unsigned int id)
Removes interactions of periodic particles when the periodic particles get deleted (see DPMBase::remo...
~InteractionHandler() override
Destructor, it destructs the InteractionHandler and all BaseInteraction it contains.
InteractionHandler()
Default constructor, it creates an empty InteractionHandler.
BaseParticle * getPeriodicFromParticle() const
Returns the 'original' particle this one's a periodic copy of.
Definition: BaseParticle.h:338
void deleteMPIInteractionDataArray(void *dataArray)
deletes an MPIInteractionDataArray
BaseInteractable * getI()
Returns a pointer to the second object involved in the interaction (often a wall or a particle)...
static Mdouble dot(const Vec3D &a, const Vec3D &b)
Calculates the dot product of two Vec3D: .
Definition: Vector.cc:76
void setTimeStamp(unsigned timeStamp)
Updates the time step of the interacting. Note, time steps used to find completed interactions...
Mdouble getTimeStamp() const
Returns an Mdouble which is the time stamp of the interaction.
Mdouble getMeanOverlap() const
The mean overlap of all interactions.
Stores information about interactions between two interactable objects; often particles but could be ...
std::string getName() const override
Returns the name of the object.
FileType
With FileType options, one is able to choose if data is to be read/written from/into no or single or ...
Definition: File.h:40
void deleteEmptyInteraction(BaseInteraction *interaction) const
Deletes an empty interaction.
BaseInteractable * getP()
Returns a pointer to first object involved in the interaction (normally a particle).
void getInteractionDetails(void *interactionData, unsigned int index, unsigned int &identificationP, unsigned int &identificationI, bool &isWallInteraction, unsigned &timeStamp)
reads the basic interaction details from an MPIInteractionDataArray
file will not be created/read
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:424
BaseInteraction * getInteraction(BaseInteractable *P, BaseInteractable *I, unsigned timeStamp)
Returns the Interaction between the BaseInteractable's P and I.
virtual void removeObject(unsigned const int index)
Removes an Object from the BaseHandler.
ParticleHandler particleHandler
An object of the class ParticleHandler, contains the pointers to all the particles created...
Definition: DPMBase.h:1329
Container to store Interaction objects.
BaseInteraction * addInteraction(BaseInteractable *P, BaseInteractable *I, unsigned timeStamp)
T * getObject(const unsigned int id)
Gets a pointer to the Object at the specified index in the BaseHandler.
Definition: BaseHandler.h:613
void setHandler(InteractionHandler *handler)
Sets the pointer to the interaction hander which is storing this interaction.
virtual BaseInteraction * getEmptyInteraction() const =0
#define UNUSED
Definition: GeneralDefine.h:39
double getLiquidBridgeVolume() const
std::vector< std::vector< BaseInteraction * > > newObjects_
void setSpecies(const BaseSpecies *species)
Set the Species of the interaction; note this can either be a Species or MixedSpecies.
Vec3D getMPISum(Vec3D &val)
Sums the values over all processors using MPI_reduce.
void addObject(BaseInteraction *I) override
Adds an Interaction to the InteractionHandler.
std::vector< BaseInteraction * > objects_
The actual list of Object pointers.
Definition: BaseHandler.h:302
Container to store the pointers to all objects that one creates in a simulation.
Definition: BaseHandler.h:50
Basic class for walls.
Definition: BaseWall.h:47
virtual unsigned int getNumberOfObjects() const
Gets the number of real Object in this BaseHandler. (i.e. no mpi or periodic particles) ...
SpeciesHandler speciesHandler
A handler to that stores the species type i.e. LinearViscoelasticSpecies, etc.
Definition: DPMBase.h:1319
std::map< unsigned, BaseWall * > wallById
void setWriteVTK(FileType f)
const std::vector< BaseInteraction * > & getInteractions() const
Returns a list of interactions which belong to this interactable.
virtual void addObject(T *object)
Adds a new Object to the BaseHandler.
Definition: BaseHandler.h:431
Defines the basic properties that a interactable object can have.
virtual BaseInteraction * getNewInteraction(BaseInteractable *P, BaseInteractable *I, unsigned timeStamp) const =0
returns new Interaction object.
unsigned int getIndSpecies() const
Returns the index of the species associated with the interactable object.
Definition: Vector.h:49
DPMBase * getDPMBase()
Gets the problem that is solved using this handler.
void readAndAddObject(std::istream &is) override
Reads an Interaction into the InteractionHandler from restart data.
BaseInteraction * getExistingInteraction(const BaseInteractable *P, const BaseInteractable *I) const
Returns the Interaction between the BaseInteractable's P and I if it exists, otherwise returns a null...
virtual void clear()
Empties the whole BaseHandler by removing all Objects and setting all other variables to 0...
virtual void deleteMPIInteractionDataArray(void *dataArray)