revision: v0.14
helpers Namespace Reference

Classes

class  KAndDisp
 return type specifically for fuctions returning k and disp at once More...
 
class  KAndDispAndKtAndDispt
 Calculates the collision time for a given stiffness, dissipation, and effective mass. More...
 

Functions

std::string lower (std::string s)
 
MERCURY_DEPRECATED KAndDispAndKtAndDispt computeDisptFromCollisionTimeAndRestitutionCoefficientAndTangentialRestitutionCoefficientAndEffectiveMass (Mdouble tc, Mdouble r, Mdouble beta, Mdouble mass)
 Set disp, k, dispt and kt such that is matches a given collision time tc and a normal and tangential restitution coefficient r, beta for a collision of effective/reduced mass m. From Deen...Kuipers2006, eq. 43 and 30. More...
 
MERCURY_DEPRECATED Mdouble getMaximumVelocity (Mdouble k, Mdouble disp, Mdouble radius, Mdouble mass)
 Calculates the maximum relative velocity allowed for a normal collision of two particles of radius r and particle mass m (for higher velocities particles could pass through each other) More...
 
unsigned int getSaveCountFromNumberOfSavesAndTimeMaxAndTimeStep (unsigned int numberOfSaves, Mdouble timeMax, Mdouble timeStep)
 Returns the correct saveCount if the total number of saves, the final time and the time step is known. More...
 
void getLineFromStringStream (std::istream &in, std::stringstream &out)
 Reads a line from one stringstream into another, and prepares the latter for reading in. More...
 
bool writeToFile (std::string filename, std::string filecontent)
 Writes a string to a file. More...
 
void writeCommandLineToFile (const std::string filename, const int argc, char *const argv[])
 Writes a string to a file. More...
 
void gnuplot (std::string command)
 Plots to a gnuplot window. More...
 
bool addToFile (std::string filename, std::string filecontent)
 Adds a string to an existing file. More...
 
bool fileExists (std::string strFilename)
 Function to check if a file exists, is used to check if a run has already need done. More...
 
bool openFile (std::fstream &file, std::string filename, std::fstream::openmode mode)
 Provides a simple interface for opening a file. More...
 
Mdouble getEffectiveMass (Mdouble mass0, Mdouble mass1)
 Calculates the effective mass of a particle pair, i.e. half the harmonic mean of two particle masses. More...
 
std::vector< doublereadArrayFromFile (std::string filename, int &n, int &m)
 
void more (std::string filename, unsigned nLines=constants::unsignedMax)
 
template<typename T >
std::string to_string (const T &n)
 
std::string to_string (Mdouble value, unsigned precision)
 
template<typename T >
bool readOptionalVariable (std::istream &is, const std::string &name, T &variable)
 Reads optional variables in the restart file. More...
 
void loadingTest (const ParticleSpecies *species, Mdouble displacement, Mdouble velocity, Mdouble radius, std::string name)
 
void normalAndTangentialLoadingTest (const ParticleSpecies *species, Mdouble displacement, Mdouble tangentialDisplacement, Mdouble velocity, Mdouble radius, std::string name)
 
void objectivenessTest (const ParticleSpecies *species, Mdouble displacement, Mdouble tangentialDisplacement, Mdouble velocity, Mdouble radius, std::string name)
 
bool compare (std::istream &is, std::string s)
 Checks if the next argument in the input stream is a certain string. More...
 
template<typename T >
readFromFile (std::string fileName, std::string varName, T value)
 
bool readFromCommandLine (int argc, char *argv[], std::string varName)
 
template<typename T >
readFromCommandLine (int argc, char *argv[], std::string varName, T value)
 
template<typename T , size_t n>
std::array< T, nreadArrayFromCommandLine (int argc, char *argv[], std::string varName, std::array< T, n > value)
 
template<typename T >
std::vector< T > readVectorFromCommandLine (int argc, char *argv[], std::string varName, size_t n, std::vector< T > values)
 
template<>
std::string readFromCommandLine< std::string > (int argc, char *argv[], std::string varName, std::string value)
 
bool removeFromCommandline (int &argc, char *argv[], std::string varName, int nArgs)
 
void check (double real, double ideal, double error, std::string errorMessage)
 
void check (Vec3D real, Vec3D ideal, double error, std::string errorMessage)
 
void check (Matrix3D real, Matrix3D ideal, double error, std::string errorMessage)
 
void check (MatrixSymmetric3D real, MatrixSymmetric3D ideal, double error, std::string errorMessage)
 
std::string getPath ()
 
Mdouble getRealTime ()
 
bool isNext (std::istream &is, const std::string name)
 
bool createDirectory (std::string)
 
Mdouble round (const Mdouble value, unsigned precision)
 
Mdouble getRayleighTime (Mdouble radius, Mdouble shearModulus, Mdouble poisson, Mdouble density)
 
std::vector< Mdoublelinspace (Mdouble a, Mdouble b, int N)
 

Function Documentation

◆ addToFile()

bool helpers::addToFile ( std::string  filename,
std::string  filecontent 
)

Adds a string to an existing file.

485 {
486  std::fstream file;
487  file.open(filename.c_str(), std::ios::app);
488  if (file.fail())
489  {
490  logger(INFO, "Error in writeToFile: file could not be opened");
491  return false;
492  }
493  file << filecontent;
494  file.close();
495  return true;
496 }

References INFO, and logger.

◆ check() [1/4]

◆ check() [2/4]

void helpers::check ( Matrix3D  real,
Matrix3D  ideal,
double  error,
std::string  errorMessage 
)
927 {
928  checkTemplate(real,ideal,error,whatIsChecked);
929 }

References checkTemplate().

◆ check() [3/4]

void helpers::check ( MatrixSymmetric3D  real,
MatrixSymmetric3D  ideal,
double  error,
std::string  errorMessage 
)
922 {
923  checkTemplate(real,ideal,error,whatIsChecked);
924 }

References checkTemplate().

◆ check() [4/4]

void helpers::check ( Vec3D  real,
Vec3D  ideal,
double  error,
std::string  errorMessage 
)
917 {
918  checkTemplate(real,ideal,error,whatIsChecked);
919 }

References checkTemplate().

◆ compare()

bool helpers::compare ( std::istream &  is,
std::string  s 
)

Checks if the next argument in the input stream is a certain string.

