revision: v0.14
DPMBase.h
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 
26 #ifndef DPMBase_H
27 #define DPMBase_H
28 
29 //so that the user doesn't have to include string/io manipulations:
30 #include <string>
31 #include <iomanip>
32 // include file for signal handler types
33 #include <csignal>
34 //The vector class contains a 3D vector class.
35 #include "Math/Vector.h"
36 //This class defines the particle handler
37 #include "ParticleHandler.h"
38 //This class defines the base particle (such that not every Driver has to include it)
39 #include "Particles/BaseParticle.h"
41 //This class defines the wall handler
42 #include "WallHandler.h"
43 //This class defines the boundary handler
44 #include "BoundaryHandler.h"
46 #include "DomainHandler.h"
47 //This class defines the interaction handler
48 #include "InteractionHandler.h"
49 //This class defines the Species handler
50 #include "SpeciesHandler.h"
51 //This class defines the cg handler
52 #include "CG/CGHandler.h"
53 //This class defines the possibleContact lists
54 #ifdef CONTACT_LIST_HGRID
55 #include "PossibleContactList.h"
56 #endif
57 //This class defines the random number generator
58 #include "Math/RNG.h"
59 #include "Domain.h"
64 #include "MercuryTime.h"
65 
66 
75 class DPMBase
76 {
81 public:
82 
88  void constructor();
89 
93  DPMBase();
94 
98  DPMBase(const DPMBase& other);
99 
103  virtual ~DPMBase();
104 
109  static void incrementRunNumberInFile();
110 
114  static int readRunNumberFromFile();
115 
120  void autoNumber();
121 
126  std::vector<int> get1DParametersFromRunNumber(int size_x) const;
127 
132  std::vector<int> get2DParametersFromRunNumber(int size_x, int size_y) const;
133 
138  std::vector<int> get3DParametersFromRunNumber(int size_x, int size_y, int size_z) const;
139 
143  int launchNewRun(const char* name, bool quick = false);
144 
145  //setters and getters
146 
150  void setRunNumber(int runNumber);
151 
155  int getRunNumber() const;
156 
160 // void sendParticlesToRoot(unsigned int &numberOfLocalParticles, unsigned int processorId, ParticleHandler &tempParticleHandler) const;
161 
165  virtual void decompose();
166 
171  void solve();
172 
177  virtual void computeOneTimeStep();
178 
183  void checkSettings();
184 
188  void forceWriteOutputFiles();
189 
194  virtual void writeOutputFiles();
195 
199  void solve(int argc, char* argv[]);
200 
211  virtual void setupInitialConditions();
212 
217  virtual void writeXBallsScript() const;
218 
222  virtual Mdouble getInfo(const BaseParticle& P) const;
223 
225 
233  virtual void writeRestartFile();
234 
235  void writeDataFile();
236 
237  void writeEneFile();
238 
239  void writeFStatFile();
240 
241  void fillDomainWithParticles(unsigned N=50);
242 
243  enum class ReadOptions : int {
244  ReadAll,
247  };
248 
255 
260  int readRestartFile(std::string fileName, ReadOptions opt = ReadOptions::ReadAll);
261 
262 // /*!
263 // * \brief Loads all MD data and plots statistics for all time steps in the .data file
264 // */
265 // void statisticsFromRestartData(const char *name);
267 
271  virtual void write(std::ostream& os, bool writeAllParticles = true) const;
272 
277  virtual void read(std::istream& is, ReadOptions opt = ReadOptions::ReadAll);
278 
282  virtual BaseWall* readUserDefinedWall(const std::string& type) const
283  { return nullptr; }
284 
289  virtual void readOld(std::istream& is);
290 
294  bool readDataFile(std::string fileName = "", unsigned int format = 0);
295 
299  bool readParAndIniFiles(std::string fileName);
300 
306  bool readNextDataFile(unsigned int format = 0);
307 
311  void readNextFStatFile();
312 
316  bool findNextExistingDataFile(Mdouble tMin, bool verbose = true);
317 
322  bool readArguments(int argc, char* argv[]);
323 
327  virtual bool readNextArgument(int& i, int argc, char* argv[]);
328 
332  virtual bool checkParticleForInteraction(const BaseParticle& P);
333 
337  virtual bool checkParticleForInteractionLocal(const BaseParticle& P);
338 
340 
342 
348 
349  //getters and setters
350 
357  File& getDataFile();
358 
365  File& getEneFile();
366 
373  File& getFStatFile();
374 
381  File& getRestartFile();
382 
389  File& getStatFile();
390 
395 
402  const File& getDataFile() const;
403 
410  const File& getEneFile() const;
411 
418  const File& getFStatFile() const;
419 
426  const File& getRestartFile() const;
427 
434  const File& getStatFile() const;
435 
439  const File& getInteractionFile() const;
440 
444  const std::string& getName() const;
445 
449  void setName(const std::string& name);
450 
454  void setName(const char* name);
455 
459  void setSaveCount(unsigned int saveCount);
460 
464  void setFileType(FileType fileType);
465 
469  void setOpenMode(std::fstream::openmode openMode);
470 
471 
472  //other member functions
473 
477  void resetFileCounter();
478 
482  void closeFiles();
483 
487  void setLastSavedTimeStep(unsigned int nextSavedTimeStep);
488 
492  Mdouble getTime() const;
493 
497  Mdouble getNextTime() const;
498 
503  unsigned int getNumberOfTimeSteps() const;
504 
508  void setTime(Mdouble time);
509 
513  void setTimeMax(Mdouble newTMax);
514 
518  Mdouble getTimeMax() const;
519 
523  void setLogarithmicSaveCount(Mdouble logarithmicSaveCountBase);
524 
528  void setNToWrite(int nToWrite);
529 
533  int getNToWrite() const;
534 
535 #ifdef CONTACT_LIST_HGRID
536 
539  PossibleContactList& getPossibleContactList();
540 #endif
541 
552  void setRotation(bool rotation)
553  { rotation_ = rotation; }
554 
559  bool getRotation() const
560  { return rotation_; }
561 
565  void setWallsWriteVTK(FileType writeWallsVTK);
566 
570  void setWallsWriteVTK(bool);
571 
575  void setInteractionsWriteVTK(bool);
576 
580  void setParticlesWriteVTK(bool writeParticlesVTK);
581 
582  void setSuperquadricParticlesWriteVTK(bool writeSuperquadricParticlesVTK);
583 
587  FileType getWallsWriteVTK() const;
588 
592  bool getParticlesWriteVTK() const;
593 
595 
600  Mdouble getXMin() const
601  { return min_.x(); }
602 
607  Mdouble getXMax() const
608  { return max_.x(); }
609 
613  Mdouble getYMin() const
614  { return min_.y(); }
615 
619  Mdouble getYMax() const
620  { return max_.y(); }
621 
625  Mdouble getZMin() const
626  { return min_.z(); }
627 
631  Mdouble getZMax() const
632  { return max_.z(); }
633 
634  /*
635  * \brief Returns the minimum coordinates of the problem domain.
636  */
637  Vec3D getMin() const
638  { return min_; }
639 
640  /*
641  * \brief Returns the maximum coordinates of the problem domain.
642  */
643  Vec3D getMax() const
644  { return max_; }
645 
650  void setXMin(Mdouble newXMin);
651 
656  void setYMin(Mdouble newYMin);
657 
662  void setZMin(Mdouble newZMin);
663 
668  void setXMax(Mdouble newXMax);
669 
674  void setYMax(Mdouble newYMax);
675 
680  void setZMax(Mdouble newZMax);
681 
685  void setMax(const Vec3D& max);
686 
690  void setMax(Mdouble, Mdouble, Mdouble);
691 
695  void setDomain(const Vec3D& min, const Vec3D& max);
696 
700  void setMin(const Vec3D& min);
701 
705  void setMin(Mdouble, Mdouble, Mdouble);
706 
707 
711  void setTimeStep(Mdouble newDt);
712 
716  Mdouble getTimeStep() const;
717 
718  /* Sets the number of omp threads */
719  void setNumberOfOMPThreads(int numberOfOMPThreads);
720 
721  /* Returns the number of omp threads */
722  int getNumberOfOMPThreads() const;
723 
727  void setXBallsColourMode(int newCMode);
728 
732  int getXBallsColourMode() const;
733 
737  void setXBallsVectorScale(double newVScale);
738 
742  double getXBallsVectorScale() const;
743 
747  void setXBallsAdditionalArguments(std::string newXBArgs);
748 
752  std::string getXBallsAdditionalArguments() const;
753 
758  void setXBallsScale(Mdouble newScale);
759 
763  double getXBallsScale() const;
764 
768  void setGravity(Vec3D newGravity);
769 
773  Vec3D getGravity() const;
774 
778  void setDimension(unsigned int newDim);
779 
783  void setSystemDimensions(unsigned int newDim);
784 
788  unsigned int getSystemDimensions() const;
789 
793  void setParticleDimensions(unsigned int particleDimensions);
794 
798  unsigned int getParticleDimensions() const;
799 
804  std::string getRestartVersion() const;
805 
809  void setRestartVersion(std::string newRV);
810 
814  bool getRestarted() const;
815 
819  void setRestarted(bool newRestartedFlag);
820 
824  bool getAppend() const;
825 
829  void setAppend(bool newAppendFlag);
830 
834  Mdouble getElasticEnergy() const;
835 
839  Mdouble getKineticEnergy() const;
840 
845 
850 
851  Mdouble getTotalEnergy() const;
852 
856  Mdouble getTotalMass() const;
857 
861  Vec3D getCentreOfMass() const;
862 
866  Vec3D getTotalMomentum() const;
867 
871  double getCPUTime() { return clock_.getCPUTime(); }
872 
876  double getWallTime() { return clock_.getWallTime(); }
877 
878 
882  static bool areInContact(const BaseParticle* pI, const BaseParticle* pJ);
883 
885 
888  virtual void hGridInsertParticle(BaseParticle* obj UNUSED);
889 
893  virtual void hGridUpdateParticle(BaseParticle* obj UNUSED);
894 
898  virtual void hGridRemoveParticle(BaseParticle* obj UNUSED);
899 
903  virtual void hGridUpdateMove(BaseParticle*, Mdouble);
904 
908  bool mpiIsInCommunicationZone(BaseParticle* particle);
909 
914 
919 
923  void updateGhostGrid(BaseParticle* P);
924 
929  virtual void gatherContactStatistics(unsigned int index1, int index2, Vec3D Contact, Mdouble delta, Mdouble ctheta,
930  Mdouble fdotn, Mdouble fdott, Vec3D P1_P2_normal_, Vec3D P1_P2_tangential);
931 
935  void setNumberOfDomains(std::vector<unsigned> direction);
936 
937  enum class DomainSplit {X, Y, Z, XY, XZ, YZ, XYZ};
938 
944  void splitDomain(DomainSplit domainSplit);
945 
949  std::vector<unsigned> getNumberOfDomains();
950 
955 
956  void removeOldFiles() const;
957 
961  virtual void hGridGetInteractingParticleList(BaseParticle* obj, std::vector<BaseParticle*>& list)
962  {};
963 
964  virtual void computeWallForces(BaseWall* w);
965 
969  virtual bool getHGridUpdateEachTimeStep() const;
970 
972  void setMeanVelocity(Vec3D V_mean_goal);
973 
975  void setMeanVelocityAndKineticEnergy(Vec3D V_mean_goal, Mdouble Ek_goal);
976 
978  Mdouble getTotalVolume() const;
979 
981  Matrix3D getKineticStress() const;
982 
984  Matrix3D getStaticStress() const;
985 
987  Matrix3D getTotalStress() const;
988 
992  virtual void handleParticleRemoval(unsigned int id);
993 
997  virtual void handleParticleAddition(unsigned int id, BaseParticle* p);
998 
999  //functions that should only be used in the class definitions
1000 protected:
1001 
1006  virtual void computeAllForces();
1007 
1012  virtual void computeInternalForces(BaseParticle*);
1013 
1019 
1023  virtual void computeExternalForces(BaseParticle*);
1024 
1029 
1034  virtual void actionsOnRestart();
1035 
1040  virtual void actionsBeforeTimeLoop();
1041 
1046  virtual void hGridActionsBeforeTimeLoop();
1047 
1052  virtual void hGridActionsBeforeTimeStep();
1053 
1058  virtual void actionsBeforeTimeStep();
1059 
1068  virtual void computeAdditionalForces() {}
1069 
1074  virtual void actionsAfterSolve();
1075 
1080  virtual void actionsAfterTimeStep();
1081 
1082  void writeVTKFiles() const;
1083 
1088  virtual void outputXBallsData(std::ostream& os) const;
1089 
1094  virtual void outputXBallsDataParticle(unsigned int i, unsigned int format, std::ostream& os) const;
1095 
1099  virtual void writeEneHeader(std::ostream& os) const;
1100 
1104  virtual void writeFstatHeader(std::ostream& os) const;
1105 
1109  virtual void writeEneTimeStep(std::ostream& os) const;
1110 
1111  // Functions for statistics
1115  virtual void initialiseStatistics();
1116 
1120  virtual void outputStatistics();
1121 
1125  void gatherContactStatistics();
1126 
1130  virtual void processStatistics(bool);
1131 
1135  virtual void finishStatistics();
1136 
1143  virtual void integrateBeforeForceComputation();
1144 
1148  virtual void integrateAfterForceComputation();
1149 
1155  virtual void checkInteractionWithBoundaries();
1156 
1160  virtual void hGridActionsBeforeIntegration();
1161 
1165  virtual void hGridActionsAfterIntegration();
1166 
1170  void setFixedParticles(unsigned int n);
1171 
1175  virtual void printTime() const;
1176 
1180  virtual bool continueSolve() const;
1181 
1186  void outputInteractionDetails() const;
1187 
1191  bool isTimeEqualTo(Mdouble time) const;
1192 
1197 
1204 
1211 
1212  void deleteGhostParticles(std::set<BaseParticle*>& particlesToBeDeleted);
1213 
1214  void synchroniseParticle(BaseParticle*, unsigned fromProcessor = 0);
1215 
1220 
1221 
1225  static void signalHandler(int signal);
1226 
1230  void setSoftStop();
1231 
1232 private:
1233 
1237  static volatile sig_atomic_t continueFlag_;
1238 
1239  /*The number of openmp (symmetric multiprocessing threads)*/
1241 
1245  unsigned int systemDimensions_;
1246 
1250  unsigned int particleDimensions_;
1251 
1256 
1260  std::vector<unsigned> numberOfDomains_;
1261 
1267 
1272 
1276  unsigned int numberOfTimeSteps_;
1277 
1282 
1287 
1293  std::string restartVersion_;
1294 
1299 
1304  bool append_;
1305 
1312 
1317 
1322 
1324 
1326 
1328 
1330 
1332 
1333  //This is the private data that is only used by the xballs output
1334 
1340 
1345 
1350 
1355 
1360 
1364  std::string name_;
1365 
1366  // defines a Macro for creating an instance of class PossibleContactList. See PossibleContactList.h
1367 #ifdef CONTACT_LIST_HGRID
1368  PossibleContactList possibleContactList;
1369 #endif
1370 
1375 
1380 
1381 public:
1386 
1391 
1396 
1401 
1406 
1411 
1416 
1421 
1426 
1427 
1432 
1437 
1442 
1447 
1452 
1457 
1463 
1468 
1470 };
1471 
1472 #endif
BaseHandler::setDPMBase
void setDPMBase(DPMBase *DPMBase)
Sets the problem that is solved using this handler.
Definition: BaseHandler.h:718
DPMBase::setMax
void setMax(const Vec3D &max)
Sets the maximum coordinates of the problem domain.
Definition: DPMBase.cc:1073
ParticleHandler::getNumberOfRealObjects
unsigned int getNumberOfRealObjects() const
Returns the number of real objects (on all processors)
Definition: ParticleHandler.cc:1294
DPMBase::xBallsScale_
Mdouble xBallsScale_
sets the xballs argument scale (see XBalls/xballs.txt)
Definition: DPMBase.h:1349
DPMBase::processStatistics
virtual void processStatistics(bool)
Definition: DPMBase.cc:1909
LinearViscoelasticSlidingFrictionSpecies.h
DPMBase::getInteractionFile
File & getInteractionFile()
Return a reference to the file InteractionsFile.
Definition: DPMBase.cc:343
DomainHandler::addNewParticles
void addNewParticles()
Definition: DomainHandler.cc:429
BaseParticle::getInteractionWith
BaseInteraction * getInteractionWith(BaseParticle *P, unsigned timeStamp, InteractionHandler *interactionHandler) override
Checks if particle is in interaction with given particle P, and if so, returns vector of pointer to t...
Definition: BaseParticle.cc:678
BoundaryHandler
Container to store pointers to all BaseBoundary objects.
Definition: BoundaryHandler.h:39
SuperQuadricParticleVtkWriter.h
mathsFunc::square
T square(const T val)
squares a number
Definition: ExtendedMath.h:106
DPMBase::rotation_
bool rotation_
A flag to turn on/off particle rotation. true will enable particle rotation. false will disable parti...
Definition: DPMBase.h:1311
DPMBase::setName
void setName(const std::string &name)
Allows to set the name of all the files (ene, data, fstat, restart, stat)
Definition: DPMBase.cc:420
SlidingFrictionSpecies.h
File::getName
const std::string & getName() const
Allows to access the file name, e.g., "problem.data".
Definition: File.cc:165
File::close
void close()
Closes the file by calling fstream_.close()
Definition: File.cc:407
InteractionHandler::addNewObjectsOMP
void addNewObjectsOMP()
Definition: InteractionHandler.cc:132
BaseInteractable::getSpecies
const ParticleSpecies * getSpecies() const
Returns a pointer to the species of this BaseInteractable.
Definition: BaseInteractable.h:108
DPMBase::computeExternalForces
virtual void computeExternalForces(BaseParticle *)
Computes the external forces, such as gravity, acting on particles.
Definition: DPMBase.cc:3074
DPMBase::setRunNumber
void setRunNumber(int runNumber)
This sets the counter/Run number, overriding the defaults.
Definition: DPMBase.cc:603
File::setFileType
void setFileType(FileType fileType)
Sets the type of file needed to write into or read from. File::fileType_.
Definition: File.cc:215
DPMBase::readRestartFile
bool readRestartFile(ReadOptions opt=ReadOptions::ReadAll)
Reads all the particle data corresponding to a given, existing . restart file (for more details regar...
Definition: DPMBase.cc:2921
DPMBase::setTimeStep
void setTimeStep(Mdouble newDt)
Sets a new value for the simulation time step.
Definition: DPMBase.cc:1225
File::openWriteNoAppend
bool openWriteNoAppend(unsigned)
Definition: File.cc:398
DPMBase::getParticleDimensions
unsigned int getParticleDimensions() const
Returns the particle dimensionality.
Definition: DPMBase.cc:1458
DPMBase::nToWrite_
int nToWrite_
number of elements to write to a screen
Definition: DPMBase.h:1379
DPMBase::append_
bool append_
A flag to determine if the file has to be appended or not. See DPMBase::Solve() for example.
Definition: DPMBase.h:1304
ParticleVtkWriter
Definition: ParticleVtkWriter.h:34
DPMBase::setMin
void setMin(const Vec3D &min)
Sets the minimum coordinates of the problem domain.
Definition: DPMBase.cc:1109
to_string_padded
std::string to_string_padded(unsigned int value)
Pads the number This function tries to pad the number to 4 digits, which is used when you create mult...
Definition: File.cc:44
DPMBase::computeOneTimeStep
virtual void computeOneTimeStep()
Performs everything needed for one time step, used in the time-loop of solve().
Definition: DPMBase.cc:4139
DPMBase::readRunNumberFromFile
static int readRunNumberFromFile()
Read the run number or the counter from the counter file (COUNTER_DONOTDEL)
Definition: DPMBase.cc:550
ParticleHandler::getLargestInteractionRadius
Mdouble getLargestInteractionRadius() const
Returns the largest interaction radius.
Definition: ParticleHandler.cc:768
DPMBase::DomainSplit::YZ
@ YZ
DPMBase::getSuperquadricParticlesWriteVTK
bool getSuperquadricParticlesWriteVTK() const
Definition: DPMBase.cc:984
DPMBase::getTimeStep
Mdouble getTimeStep() const
Returns the simulation time step.
Definition: DPMBase.cc:1241
DPMBase::getVtkWriter
ParticleVtkWriter * getVtkWriter() const
Definition: DPMBase.cc:5130
CGHandler
Container that stores all CG objects.
Definition: CGHandler.h:65
LinearViscoelasticSlidingFrictionSpecies
Species< LinearViscoelasticNormalSpecies, SlidingFrictionSpecies > LinearViscoelasticSlidingFrictionSpecies
Definition: LinearViscoelasticSlidingFrictionSpecies.h:34
File::saveCurrentTimeStepNoFileTypeCheck
bool saveCurrentTimeStepNoFileTypeCheck(unsigned int ntimeSteps)
Definition: File.cc:317
Vector.h
DPMBase::hGridActionsBeforeTimeLoop
virtual void hGridActionsBeforeTimeLoop()
A virtual function that allows one to carry out hGrid operations before the start of the time loop.
Definition: DPMBase.cc:1667
PeriodicBoundaryHandler
Container to store pointers to all BasePeriodicBoundary objects.
Definition: PeriodicBoundaryHandler.h:46
File::setlogarithmicSaveCount
void setlogarithmicSaveCount(const Mdouble logarithmicSaveCountBase)
the function to set the user input base of logarithmic saving count
Definition: File.cc:283
DPMBase::getCPUTime
double getCPUTime()
Definition: DPMBase.h:871
PeriodicBoundary.h
DPMBase::mpiIsInCommunicationZone
bool mpiIsInCommunicationZone(BaseParticle *particle)
Checks if the position of the particle is in an mpi communication zone or not.
Definition: DPMBase.cc:1755
DPMBase::writePythonFileForVTKVisualisation
void writePythonFileForVTKVisualisation() const
Definition: DPMBase.cc:2155
DPMBase::setXBallsColourMode
void setXBallsColourMode(int newCMode)
Set the xballs output mode.
Definition: DPMBase.cc:1291
DomainHandler.h
ParticleHandler::getMomentum
Vec3D getMomentum() const
Definition: ParticleHandler.cc:662
PossibleContactList.h
Flusher::FLUSH
@ FLUSH
DPMBase::actionsBeforeTimeLoop
virtual void actionsBeforeTimeLoop()
A virtual function. Allows one to carry out any operations before the start of the time loop.
Definition: DPMBase.cc:1660
DPMBase::clock_
Time clock_
record when the simulation started
Definition: DPMBase.h:1467
InteractionHandler
Container to store Interaction objects.
Definition: InteractionHandler.h:44
ParticleHandler::actionsAfterTimeStep
void actionsAfterTimeStep()
Definition: ParticleHandler.cc:1362
DPMBase::readSpeciesFromDataFile_
bool readSpeciesFromDataFile_
Determines if the last column of the data file is interpreted as the info parameter during restart.
Definition: DPMBase.h:1374
DPMBase::setParticlesWriteVTK
void setParticlesWriteVTK(bool writeParticlesVTK)
Sets whether particles are written in a VTK file.
Definition: DPMBase.cc:934
SpeciesHandler::clear
void clear() override
Empties the whole BaseHandler by removing all Objects and setting all other variables to 0.
Definition: SpeciesHandler.h:54
BoundaryHandler::setWriteVTK
void setWriteVTK(bool writeVTK)
Definition: BoundaryHandler.h:88
DPMBase::readDataFile
bool readDataFile(std::string fileName="", unsigned int format=0)
This allows particle data to be reloaded from data files.
Definition: DPMBase.cc:2346
Quaternion::setAngleZ
void setAngleZ(Mdouble psi)
Converts the rotation angle in the XY plane into a quaternion (for Mercury2D). See Wikipedia for deta...
Definition: Quaternion.cc:492
BaseParticle::fixParticle
void fixParticle()
Fix Particle function. It fixes a Particle by setting its inverse mass and inertia and velocities to ...
Definition: BaseParticle.cc:165
DPMBase::DomainSplit::XY
@ XY
DPMBase::actionsAfterTimeStep
virtual void actionsAfterTimeStep()
A virtual function which allows to define operations to be executed after time step.
Definition: DPMBase.cc:1869
DPMBase::getRestarted
bool getRestarted() const
Returns the flag denoting if the simulation was restarted or not.
Definition: DPMBase.cc:1484
PeriodicBoundary::shiftPosition
virtual void shiftPosition(BaseParticle *p) const override
shifts the particle
Definition: PeriodicBoundary.cc:219
DPMBase::incrementRunNumberInFile
static void incrementRunNumberInFile()
Increment the run Number (counter value) stored in the file_counter (COUNTER_DONOTDEL) by 1 and store...
Definition: DPMBase.cc:625
DPMBase::getStaticStress
Matrix3D getStaticStress() const
Calculate the static stress tensor in the system averaged over the whole volume.
Definition: DPMBase.cc:5237
BaseParticle::getPeriodicFromParticle
BaseParticle * getPeriodicFromParticle() const
Returns the 'original' particle this one's a periodic copy of.
Definition: BaseParticle.h:341
DPMBase::getXBallsAdditionalArguments
std::string getXBallsAdditionalArguments() const
Returns the additional arguments for xballs.
Definition: DPMBase.cc:1346
InteractionHandler::resetNewObjectsOMP
void resetNewObjectsOMP()
Definition: InteractionHandler.cc:125
DPMBase::restartVersion_
std::string restartVersion_
Previous versions of MercuryDPM had a different restart file format, the below member variable allows...
Definition: DPMBase.h:1293
DPMBase::setParticleDimensions
void setParticleDimensions(unsigned int particleDimensions)
Sets the particle dimensionality.
Definition: DPMBase.cc:1439
DPMBase::setLogarithmicSaveCount
void setLogarithmicSaveCount(Mdouble logarithmicSaveCountBase)
Sets File::logarithmicSaveCount_ for all files (ene, data, fstat, restart, stat)
Definition: DPMBase.cc:5283
ParticleHandler::getRotationalEnergy
Mdouble getRotationalEnergy() const
Definition: ParticleHandler.cc:582
DPMBase::setSoftStop
void setSoftStop()
function for setting sigaction constructor.
Definition: DPMBase.cc:5003
DPMBase::hGridUpdateParticle
virtual void hGridUpdateParticle(BaseParticle *obj UNUSED)
Definition: DPMBase.cc:1695
DPMBase::signalHandler
static void signalHandler(int signal)
signal handler function.
Definition: DPMBase.cc:4977
DPMBase::read
virtual void read(std::istream &is, ReadOptions opt=ReadOptions::ReadAll)
Reads all data from a restart file, e.g. domain data and particle data.
Definition: DPMBase.cc:3531
DPMBase::hGridActionsBeforeIntegration
virtual void hGridActionsBeforeIntegration()
This function has to be called before integrateBeforeForceComputation.
Definition: DPMBase.cc:1931
DPMBase::writeWallsVTK_
FileType writeWallsVTK_
A flag to turn on/off the vtk writer for walls.
Definition: DPMBase.h:1316
Domain::isInCommunicationZone
bool isInCommunicationZone(BaseParticle *particle)
Check if the particle is in the communication zone of the current domain.
Definition: Domain.cc:441
YAXIS
@ YAXIS
Definition: GeneralDefine.h:78
DPMBase::launchNewRun
int launchNewRun(const char *name, bool quick=false)
This launches a code from within this code. Please pass the name of the code to run.
Definition: DPMBase.cc:772
DPMBase::getNextTime
Mdouble getNextTime() const
Returns the current simulation time.
Definition: DPMBase.cc:813
DPMBase::getXMax
Mdouble getXMax() const
If the length of the problem domain in x-direction is XMax - XMin, then getXMax() returns XMax.
Definition: DPMBase.h:607
BoundaryHandler::readAndAddObject
void readAndAddObject(std::istream &is) final
Reads BaseBoundary into the BoundaryHandler from restart data.
Definition: BoundaryHandler.cc:209
BaseInteractable::setPosition
void setPosition(const Vec3D &position)
Sets the position of this BaseInteractable.
Definition: BaseInteractable.h:239
Time::tic
void tic()
This is like a start button of a stopwatch. Assigns the variable start with the current number of clo...
Definition: MercuryTime.h:55
DPMBase::getParticlesWriteVTK
bool getParticlesWriteVTK() const
Returns whether particles are written in a VTK file.
Definition: DPMBase.cc:976
DPMBase::getStatFile
MERCURY_DEPRECATED File & getStatFile()
The non const version. Allows to edit the File::statFile.
Definition: DPMBase.cc:335
BaseInteractable::addTorque
void addTorque(const Vec3D &addTorque)
Adds an amount to the torque on this BaseInteractable.
Definition: BaseInteractable.cc:132
DPMBase::readOld
virtual void readOld(std::istream &is)
Reads all data from a restart file, e.g. domain data and particle data; old version.
Definition: DPMBase.cc:3774
BaseWall
Basic class for walls.
Definition: BaseWall.h:48
logger
Logger< MERCURY_LOGLEVEL > logger("MercuryKernel")
Definition of different loggers with certain modules. A user can define its own custom logger here.
PeriodicBoundaryHandler::addNewParticle
void addNewParticle(BaseParticle *particle)
Adds a new particle to the periodic list.
Definition: PeriodicBoundaryHandler.cc:333
DPMBase::hGridRemoveParticle
virtual void hGridRemoveParticle(BaseParticle *obj UNUSED)
Definition: DPMBase.cc:1702
DPMBase::checkInteractionWithBoundaries
virtual void checkInteractionWithBoundaries()
There are a range of boundaries one could implement depending on ones' problem. This methods checks f...
Definition: DPMBase.cc:3177
ParticleVtkWriter.h
DPMBase::random
RNG random
This is a random generator, often used for setting up the initial conditions etc.....
Definition: DPMBase.h:1390
PossibleContact::getP2
BaseParticle * getP2()
Gets a pointer to the second BaseParticle in this PossibleContact.
Definition: PossibleContact.h:139
RNG::read
void read(std::istream &is)
Definition: RNG.cc:58
DPMBase::cgHandler
CGHandler cgHandler
Object of the class cgHandler.
Definition: DPMBase.h:1431
DPMBase::numberOfTimeSteps_
unsigned int numberOfTimeSteps_
Stores the number of time steps.
Definition: DPMBase.h:1276
DPMBase::getInfo
virtual Mdouble getInfo(const BaseParticle &P) const
A virtual function that returns some user-specified information about a particle.
Definition: DPMBase.cc:1633
ParticleHandler::getNumberOfObjects
unsigned int getNumberOfObjects() const override
Returns the number of objects in the container. In parallel code this practice is forbidden to avoid ...
Definition: ParticleHandler.cc:1315
Vec3D::X
Mdouble X
the vector components
Definition: Vector.h:65
ParticleHandler::getMassTimesPosition
Vec3D getMassTimesPosition() const
Definition: ParticleHandler.cc:632
InteractionHandler::addInteraction
BaseInteraction * addInteraction(BaseInteractable *P, BaseInteractable *I, unsigned timeStamp)
Definition: InteractionHandler.cc:117
DPMBase::periodicBoundaryHandler
PeriodicBoundaryHandler periodicBoundaryHandler
Internal handler that deals with periodic boundaries, especially in a parallel build.
Definition: DPMBase.h:1415
DPMBase::setMeanVelocity
void setMeanVelocity(Vec3D V_mean_goal)
This function will help you set a fixed kinetic energy and mean velocity in your system.
Definition: DPMBase.cc:5141
Vec3D::dot
static Mdouble dot(const Vec3D &a, const Vec3D &b)
Calculates the dot product of two Vec3D: .
Definition: Vector.cc:76
File::getCounter
unsigned int getCounter() const
In case of multiple files, File::getCounter() returns the the number (FILE::Counter_) of the next fil...
Definition: File.cc:223
SpeciesHandler.h
DPMBase::synchroniseParticle
void synchroniseParticle(BaseParticle *, unsigned fromProcessor=0)
Definition: DPMBase.cc:4946
DPMBase::DPMBase
DPMBase()
Constructor that calls the "void constructor()".
Definition: DPMBase.cc:193
PeriodicBoundaryHandler::updateStatus
void updateStatus(std::set< BaseParticle * > &ghostParticlesToBeDeleted)
Updates the positions/velocity of ghost particles and accordingly the status of these particles.
Definition: PeriodicBoundaryHandler.cc:151
WallHandler::actionsAfterParticleGhostUpdate
void actionsAfterParticleGhostUpdate()
Calls the method actionsAfterParticleGhostUpdate of every wall in the handler.
Definition: WallHandler.cc:479
DPMBase::outputStatistics
virtual void outputStatistics()
Definition: DPMBase.cc:1884
DPMBase::setYMax
void setYMax(Mdouble newYMax)
Sets the value of YMax, the upper bound of the problem domain in the y-direction.
Definition: DPMBase.cc:1182
DPMBase::constructor
void constructor()
A function which initialises the member variables to default values, so that the problem can be solve...
Definition: DPMBase.cc:206
DPMBase::wallVTKWriter_
WallVTKWriter wallVTKWriter_
Definition: DPMBase.h:1327
initialiseMPI
void initialiseMPI()
Inialises the MPI library.
Definition: MpiContainer.cc:137
BaseHandler::read
void read(std::istream &is)
Reads all objects from restart data.
Definition: BaseHandler.h:543
BaseInteraction::setFStatData
void setFStatData(std::fstream &fstat, BaseParticle *P, BaseWall *I)
Definition: BaseInteraction.cc:947
copyDataFromMPIParticleToParticle
void copyDataFromMPIParticleToParticle(MPIParticle *bP, BaseParticle *p, ParticleHandler *particleHandler)
Copies data from an MPIParticle class to a BaseParticle and sets the particleHandler and species.
Definition: MpiDataClass.cc:105
getSVNURL
const std::string getSVNURL()
Definition: CMakeDefinitions.cc:46
DPMBase::ReadOptions
ReadOptions
Definition: DPMBase.h:243
DPMBase::splitDomain
void splitDomain(DomainSplit domainSplit)
Definition: DPMBase.cc:5059
DomainHandler::initialise
void initialise()
Definition: DomainHandler.cc:434
BaseHandler::getSize
unsigned int getSize() const
Gets the size of the particleHandler (including mpi and periodic particles)
Definition: BaseHandler.h:655
MPISphericalParticle::copyDataFromParticleToMPIParticle
void copyDataFromParticleToMPIParticle(BaseParticle *p)
Definition: MpiDataClass.cc:131
DPMBase::writeXBallsScript
virtual void writeXBallsScript() const
This writes a script which can be used to load the xballs problem to display the data just generated.
File::setSaveCount
void setSaveCount(unsigned int saveCount)
Sets File::saveCount_.
Definition: File.cc:273
BaseParticle::integrateBeforeForceComputation
void integrateBeforeForceComputation(double time, double timeStep)
First step of Velocity Verlet integration.
Definition: BaseParticle.cc:708
Domain::cleanCommunicationLists
void cleanCommunicationLists()
Removes nullptrs from boundaryParticleList_ and boundaryParticleListNeighbour_.
Definition: Domain.cc:1742
DPMBase::restarted_
bool restarted_
A bool to check if the simulation was restarted or not, ie. if setupInitialConditionsShould be run an...
Definition: DPMBase.h:1298
DomainHandler
Container to store all Domain.
Definition: DomainHandler.h:47
ParticleHandler.h
DPMBase::checkParticleForInteractionLocal
virtual bool checkParticleForInteractionLocal(const BaseParticle &P)
Checks if a particle P has any interaction with walls or other particles in the local domain.
Definition: DPMBase.cc:4738
DPMBase::getNumberOfOMPThreads
int getNumberOfOMPThreads() const
Definition: DPMBase.cc:1277
DPMBase::outputInteractionDetails
void outputInteractionDetails() const
Displays the interaction details corresponding to the pointer objects in the interaction handler.
Definition: DPMBase.cc:5019
DPMBase::writeEneTimeStep
virtual void writeEneTimeStep(std::ostream &os) const
Write the global kinetic, potential energy, etc. in the system.
Definition: DPMBase.cc:2095
DPMBase::integrateBeforeForceComputation
virtual void integrateBeforeForceComputation()
Update particles' and walls' positions and velocities before force computation.
Definition: DPMBase.cc:3133
DPMBase::computeForcesDueToWalls
void computeForcesDueToWalls(BaseParticle *, BaseWall *)
Computes the forces on the particles due to the walls (normals are outward normals)
Definition: DPMBase.cc:3094
SphericalParticleVtkWriter
Definition: SphericalParticleVtkWriter.h:35
FileType::NO_FILE
@ NO_FILE
file will not be created/read
DPMBase::restartFile
File restartFile
An instance of class File to handle in- and output into a .restart file.
Definition: DPMBase.h:1451
DPMBase::readParAndIniFiles
bool readParAndIniFiles(std::string fileName)
Allows the user to read par.ini files (useful to read files produced by the MDCLR simulation code - e...
Definition: DPMBase.cc:2382
DPMBase::writeRestartFile
virtual void writeRestartFile()
Stores all the particle data for current save time step to a "restart" file, which is a file simply i...
Definition: DPMBase.cc:2854
PARTICLE
@ PARTICLE
Definition: MpiContainer.h:67
DPMBase::numberOfOMPThreads_
int numberOfOMPThreads_
Definition: DPMBase.h:1240
MPIContainer::getProcessorID
std::size_t getProcessorID()
Reduces a scalar on all processors to one scalar on a target processor.
Definition: MpiContainer.cc:113
DPMBase::closeFiles
void closeFiles()
Closes all files (ene, data, fstat, restart, stat) that were opened to read or write.
Definition: DPMBase.cc:500
DPMBase::get3DParametersFromRunNumber
std::vector< int > get3DParametersFromRunNumber(int size_x, int size_y, int size_z) const
This turns a counter into 3 indices, which is a useful feature for performing a 3D parameter study....
Definition: DPMBase.cc:732
DPMBase::handleParticleRemoval
virtual void handleParticleRemoval(unsigned int id)
Handles the removal of particles from the particleHandler.
Definition: DPMBase.cc:5298
INFO
LL< Log::INFO > INFO
Info log level.
Definition: Logger.cc:55
DPMBase::setWallsWriteVTK
void setWallsWriteVTK(FileType writeWallsVTK)
Sets whether walls are written into a VTK file.
Definition: DPMBase.cc:908
Time::getCPUTime
Mdouble getCPUTime()
Definition: MercuryTime.h:76
BoundaryVTKWriter
Definition: BoundaryVTKWriter.h:34
MpiContainer.h
ParticleSpecies
Definition: ParticleSpecies.h:37
BaseHandler::setStorageCapacity
void setStorageCapacity(const unsigned int N)
Sets the storage capacity of this BaseHandler.
Definition: BaseHandler.h:669
ParticleHandler::removeGhostObject
void removeGhostObject(unsigned int index)
Removes a BaseParticle from the ParticleHandler without a global check, this is only to be done for m...
Definition: ParticleHandler.cc:418
BaseBoundary
Definition: BaseBoundary.h:49
CubeInsertionBoundary
It's an insertion boundary which has cuboidal shape (yes, 'CuboidalInsertionBoundary' would have been...
Definition: CubeInsertionBoundary.h:42
DPMBase::paoloParticleHandler
ParticleHandler paoloParticleHandler
Fake particleHandler created by Paolo needed temporary by just Paolo.
Definition: DPMBase.h:1400
DPMBase::getDataFile
MERCURY_DEPRECATED File & getDataFile()
The non const version. Allows one to edit the File::dataFile.
Definition: DPMBase.cc:303
DPMBase::hGridUpdateMove
virtual void hGridUpdateMove(BaseParticle *, Mdouble)
Definition: DPMBase.cc:1924
InsertionBoundary::insertParticles
void insertParticles(DPMBase *md)
Fill a certain domain with particles.
Definition: InsertionBoundary.cc:328
DPMBase::setFixedParticles
void setFixedParticles(unsigned int n)
Sets a number, n, of particles in the particleHandler as "fixed particles".
Definition: DPMBase.cc:1951
BaseWall::checkInteractions
virtual void checkInteractions(InteractionHandler *interactionHandler, unsigned int timeStamp)
Check if all interactions are valid.
Definition: BaseWall.h:219
BaseParticle::setRadius
virtual void setRadius(Mdouble radius)
Sets the particle's radius_ (and adjusts the mass_ accordingly, based on the particle's species)
Definition: BaseParticle.cc:542
DPMBase::writeOutputFiles
virtual void writeOutputFiles()
Writes simulation data to all the main Mercury files: .data, .ene, .fstat, .xballs and ....
Definition: DPMBase.cc:3881
Vec3D
Definition: Vector.h:50
DPMBase::getWallsWriteVTK
FileType getWallsWriteVTK() const
Returns whether walls are written in a VTK file.
Definition: DPMBase.cc:965
DPMBase::solve
void solve()
The work horse of the code.
Definition: DPMBase.cc:4003
DPMBase::interactionFile
File interactionFile
File class to handle in- and output into .interactions file. This file hold information about interac...
Definition: DPMBase.h:1462
BaseInteractable::setOrientation
void setOrientation(const Quaternion &orientation)
Sets the orientation of this BaseInteractable.
Definition: BaseInteractable.h:260
WallVTKWriter.h
DPMBase::computeAdditionalForces
virtual void computeAdditionalForces()
A virtual function which allows to define operations to be executed prior to the OMP force collect.
Definition: DPMBase.h:1068
DPMBase::fStatFile
File fStatFile
An instance of class File to handle in- and output into a .fstat file.
Definition: DPMBase.h:1441
DPMBase::numberOfDomains_
std::vector< unsigned > numberOfDomains_
Vector containing the number of domains in x-,y- and z-direction, required for parallel computations.
Definition: DPMBase.h:1260
Interaction.h
verbose
bool verbose
Definition: statXZ.cpp:28
DPMBase::handleParticleAddition
virtual void handleParticleAddition(unsigned int id, BaseParticle *p)
Definition: DPMBase.cc:5313
CGHandler::evaluate
void evaluate()
Contains the code executed at each time step.
Definition: CGHandler.cc:97
DPMBase::writeFstatHeader
virtual void writeFstatHeader(std::ostream &os) const
Writes a header with a certain format for FStat file.
Definition: DPMBase.cc:2034
DPMBase::DomainSplit::Z
@ Z
DPMBase::outputXBallsDataParticle
virtual void outputXBallsDataParticle(unsigned int i, unsigned int format, std::ostream &os) const
This function writes out the particle locations into an output stream in a format the XBalls program ...
BaseInteraction
Stores information about interactions between two interactable objects; often particles but could be ...
Definition: BaseInteraction.h:60
DPMBase::getYMin
Mdouble getYMin() const
If the length of the problem domain in y-direction is YMax - YMin, then getYMin() returns YMin.
Definition: DPMBase.h:613
File
Definition: File.h:81
DPMBase::getCurrentDomain
Domain * getCurrentDomain()
Function that returns a pointer to the domain corresponding to the processor.
Definition: DPMBase.cc:5125
WallHandler
Container to store all BaseWall.
Definition: WallHandler.h:43
DPMBase::vtkWriter_
ParticleVtkWriter * vtkWriter_
Definition: DPMBase.h:1325
DPMBase::setZMax
void setZMax(Mdouble newZMax)
Sets the value of ZMax, the upper bound of the problem domain in the z-direction.
Definition: DPMBase.cc:1208
File::open
bool open()
Checks if the file stream fstream_ has any issues while opening. Alongside, it also increments the ne...
Definition: File.cc:347
DPMBase::xBallsColourMode_
int xBallsColourMode_
XBalls is a package to view the particle data. As an alternative MercuryDPM also supports ParaView....
Definition: DPMBase.h:1339
Mdouble
double Mdouble
Definition: GeneralDefine.h:34
WallHandler::readAndAddObject
void readAndAddObject(std::istream &is) final
Create a new wall in the WallHandler, based on the information provided in a restart file.
Definition: WallHandler.cc:281
DPMBase::setNToWrite
void setNToWrite(int nToWrite)
set the number of elements to write to the screen
Definition: DPMBase.cc:850
PeriodicBoundaryHandler::cleanCommunicationLists
void cleanCommunicationLists()
Definition: PeriodicBoundaryHandler.cc:1775
DPMBase::boundaryHandler
BoundaryHandler boundaryHandler
An object of the class BoundaryHandler which concerns insertion and deletion of particles into or fro...
Definition: DPMBase.h:1410
DPMBase::setXBallsScale
void setXBallsScale(Mdouble newScale)
Sets the scale of the view (either normal, zoom in or zoom out) to display in xballs....
Definition: DPMBase.cc:1354
DPMBase::runNumber_
int runNumber_
This stores the run number for saving.
Definition: DPMBase.h:1359
DPMBase::fillDomainWithParticles
void fillDomainWithParticles(unsigned N=50)
Definition: DPMBase.cc:2899
InteractionHandler::actionsAfterTimeStep
void actionsAfterTimeStep()
Definition: InteractionHandler.cc:398
DPMBase::autoNumber
void autoNumber()
The autoNumber() function calls three functions: setRunNumber(), readRunNumberFromFile() and incremen...
Definition: DPMBase.cc:536
BaseInteractable::sumForceTorqueOMP
void sumForceTorqueOMP()
Definition: BaseInteractable.cc:162
Matrix3D
Implementation of a 3D matrix.
Definition: Matrix.h:38
DomainHandler::setInteractionDistance
void setInteractionDistance(Mdouble interactionDistance)
Sets the interaction distance of the domain handler.
Definition: DomainHandler.cc:299
NUMBER_OF_PROCESSORS
#define NUMBER_OF_PROCESSORS
For the MPI communication routines this quantity is often required. defining this macro makes the cod...
Definition: GeneralDefine.h:62
CGHandler.h
BaseParticle::getRadius
Mdouble getRadius() const
Returns the particle's radius.
Definition: BaseParticle.h:348
RNG.h
File::setOpenMode
void setOpenMode(std::fstream::openmode openMode)
Allows the user to Sets File::openMode_.
Definition: File.cc:247
DPMBase::name_
std::string name_
the name of the problem, used, e.g., for the output files
Definition: DPMBase.h:1364
DPMBase::min_
Vec3D min_
These vectors are used for the XBalls domain, and occasionally people use it to add walls.
Definition: DPMBase.h:1265
BaseParticle::unfix
void unfix()
Unfix Particle function, which required a reference to the Species vector. It unfixes a Particle by c...
Definition: BaseParticle.cc:316
DPMBase::writeEneFile
void writeEneFile()
Definition: DPMBase.cc:2873
DPMBase::getXBallsScale
double getXBallsScale() const
Returns the scale of the view in xballs.
Definition: DPMBase.cc:1363
ParticleHandler::getLargestParticle
BaseParticle * getLargestParticle() const
Returns the pointer of the largest particle in the particle handler. When mercury is running in paral...
Definition: ParticleHandler.cc:530
InteractionHandler::write
void write(std::ostream &os) const
Writes the InteractionHandler to an output stream, for example a restart file.
Definition: InteractionHandler.cc:436
DPMBase::actionsOnRestart
virtual void actionsOnRestart()
A virtual function where the users can add extra code which is executed only when the code is restart...
Definition: DPMBase.cc:1674
DPMBase::checkSettings
void checkSettings()
Checks if the essentials are set properly to go ahead with solving the problem.
Definition: DPMBase.cc:3831
BaseInteractable::setVelocity
void setVelocity(const Vec3D &velocity)
set the velocity of the BaseInteractable.
Definition: BaseInteractable.cc:350
Domain::containsParticle
bool containsParticle(BaseParticle *particle, Mdouble offset=0.0)
Check to see if a given particle is within the current domain.
Definition: Domain.cc:400
BaseVTKWriter::setFileCounter
void setFileCounter(unsigned fileCounter)
Definition: BaseVTKWriter.h:61
DPMBaseXBalls.icc
BoundaryHandler::boundaryActionsBeforeTimeLoop
void boundaryActionsBeforeTimeLoop()
Definition: BoundaryHandler.cc:261
BoundaryHandler.h
BaseParticle::copy
virtual BaseParticle * copy() const =0
Particle copy method. It calls to copy constructor of this Particle, useful for polymorphism.
File::saveCurrentTimeStep
bool saveCurrentTimeStep(unsigned int ntimeSteps)
determined if this time step has to be written; if so, opens the output file
Definition: File.cc:312
DPMBase::mpiInsertParticleCheck
bool mpiInsertParticleCheck(BaseParticle *P)
Function that checks if the mpi particle should really be inserted by the current domain.
Definition: DPMBase.cc:1723
SpeciesHandler
Container to store all ParticleSpecies.
Definition: SpeciesHandler.h:37
RNG::randomise
void randomise()
sets the random variables such that they differ for each run
Definition: RNG.cc:97
BaseParticle::isMPIParticle
bool isMPIParticle() const
Indicates if this particle is a ghost in the MPI domain.
Definition: BaseParticle.cc:175
RNG
This is a class that generates random numbers i.e. named the Random Number Generator (RNG).
Definition: RNG.h:53
double
DPMBase::getMax
Vec3D getMax() const
Definition: DPMBase.h:643
BaseInteractable::setAngularVelocity
void setAngularVelocity(const Vec3D &angularVelocity)
set the angular velocity of the BaseInteractble.
Definition: BaseInteractable.cc:360
DPMBase::systemDimensions_
unsigned int systemDimensions_
The dimensions of the simulation i.e. 2D or 3D.
Definition: DPMBase.h:1245
VERBOSE
LL< Log::VERBOSE > VERBOSE
Verbose information.
Definition: Logger.cc:57
ERROR
LL< Log::ERROR > ERROR
Error log level.
Definition: Logger.cc:53
DPMBase::write
virtual void write(std::ostream &os, bool writeAllParticles=true) const
Definition: DPMBase.cc:3409
DPMBase::readUserDefinedWall
virtual BaseWall * readUserDefinedWall(const std::string &type) const
Allows you to read in a wall defined in a Driver directory; see USER/Luca/ScrewFiller.
Definition: DPMBase.h:282
MPIContainer::gather
void gather(T &send_t, T *receive_t)
Gathers a scaler from all processors to a vector of scalars on the root.
Definition: MpiContainer.h:428
DPMBase::getKineticEnergy
Mdouble getKineticEnergy() const
Returns the global kinetic energy stored in the system.
Definition: DPMBase.cc:1535
DPMBase::readNextFStatFile
void readNextFStatFile()
Reads the next fstat file.
Definition: DPMBase.cc:2779
Logger.h
MERCURY_DEPRECATED
#define MERCURY_DEPRECATED
Definition: GeneralDefine.h:37
DPMBase::getFStatFile
MERCURY_DEPRECATED File & getFStatFile()
The non const version. Allows to edit the File::fStatFile.
Definition: DPMBase.cc:319
BaseParticle::setSpecies
void setSpecies(const ParticleSpecies *species)
Definition: BaseParticle.cc:804
CGHandler::finish
void finish()
Contains the code executed after the last time step.
Definition: CGHandler.cc:113
Vec3D::z
Mdouble & z()
RW reference to Z.
Definition: Vector.h:368
SpeciesHandler::addObject
void addObject(ParticleSpecies *S) override
Adds a new ParticleSpecies to the SpeciesHandler.
Definition: SpeciesHandler.cc:797
DPMBase::setDomain
void setDomain(const Vec3D &min, const Vec3D &max)
Sets the minimum coordinates of the problem domain.
Definition: DPMBase.cc:1089
CMakeDefinitions.h
DomainHandler::getCurrentDomain
Domain * getCurrentDomain()
Gets the domain assigned to the processor.
Definition: DomainHandler.cc:250
BoundaryVTKWriter::writeVTK
void writeVTK() const override
writes a vtk file
Definition: BoundaryVTKWriter.cc:29
InteractionHandler::read
void read(std::istream &is)
Definition: InteractionHandler.cc:457
Contact
Definition: SpeedTestParticleInteractions.cpp:36
XAXIS
@ XAXIS
Definition: GeneralDefine.h:78
File::setName
void setName(const std::string &name)
Sets the file name, e.g. "Name.data".
Definition: File.cc:198
Time::getWallTime
Mdouble getWallTime()
Definition: MercuryTime.h:84
DPMBase::getTime
Mdouble getTime() const
Returns the current simulation time.
Definition: DPMBase.cc:805
Vec3D::x
Mdouble & x()
RW reference to X.
Definition: Vector.h:344
File::getSaveCount
unsigned int getSaveCount() const
Gets File::saveCount_.
Definition: File.cc:255
DPMBase::timeMax_
Mdouble timeMax_
Stores the duration of the simulation.
Definition: DPMBase.h:1286
DPMBase::setTime
void setTime(Mdouble time)
Sets a new value for the current simulation time.
Definition: DPMBase.cc:833
PeriodicBoundaryHandler::initialise
void initialise()
Initialises the communication list vectors as they can not be determined on compile time.
Definition: PeriodicBoundaryHandler.cc:1790
InteractionHandler::getWriteVTK
FileType getWriteVTK() const
Definition: InteractionHandler.cc:558
DPMBase::setRestartVersion
void setRestartVersion(std::string newRV)
Sets restart_version.
Definition: DPMBase.cc:1475
DPMBase::DomainSplit
DomainSplit
Definition: DPMBase.h:937
WallVTKWriter::writeVTK
void writeVTK() const override
Definition: WallVTKWriter.cc:52
DPMBase::setXMax
void setXMax(Mdouble newXMax)
Sets the value of XMax, the upper bound of the problem domain in the x-direction.
Definition: DPMBase.cc:1156
WARN
LL< Log::WARN > WARN
Warning log level.
Definition: Logger.cc:54
PeriodicBoundaryHandler::flushParticles
void flushParticles(std::set< BaseParticle * > &particlesToBeFlushed)
Removes particles from the periodiocParticleList_ and periociGhostList_.
Definition: PeriodicBoundaryHandler.cc:1708
InfiniteWall.h
ParticleHandler::clear
void clear() override
Empties the whole ParticleHandler by removing all BaseParticle.
Definition: ParticleHandler.cc:973
DPMBase::setOpenMode
void setOpenMode(std::fstream::openmode openMode)
Sets File::openMode_ for all files (ene, data, fstat, restart, stat)
Definition: DPMBase.cc:487
PossibleContact::getNext
PossibleContact * getNext()
Gets the next PossibleContact in the general linked list of PossibleContact.
Definition: PossibleContact.h:167
Log::FATAL
@ FATAL
DPMBase::writeEneHeader
virtual void writeEneHeader(std::ostream &os) const
Writes a header with a certain format for ENE file.
Definition: DPMBase.cc:2005
DPMBase::speciesHandler
SpeciesHandler speciesHandler
A handler to that stores the species type i.e. LinearViscoelasticSpecies, etc.
Definition: DPMBase.h:1385
DPMBase::xBallsVectorScale_
Mdouble xBallsVectorScale_
sets the xballs argument vscale (see XBalls/xballs.txt)
Definition: DPMBase.h:1344
PeriodicBoundary::getDistance
Mdouble getDistance(const BaseParticle &p) const override
Returns the distance of the edge to the particle.
Definition: PeriodicBoundary.cc:197
helpers::round
Mdouble round(const Mdouble value, unsigned precision)
Definition: Helpers.cc:598
DPMBase::domainHandler
DomainHandler domainHandler
An object of the class DomainHandler which deals with parallel code.
Definition: DPMBase.h:1420
DPMBase::writeVTKFiles
void writeVTKFiles() const
Definition: DPMBase.cc:2119
DPMBase::getKineticStress
Matrix3D getKineticStress() const
Calculate the kinetic stress tensor in the system averaged over the whole volume.
Definition: DPMBase.cc:5215
DPMBase::writeFStatFile
void writeFStatFile()
Definition: DPMBase.cc:2883
ParticleHandler::getNumberOfRealObjectsLocal
unsigned int getNumberOfRealObjectsLocal() const
Returns the number of real objects on a local domain. MPI particles and periodic particles are neglec...
Definition: ParticleHandler.cc:1273
DPMBase::getElasticEnergy
Mdouble getElasticEnergy() const
Returns the global elastic energy within the system.
Definition: DPMBase.cc:1521
BaseParticle::setMPIParticle
void setMPIParticle(bool flag)
Flags the mpi particle status.
Definition: BaseParticle.cc:185
BoundaryHandler::getWriteVTK
bool getWriteVTK() const
Definition: BoundaryHandler.h:90
DPMBase::isTimeEqualTo
bool isTimeEqualTo(Mdouble time) const
Checks whether the input variable "time" is the current time in the simulation.
Definition: DPMBase.cc:5039
DPMBase::setInteractionsWriteVTK
void setInteractionsWriteVTK(bool)
Sets whether interactions are written into a VTK file.
Definition: DPMBase.cc:924
DPMBase.h
DPMBase::readNextArgument
virtual bool readNextArgument(int &i, int argc, char *argv[])
Interprets the i^th command-line argument.
Definition: DPMBase.cc:4350
ParticleHandler::readAndAddObject
void readAndAddObject(std::istream &is) override
Definition: ParticleHandler.cc:1084
BaseInteractable::addForce
void addForce(const Vec3D &addForce)
Adds an amount to the force on this BaseInteractable.
Definition: BaseInteractable.cc:116
BaseInteractable::resetForceTorque
void resetForceTorque(int numberOfOMPthreads)
Definition: BaseInteractable.cc:141
SpeciesHandler::write
virtual void write(std::ostream &os) const
Write all the species and mixed species to an output stream.
Definition: SpeciesHandler.cc:836
DPMBase::continueSolve
virtual bool continueSolve() const
A virtual function for deciding whether to continue the simulation, based on a user-specified criteri...
Definition: DPMBase.cc:1981
SphericalParticle
A spherical particle is the most simple particle used in MercuryDPM.
Definition: SphericalParticle.h:37
Domain
The simulation can be subdivided into Domain's used in parallel code.
Definition: Domain.h:64
DPMBase::wallHandler
WallHandler wallHandler
An object of the class WallHandler. Contains pointers to all the walls created.
Definition: DPMBase.h:1405
MPIContainer::initialiseMercuryMPITypes
void initialiseMercuryMPITypes(const SpeciesHandler &speciesHandler)
Creates the MPI types required for communication of Mercury data through the MPI interface.
Definition: MpiContainer.cc:74
ZAXIS
@ ZAXIS
Definition: GeneralDefine.h:78
getVersion
const std::string getVersion()
Definition: CMakeDefinitions.cc:51
BaseInteractable::integrateAfterForceComputation
void integrateAfterForceComputation(double time, double timeStep)
This is part of the integration routine for objects with infinite mass.
Definition: BaseInteractable.cc:611
DPMBase::getTotalEnergy
Mdouble getTotalEnergy() const
Definition: DPMBase.cc:1582
BaseHandler::getNumberOfObjects
virtual unsigned int getNumberOfObjects() const
Gets the number of real Object in this BaseHandler. (i.e. no mpi or periodic particles)
Definition: BaseHandler.h:648
Domain::addParticle
void addParticle(BaseParticle *particle)
Initialises a single particle which is added during the simulation.
Definition: Domain.cc:1610
DPMBase::time_
Mdouble time_
Stores the current simulation time.
Definition: DPMBase.h:1271
File::setLastSavedTimeStep
void setLastSavedTimeStep(unsigned int lastSavedTimeStep)
Sets File::nextSavedTimeStep_.
Definition: File.cc:302
MPIContainer::sync
void sync()
Process all pending asynchronous communication requests before continuing.
Definition: MpiContainer.h:152
DPMBase::setDimension
void setDimension(unsigned int newDim)
Sets both the system dimensions and the particle dimensionality.
Definition: DPMBase.cc:1394
DPMBase::getNToWrite
int getNToWrite() const
get the number of elements to write to the
Definition: DPMBase.cc:861
MPIContainer
This class contains all information and functions required for communication between processors.
Definition: MpiContainer.h:130
ParticleHandler::getSmallestParticleLocal
BaseParticle * getSmallestParticleLocal() const
Gets a pointer to the smallest BaseParticle (by interactionRadius) in this ParticleHandler of the loc...
Definition: ParticleHandler.cc:498
DPMBase::getXBallsVectorScale
double getXBallsVectorScale() const
Returns the scale of vectors used in xballs.
Definition: DPMBase.cc:1321
DPMBase::setTimeMax
void setTimeMax(Mdouble newTMax)
Sets a new value for the maximum simulation duration.
Definition: DPMBase.cc:870
BaseParticle::integrateAfterForceComputation
void integrateAfterForceComputation(double time, double timeStep)
Second step of Velocity Verlet integration.
Definition: BaseParticle.cc:748
DPMBase::updateGhostGrid
void updateGhostGrid(BaseParticle *P)
Checks if the Domain/periodic interaction distance needs to be updated and updates it accordingly.
Definition: DPMBase.cc:1828
Vec3D::Y
Mdouble Y
Definition: Vector.h:65
DPMBase::getHGridUpdateEachTimeStep
virtual bool getHGridUpdateEachTimeStep() const
Definition: DPMBase.cc:1709
DPMBase::getCentreOfMass
Vec3D getCentreOfMass() const
JMFT: Return the centre of mass of the system, excluding fixed particles.
Definition: DPMBase.cc:1605
DPMBase::setGravity
void setGravity(Vec3D newGravity)
Sets a new value for the gravitational acceleration.
Definition: DPMBase.cc:1374
MPIParticle
Data class to send a particle over MPI.
Definition: MpiDataClass.h:81
DPMBase::getRotationalEnergy
Mdouble getRotationalEnergy() const
JMFT Returns the global rotational energy stored in the system.
Definition: DPMBase.cc:1568
DPMBase::resetFileCounter
void resetFileCounter()
Resets the file counter for each file i.e. for ene, data, fstat, restart, stat)
Definition: DPMBase.cc:469
NEVER
const unsigned NEVER
Definition: File.h:35
Domain::flushParticles
void flushParticles(std::set< BaseParticle * > &toBeDeletedList)
Particles that are going to be deleted from the simulation are flushed out of the communication bound...
Definition: Domain.cc:1698
ParticleHandler::getCentreOfMass
Vec3D getCentreOfMass() const
Definition: ParticleHandler.cc:650
InfiniteWall
A infinite wall fills the half-space {point: (position_-point)*normal_<=0}.
Definition: InfiniteWall.h:48
FATAL
LL< Log::FATAL > FATAL
Definition of the different loglevels by its wrapper class LL. These are used as tags in template met...
Definition: Logger.cc:52
DPMBase::getGravitationalEnergy
Mdouble getGravitationalEnergy() const
Returns the global gravitational potential energy stored in the system.
Definition: DPMBase.cc:1552
helpers::getLineFromStringStream
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:423
BaseWall.h
DomainHandler::createMesh
void createMesh(std::vector< Mdouble > &simulationMin, std::vector< Mdouble > &simulationMax, std::vector< unsigned > &numberOfDomains, bool open)
Creates a Cartesian square mesh in 3D.
Definition: DomainHandler.cc:93
DPMBase::getName
const std::string & getName() const
Returns the name of the file. Does not allow to change it though.
Definition: DPMBase.cc:397
DPMBase::getRunNumber
int getRunNumber() const
This returns the current value of the counter (runNumber_)
Definition: DPMBase.cc:614
DPMBase::writeDataFile
void writeDataFile()
Definition: DPMBase.cc:2864
constants::i
const std::complex< Mdouble > i
Definition: ExtendedMath.h:51
DPMBase::DomainSplit::X
@ X
DPMBase::setRotation
void setRotation(bool rotation)
Sets whether particle rotation is enabled or disabled.
Definition: DPMBase.h:552
Domain.h
BaseWall::setSpecies
void setSpecies(const ParticleSpecies *species)
Defines the species of the current wall.
Definition: BaseWall.cc:171
DPMBase::getZMin
Mdouble getZMin() const
If the length of the problem domain in z-direction is ZMax - ZMin, then getZMin() returns ZMin.
Definition: DPMBase.h:625
Quaternion
This class contains the 4 components of a quaternion and the standard operators and functions needed ...
Definition: Quaternion.h:63
DPMBase::setYMin
void setYMin(Mdouble newYMin)
Sets the value of YMin, the lower bound of the problem domain in the y-direction.
Definition: DPMBase.cc:1025
helpers::fileExists
bool fileExists(std::string strFilename)
Function to check if a file exists, is used to check if a run has already need done.
Definition: Helpers.cc:502
DPMBase::statFile
File statFile
An instance of class File to handle in- and output into a .stat file.
Definition: DPMBase.h:1456
DPMBase::setFileType
void setFileType(FileType fileType)
Sets File::fileType_ for all files (ene, data, fstat, restart, stat)
Definition: DPMBase.cc:457
DPMBase::findNextExistingDataFile
bool findNextExistingDataFile(Mdouble tMin, bool verbose=true)
Finds and opens the next data file, if such a file exists.
Definition: DPMBase.cc:2563
DPMBase::decompose
virtual void decompose()
Sends particles from processorId to the root processor.
Definition: DPMBase.cc:3920
DPMBase::ReadOptions::ReadAll
@ ReadAll
BaseHandler::getLastObject
T * getLastObject()
Gets a pointer to the last Object in this BaseHandler.
Definition: BaseHandler.h:634
DPMBase::setNumberOfDomains
void setNumberOfDomains(std::vector< unsigned > direction)
Sets the number of domains in x-,y- and z-direction. Required for parallel computations.
Definition: DPMBase.cc:5049
DPMBase::get1DParametersFromRunNumber
std::vector< int > get1DParametersFromRunNumber(int size_x) const
This turns a counter into 1 index, which is a useful feature for performing 1D parameter study....
Definition: DPMBase.cc:667
DPMBase::setupInitialConditions
virtual void setupInitialConditions()
This function allows to set the initial conditions for our problem to be solved, by default particle ...
Definition: DPMBase.cc:1989
DPMBase::getNumberOfDomains
std::vector< unsigned > getNumberOfDomains()
returns the number of domains
Definition: DPMBase.cc:5116
DPMBase::initialiseStatistics
virtual void initialiseStatistics()
Definition: DPMBase.cc:1876
DPMBase::setSuperquadricParticlesWriteVTK
void setSuperquadricParticlesWriteVTK(bool writeSuperquadricParticlesVTK)
Definition: DPMBase.cc:948
DPMBase::actionsAfterSolve
virtual void actionsAfterSolve()
A virtual function which allows to define operations to be executed after the solve().
Definition: DPMBase.cc:1862
BaseInteractable::setOrientationViaEuler
void setOrientationViaEuler(Vec3D eulerAngle)
Sets the orientation of this BaseInteractable by defining the euler angles.
Definition: BaseInteractable.cc:204
BaseHandler::getObject
T * getObject(const unsigned int id)
Gets a pointer to the Object at the specified index in the BaseHandler.
Definition: BaseHandler.h:613
DPMBase::checkParticleForInteractionLocalPeriodic
bool checkParticleForInteractionLocalPeriodic(const BaseParticle &P)
Definition: DPMBase.cc:4688
File::openWrite
bool openWrite(unsigned)
First sets openmode to write (and append in some cases), then calls open().
Definition: File.cc:381
RNG::write
void write(std::ostream &os) const
Definition: RNG.cc:77
BaseObject::getId
unsigned int getId() const
Returns the unique identifier of any particular object.
Definition: BaseObject.h:125
Vec3D::y
Mdouble & y()
RW reference to Y.
Definition: Vector.h:356
SuperQuadricParticleVtkWriter
Definition: SuperQuadricParticleVtkWriter.h:34
DPMBase::forceWriteOutputFiles
void forceWriteOutputFiles()
Writes output files immediately, even if the current time step was not meant to be written....
Definition: DPMBase.cc:3860
DPMBase::setNumberOfOMPThreads
void setNumberOfOMPThreads(int numberOfOMPThreads)
Definition: DPMBase.cc:1248
DPMBase::eneFile
File eneFile
An instance of class File to handle in- and output into a .ene file.
Definition: DPMBase.h:1446
File::setCounter
void setCounter(unsigned int counter)
Allows the user to set the file counter according to his need. Sets File::counter_.
Definition: File.cc:231
BaseParticle
Definition: BaseParticle.h:54
Vec3D::getLengthSquared
static Mdouble getLengthSquared(const Vec3D &a)
Calculates the squared length of a Vec3D: .
Definition: Vector.h:316
DPMBase::interactionVTKWriter_
InteractionVTKWriter interactionVTKWriter_
Definition: DPMBase.h:1329
BaseInteractable::getInteractions
const std::vector< BaseInteraction * > & getInteractions() const
Returns a list of interactions which belong to this interactable.
Definition: BaseInteractable.h:277
UNUSED
#define UNUSED
Definition: GeneralDefine.h:39
DPMBase::getMin
Vec3D getMin() const
Definition: DPMBase.h:637
quick
bool quick
Definition: LongPeriodicChute.cpp:36
DPMBase::DomainSplit::XYZ
@ XYZ
DPMBase::getSystemDimensions
unsigned int getSystemDimensions() const
Returns the system dimensionality.
Definition: DPMBase.cc:1421
BaseHandler::copyAndAddObject
std::enable_if<!std::is_pointer< U >::value, U * >::type copyAndAddObject(const U &object)
Creates a copy of a Object and adds it to the BaseHandler.
Definition: BaseHandler.h:379
DPMBase::getWallTime
double getWallTime()
Definition: DPMBase.h:876
DPMBase::computeAllForces
virtual void computeAllForces()
Computes all the forces acting on the particles using the BaseInteractable::setForce() and BaseIntera...
Definition: DPMBase.cc:3302
BaseParticle::getMaxInteractionRadius
Mdouble getMaxInteractionRadius() const
Returns the particle's interaction radius, which might be different from radius_ (e....
Definition: BaseParticle.h:362
ParticleHandler::getKineticEnergy
Mdouble getKineticEnergy() const
Definition: ParticleHandler.cc:553
n
const unsigned n
Definition: CG3DPackingUnitTest.cpp:32
Time
Allows for timing the algorithms; accurate up to 0.01 sec.
Definition: MercuryTime.h:44
DPMBase::timeStep_
Mdouble timeStep_
Stores the simulation time step.
Definition: DPMBase.h:1281
DPMBase::integrateAfterForceComputation
virtual void integrateAfterForceComputation()
Update particles' and walls' positions and velocities after force computation.
Definition: DPMBase.cc:3219
DPMBase::getAppend
bool getAppend() const
Returns whether the "append" option is on or off.
Definition: DPMBase.cc:1501
PossibleContact::getP1
BaseParticle * getP1()
Gets a pointer to the first BaseParticle in this PossibleContact.
Definition: PossibleContact.h:130
Vec3D::cross
static Vec3D cross(const Vec3D &a, const Vec3D &b)
Calculates the cross product of two Vec3D: .
Definition: Vector.cc:163
DPMBase::computeWallForces
virtual void computeWallForces(BaseWall *w)
Definition: DPMBase.cc:5263
RNG::setRandomSeed
void setRandomSeed(unsigned long int new_seed)
This is the seed for the random number generator (note the call to seed_LFG is only required really i...
Definition: RNG.cc:52
PeriodicBoundaryHandler::addNewParticles
void addNewParticles()
Adds new particles to the periodic particle lists.
Definition: PeriodicBoundaryHandler.cc:300
InteractionVTKWriter.h
DPMBase::writeParticlesVTK_
bool writeParticlesVTK_
A flag to turn on/off the vtk writer for particles.
Definition: DPMBase.h:1321
DomainHandler::setCurrentDomainIndex
void setCurrentDomainIndex(unsigned int index)
This sets a domain to the processor.
Definition: DomainHandler.cc:240
BaseHandler::begin
const std::vector< T * >::const_iterator begin() const
Gets the begin of the const_iterator over all Object in this BaseHandler.
Definition: BaseHandler.h:690
DPMBase::getXBallsColourMode
int getXBallsColourMode() const
Get the xballs colour mode (CMode).
Definition: DPMBase.cc:1301
Quaternion::setEuler
void setEuler(const Vec3D &e)
Convert Euler angles to a quaternion. See Wikipedia for details.
Definition: Quaternion.cc:473
DPMBase::ReadOptions::ReadNoParticlesAndInteractions
@ ReadNoParticlesAndInteractions
DPMBase::getRestartVersion
std::string getRestartVersion() const
This is to take into account for different Mercury versions. Returns the version of the restart file.
Definition: DPMBase.cc:1466
DPMBase::getTotalMass
Mdouble getTotalMass() const
JMFT: Return the total mass of the system, excluding fixed particles.
Definition: DPMBase.cc:1589
getSVNRevision
const int getSVNRevision()
Definition: CMakeDefinitions.cc:41
MpiDataClass.h
DPMBase::getEneFile
MERCURY_DEPRECATED File & getEneFile()
The non const version. Allows to edit the File::eneFile.
Definition: DPMBase.cc:311
Time::toc
Mdouble toc()
This is like a stop button of a stopwatch. Assigns the variable finish to the current value of ticks ...
Definition: MercuryTime.h:66
RNG::getRandomNumber
Mdouble getRandomNumber()
This is a random generating routine can be used for initial positions.
Definition: RNG.cc:142
DPMBase::checkAndDuplicatePeriodicParticles
void checkAndDuplicatePeriodicParticles()
For simulations using periodic boundaries, checks and adds particles when necessary into the particle...
Definition: DPMBase.cc:4865
DPMBase::setXBallsAdditionalArguments
void setXBallsAdditionalArguments(std::string newXBArgs)
Set the additional arguments for xballs.
Definition: DPMBase.cc:1338
DPMBase::get2DParametersFromRunNumber
std::vector< int > get2DParametersFromRunNumber(int size_x, int size_y) const
This turns a counter into 2 indices which is a very useful feature for performing a 2D study....
Definition: DPMBase.cc:695
SphericalParticle.h
CGHandler::initialise
void initialise()
Contains the code executed before the first time step.
Definition: CGHandler.cc:90
DPMBase::finishStatistics
virtual void finishStatistics()
Definition: DPMBase.cc:1916
DPMBase::setRestarted
void setRestarted(bool newRestartedFlag)
Allows to set the flag stating if the simulation is to be restarted or not.
Definition: DPMBase.cc:1492
helpers::writeToFile
bool writeToFile(std::string filename, std::string filecontent)
Writes a string to a file.
Definition: Helpers.cc:445
PeriodicBoundary
Defines a pair of periodic walls. Inherits from BaseBoundary.
Definition: PeriodicBoundary.h:41
DPMBase::getGravity
Vec3D getGravity() const
Returns the gravitational acceleration.
Definition: DPMBase.cc:1382
DPMBase::getYMax
Mdouble getYMax() const
If the length of the problem domain in y-direction is YMax - YMin, then getYMax() returns XMax.
Definition: DPMBase.h:619
DPMBase::performGhostVelocityUpdate
void performGhostVelocityUpdate()
updates the final time-step velocity of the ghost particles
Definition: DPMBase.cc:4964
PossibleContactList
Manages the linked list of PossibleContact.
Definition: PossibleContactList.h:44
BaseHandler::clear
virtual void clear()
Empties the whole BaseHandler by removing all Objects and setting all other variables to 0.
Definition: BaseHandler.h:528
constants::NaN
const Mdouble NaN
Definition: GeneralDefine.h:43
DPMBase::getZMax
Mdouble getZMax() const
If the length of the problem domain in z-direction is ZMax - ZMin, then getZMax() returns ZMax.
Definition: DPMBase.h:631
DPMBase::DomainSplit::Y
@ Y
SphericalParticleVtkWriter.h
DPMBase::outputXBallsData
virtual void outputXBallsData(std::ostream &os) const
This function writes the location of the walls and particles in a format the XBalls program can read....
Definition: DPMBase.cc:2268
DPMBase::setXMin
void setXMin(Mdouble newXMin)
Sets the value of XMin, the lower bound of the problem domain in the x-direction.
Definition: DPMBase.cc:1001
ParticleHandler::getMass
Mdouble getMass() const
Definition: ParticleHandler.cc:607
PeriodicBoundary::set
void set(Vec3D normal, Mdouble distanceLeft, Mdouble distanceRight)
Defines a PeriodicBoundary by its normal and positions.
Definition: PeriodicBoundary.cc:84
MPIContainer::broadcast
std::enable_if< std::is_scalar< T >::value, void >::type broadcast(T &t, int fromProcessor=0)
Broadcasts a scalar from the root to all other processors.
Definition: MpiContainer.h:441
ParticleHandler
Container to store all BaseParticle.
Definition: ParticleHandler.h:48
Vec3D::getDistanceSquared
static Mdouble getDistanceSquared(const Vec3D &a, const Vec3D &b)
Calculates the squared distance between two Vec3D: .
Definition: Vector.h:295
PeriodicBoundaryHandler::setInteractionDistance
void setInteractionDistance(Mdouble interactionDistance)
Sets the interaction distance.
Definition: PeriodicBoundaryHandler.cc:128
DPMBase::computeInternalForces
virtual void computeInternalForces(BaseParticle *)
Computes the internal forces on particle i (internal in the sense that the sum over all these forces ...
Definition: DPMBase.cc:3393
DPMBase::writeSuperquadricParticlesVTK_
bool writeSuperquadricParticlesVTK_
Definition: DPMBase.h:1323
DPMBase::printTime
virtual void printTime() const
Displays the current simulation time and the maximum simulation duration.
Definition: DPMBase.cc:1961
ParticleHandler::getLargestParticleLocal
BaseParticle * getLargestParticleLocal() const
Gets a pointer to the largest BaseParticle (by interactionRadius) in the ParticleHandler of the local...
Definition: ParticleHandler.cc:522
DPMBase::hGridGetInteractingParticleList
virtual void hGridGetInteractingParticleList(BaseParticle *obj, std::vector< BaseParticle * > &list)
Creates a list of neighbour particles obtained from the hgrid.
Definition: DPMBase.h:961
BaseParticle::isFixed
bool isFixed() const override
Is fixed Particle function. It returns whether a Particle is fixed or not, by checking its inverse Ma...
Definition: BaseParticle.h:93
PROCESSOR_ID
#define PROCESSOR_ID
Definition: GeneralDefine.h:63
BoundaryVTKWriter.h
ParticleHandler::removeObject
void removeObject(unsigned int index) override
Removes a BaseParticle from the ParticleHandler.
Definition: ParticleHandler.cc:390
DomainHandler::getInteractionDistance
Mdouble getInteractionDistance()
Gets the interaction distance of the domain handler.
Definition: DomainHandler.cc:324
DPMBase::getRotation
bool getRotation() const
Indicates whether particle rotation is enabled or disabled.
Definition: DPMBase.h:559
DPMBase::particleHandler
ParticleHandler particleHandler
An object of the class ParticleHandler, contains the pointers to all the particles created.
Definition: DPMBase.h:1395
DPMBase::readSpeciesFromDataFile
void readSpeciesFromDataFile(bool read=true)
Definition: DPMBase.h:341
DPMBase::importParticlesAs
void importParticlesAs(ParticleHandler &particleHandler, InteractionHandler &interactionHandler, const ParticleSpecies *species)
Copies particles, interactions assigning species from a local simulation to a global one....
Definition: DPMBase.cc:4788
Vec3D::Z
Mdouble Z
Definition: Vector.h:65
DPMBase::gravity_
Vec3D gravity_
Gravity vector.
Definition: DPMBase.h:1255
DPMBase::performGhostParticleUpdate
void performGhostParticleUpdate()
When the Verlet scheme updates the positions and velocities of particles, ghost particles will need a...
Definition: DPMBase.cc:4899
BaseParticle::isPeriodicGhostParticle
bool isPeriodicGhostParticle() const
Indicates if this particle is a ghost in the periodic boundary.
Definition: BaseParticle.cc:291
DPMBase::removeDuplicatePeriodicParticles
void removeDuplicatePeriodicParticles()
Removes periodic duplicate Particles.
Definition: DPMBase.cc:4822
DPMBase::setSaveCount
void setSaveCount(unsigned int saveCount)
Sets File::saveCount_ for all files (ene, data, fstat, restart, stat)
Definition: DPMBase.cc:406
FileType
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:41
helpers::getPath
std::string getPath()
Definition: Helpers.cc:931
BaseVTKWriter::writeVTK
virtual void writeVTK() const =0
CubeInsertionBoundary.h
DPMBase::particleDimensions_
unsigned int particleDimensions_
determines if 2D or 3D particle volume is used for mass calculations
Definition: DPMBase.h:1250
BaseParticle::isInContactWith
virtual bool isInContactWith(const BaseParticle *P) const
Get whether or not this particle is in contact with the given particle.
Definition: BaseParticle.cc:851
DPMBase::insertGhostParticle
void insertGhostParticle(BaseParticle *P)
This function inserts a particle in the mpi communication boundaries.
Definition: DPMBase.cc:1802
DPMBase::ReadOptions::ReadNoInteractions
@ ReadNoInteractions
DPMBase::setXBallsVectorScale
void setXBallsVectorScale(double newVScale)
Set the scale of vectors in xballs.
Definition: DPMBase.cc:1311
DPMBase::xBallsAdditionalArguments_
std::string xBallsAdditionalArguments_
A string of additional arguments for xballs can be specified (see XBalls/xballs.txt)....
Definition: DPMBase.h:1354
DomainHandler::updateStatus
void updateStatus(std::set< BaseParticle * > &particlesToBeDeleted)
Definition: DomainHandler.cc:419
DPMBase::removeOldFiles
void removeOldFiles() const
Definition: DPMBase.cc:4256
BaseInteractable::integrateBeforeForceComputation
void integrateBeforeForceComputation(double time, double timeStep)
This is part of integrate routine for objects with infinite mass.
Definition: BaseInteractable.cc:538
WallHandler.h
DPMBase::DomainSplit::XZ
@ XZ
DPMBase::gatherContactStatistics
void gatherContactStatistics()
Definition: DPMBase.cc:1889
BaseParticle::getSumOfInteractionRadii
Mdouble getSumOfInteractionRadii(const BaseParticle *particle) const
returns the sum of the radii plus the interactionDistance
Definition: BaseParticle.h:379
BaseParticle::getMass
Mdouble getMass() const
Returns the particle's mass.
Definition: BaseParticle.h:322
InfiniteWall::set
void set(Vec3D normal, Vec3D point)
Defines a standard wall, given an outward normal vector s.t. normal*x=normal*point for all x of the w...
Definition: InfiniteWall.cc:118
File::getFstream
std::fstream & getFstream()
Allows to access the member variable File::fstream_.
Definition: File.cc:153
WallVTKWriter
Definition: WallVTKWriter.h:34
helpers::compare
bool compare(std::istream &is, std::string s)
Checks if the next argument in the input stream is a certain string.
Definition: Helpers.cc:858
DPMBase::dataFile
File dataFile
An instance of class File to handle in- and output into a .data file.
Definition: DPMBase.h:1436
InteractionVTKWriter::writeVTK
void writeVTK() const override
writes a vtk file
Definition: InteractionVTKWriter.cc:29
InteractionHandler.h
BaseHandler::end
const std::vector< T * >::const_iterator end() const
Gets the end of the const_iterator over all BaseBoundary in this BaseHandler.
Definition: BaseHandler.h:704
BaseParticle.h
ParticleHandler::write
void write(std::ostream &os) const
units::name
std::string name
Definition: MercuryProb.h:48
InteractionHandler::eraseOldInteractions
void eraseOldInteractions(unsigned)
erases interactions which have an old timestamp.
Definition: InteractionHandler.cc:380
BaseParticle::setPeriodicGhostParticle
void setPeriodicGhostParticle(bool flag)
Flags the status of the particle to be a ghost in periodic boundary or not.
Definition: BaseParticle.cc:296
DPMBase::hGridActionsAfterIntegration
virtual void hGridActionsAfterIntegration()
This function has to be called after integrateBeforeForceComputation.
Definition: DPMBase.cc:1938
DPMBase::~DPMBase
virtual ~DPMBase()
virtual destructor
Definition: DPMBase.cc:295
PeriodicBoundaryHandler.h
DPMBase::getNumberOfTimeSteps
unsigned int getNumberOfTimeSteps() const
Returns the current counter of time-steps, i.e. the number of time-steps that the simulation has unde...
Definition: DPMBase.cc:821
DPMBase::areInContact
static bool areInContact(const BaseParticle *pI, const BaseParticle *pJ)
Checks if two particle are in contact or is there any positive overlap.
Definition: DPMBase.cc:1652
DPMBase::readNextDataFile
bool readNextDataFile(unsigned int format=0)
Reads the next data file with default format=0. However, one can modify the format based on whether t...
Definition: DPMBase.cc:2600
CubeInsertionBoundary::set
void set(BaseParticle *particleToCopy, unsigned int maxFailed, Vec3D posMin, Vec3D posMax, Vec3D velMin, Vec3D velMax)
Sets the properties of the InsertionBoundary for mutliple different particle types.
Definition: CubeInsertionBoundary.cc:107
DPMBase::getTotalMomentum
Vec3D getTotalMomentum() const
JMFT: Return the total momentum of the system, excluding fixed particles.
Definition: DPMBase.cc:1615
DPMBase::setZMin
void setZMin(Mdouble newZMin)
Sets the value of ZMin, the lower bound of the problem domain in the z-direction.
Definition: DPMBase.cc:1049
InteractionVTKWriter
Definition: InteractionVTKWriter.h:34
DPMBase::getTotalStress
Matrix3D getTotalStress() const
Calculate the total stress tensor in the system averaged over the whole volume.
Definition: DPMBase.cc:5257
BaseParticle::getInfo
virtual Mdouble getInfo() const
Returns some user-defined information about this object (by default, species ID).
Definition: BaseParticle.cc:352
File::getFileType
FileType getFileType() const
Gets the file type e.g. NOFILE, ONEFILE and MULTIPLE FILES. File::fileType_.
Definition: File.cc:207
DPMBase::setLastSavedTimeStep
void setLastSavedTimeStep(unsigned int nextSavedTimeStep)
Sets the next time step for all the files (ene, data, fstat, restart, stat) at which the data is to b...
Definition: DPMBase.cc:516
DPMBase::deleteGhostParticles
void deleteGhostParticles(std::set< BaseParticle * > &particlesToBeDeleted)
Definition: DPMBase.cc:4922
helpers::to_string
std::string to_string(const T &n)
Definition: Helpers.h:227
DPMBase::getXMin
Mdouble getXMin() const
If the length of the problem domain in x-direction is XMax - XMin, then getXMin() returns XMin.
Definition: DPMBase.h:600
DPMBase::actionsBeforeTimeStep
virtual void actionsBeforeTimeStep()
A virtual function which allows to define operations to be executed before the new time step.
Definition: DPMBase.cc:1855
BaseWall::getInteractionWith
BaseInteraction * getInteractionWith(BaseParticle *p, unsigned timeStamp, InteractionHandler *interactionHandler) override
Returns the interaction between this wall and a given particle, nullptr if there is no interaction.
Definition: BaseWall.cc:369
InteractionHandler::setWriteVTK
void setWriteVTK(FileType f)
Definition: InteractionHandler.cc:553
DPMBase::setMeanVelocityAndKineticEnergy
void setMeanVelocityAndKineticEnergy(Vec3D V_mean_goal, Mdouble Ek_goal)
This function will help you set a fixed kinetic energy and mean velocity in your system.
Definition: DPMBase.cc:5158
MercuryTime.h
DPMBase::hGridInsertParticle
virtual void hGridInsertParticle(BaseParticle *obj UNUSED)
Definition: DPMBase.cc:1688
logWriteAndDie
MERCURY_DEPRECATED void logWriteAndDie(const std::string &module, std::string message)
todo strcmp relies on this, should be changed to more modern version
Definition: DPMBase.cc:79
DPMBase::interactionHandler
InteractionHandler interactionHandler
An object of the class InteractionHandler.
Definition: DPMBase.h:1425
operator<<
std::ostream & operator<<(std::ostream &os, DPMBase &md)
Definition: DPMBase.cc:96
InteractionHandler::removeObjectKeepingPeriodics
void removeObjectKeepingPeriodics(unsigned int id)
Removes interactions of periodic particles when the periodic particles get deleted (see DPMBase::remo...
Definition: InteractionHandler.cc:319
DPMBase::continueFlag_
static volatile sig_atomic_t continueFlag_
Definition: DPMBase.h:1237
BaseInteractable::getPosition
const Vec3D & getPosition() const
Returns the position of this BaseInteractable.
Definition: BaseInteractable.h:218
PossibleContact
Class that describes a possible contact between two BaseParticle.
Definition: PossibleContact.h:42
DPMBase::getTimeMax
Mdouble getTimeMax() const
Returns the maximum simulation duration.
Definition: DPMBase.cc:885
DPMBase::getRestartFile
MERCURY_DEPRECATED File & getRestartFile()
The non const version. Allows to edit the File::restartFile.
Definition: DPMBase.cc:327
int
Matrix3D::dyadic
static Matrix3D dyadic(const Vec3D &a, const Vec3D &b)
Calculates the dyadic product of a two Vec3D: .
Definition: Matrix.cc:323
BaseVTKWriter::getFileCounter
unsigned getFileCounter() const
Definition: BaseVTKWriter.h:56
DPMBase::checkParticleForInteraction
virtual bool checkParticleForInteraction(const BaseParticle &P)
Checks whether a particle P has any interaction with walls or other particles.
Definition: DPMBase.cc:4639
DPMBase::computeInternalForce
virtual void computeInternalForce(BaseParticle *, BaseParticle *)
Computes the forces between two particles (internal in the sense that the sum over all these forces i...
Definition: DPMBase.cc:3019
DPMBase::hGridActionsBeforeTimeStep
virtual void hGridActionsBeforeTimeStep()
A virtual function that allows one to set or execute hGrid parameters or operations before every simu...
Definition: DPMBase.cc:1681
File::getFullName
const std::string getFullName() const
Also allows to access the file name, however with additional information which is the file counter,...
Definition: File.cc:170
DPMBase
The DPMBase header includes quite a few header files, defining all the handlers, which are essential....
Definition: DPMBase.h:76
ParticleHandler::computeAllMasses
void computeAllMasses(unsigned int indSpecies)
Computes the mass for all BaseParticle of the given species in this ParticleHandler.
Definition: ParticleHandler.cc:1205
DPMBase::setAppend
void setAppend(bool newAppendFlag)
Sets whether the "append" option is on or off.
Definition: DPMBase.cc:1513
MPIContainer::Instance
static MPIContainer & Instance()
fetch the instance to be used for communication
Definition: MpiContainer.h:134
DPMBase::max_
Vec3D max_
Definition: DPMBase.h:1266
DPMBase::boundaryVTKWriter_
BoundaryVTKWriter boundaryVTKWriter_
Definition: DPMBase.h:1331
DPMBase::readArguments
bool readArguments(int argc, char *argv[])
Can interpret main function input arguments that are passed by the driver codes.
Definition: DPMBase.cc:4225
DPMBase::getTotalVolume
Mdouble getTotalVolume() const
Get the total volume of the cuboid system.
Definition: DPMBase.cc:5206
DPMBase::setSystemDimensions
void setSystemDimensions(unsigned int newDim)
Sets the system dimensionality.
Definition: DPMBase.cc:1408