859 {
860  // Get current position
861  //check if the next line starts with 'interactionFile'; otherwise, skip interaction
862  int len = is.tellg();
863  std::string dummy;
864  is >> dummy;
865  if (dummy != s)
866  {
867  is.seekg(len, std::ios_base::beg);
868  logger(INFO, "helpers::compare: Next stream value (%) is not %", dummy, s);
869  return false;
870  }
871  return true;
872 }

References INFO, and logger.

Referenced by DPMBase::read(), File::write(), and DPMBase::write().

◆ computeDisptFromCollisionTimeAndRestitutionCoefficientAndTangentialRestitutionCoefficientAndEffectiveMass()

helpers::KAndDispAndKtAndDispt helpers::computeDisptFromCollisionTimeAndRestitutionCoefficientAndTangentialRestitutionCoefficientAndEffectiveMass ( Mdouble  tc,
Mdouble  r,
Mdouble  beta,
Mdouble  mass 
)

Set disp, k, dispt and kt such that is matches a given collision time tc and a normal and tangential restitution coefficient r, beta for a collision of effective/reduced mass m. From Deen...Kuipers2006, eq. 43 and 30.

Todo:
what should be used instead of this function?
332 {
334  ans.disp = -2.0 * mass * log(r) / tc;
335  ans.k = mass * (mathsFunc::square(constants::pi / tc) + mathsFunc::square(ans.disp / (2.0 * mass)));
336  ans.kt = 2.0 / 7.0 * ans.k * (mathsFunc::square(constants::pi) + mathsFunc::square(log(beta))) /
338  if (beta != 0.0)
339  ans.dispt = -2 * log(beta) *
340  sqrt(1.0 / 7.0 * mass * ans.kt / (mathsFunc::square(constants::pi) + mathsFunc::square(log(beta))));
341  else
342  ans.dispt = 2. * sqrt(1.0 / 7.0 * mass * ans.kt);
343  return ans;
344 }

References mathsFunc::beta(), helpers::KAndDispAndKtAndDispt::disp, helpers::KAndDispAndKtAndDispt::dispt, helpers::KAndDispAndKtAndDispt::k, helpers::KAndDispAndKtAndDispt::kt, mathsFunc::log(), constants::pi, and mathsFunc::square().

◆ createDirectory()

bool helpers::createDirectory ( std::string  path)
1035  {
1036  //see https://stackoverflow.com/questions/20358455/cross-platform-way-to-make-a-directory
1037  mode_t nMode = 0733; // UNIX style permissions
1038  int nError = 0;
1039 #if defined(_WIN32)
1040  nError = _mkdir(path.c_str()); // can be used on Windows
1041 #else
1042  nError = mkdir(path.c_str(),nMode); // can be used on non-Windows
1043 #endif
1044  if (nError != 0) {
1045  // handle your error here
1046  }
1047  return false;
1048 }

◆ fileExists()

bool helpers::fileExists ( std::string  strFilename)

Function to check if a file exists, is used to check if a run has already need done.

This is a FileExist routine, which is used to test if a run have already need preformed, allows me to plug holes in parm studies.

503 {
504  struct stat stFileInfo;
505  bool blnReturn;
506  int intStat;
507 
508  // Attempt to get the file attributes
509 
510  intStat = stat(strFilename.c_str(), &stFileInfo);
511  if (intStat == 0)
512  {
513  // We were able to get the file attributes
514  // so the file obviously exists.
515  blnReturn = true;
516  }
517  else
518  {
519  // We were not able to get the file attributes.
520  // This may mean that we don't have permission to
521  // access the folder which contains this file. If you
522  // need to do that level of checking, lookup the
523  // return values of stat which will give you
524  // more details on why stat failed.
525  blnReturn = false;
526  }
527 
528  return blnReturn;
529 }

Referenced by ParameterStudy1DDemo::actionsBeforeTimeLoop(), ParameterStudy2DDemo::actionsBeforeTimeLoop(), ParameterStudy3DDemo::actionsBeforeTimeLoop(), main(), DPMBase::readRestartFile(), FlowRule::run(), AngleOfRepose::run(), and SegregationPeriodic::setupInitialConditions().

◆ getEffectiveMass()

Mdouble helpers::getEffectiveMass ( Mdouble  mass0,
Mdouble  mass1 
)

Calculates the effective mass of a particle pair, i.e. half the harmonic mean of two particle masses.

The effective mass is an important parameter in a collision. E.g. the collision time and the restitution coefficient are functions of the effective mass.

Parameters
[in]mass0The mass of the first particle.
[in]mass1The mass of the second particle.
Returns
the effective mass of the particle pair.
556 {
557  return mass0 * mass1 / (mass0 + mass1);
558 }

Referenced by main().

◆ getLineFromStringStream()

void helpers::getLineFromStringStream ( std::istream &  in,
std::stringstream &  out 
)

Reads a line from one stringstream into another, and prepares the latter for reading in.

This function is used to avoid errors from reading in old or manually modified restart files. Instead of reading variable by variable directly from the restart stringstream, a full line is read first, from which the variables are read. Thus, if a line has the wrong number of arguments, it might affect the reading of the current line, but correctly reads the next line.

Example of usage:

std::stringstream line; std::stringstream is = restartFile.getFStream(); helpers::getLineFromStringStream(is, line); std::string dummy; line >> dummy;

Parameters
[in]inthe stringstream from which a line is read out should be initialized as std::stringstream(std::stringstream::out)
[out]outthe stringstream into which the line is read; should be initialized as std::stringstream(std::stringstream::in | std::stringstream::out)
424 {
425  std::string line_string;
426  getline(in, line_string);
427  out.str(std::move(line_string));
428  out.clear();
429 }

Referenced by main(), BaseHandler< T >::read(), InteractionHandler::read(), DPMBase::read(), BaseCluster::read(), MercuryBase::read(), SpeciesHandler::readAndAddObject(), WallHandler::readAndCreateOldObject(), DPMBase::readNextDataFile(), BoundaryHandler::readOldObject(), and SpeciesHandler::readOldObject().

◆ getMaximumVelocity()

Mdouble helpers::getMaximumVelocity ( Mdouble  k,
Mdouble  disp,
Mdouble  radius,
Mdouble  mass 
)

Calculates the maximum relative velocity allowed for a normal collision of two particles of radius r and particle mass m (for higher velocities particles could pass through each other)

Todo:
what should be used instead of this function?
347 {
348  // note: underestimate based on energy argument,
349  // Ekin = 2*1/2*m*(v/2)^2 = 1/2*k*(2*r)^2, gives
350  // return radius * sqrt(8.0*k/mass);
351 
352  // with dissipation, see S. Luding, Collisions & Contacts between two particles, eq 34
353  Mdouble w = sqrt(k / mass - mathsFunc::square(disp / (2.0 * mass)));
354  Mdouble w0 = sqrt(k / mass);
355  Mdouble DispCorrection = exp(-disp / mass / w) * asin(w / w0);
356  //std::cout << "DispCorrection" << DispCorrection << std::endl;
357  return radius * sqrt(8.0 * k / mass) / DispCorrection;
358 }

References mathsFunc::exp(), and mathsFunc::square().

Referenced by ChutePeriodic::setupInitialConditions(), Vreman::setupInitialConditions(), and SilbertPeriodic::setupInitialConditions().

◆ getPath()

std::string helpers::getPath ( )
932 {
933  char cCurrentPath[FILENAME_MAX];
934  if (not GetCurrentDir(cCurrentPath, sizeof(cCurrentPath))) {
935  logger(WARN, "Get current dir failed: %", cCurrentPath);
936  }
937  return std::string(cCurrentPath);
938 }

References GetCurrentDir, logger, and WARN.

Referenced by DPMBase::writePythonFileForVTKVisualisation(), and HorizontalMixer::writeScript().

◆ getRayleighTime()

Mdouble helpers::getRayleighTime ( Mdouble  radius,
Mdouble  shearModulus,
Mdouble  poisson,
Mdouble  density 
)
1050  {
1051  return constants::pi*radius*sqrt(density/shearModulus)/(0.1631*poisson+0.8766);
1052 }

References constants::pi.

Referenced by NozzleDemo::setupInitialConditions(), and NozzleSelfTest::setupInitialConditions().

◆ getRealTime()

Mdouble helpers::getRealTime ( )
941 {
942  // record start time
943  static auto start = std::chrono::steady_clock::now();
944  auto end = std::chrono::steady_clock::now();
945  std::chrono::duration<double> diff = end - start;
946  return diff.count();
947 }

◆ getSaveCountFromNumberOfSavesAndTimeMaxAndTimeStep()

unsigned int helpers::getSaveCountFromNumberOfSavesAndTimeMaxAndTimeStep ( unsigned int  numberOfSaves,
Mdouble  timeMax,
Mdouble  timeStep 
)

Returns the correct saveCount if the total number of saves, the final time and the time step is known.

MercuryDPM uses the DPMBase::setSaveCount to determine how often output is written. However, if the user wants to set the total number of saves instead of the saveCount, he can use this function to calculate the correct saveCount, assuming that the final time and the mean time step is known.

Example of use:

setSaveCount(helpers::getSaveCountFromNumberOfSavesAndTimeMaxAndTimeStep(numberOfSaves, getTimeMax(), getTimeStep()));

Parameters
[in]numberOfSavesthe total number of output files the user wants at the end of the simulation.
[in]timeMaxthe final time of the simulation
[in]timeStepthe mean time step used during the simulation
Returns
the saveCount value that should be used to get the desired number of saves.
376 {
377  if (numberOfSaves > 0 && timeMax > 0 && timeStep > 0)
378  {
379  return static_cast<unsigned int>(ceil(
380  (timeMax + timeStep) / timeStep / static_cast<double>(numberOfSaves - 1)));
381  }
382  else
383  {
384  logger(ERROR,
385  "[Helpers::getSaveCountFromNumberOfSavesAndTimeMaxAndTimeStep()] numberOfSaves: %, timeMax: %, "
386  "timestep: %\n Arguments need to be positive",
387  numberOfSaves, timeMax, timeStep);
388  }
389 }

References ERROR, and logger.

Referenced by ChutePeriodicDemo::ChutePeriodicDemo(), ChuteWithWedge::ChuteWithWedge(), main(), ChutePeriodic::setup(), AxisymmetricHopper::setupInitialConditions(), MovingIntersectionOfWallsUnitTest_Basic::setupInitialConditions(), MovingWall::setupInitialConditions(), SeparateFilesSelfTest::setupInitialConditions(), and Slide::Slide().

◆ gnuplot()

void helpers::gnuplot ( std::string  command)

Plots to a gnuplot window.

470 {
471 #ifdef __CYGWIN__
472  logger(WARN, "[helpers::gnuplot] is not supported on Cygwin");
473 #else
474 #ifdef WINDOWS
475  logger(WARN, "[helpers::gnuplot] is not supported on Windows");
476 #else
477  FILE* pipe = popen("gnuplot -persist", "w");
478  fprintf(pipe, "%s", command.c_str());
479  fflush(pipe);
480 #endif
481 #endif
482 }

References logger, and WARN.

◆ isNext()

bool helpers::isNext ( std::istream &  is,
const std::string  name 
)

reads next value in stream as a string and compares it with name. If name is equal to string, the function outputs true. If name is not equal to string, the function undoes the read by setting seekg, and outputs false.

Parameters
is
name
Returns
957  {
958  std::string dummy;
959  auto pos = is.tellg();
960  is >> dummy;
961  if (dummy != name) {
962  is.seekg(pos);
963  return false;
964  } else {
965  return true;
966  }
967 }

References units::name.

Referenced by BaseWall::read().

◆ linspace()

std::vector< Mdouble > helpers::linspace ( Mdouble  a,
Mdouble  b,
int  N 
)
888 {
889  Mdouble dx = (Max - Min) / static_cast<Mdouble>(numberOfBins - 1);
890  Mdouble val;
891  std::vector<Mdouble> linearVector(numberOfBins);
892  typename std::vector<Mdouble>::iterator x;
893  for (x = linearVector.begin(), val = Min; x != linearVector.end(); ++x, val += dx)
894  {
895  *x = val;
896  }
897  // ensure that last value is equal to Max.
898  linearVector.pop_back();
899  linearVector.push_back(Max);
900  return linearVector;
901 }

Referenced by PSD::setDistributionNormal(), and PSD::setDistributionUniform().

◆ loadingTest()

void helpers::loadingTest ( const ParticleSpecies species,
Mdouble  displacement,
Mdouble  velocity,
Mdouble  radius,
std::string  name 
)

Used to test the loading/unloading properties of a contact law

Creates a DPMBase with a particles of unit size and a flat wall and loads/unloads the particle-wall contact

Parameters
[in]speciesparticle species specifying the contact law
[in]displacementpeak displacement before unloading
[in]velocityloading/unloading velocity
620 {
621  class LoadingTest : public DPMBase
622  {
623  const ParticleSpecies* species;
624  Mdouble displacement;
625  Mdouble velocity;
626  Mdouble radius;
627  public:
628  //public variables
629  LoadingTest(const ParticleSpecies* species, Mdouble displacement, Mdouble velocity, Mdouble radius)
630  : species(species), displacement(displacement), velocity(velocity), radius(radius)
631  {}
632 
633  void setupInitialConditions() override
634  {
635  //setName("LoadingTest"+species->getName());
636  setTimeMax(2.0 * displacement / velocity);
637  setTimeStep(2e-3 * getTimeMax());
638  setSaveCount(1);
639  setFileType(FileType::NO_FILE);
640  fStatFile.setFileType(FileType::ONE_FILE);
641 
642  setMax({radius, radius, radius + radius});
643  setMin({-radius, -radius, 0});
646 
647  speciesHandler.copyAndAddObject(*species);
648 
650  p.setSpecies(speciesHandler.getObject(0));
651  p.setRadius(radius);
652  p.setPosition({0, 0, radius});
653  particleHandler.copyAndAddObject(p);
654 
655  InfiniteWall w;
656  w.setSpecies(speciesHandler.getObject(0));
657  w.set(Vec3D(0, 0, -1), Vec3D(0.0, 0.0, 0.0));
658  wallHandler.copyAndAddObject(w);
659  }
660 
661  void actionsBeforeTimeStep() override
662  {
663  BaseParticle* p = particleHandler.getLastObject();
664  logger.assert_debug(p,"Empty particle handler");
665  p->setAngularVelocity({0, 0, 0});
666 
667  //Moving particle normally into surface
668  if (getTime() <= displacement / velocity)
669  {
670  p->setVelocity({0, 0, velocity});
671  p->setPosition({0, 0, radius - velocity * getTime()});
672  }
673  else
674  {
675  p->setVelocity({0, 0, -velocity});
676  p->setPosition({0, 0, radius - displacement - displacement + velocity * getTime()});
677  }
678  }
679  } test(species, displacement, velocity, radius);
680  test.setName(name);
681  test.solve();
682  writeToFile(test.getName() + ".gnu", "plot '" + test.getName() + ".fstat' u 7:9 w lp");
683  logger(INFO, "finished loading test: run 'gnuplot %.gnu' to view output", test.getName());
684 }

References DPMBase::actionsBeforeTimeStep(), DPMBase::getTime(), DPMBase::getTimeMax(), INFO, logger, units::name, InfiniteWall::set(), BaseInteractable::setAngularVelocity(), DPMBase::setFileType(), DPMBase::setMax(), DPMBase::setMin(), DPMBase::setParticleDimensions(), BaseInteractable::setPosition(), BaseParticle::setRadius(), DPMBase::setSaveCount(), BaseParticle::setSpecies(), BaseWall::setSpecies(), DPMBase::setSystemDimensions(), DPMBase::setTimeMax(), DPMBase::setTimeStep(), DPMBase::setupInitialConditions(), BaseInteractable::setVelocity(), and writeToFile().

Referenced by main().

◆ lower()

std::string helpers::lower ( std::string  s)
55  {
56  std::transform(s.begin(), s.end(), s.begin(), ::tolower);
57  return s;
58 }

Referenced by WallHandler::readTriangleWall().

◆ more()

void helpers::more ( std::string  filename,
unsigned  nLines = constants::unsignedMax 
)
581 {
582  if (nLines != constants::unsignedMax)
583  logger(INFO, "First % lines of %:\n", Flusher::NO_FLUSH, nLines, filename);
584  std::fstream file;
585  file.open(filename.c_str(), std::ios::in);
586  if (file.fail())
587  logger(ERROR, "Error in readArrayFromFile: file could not be opened");
588  std::string line;
589  for (unsigned i = 0; i < nLines; i++)
590  {
591  if (file.eof()) break;
592  getline(file, line);
593  logger(INFO, " %\n", line);
594  }
595  file.close();
596 }

References ERROR, FLUSH, constants::i, INFO, logger, and constants::unsignedMax.

Referenced by commandLineCG(), and main().

◆ normalAndTangentialLoadingTest()

void helpers::normalAndTangentialLoadingTest ( const ParticleSpecies species,
Mdouble  displacement,
Mdouble  tangentialDisplacement,
Mdouble  velocity,
Mdouble  radius,
std::string  name 
)

Used to test the tangential loading/unloading/reloading properties of a frictional contact law

Creates a DPMBase with a particles of unit size and a flat wall and loads/unloads/reloads the particle-wall contact in tangential direction

Parameters
[in]speciesparticle species specifying the contact law
[in]displacementpeak displacement before unloading
[in]velocityloading/unloading velocity
695 {
696  class LoadingTest : public DPMBase
697  {
698  const ParticleSpecies* species;
699  Mdouble displacement;
700  Mdouble tangentialDisplacement;
701  Mdouble velocity;
702  Mdouble radius;
703  public:
704  //public variables
705  LoadingTest(const ParticleSpecies* species, Mdouble displacement, Mdouble tangentialDisplacement,
706  Mdouble velocity, Mdouble radius)
707  : species(species), displacement(displacement), tangentialDisplacement(tangentialDisplacement),
708  velocity(velocity), radius(radius)
709  {}
710 
711  void setupInitialConditions() override
712  {
713  //setName("TangentialLoadingTest"+species->getName());
714  setTimeMax(4.0 * tangentialDisplacement / velocity);
715  setTimeStep(4e-4 * getTimeMax());
716  setSaveCount(1);
717  setFileType(FileType::NO_FILE);
718  fStatFile.setFileType(FileType::ONE_FILE);
719 
720  setMax({radius, radius, radius + radius});
721  setMin({-radius, -radius, 0});
724 
725  speciesHandler.copyAndAddObject(*species);
726 
728  p.setSpecies(speciesHandler.getObject(0));
729  p.setRadius(radius);
730  p.setPosition({0, 0, radius - displacement});
731  particleHandler.copyAndAddObject(p);
732 
733  InfiniteWall w;
734  w.setSpecies(speciesHandler.getObject(0));
735  w.set(Vec3D(0, 0, -1), Vec3D(0.0, 0.0, 0.0));
736  wallHandler.copyAndAddObject(w);
737  }
738 
739  void actionsBeforeTimeStep() override
740  {
741  BaseParticle* p = particleHandler.getLastObject();
742  logger.assert_debug(p,"Empty particle handler");
743  p->setAngularVelocity({0, 0, 0});
744 
745  //Moving particle cyclically right and left between +-tangentialDisplacement
746  bool moveRight = static_cast<int>(getTime() / (2.0*tangentialDisplacement / velocity) +0.5)%2==0;
747  if (moveRight)
748  {
749  p->setVelocity({-velocity, 0, 0});
750  p->setPosition({tangentialDisplacement - velocity * getTime(), 0, radius - displacement});
751  }
752  else
753  {
754  p->setVelocity({velocity, 0, 0});
755  p->setPosition({-2*tangentialDisplacement + velocity * getTime(), 0, radius - displacement});
756  }
757  }
758 
759  } test(species, displacement, tangentialDisplacement, velocity, radius);
760  test.setName(name);
761  test.solve();
762  writeToFile(test.getName() + ".gnu", "plot '" + test.getName() + ".fstat' u 8:($10*$14) w lp");
763  logger(INFO, "finished tangential loading test: run 'gnuplot %.gnu' to view output", test.getName());
764 }

References DPMBase::actionsBeforeTimeStep(), DPMBase::getTime(), DPMBase::getTimeMax(), INFO, logger, units::name, InfiniteWall::set(), BaseInteractable::setAngularVelocity(), DPMBase::setFileType(), DPMBase::setMax(), DPMBase::setMin(), DPMBase::setParticleDimensions(), BaseInteractable::setPosition(), BaseParticle::setRadius(), DPMBase::setSaveCount(), BaseParticle::setSpecies(), BaseWall::setSpecies(), DPMBase::setSystemDimensions(), DPMBase::setTimeMax(), DPMBase::setTimeStep(), DPMBase::setupInitialConditions(), BaseInteractable::setVelocity(), and writeToFile().

Referenced by main().

◆ objectivenessTest()

void helpers::objectivenessTest ( const ParticleSpecies species,
Mdouble  displacement,
Mdouble  tangentialDisplacement,
Mdouble  velocity,
Mdouble  radius,
std::string  name 
)

Used to test the tangential loading/unloading/reloading properties of a frictional contact law

Creates a DPMBase with a particles of unit size and a flat wall, loads the particle-wall contact in normal and tangential direction, then rotates.

Parameters
[in]speciesparticle species specifying the contact law
[in]displacementpeak displacement before unloading
[in]velocityloading/unloading velocity
774 {
775  class ObjectivenessTest : public DPMBase
776  {
777  const ParticleSpecies* species;
778  Mdouble displacement;
779  Mdouble tangentialDisplacement;
780  Mdouble velocity;
781  Mdouble radius;
782  public:
783  //public variables
784  ObjectivenessTest(const ParticleSpecies* species, Mdouble displacement, Mdouble tangentialDisplacement,
785  Mdouble velocity, Mdouble radius)
786  : species(species), displacement(displacement), tangentialDisplacement(tangentialDisplacement),
787  velocity(velocity), radius(radius)
788  {}
789 
790  void setupInitialConditions() override
791  {
792  //setName("ObjectivenessTest"+species->getName());
793  setTimeMax((tangentialDisplacement + 0.5 * constants::pi * radius) / velocity);
794  setTimeStep(1e-4 * getTimeMax());
795  setSaveCount(20);
796  setFileType(FileType::NO_FILE);
797  dataFile.setFileType(FileType::ONE_FILE);
798  fStatFile.setFileType(FileType::ONE_FILE);
799 
800  setMax(radius * Vec3D(2, 2, 2));
801  setMin(radius * Vec3D(-2, -2, -2));
804 
805  speciesHandler.copyAndAddObject(*species);
806 
808  p.setSpecies(speciesHandler.getObject(0));
809  p.setRadius(radius);
810  p.setPosition({0, radius - displacement, 0});
811  particleHandler.copyAndAddObject(p);
812  p.setPosition({0, -radius + displacement, 0});
813  particleHandler.copyAndAddObject(p);
814  }
815 
816  void actionsBeforeTimeStep() override
817  {
818  BaseParticle* p = particleHandler.getObject(0);
819  BaseParticle* q = particleHandler.getLastObject();
820  logger.assert_debug(p,"Empty particle handler");
821  logger.assert_debug(q,"Empty particle handler");
822 
823  //Moving particle normally into surface
824  if (getTime() <= tangentialDisplacement / velocity)
825  {
826  p->setAngularVelocity({0, 0, 0});
827  p->setVelocity({velocity, 0, 0});
828  p->setPosition({-tangentialDisplacement + velocity * getTime(), radius - displacement, 0});
829  q->setAngularVelocity({0, 0, 0});
830  q->setVelocity({-velocity, 0, 0});
831  q->setPosition({tangentialDisplacement - velocity * getTime(), -radius + displacement, 0});
832  }
833  else
834  {
835  Mdouble angle = velocity / (radius - displacement) * (getTime() - tangentialDisplacement / velocity);
836  Mdouble s = sin(angle);
837  Mdouble c = cos(angle);
838  p->setAngularVelocity(velocity / (radius - displacement) * Vec3D(0, 0, -1));
839  //p->setAngularVelocity(Vec3D(0,0,0));
840  p->setOrientation({1, 0, 0, 0});
841  p->setVelocity(velocity * Vec3D(c, -s, 0));
842  //p->setVelocity(Vec3D(0,0,0));
843  p->setPosition((radius - displacement) * Vec3D(s, c, 0));
845  q->setOrientation(-p->getOrientation());
846  q->setVelocity(-p->getVelocity());
847  q->setPosition(-p->getPosition());
848  }
849  }
850 
851  } test(species, displacement, tangentialDisplacement, velocity, radius);
852  test.setName(name);
853  test.solve();
854  writeToFile(test.getName() + ".gnu", "set size ratio -1; plot '" + test.getName() + ".fstat' u 14:15 every 2 w lp");
855  logger(INFO, "finished objectiveness test: run 'gnuplot %.gnu' to view output", test.getName());
856 }

References DPMBase::actionsBeforeTimeStep(), mathsFunc::cos(), BaseInteractable::getAngularVelocity(), BaseInteractable::getOrientation(), BaseInteractable::getPosition(), DPMBase::getTime(), DPMBase::getTimeMax(), BaseInteractable::getVelocity(), INFO, logger, units::name, constants::pi, BaseInteractable::setAngularVelocity(), DPMBase::setFileType(), DPMBase::setMax(), DPMBase::setMin(), BaseInteractable::setOrientation(), DPMBase::setParticleDimensions(), BaseInteractable::setPosition(), BaseParticle::setRadius(), DPMBase::setSaveCount(), BaseParticle::setSpecies(), DPMBase::setSystemDimensions(), DPMBase::setTimeMax(), DPMBase::setTimeStep(), DPMBase::setupInitialConditions(), BaseInteractable::setVelocity(), mathsFunc::sin(), and writeToFile().

Referenced by main().

◆ openFile()

bool helpers::openFile ( std::fstream &  file,
std::string  filename,
std::fstream::openmode  mode 
)

Provides a simple interface for opening a file.

Provides a simple interface for opening a file, in order to avoid that the user has to learn the syntax for opening a file.

Parameters
[out]fileThe std::fstream object that the user can write to.
[in]filenameThe name of the file.
[in]modeThe openmode of the file, typically std::fstream::out or std::fstream::in.
Returns
true is the file was successfully opened, false else.
540 {
541  file.open(filename.c_str(), mode);
542  if (file.fail())
543  return false;
544  else
545  return true;
546 }

Referenced by main().

◆ readArrayFromCommandLine()

template<typename T , size_t n>
std::array<T,n> helpers::readArrayFromCommandLine ( int  argc,
char argv[],
std::string  varName,
std::array< T, n value 
)
368 {
369  for (unsigned i=0; i<argc-1; ++i) {
370  if (varName == argv[i]) {
371  unsigned j = i+1;
372  std::stringstream out;
373  for (auto& v : value) {
374  v = atof(argv[j]);
375  out << v << ' ';
376  ++j;
377  }
378  logger(INFO, "readFromCommandLine: % set to % ", varName.substr(1), out.str());
379  return value;
380  }
381  }
382  //if the variable is not found
383  std::stringstream out;
384  for (auto& v : value) out << v << ' ';
385  logger(INFO, "readFromCommandLine: % set to default value % ", varName.substr(1), out.str());
386  return value;
387 }

References constants::i, INFO, and logger.

◆ readArrayFromFile()

std::vector< double > helpers::readArrayFromFile ( std::string  filename,
int n,
int m 
)
561 {
562  std::fstream file;
563  file.open(filename.c_str(), std::ios::in);
564  if (file.fail())
565  {
566  logger(ERROR, "Error in readArrayFromFile: file could not be opened");
567  }
568  file >> n >> m;
569  std::vector<double> v;
570  Mdouble val;
571  for (int i = 0; i < n * m; i++)
572  {
573  file >> val;
574  v.push_back(val);
575  }
576  file.close();
577  return v;
578 }

References ERROR, constants::i, logger, and n.

◆ readFromCommandLine() [1/2]

bool helpers::readFromCommandLine ( int  argc,
char argv[],
std::string  varName 
)

Returns true if command line arguments contain varName, false else Usage example: if (readFromCommandLine(argc, argv, '-verbose')) ...

Parameters
argcpass through number of command line arguments
argvpass through values of command line arguments
varNamename of command line arguments that is required to return true
Returns
true or false
875 {
876  for (unsigned i=0; i<argc; ++i) {
877  if (varName == argv[i]) {
878  logger(INFO, "readFromCommandLine: % set to true",varName.substr(1));
879  return true;
880  }
881  }
882  //if the variable is not found
883  logger(INFO, "readFromCommandLine: % set to default value false",varName.substr(1));
884  return false;
885 }

References constants::i, INFO, and logger.

Referenced by main().

◆ readFromCommandLine() [2/2]

template<typename T >
T helpers::readFromCommandLine ( int  argc,
char argv[],
std::string  varName,
value 
)
353 {
354  for (unsigned i=0; i<argc-1; ++i) {
355  if (varName == argv[i]) {
356  value = atof(argv[i+1]);
357  logger(INFO, "readFromCommandLine: % set to % ", varName.substr(1), value);
358  return value;
359  }
360  }
361  //if the variable is not found
362  logger(INFO, "readFromCommandLine: % set to default value % ", varName.substr(1), value);
363  return value;
364 }

References constants::i, INFO, and logger.

◆ readFromCommandLine< std::string >()

template<>
std::string helpers::readFromCommandLine< std::string > ( int  argc,
char argv[],
std::string  varName,
std::string  value 
)

◆ readFromFile()

template<typename T >
T helpers::readFromFile ( std::string  fileName,
std::string  varName,
value 
)

Allows a quick read-in from a parameter file.

For example, the following code reads in the variable time from the file in:

double time = readFromFile("in","time",24);

The in file needs to contain the string time, followed by the value, e.g.

time 20

If the file cannot be opened, or the parameter is not found, the variable is set to the default value specified by the third argument.

Parameters
fileNamename of input
varNamevariable name as it appears in the input file
valuedefault value (used if the parameter could not be read)
Returns
value of variable
312 {
313  //open filestream
314  std::ifstream is(fileName.c_str(), std::ios::in);
315  if (is.fail())
316  {
317  logger(INFO, "readFromFile: file % could not be opened, variable % set to default value %",
318  fileName, varName, value);
319  return value;
320  }
321 
322  //read in variables, until the right one is fount
323  std::string s;
324  while (!is.eof())
325  {
326  is >> s;
327  if (!s.compare(varName))
328  {
329  is >> value;
330  logger(INFO, "readFromFile: variable % set to % ", varName, value);
331  return value;
332  }
333  }
334 
335  //if the right variable is never found
336  logger(WARN, "readFromFile: variable % not set in file %, using default value % ", varName, fileName, value);
337  return value;
338 }

References INFO, logger, and WARN.

Referenced by main().

◆ readOptionalVariable()

template<typename T >
bool helpers::readOptionalVariable ( std::istream &  is,
const std::string &  name,
T &  variable 
)

Reads optional variables in the restart file.

A variable is stored in the restart file by storing the variables name and value, e.g. " name value" If a variable is always written to the restart file, it can be read-in like this: is >> dummy >> variable; If a variable is optional, the variable name has to be checked, and should be read in like this: readOptionalVariable(is,"name",variable);

Todo:
readOptionalVariable should check the full variable name, not just the next value. However, I don't know how to put the location in the ifstream back.
248 {
251  const auto pos = is.tellg();
252  std::string dummy;
253  is >> dummy;
254  if (dummy == name)
255  {
256  is >> variable;
257  return true;
258  } else {
259  is.seekg(pos);
260  return false;
261  }
262 }

References units::name.

Referenced by LiquidMigrationWilletSpecies::read(), InsertionBoundary::read(), BaseInteraction::read(), BaseSpecies::read(), InfiniteWall::read(), and BaseObject::read().

◆ readVectorFromCommandLine()

template<typename T >
std::vector<T> helpers::readVectorFromCommandLine ( int  argc,
char argv[],
std::string  varName,
size_t  n,
std::vector< T >  values 
)
391 {
392  for (unsigned i=0; i<argc-1; ++i) {
393  if (varName == argv[i]) {
394  // read until the next argument starts
395  values.resize(0);
396  std::stringstream out;
397  for (int j = i+1; j<argc and argv[j][0]!='-'; ++j) {
398  values.push_back(atof(argv[j]));
399  out << values.back() << ' ';
400  }
401  logger(INFO, "readFromCommandLine: % set to % ", varName.substr(1), out.str());
402  return values;
403  }
404  }
405  //if the variable is not found
406  std::stringstream out;
407  for (auto& v: values) out << v << ' ';
408  logger(INFO, "readFromCommandLine: % set to default value % ", varName.substr(1), out.str());
409  return values;
410 }

References constants::i, INFO, and logger.

◆ removeFromCommandline()

bool helpers::removeFromCommandline ( int argc,
char argv[],
std::string  varName,
int  nArgs 
)
Parameters
[in]argc
[in]argv
[in]varNameThe name of the commandline argument to be removed
[in]nArgsThe number of additional command line arguments to remove
Returns
boolean True if the argument was found and removed

This function does hide the desired argument from the supplied argv argc combination. It does it by moving the specified arguent to the end of the supplied argv and reducing argc by the coorect number. Other pieces of code that rely on argc should therefor no longer see the hidden argument.

If used in comination with readFromCommandLine(...), this function allows handling of arguments that are not seen by solve(), even if commandline arguments are passed.

1000 {
1001  unsigned int i, j;
1002 
1003  for (i=0; i<argc; ++i) {
1004  if (varName == argv[i])
1005  {
1006  char *tmp[nArgs+1];
1007 
1008  for (j=i; j<argc-1-nArgs; j++)
1009  {
1010  // Store the pointers of the handled argument
1011  if ( j < i + 1 + nArgs)
1012  {
1013  tmp[j-i] = argv[j];
1014  }
1015 
1016  // Move the pointers after the argument to the front
1017  argv[j] = argv[j + nArgs + 1];
1018 
1019  }
1020 
1021  // Move the stored argument to the end
1022  for (j=argc-1-nArgs; j< argc; j++)
1023  {
1024  argv[j] = tmp[j + 1 + nArgs - argc];
1025  }
1026  argc -= nArgs + 1;
1027 
1028  return true;
1029  }
1030  }
1031 
1032  return false;
1033 }

References constants::i.

Referenced by main().

◆ round()

◆ to_string() [1/2]

◆ to_string() [2/2]

std::string helpers::to_string ( Mdouble  value,
unsigned  precision 
)
606 {
607  std::ostringstream stm;
608  stm << round(value,precision);
609  return stm.str();
610 }

References round().

◆ writeCommandLineToFile()

void helpers::writeCommandLineToFile ( const std::string  filename,
const int  argc,
char *const  argv[] 
)

Writes a string to a file.

460 {
461  std::stringstream ss;
462  for (int i=0; i<argc; ++i) {
463  ss << argv[i] << ' ';
464  }
465  writeToFile(filename,ss.str());
466 }

References constants::i, and writeToFile().

◆ writeToFile()

bool helpers::writeToFile ( std::string  filename,
std::string  filecontent 
)

Writes a string to a file.

Provides a simple interface for writing a string to a file. This function is mainly used to create ini or restart file that the code later reads back in.

Example of usage:

helpers::writeToFile("RestartUnitTest.ini", "1 0 0 0 0 1 1 1\n" "0.5 0.5 0 0 0 0.5 0 0 0 0 0 0 0 0\n");

Parameters
[in]filenamethe name of the file
[in]filecontentthe content
Returns
true on success.
446 {
447  std::fstream file;
448  file.open(filename.c_str(), std::ios::out);
449  if (file.fail())
450  {
451  logger(WARN, "Error in writeToFile: file could not be opened");
452  return false;
453  }
454  file << filecontent;
455  file.close();
456  return true;
457 }

References logger, and WARN.

Referenced by statistics_while_running< T >::finishStatistics(), loadingTest(), main(), normalAndTangentialLoadingTest(), objectivenessTest(), LevelSetWall::setShapeFourSided(), GetDistanceAndNormalForIntersectionOfWalls::setupInitialConditions(), GetDistanceAndNormalForScrew::setupInitialConditions(), GetDistanceAndNormalForTriangleWall::setupInitialConditions(), GetDistanceAndNormalForTriangleWalls::setupInitialConditions(), SingleParticle< SpeciesType >::setupInitialConditions(), SinterPair::SinterPair(), writeCommandLineToFile(), DPMBase::writePythonFileForVTKVisualisation(), HorizontalMixer::writeScript(), LevelSetWall::writeToFile(), and Screw::writeVTK().

helpers::KAndDispAndKtAndDispt::k
Mdouble k
Definition: Helpers.h:151
DPMBase::setMax
void setMax(const Vec3D &max)
Sets the maximum coordinates of the problem domain.
Definition: DPMBase.cc:1073
BaseInteractable::getAngularVelocity
virtual const Vec3D & getAngularVelocity() const
Returns the angular velocity of this interactable.
Definition: BaseInteractable.cc:341
mathsFunc::square
T square(const T val)
squares a number
Definition: ExtendedMath.h:106
constants::unsignedMax
const unsigned unsignedMax
Definition: GeneralDefine.h:46
DPMBase::setTimeStep
void setTimeStep(Mdouble newDt)
Sets a new value for the simulation time step.
Definition: DPMBase.cc:1225
DPMBase::setMin
void setMin(const Vec3D &min)
Sets the minimum coordinates of the problem domain.
Definition: DPMBase.cc:1109
constants::pi
const Mdouble pi
Definition: ExtendedMath.h:45
Flusher::FLUSH
@ FLUSH
helpers::KAndDispAndKtAndDispt::kt
Mdouble kt
Definition: Helpers.h:153
DPMBase::setParticleDimensions
void setParticleDimensions(unsigned int particleDimensions)
Sets the particle dimensionality.
Definition: DPMBase.cc:1439
BaseInteractable::setPosition
void setPosition(const Vec3D &position)
Sets the position of this BaseInteractable.
Definition: BaseInteractable.h:239
logger
Logger< MERCURY_LOGLEVEL > logger("MercuryKernel")
Definition of different loggers with certain modules. A user can define its own custom logger here.
mathsFunc::exp
Mdouble exp(Mdouble Exponent)
Definition: ExtendedMath.cc:84
BaseInteractable::getOrientation
const Quaternion & getOrientation() const
Returns the orientation of this BaseInteractable.
Definition: BaseInteractable.h:230
INFO
LL< Log::INFO > INFO
Info log level.
Definition: Logger.cc:55
ParticleSpecies
Definition: ParticleSpecies.h:37
mathsFunc::log
Mdouble log(Mdouble Power)
Definition: ExtendedMath.cc:104
BaseInteractable::getVelocity
virtual const Vec3D & getVelocity() const
Returns the velocity of this interactable.
Definition: BaseInteractable.cc:329
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
Vec3D
Definition: Vector.h:50
BaseInteractable::setOrientation
void setOrientation(const Quaternion &orientation)
Sets the orientation of this BaseInteractable.
Definition: BaseInteractable.h:260
Mdouble
double Mdouble
Definition: GeneralDefine.h:34
checkTemplate
void checkTemplate(T real, T ideal, double error, std::string whatIsChecked)
Definition: Helpers.cc:904
BaseInteractable::setVelocity
void setVelocity(const Vec3D &velocity)
set the velocity of the BaseInteractable.
Definition: BaseInteractable.cc:350
BaseInteractable::setAngularVelocity
void setAngularVelocity(const Vec3D &angularVelocity)
set the angular velocity of the BaseInteractble.
Definition: BaseInteractable.cc:360
helpers::KAndDispAndKtAndDispt::dispt
Mdouble dispt
Definition: Helpers.h:154
ERROR
LL< Log::ERROR > ERROR
Error log level.
Definition: Logger.cc:53
mathsFunc::sin
Mdouble sin(Mdouble x)
Definition: ExtendedMath.cc:44
BaseParticle::setSpecies
void setSpecies(const ParticleSpecies *species)
Definition: BaseParticle.cc:804
DPMBase::getTime
Mdouble getTime() const
Returns the current simulation time.
Definition: DPMBase.cc:805
WARN
LL< Log::WARN > WARN
Warning log level.
Definition: Logger.cc:54
helpers::round
Mdouble round(const Mdouble value, unsigned precision)
Definition: Helpers.cc:598
GetCurrentDir
#define GetCurrentDir
Definition: Helpers.cc:51
SphericalParticle
A spherical particle is the most simple particle used in MercuryDPM.
Definition: SphericalParticle.h:37
helpers::KAndDispAndKtAndDispt::disp
Mdouble disp
Definition: Helpers.h:152
DPMBase::setTimeMax
void setTimeMax(Mdouble newTMax)
Sets a new value for the maximum simulation duration.
Definition: DPMBase.cc:870
InfiniteWall
A infinite wall fills the half-space {point: (position_-point)*normal_<=0}.
Definition: InfiniteWall.h:48
constants::i
const std::complex< Mdouble > i
Definition: ExtendedMath.h:51
BaseWall::setSpecies
void setSpecies(const ParticleSpecies *species)
Defines the species of the current wall.
Definition: BaseWall.cc:171
DPMBase::setFileType
void setFileType(FileType fileType)
Sets File::fileType_ for all files (ene, data, fstat, restart, stat)
Definition: DPMBase.cc:457
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
BaseParticle
Definition: BaseParticle.h:54
helpers::KAndDispAndKtAndDispt
Calculates the collision time for a given stiffness, dissipation, and effective mass.
Definition: Helpers.h:149
n
const unsigned n
Definition: CG3DPackingUnitTest.cpp:32
helpers::writeToFile
bool writeToFile(std::string filename, std::string filecontent)
Writes a string to a file.
Definition: Helpers.cc:445
mathsFunc::beta
Mdouble beta(Mdouble z, Mdouble w)
This is the beta function, returns the approximation based on cmath's implementation of ln(gamma)
Definition: ExtendedMath.cc:164
DPMBase::setSaveCount
void setSaveCount(unsigned int saveCount)
Sets File::saveCount_ for all files (ene, data, fstat, restart, stat)
Definition: DPMBase.cc:406
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
units::name
std::string name
Definition: MercuryProb.h:48
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
mathsFunc::cos
Mdouble cos(Mdouble x)
Definition: ExtendedMath.cc:64
BaseInteractable::getPosition
const Vec3D & getPosition() const
Returns the position of this BaseInteractable.
Definition: BaseInteractable.h:218
DPMBase::getTimeMax
Mdouble getTimeMax() const
Returns the maximum simulation duration.
Definition: DPMBase.cc:885
DPMBase
The DPMBase header includes quite a few header files, defining all the handlers, which are essential....
Definition: DPMBase.h:76
helpers::readFromFile
T readFromFile(std::string fileName, std::string varName, T value)
Definition: Helpers.h:311
DPMBase::setSystemDimensions
void setSystemDimensions(unsigned int newDim)
Sets the system dimensionality.
Definition: DPMBase.cc:1408