MercuryDPM  Trunk
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Helpers.cc
Go to the documentation of this file.
1 //Copyright (c) 2013-2020, The MercuryDPM Developers Team. All rights reserved.
2 //For the list of developers, see <http://www.MercuryDPM.org/Team>.
3 //
4 //Redistribution and use in source and binary forms, with or without
5 //modification, are permitted provided that the following conditions are met:
6 // * Redistributions of source code must retain the above copyright
7 // notice, this list of conditions and the following disclaimer.
8 // * Redistributions in binary form must reproduce the above copyright
9 // notice, this list of conditions and the following disclaimer in the
10 // documentation and/or other materials provided with the distribution.
11 // * Neither the name MercuryDPM nor the
12 // names of its contributors may be used to endorse or promote products
13 // derived from this software without specific prior written permission.
14 //
15 //THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND
16 //ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
17 //WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
18 //DISCLAIMED. IN NO EVENT SHALL THE MERCURYDPM DEVELOPERS TEAM BE LIABLE FOR ANY
19 //DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
20 //(INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
21 //LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
22 //ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
23 //(INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
24 //SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
25 
26 #include "Helpers.h"
27 #include <cmath>
28 #include <sys/stat.h>
29 #include <cstdio>
30 #include <sstream>
31 #include <string>
32 #include <Logger.h>
33 #include <DPMBase.h>
34 #include <Walls/InfiniteWall.h>
35 #include <Particles/BaseParticle.h>
36 #include <Species/BaseSpecies.h>
37 //#include <cassert>
38 #include "Math/ExtendedMath.h"
39 #include <numeric>
40 #include <chrono>
41 #include <sys/types.h> // required for stat.h
42 
43 #ifdef WINDOWS
44 #include <direct.h>
45 #define GetCurrentDir _getcwd
46 #else
47 
48 #include <unistd.h>
50 
51 #define GetCurrentDir getcwd
52 #endif
53 
54 
55 std::string helpers::lower(std::string s) {
56  std::transform(s.begin(), s.end(), s.begin(), ::tolower);
57  return s;
58 }
59 
60 //helpers::KAndDisp helpers::computeKAndDispFromCollisionTimeAndRestitutionCoefficientAndEffectiveMass(Mdouble tc, Mdouble r, Mdouble mass)
61 //{
62 // helpers::KAndDisp ans;
63 // ans.k = computeKFromCollisionTimeAndRestitutionCoefficientAndEffectiveMass(tc, r, mass);
64 // ans.disp = computeDispFromCollisionTimeAndRestitutionCoefficientAndEffectiveMass(tc, r, mass);
65 // return ans;
66 //}
67 //
68 //Mdouble helpers::computeCollisionTimeFromKAndDispAndEffectiveMass(Mdouble k, Mdouble disp, Mdouble mass)
69 //{
70 // if (k <= 0)
71 // {
72 // std::cerr << "Error in helpers::computeCollisionTimeFromKandDispAndEffectiveMass(Mdouble k, Mdouble disp, Mdouble mass) stiffness is not set or has an unexpected value, (stiffness=" << k << ")" << std::endl;
73 // exit(-1);
74 // }
75 // if (disp < 0)
76 // {
77 // std::cerr << "Error in helpers::computeCollisionTimeFromKandDispAndEffectiveMass(Mdouble k, Mdouble disp, Mdouble mass) dissipation is not set or has an unexpected value, (dissipation=" << disp << ")" << std::endl;
78 // exit(-1);
79 // }
80 // if (mass <= 0)
81 // {
82 // std::cerr << "Error in helpers::computeCollisionTimeFromKandDispAndEffectiveMass(Mdouble k, Mdouble disp, Mdouble mass) mass is not set or has an unexpected value (mass=" << mass << ")" << std::endl;
83 // exit(-1);
84 // }
85 // if (4.0 * k / mass - mathsFunc::square(disp / mass) < 0)
86 // {
87 // std::cerr << "Error in helpers::computeCollisionTimeFromKandDispAndEffectiveMass(Mdouble k, Mdouble disp, Mdouble mass) values for stiffness, dissipation and mass lead to an overdamped system (stiffness=" << k << " dissipation=" << disp << " mass=" << mass << ")" << std::endl;
88 // exit(-1);
89 // }
90 // return 2.0 * constants::pi / std::sqrt(4.0 * k / mass - mathsFunc::square(disp / mass));
91 //}
92 //
93 //Mdouble helpers::computeRestitutionCoefficientFromKAndDispAndEffectiveMass(Mdouble k, Mdouble disp, Mdouble mass)
94 //{
95 // if (k <= 0)
96 // {
97 // std::cerr << "Error in helpers::computeRestitutionCoefficientFromKandDispAndEffectiveMass(Mdouble k, Mdouble disp, Mdouble mass) stiffness is not set or has an unexpected value, (stiffness=" << k << ")" << std::endl;
98 // exit(-1);
99 // }
100 // if (disp < 0)
101 // {
102 // std::cerr << "Error in helpers::computeRestitutionCoefficientFromKandDispAndEffectiveMass(Mdouble k, Mdouble disp, Mdouble mass) dissipation is not set or has an unexpected value, (dissipation=" << disp << ")" << std::endl;
103 // exit(-1);
104 // }
105 // if (mass <= 0)
106 // {
107 // std::cerr << "Error in helpers::computeRestitutionCoefficientFromKandDispAndEffectiveMass(Mdouble k, Mdouble disp, Mdouble mass) mass is not set or has an unexpected value (mass=" << mass << ")" << std::endl;
108 // exit(-1);
109 // }
110 // if (4.0 * mass * k - mathsFunc::square(disp) < 0)
111 // {
112 // std::cerr << "Error in helpers::computeRestitutionCoefficientFromKandDispAndEffectiveMass(Mdouble k, Mdouble disp, Mdouble mass) values for stiffness, dissipation and mass lead to an overdamped system (stiffness=" << k << " dissipation=" << disp << " mass=" << mass << ")" << std::endl;
113 // exit(-1);
114 // }
115 // return std::exp(-disp * constants::pi / std::sqrt(4.0 * mass * k - mathsFunc::square(disp)));
116 //}
117 //
118 //Mdouble helpers::computeDispFromKAndRestitutionCoefficientAndEffectiveMass(Mdouble k, Mdouble r, Mdouble mass)
119 //{
120 // if (k <= 0)
121 // {
122 // std::cerr << "Error in helpers::computeDispFromKAndRestitutionCoefficientAndEffectiveMass(Mdouble k, Mdouble r, Mdouble mass) stiffness is not set or has an unexpected value, (stiffness=" << k << ")" << std::endl;
123 // exit(-1);
124 // }
125 // if (r < 0)
126 // {
127 // std::cerr << "Error in helpers::computeDispFromKAndRestitutionCoefficientAndEffectiveMass(Mdouble k, Mdouble r, Mdouble mass) restitution coefficient is not set or has an unexpected value, (restitution coefficient=" << r << ")" << std::endl;
128 // exit(-1);
129 // }
130 // if (mass <= 0)
131 // {
132 // std::cerr << "Error in helpers::computeDispFromKAndRestitutionCoefficientAndEffectiveMass(Mdouble k, Mdouble r, Mdouble mass) mass is not set or has an unexpected value (mass=" << mass << ")" << std::endl;
133 // exit(-1);
134 // }
135 // return -2.0 * sqrt(mass * k / (constants::sqr_pi + mathsFunc::square(std::log(r)))) * std::log(r);
136 //}
137 //
138 //Mdouble helpers::computeCollisionTimeFromKAndRestitutionCoefficientAndEffectiveMass(Mdouble k, Mdouble r, Mdouble mass)
139 //{
140 // if (k <= 0)
141 // {
142 // std::cerr << "Error in helpers::computeCollisionTimeFromKAndRestitutionCoefficientAndEffectiveMass(Mdouble k, Mdouble r, Mdouble mass) stiffness is not set or has an unexpected value, (stiffness=" << k << ")" << std::endl;
143 // exit(-1);
144 // }
145 // if (r < 0)
146 // {
147 // std::cerr << "Error in helpers::computeCollisionTimeFromKAndRestitutionCoefficientAndEffectiveMass(Mdouble k, Mdouble r, Mdouble mass) restitution coefficient is not set or has an unexpected value, (restitution coefficient=" << r << ")" << std::endl;
148 // exit(-1);
149 // }
150 // if (mass <= 0)
151 // {
152 // std::cerr << "Error in helpers::computeCollisionTimeFromKAndRestitutionCoefficientAndEffectiveMass(Mdouble k, Mdouble r, Mdouble mass) mass is not set or has an unexpected value (mass=" << mass << ")" << std::endl;
153 // exit(-1);
154 // }
155 // return sqrt(mass / k * (constants::sqr_pi + mathsFunc::square(std::log(r))));
156 //}
157 //
158 //Mdouble helpers::computeDispFromKAndCollisionTimeAndEffectiveMass(Mdouble k, Mdouble tc, Mdouble mass)
159 //{
160 // if (k <= 0)
161 // {
162 // std::cerr << "Error in helpers::computeDispFromKAndCollisionTimeAndEffectiveMass(Mdouble k, Mdouble tc, Mdouble mass) stiffness is not set or has an unexpected value, (stiffness=" << k << ")" << std::endl;
163 // exit(-1);
164 // }
165 // if (tc <= 0)
166 // {
167 // std::cerr << "Error in helpers::computeDispFromKAndCollisionTimeAndEffectiveMass(Mdouble k, Mdouble tc, Mdouble mass) collision time is not set or has an unexpected value, (collision time=" << tc << ")" << std::endl;
168 // exit(-1);
169 // }
170 // if (mass <= 0)
171 // {
172 // std::cerr << "Error in helpers::computeDispFromKAndCollisionTimeAndEffectiveMass(Mdouble k, Mdouble tc, Mdouble mass) mass is not set or has an unexpected value (mass=" << mass << ")" << std::endl;
173 // exit(-1);
174 // }
175 // if (mass * k - constants::sqr_pi * mathsFunc::square(mass / tc) < 0)
176 // {
177 // std::cerr << "Error in helpers::computeDispFromKAndCollisionTimeAndEffectiveMass(Mdouble k, Mdouble disp, Mdouble mass) values for stiffness, collision time and mass lead to an overdamped system (stiffness=" << k << " collision time=" << tc << " mass=" << mass << ")" << std::endl;
178 // exit(-1);
179 // }
180 // return 2.0 * std::sqrt(mass * k - constants::sqr_pi * mathsFunc::square(mass / tc));
181 //}
182 //
183 //Mdouble helpers::computeRestitutionCoefficientFromKAndCollisionTimeAndEffectiveMass(Mdouble k, Mdouble tc, Mdouble mass)
184 //{
185 // if (k <= 0)
186 // {
187 // std::cerr << "Error in helpers::computeRestitutionCoefficientFromKAndCollisionTimeAndEffectiveMass(Mdouble k, Mdouble tc, Mdouble mass) stiffness is not set or has an unexpected value, (stiffness=" << k << ")" << std::endl;
188 // exit(-1);
189 // }
190 // if (tc <= 0)
191 // {
192 // std::cerr << "Error in helpers::computeRestitutionCoefficientFromKAndCollisionTimeAndEffectiveMass(Mdouble k, Mdouble tc, Mdouble mass) collision time is not set or has an unexpected value, (collision time=" << tc << ")" << std::endl;
193 // exit(-1);
194 // }
195 // if (mass <= 0)
196 // {
197 // std::cerr << "Error in helpers::computeRestitutionCoefficientFromKAndCollisionTimeAndEffectiveMass(Mdouble k, Mdouble tc, Mdouble mass) mass is not set or has an unexpected value (mass=" << mass << ")" << std::endl;
198 // exit(-1);
199 // }
200 // if (k / mass * mathsFunc::square(tc) - constants::sqr_pi < 0)
201 // {
202 // std::cerr << "Error in helpers::computeRestitutionCoefficientFromKAndCollisionTimeAndEffectiveMass(Mdouble k, Mdouble disp, Mdouble mass) values for stiffness, collision time and mass lead to an overdamped system (stiffness=" << k << " collision time=" << tc << " mass=" << mass << ")" << std::endl;
203 // exit(-1);
204 // }
205 // return std::exp(-std::sqrt(k / mass * mathsFunc::square(tc) - constants::sqr_pi));
206 //}
207 //
208 //Mdouble helpers::computeDispFromCollisionTimeAndRestitutionCoefficientAndEffectiveMass(Mdouble tc, Mdouble r, Mdouble mass)
209 //{
210 // if (tc <= 0)
211 // {
212 // std::cerr << "Error in helpers::computeDispFromCollisionTimeAndRestitutionCoefficientAndEffectiveMass(Mdouble tc, Mdouble r, Mdouble mass) collision time is not set or has an unexpected value, (collision time=" << tc << ")" << std::endl;
213 // exit(-1);
214 // }
215 // if (r < 0)
216 // {
217 // std::cerr << "Error in helpers::computeDispFromCollisionTimeAndRestitutionCoefficientAndEffectiveMass(Mdouble tc, Mdouble r, Mdouble mass) restitution coefficient not set or has an unexpected value, (restitution coefficient=" << r << ")" << std::endl;
218 // exit(-1);
219 // }
220 // if (mass <= 0)
221 // {
222 // std::cerr << "Error in helpers::computeDispFromCollisionTimeAndRestitutionCoefficientAndEffectiveMass(Mdouble tc, Mdouble r, Mdouble mass) mass is not set or has an unexpected value (mass=" << mass << ")" << std::endl;
223 // exit(-1);
224 // }
225 // return -2.0 * mass * std::log(r) / tc;
226 //}
227 //
228 //Mdouble helpers::computeKFromCollisionTimeAndRestitutionCoefficientAndEffectiveMass(Mdouble tc, Mdouble r, Mdouble mass)
229 //{
230 // if (tc <= 0)
231 // {
232 // std::cerr << "Error in helpers::computeKFromCollisionTimeAndRestitutionCoefficientAndEffectiveMass(Mdouble tc, Mdouble r, Mdouble mass) collision time is not set or has an unexpected value, (collision time=" << tc << ")" << std::endl;
233 // exit(-1);
234 // }
235 // if (r < 0)
236 // {
237 // std::cerr << "Error in helpers::computeKFromCollisionTimeAndRestitutionCoefficientAndEffectiveMass(Mdouble tc, Mdouble r, Mdouble mass) restitution coefficient not set or has an unexpected value, (restitution coefficient=" << r << ")" << std::endl;
238 // exit(-1);
239 // }
240 // if (mass <= 0)
241 // {
242 // std::cerr << "Error in helpers::computeKFromCollisionTimeAndRestitutionCoefficientAndEffectiveMass(Mdouble tc, Mdouble r, Mdouble mass) mass is not set or has an unexpected value (mass=" << mass << ")" << std::endl;
243 // exit(-1);
244 // }
245 // return mass * (mathsFunc::square(constants::pi / tc) + mathsFunc::square(-std::log(r) / tc));
246 //}
247 //
248 //Mdouble helpers::computeKFromCollisionTimeAndDispAndEffectiveMass(Mdouble tc, Mdouble disp, Mdouble mass)
249 //{
250 // if (tc <= 0)
251 // {
252 // std::cerr << "Error in helpers::computeKFromCollisionTimeAndDispAndEffectiveMass(Mdouble tc, Mdouble disp, Mdouble mass) collision time is not set or has an unexpected value, (collision time=" << tc << ")" << std::endl;
253 // exit(-1);
254 // }
255 // if (disp < 0)
256 // {
257 // std::cerr << "Error in helpers::computeKFromCollisionTimeAndDispAndEffectiveMass(Mdouble tc, Mdouble disp, Mdouble mass) dissipation not set or has an unexpected value, (dissipation=" << disp << ")" << std::endl;
258 // exit(-1);
259 // }
260 // if (mass <= 0)
261 // {
262 // std::cerr << "Error in helpers::computeKFromCollisionTimeAndDispAndEffectiveMass(Mdouble tc, Mdouble disp, Mdouble mass) mass is not set or has an unexpected value (mass=" << mass << ")" << std::endl;
263 // exit(-1);
264 // }
265 // return 0.25 * mathsFunc::square(disp) / mass + constants::sqr_pi * mass / mathsFunc::square(tc);
266 //}
267 //
268 //Mdouble helpers::computeRestitutionCoefficientFromCollisionTimeAndDispAndEffectiveMass(Mdouble tc, Mdouble disp, Mdouble mass)
269 //{
270 // if (tc <= 0)
271 // {
272 // std::cerr << "Error in helpers::computeRestitutionCoefficientFromCollisionTimeAndDispAndEffectiveMass(Mdouble tc, Mdouble disp, Mdouble mass) collision time is not set or has an unexpected value, (collision time=" << tc << ")" << std::endl;
273 // exit(-1);
274 // }
275 // if (disp < 0)
276 // {
277 // std::cerr << "Error in helpers::computeRestitutionCoefficientFromCollisionTimeAndDispAndEffectiveMass(Mdouble tc, Mdouble disp, Mdouble mass) dissipation not set or has an unexpected value, (dissipation=" << disp << ")" << std::endl;
278 // exit(-1);
279 // }
280 // if (mass <= 0)
281 // {
282 // std::cerr << "Error in helpers::computeRestitutionCoefficientFromCollisionTimeAndDispAndEffectiveMass(Mdouble tc, Mdouble disp, Mdouble mass) mass is not set or has an unexpected value (mass=" << mass << ")" << std::endl;
283 // exit(-1);
284 // }
285 // return std::exp(-0.5 * disp * tc / mass);
286 //}
287 //
288 //Mdouble helpers::computeKFromDispAndRestitutionCoefficientAndEffectiveMass(Mdouble disp, Mdouble r, Mdouble mass)
289 //{
290 // if (disp <= 0)
291 // {
292 // std::cerr << "Error in helpers::computeKFromDispAndRestitutionCoefficientAndEffectiveMass(Mdouble disp, Mdouble r, Mdouble mass) dissipation is not set or has an unexpected value, (dissipation time=" << disp << ")" << std::endl;
293 // exit(-1);
294 // }
295 // if (r < 0)
296 // {
297 // std::cerr << "Error in helpers::computeKFromDispAndRestitutionCoefficientAndEffectiveMass(Mdouble disp, Mdouble r, Mdouble mass) restitution coefficient not set or has an unexpected value, (restitution coefficient=" << r << ")" << std::endl;
298 // exit(-1);
299 // }
300 // if (mass <= 0)
301 // {
302 // std::cerr << "Error in helpers::computeKFromDispAndRestitutionCoefficientAndEffectiveMass(Mdouble disp, Mdouble r, Mdouble mass) mass is not set or has an unexpected value (mass=" << mass << ")" << std::endl;
303 // exit(-1);
304 // }
305 // return 0.25 * mathsFunc::square(disp)*(constants::sqr_pi / (mathsFunc::square(std::log(r))) + 1.0) / mass;
306 //}
307 //
308 //Mdouble helpers::computeCollisionTimeFromDispAndRestitutionCoefficientAndEffectiveMass(Mdouble disp, Mdouble r, Mdouble mass)
309 //{
310 // if (disp <= 0)
311 // {
312 // std::cerr << "Error in helpers::computeCollisionTimeFromDispAndRestitutionCoefficientAndEffectiveMass(Mdouble disp, Mdouble r, Mdouble mass) dissipation is not set or has an unexpected value, (dissipation time=" << disp << ")" << std::endl;
313 // exit(-1);
314 // }
315 // if (r < 0)
316 // {
317 // std::cerr << "Error in helpers::computeCollisionTimeFromDispAndRestitutionCoefficientAndEffectiveMass(Mdouble disp, Mdouble r, Mdouble mass) restitution coefficient not set or has an unexpected value, (restitution coefficient=" << r << ")" << std::endl;
318 // exit(-1);
319 // }
320 // if (mass <= 0)
321 // {
322 // std::cerr << "Error in helpers::computeCollisionTimeFromDispAndRestitutionCoefficientAndEffectiveMass(Mdouble disp, Mdouble r, Mdouble mass) mass is not set or has an unexpected value (mass=" << mass << ")" << std::endl;
323 // exit(-1);
324 // }
325 // return -2.0 * mass * std::log(r) / disp;
326 //}
328 
331  Mdouble tc, Mdouble r, Mdouble beta, Mdouble mass)
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 }
345 
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 }
359 
374 unsigned int helpers::getSaveCountFromNumberOfSavesAndTimeMaxAndTimeStep(unsigned int numberOfSaves, Mdouble timeMax,
375  Mdouble timeStep)
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: %, timestep: %",
386  numberOfSaves, timeMax, timeStep);
387  logger(ERROR, " Arguments need to be positive");
388  exit(-1);
389  }
390 }
391 
392 //seems to be unused, consider taking out \author weinhartt
393 //unsigned int helpers::getSaveCountFromNumberOfSavesPerTimeUnitAndTimeStep(unsigned int numberOfSaves, Mdouble time step)
394 //{
395 // if (numberOfSaves > 1 && time step > 0)
396 // {
397 // return static_cast<unsigned int>(ceil(1.0 / time step / static_cast<double>(numberOfSaves - 1)));
398 // }
399 // else
400 // {
401 // std::cerr << "Error in getSaveCountFromNumberOfSavesPerTimeUnitAndTimeStep (" << numberOfSaves << "," << time step << ")" << std::endl;
402 // std::cerr << "Arguments need to be positive" << std::endl;
403 // exit(-1);
404 // }
405 //}
406 
424 void helpers::getLineFromStringStream(std::istream& in, std::stringstream& out)
425 {
426  std::string line_string;
427  getline(in, line_string);
428  out.str(std::move(line_string));
429  out.clear();
430 }
431 
446 bool helpers::writeToFile(std::string filename, std::string filecontent)
447 {
448  std::fstream file;
449  file.open(filename.c_str(), std::ios::out);
450  if (file.fail())
451  {
452  logger(WARN, "Error in writeToFile: file could not be opened");
453  return false;
454  }
455  file << filecontent;
456  file.close();
457  return true;
458 }
459 
460 void helpers::writeCommandLineToFile(const std::string filename, const int argc, char * const argv[])
461 {
462  std::stringstream ss;
463  for (int i=0; i<argc; ++i) {
464  ss << argv[i] << ' ';
465  }
466  writeToFile(filename,ss.str());
467 }
468 
469 
470 void helpers::gnuplot(std::string command)
471 {
472 #ifdef __CYGWIN__
473  logger(WARN, "[helpers::gnuplot] is not supported on Cygwin");
474 #else
475 #ifdef WINDOWS
476  logger(WARN, "[helpers::gnuplot] is not supported on Windows");
477 #else
478  FILE* pipe = popen("gnuplot -persist", "w");
479  fprintf(pipe, "%s", command.c_str());
480  fflush(pipe);
481 #endif
482 #endif
483 }
484 
485 bool helpers::addToFile(std::string filename, std::string filecontent)
486 {
487  std::fstream file;
488  file.open(filename.c_str(), std::ios::app);
489  if (file.fail())
490  {
491  std::cerr << "Error in writeToFile: file could not be opened" << std::endl;
492  return false;
493  }
494  file << filecontent;
495  file.close();
496  return true;
497 }
498 
503 bool helpers::fileExists(std::string strFilename)
504 {
505  struct stat stFileInfo;
506  bool blnReturn;
507  int intStat;
508 
509  // Attempt to get the file attributes
510 
511  intStat = stat(strFilename.c_str(), &stFileInfo);
512  if (intStat == 0)
513  {
514  // We were able to get the file attributes
515  // so the file obviously exists.
516  blnReturn = true;
517  }
518  else
519  {
520  // We were not able to get the file attributes.
521  // This may mean that we don't have permission to
522  // access the folder which contains this file. If you
523  // need to do that level of checking, lookup the
524  // return values of stat which will give you
525  // more details on why stat failed.
526  blnReturn = false;
527  }
528 
529  return blnReturn;
530 }
531 
540 bool helpers::openFile(std::fstream& file, std::string filename, std::fstream::openmode mode)
541 {
542  file.open(filename.c_str(), mode);
543  if (file.fail())
544  return false;
545  else
546  return true;
547 }
548 
557 {
558  return mass0 * mass1 / (mass0 + mass1);
559 }
560 
561 std::vector<double> helpers::readArrayFromFile(std::string filename, int& n, int& m)
562 {
563  std::fstream file;
564  file.open(filename.c_str(), std::ios::in);
565  if (file.fail())
566  {
567  std::cerr << "Error in readArrayFromFile: file could not be opened" << std::endl;
568  exit(-1);
569  }
570  file >> n >> m;
571  std::vector<double> v;
572  Mdouble val;
573  for (int i = 0; i < n * m; i++)
574  {
575  file >> val;
576  v.push_back(val);
577  }
578  file.close();
579  return v;
580 }
581 
582 void helpers::more(std::string filename, unsigned nLines)
583 {
584  if (nLines != constants::unsignedMax)
585  std::cout << "First " << nLines << " lines of " << filename << ":\n";
586  std::fstream file;
587  file.open(filename.c_str(), std::ios::in);
588  if (file.fail())
589  logger(ERROR, "Error in readArrayFromFile: file could not be opened");
590  std::string line;
591  for (unsigned i = 0; i < nLines; i++)
592  {
593  if (file.eof()) break;
594  getline(file, line);
595  std::cout << " " << line << '\n';
596  }
597  file.close();
598 }
599 
600 Mdouble helpers::round(const Mdouble value, unsigned precision)
601 {
602  const Mdouble logValue = log10(value);
603  const int factor = std::pow(10, precision - std::ceil(logValue));
604  return std::round(value * factor) / factor;
605 }
606 
607 std::string helpers::to_string(const Mdouble value, unsigned precision)
608 {
609  std::ostringstream stm;
610  stm << round(value,precision);
611  return stm.str();
612 }
613 
620 void helpers::loadingTest(const ParticleSpecies* species, Mdouble displacement, Mdouble velocity, Mdouble radius,
621  std::string name)
622 {
623  class LoadingTest : public DPMBase
624  {
625  const ParticleSpecies* species;
626  Mdouble displacement;
627  Mdouble velocity;
628  Mdouble radius;
629  public:
630  //public variables
631  LoadingTest(const ParticleSpecies* species, Mdouble displacement, Mdouble velocity, Mdouble radius)
632  : species(species), displacement(displacement), velocity(velocity), radius(radius)
633  {}
634 
635  void setupInitialConditions() override
636  {
637  //setName("LoadingTest"+species->getName());
638  setTimeMax(2.0 * displacement / velocity);
639  setTimeStep(2e-3 * getTimeMax());
640  setSaveCount(1);
642  fStatFile.setFileType(FileType::ONE_FILE);
643 
644  setMax({radius, radius, radius + radius});
645  setMin({-radius, -radius, 0});
648 
649  speciesHandler.copyAndAddObject(*species);
650 
652  p.setSpecies(speciesHandler.getObject(0));
653  p.setRadius(radius);
654  p.setPosition({0, 0, radius});
655  particleHandler.copyAndAddObject(p);
656 
657  InfiniteWall w;
658  w.setSpecies(speciesHandler.getObject(0));
659  w.set(Vec3D(0, 0, -1), Vec3D(0.0, 0.0, 0.0));
660  wallHandler.copyAndAddObject(w);
661  }
662 
663  void actionsBeforeTimeStep() override
664  {
665  BaseParticle* p = particleHandler.getLastObject();
666  logger.assert(p,"Empty particle handler");
667  p->setAngularVelocity({0, 0, 0});
668 
669  //Moving particle normally into surface
670  if (getTime() <= displacement / velocity)
671  {
672  p->setVelocity({0, 0, velocity});
673  p->setPosition({0, 0, radius - velocity * getTime()});
674  }
675  else
676  {
677  p->setVelocity({0, 0, -velocity});
678  p->setPosition({0, 0, radius - displacement - displacement + velocity * getTime()});
679  }
680  }
681  } test(species, displacement, velocity, radius);
682  test.setName(name);
683  test.solve();
684  writeToFile(test.getName() + ".gnu", "plot '" + test.getName() + ".fstat' u 7:9 w lp");
685  logger(INFO, "finished loading test: run 'gnuplot %.gnu' to view output", test.getName());
686 }
687 
695  Mdouble tangentialDisplacement, Mdouble velocity, Mdouble radius,
696  std::string name)
697 {
698  class LoadingTest : public DPMBase
699  {
700  const ParticleSpecies* species;
701  Mdouble displacement;
702  Mdouble tangentialDisplacement;
703  Mdouble velocity;
704  Mdouble radius;
705  public:
706  //public variables
707  LoadingTest(const ParticleSpecies* species, Mdouble displacement, Mdouble tangentialDisplacement,
708  Mdouble velocity, Mdouble radius)
709  : species(species), displacement(displacement), tangentialDisplacement(tangentialDisplacement),
710  velocity(velocity), radius(radius)
711  {}
712 
713  void setupInitialConditions() override
714  {
715  //setName("TangentialLoadingTest"+species->getName());
716  setTimeMax(4.0 * tangentialDisplacement / velocity);
717  setTimeStep(4e-4 * getTimeMax());
718  setSaveCount(1);
720  fStatFile.setFileType(FileType::ONE_FILE);
721 
722  setMax({radius, radius, radius + radius});
723  setMin({-radius, -radius, 0});
726 
727  speciesHandler.copyAndAddObject(*species);
728 
730  p.setSpecies(speciesHandler.getObject(0));
731  p.setRadius(radius);
732  p.setPosition({0, 0, radius - displacement});
733  particleHandler.copyAndAddObject(p);
734 
735  InfiniteWall w;
736  w.setSpecies(speciesHandler.getObject(0));
737  w.set(Vec3D(0, 0, -1), Vec3D(0.0, 0.0, 0.0));
738  wallHandler.copyAndAddObject(w);
739  }
740 
741  void actionsBeforeTimeStep() override
742  {
743  BaseParticle* p = particleHandler.getLastObject();
744  logger.assert(p,"Empty particle handler");
745  p->setAngularVelocity({0, 0, 0});
746 
747  //Moving particle cyclically right and left between +-tangentialDisplacement
748  bool moveRight = static_cast<int>(getTime() / (2.0*tangentialDisplacement / velocity) +0.5)%2==0;
749  if (moveRight)
750  {
751  p->setVelocity({-velocity, 0, 0});
752  p->setPosition({tangentialDisplacement - velocity * getTime(), 0, radius - displacement});
753  }
754  else
755  {
756  p->setVelocity({velocity, 0, 0});
757  p->setPosition({-2*tangentialDisplacement + velocity * getTime(), 0, radius - displacement});
758  }
759  }
760 
761  } test(species, displacement, tangentialDisplacement, velocity, radius);
762  test.setName(name);
763  test.solve();
764  writeToFile(test.getName() + ".gnu", "plot '" + test.getName() + ".fstat' u 8:($10*$14) w lp");
765  logger(INFO, "finished tangential loading test: run 'gnuplot %.gnu' to view output", test.getName());
766 }
767 
774 void helpers::objectivenessTest(const ParticleSpecies* species, Mdouble displacement, Mdouble tangentialDisplacement,
775  Mdouble velocity, Mdouble radius, std::string name)
776 {
777  class ObjectivenessTest : public DPMBase
778  {
779  const ParticleSpecies* species;
780  Mdouble displacement;
781  Mdouble tangentialDisplacement;
782  Mdouble velocity;
783  Mdouble radius;
784  public:
785  //public variables
786  ObjectivenessTest(const ParticleSpecies* species, Mdouble displacement, Mdouble tangentialDisplacement,
787  Mdouble velocity, Mdouble radius)
788  : species(species), displacement(displacement), tangentialDisplacement(tangentialDisplacement),
789  velocity(velocity), radius(radius)
790  {}
791 
792  void setupInitialConditions() override
793  {
794  //setName("ObjectivenessTest"+species->getName());
795  setTimeMax((tangentialDisplacement + 0.5 * constants::pi * radius) / velocity);
796  setTimeStep(1e-4 * getTimeMax());
797  setSaveCount(20);
799  dataFile.setFileType(FileType::ONE_FILE);
800  fStatFile.setFileType(FileType::ONE_FILE);
801 
802  setMax(radius * Vec3D(2, 2, 2));
803  setMin(radius * Vec3D(-2, -2, -2));
806 
807  speciesHandler.copyAndAddObject(*species);
808 
810  p.setSpecies(speciesHandler.getObject(0));
811  p.setRadius(radius);
812  p.setPosition({0, radius - displacement, 0});
813  particleHandler.copyAndAddObject(p);
814  p.setPosition({0, -radius + displacement, 0});
815  particleHandler.copyAndAddObject(p);
816  }
817 
818  void actionsBeforeTimeStep() override
819  {
820  BaseParticle* p = particleHandler.getObject(0);
821  BaseParticle* q = particleHandler.getLastObject();
822  logger.assert(p,"Empty particle handler");
823  logger.assert(q,"Empty particle handler");
824 
825  //Moving particle normally into surface
826  if (getTime() <= tangentialDisplacement / velocity)
827  {
828  p->setAngularVelocity({0, 0, 0});
829  p->setVelocity({velocity, 0, 0});
830  p->setPosition({-tangentialDisplacement + velocity * getTime(), radius - displacement, 0});
831  q->setAngularVelocity({0, 0, 0});
832  q->setVelocity({-velocity, 0, 0});
833  q->setPosition({tangentialDisplacement - velocity * getTime(), -radius + displacement, 0});
834  }
835  else
836  {
837  Mdouble angle = velocity / (radius - displacement) * (getTime() - tangentialDisplacement / velocity);
838  Mdouble s = sin(angle);
839  Mdouble c = cos(angle);
840  p->setAngularVelocity(velocity / (radius - displacement) * Vec3D(0, 0, -1));
841  //p->setAngularVelocity(Vec3D(0,0,0));
842  p->setOrientation({1, 0, 0, 0});
843  p->setVelocity(velocity * Vec3D(c, -s, 0));
844  //p->setVelocity(Vec3D(0,0,0));
845  p->setPosition((radius - displacement) * Vec3D(s, c, 0));
847  q->setOrientation(-p->getOrientation());
848  q->setVelocity(-p->getVelocity());
849  q->setPosition(-p->getPosition());
850  }
851  }
852 
853  } test(species, displacement, tangentialDisplacement, velocity, radius);
854  test.setName(name);
855  test.solve();
856  writeToFile(test.getName() + ".gnu", "set size ratio -1; plot '" + test.getName() + ".fstat' u 14:15 every 2 w lp");
857  logger(INFO, "finished objectiveness test: run 'gnuplot %.gnu' to view output", test.getName());
858 }
859 
860 bool helpers::compare(std::istream& is, std::string s)
861 {
862  // Get current position
863  //check if the next line starts with 'interactionFile'; otherwise, skip interaction
864  int len = is.tellg();
865  std::string dummy;
866  is >> dummy;
867  if (dummy != s)
868  {
869  is.seekg(len, std::ios_base::beg);
870  return false;
871  logger(INFO, "helpers::compare: Next stream value (%) is not %", dummy, s);
872  }
873  return true;
874 }
875 
876 bool helpers::readFromCommandLine(int argc, char *argv[], std::string varName)
877 {
878  for (unsigned i=0; i<argc; ++i) {
879  if (varName == argv[i]) {
880  logger(INFO, "readFromCommandLine: % set to true",varName.substr(1));
881  return true;
882  }
883  }
884  //if the variable is not found
885  logger(INFO, "readFromCommandLine: % set to default value false",varName.substr(1));
886  return false;
887 }
888 
889 
890 template<class T>
891 void checkTemplate(T real, T ideal, double error, std::string whatIsChecked)
892 {
893  logger.assert_always(mathsFunc::isEqual(real, ideal, error),
894  whatIsChecked + ": % (should be %) ", real, ideal);
895  logger(INFO, whatIsChecked + ": % (correct)", real);
896 }
897 
898 void helpers::check(double real, double ideal, double error, std::string whatIsChecked)
899 {
900  checkTemplate(real,ideal,error,whatIsChecked);
901 }
902 
903 void helpers::check(Vec3D real, Vec3D ideal, double error, std::string whatIsChecked)
904 {
905  checkTemplate(real,ideal,error,whatIsChecked);
906 }
907 
908 void helpers::check(MatrixSymmetric3D real, MatrixSymmetric3D ideal, double error, std::string whatIsChecked)
909 {
910  checkTemplate(real,ideal,error,whatIsChecked);
911 }
912 
913 void helpers::check(Matrix3D real, Matrix3D ideal, double error, std::string whatIsChecked)
914 {
915  checkTemplate(real,ideal,error,whatIsChecked);
916 }
917 
918 std::string helpers::getPath()
919 {
920  char cCurrentPath[FILENAME_MAX];
921  if (not GetCurrentDir(cCurrentPath, sizeof(cCurrentPath))) {
922  logger(WARN, "Get current dir failed: %", cCurrentPath);
923  }
924  return std::string(cCurrentPath);
925 }
926 
928 {
929  // record start time
930  static auto start = std::chrono::steady_clock::now();
931  auto end = std::chrono::steady_clock::now();
932  std::chrono::duration<double> diff = end - start;
933  return diff.count();
934 }
935 
944 bool helpers::isNext(std::istream& is, const std::string name) {
945  std::string dummy;
946  auto pos = is.tellg();
947  is >> dummy;
948  if (dummy != name) {
949  is.seekg(pos);
950  return false;
951  } else {
952  return true;
953  }
954 }
955 
956 template<>
957 std::string helpers::readFromCommandLine<std::string>(int argc, char *argv[], std::string varName, std::string value)
958 {
959  for (unsigned i=0; i<argc-1; ++i) {
960  if (varName == argv[i]) {
961  value = argv[i+1];
962  logger(INFO, "readFromCommandLine: % set to % ", varName.substr(1), value);
963  return value;
964  }
965  }
966  //if the variable is not found
967  logger(INFO, "readFromCommandLine: % set to default value % ", varName.substr(1), value);
968  return value;
969 }
970 
971 bool helpers::createDirectory(std::string path) {
972  //see https://stackoverflow.com/questions/20358455/cross-platform-way-to-make-a-directory
973  mode_t nMode = 0733; // UNIX style permissions
974  int nError = 0;
975 #if defined(_WIN32)
976  nError = _mkdir(path.c_str()); // can be used on Windows
977 #else
978  nError = mkdir(path.c_str(),nMode); // can be used on non-Windows
979 #endif
980  if (nError != 0) {
981  // handle your error here
982  }
983  return false;
984 }
985 
986 Mdouble helpers::getRayleighTime(Mdouble radius, Mdouble shearModulus, Mdouble poisson, Mdouble density) {
987  return constants::pi*radius*sqrt(density/shearModulus)/(0.1631*poisson+0.8766);
988 }
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 ...
Definition: Helpers.cc:330
bool writeToFile(std::string filename, std::string filecontent)
Writes a string to a file.
Definition: Helpers.cc:446
const Vec3D & getPosition() const
Returns the position of this BaseInteractable.
The DPMBase header includes quite a few header files, defining all the handlers, which are essential...
Definition: DPMBase.h:72
std::string lower(std::string s)
Definition: Helpers.cc:55
A basic particle.
void setVelocity(const Vec3D &velocity)
set the velocity of the BaseInteractable.
void setTimeMax(Mdouble newTMax)
Sets a new value for the maximum simulation duration.
Definition: DPMBase.cc:840
Mdouble getEffectiveMass(Mdouble mass0, Mdouble mass1)
Calculates the effective mass of a particle pair, i.e. half the harmonic mean of two particle masses...
Definition: Helpers.cc:556
Logger< MERCURY_LOGLEVEL > logger("MercuryKernel")
double Mdouble
Definition: GeneralDefine.h:34
virtual void setRadius(Mdouble radius)
Sets the particle's radius_ (and adjusts the mass_ accordingly, based on the particle's species) ...
std::string getPath()
Definition: Helpers.cc:918
void setMax(const Vec3D &max)
Sets the maximum coordinates of the problem domain.
Definition: DPMBase.cc:1043
LL< Log::INFO > INFO
Info log level.
Definition: Logger.cc:55
Mdouble getRayleighTime(Mdouble radius, Mdouble shearModulus, Mdouble poisson, Mdouble density)
Definition: Helpers.cc:986
Mdouble exp(Mdouble Exponent)
Definition: ExtendedMath.cc:84
Mdouble beta(Mdouble z, Mdouble w)
This is the beta function, returns the approximation based on cmath's implementation of ln(gamma) ...
void check(double real, double ideal, double error, std::string errorMessage)
Definition: Helpers.cc:898
const std::complex< Mdouble > i
Definition: ExtendedMath.h:50
void setParticleDimensions(unsigned int particleDimensions)
Sets the particle dimensionality.
Definition: DPMBase.cc:1408
virtual const Vec3D & getAngularVelocity() const
Returns the angular velocity of this interactable.
void setSpecies(const ParticleSpecies *species)
Set disp and k such that is matches a given collision time tc and restitution coefficient r for a col...
Definition: Helpers.h:148
void 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...
Definition: Helpers.cc:694
void setSystemDimensions(unsigned int newDim)
Sets the system dimensionality.
Definition: DPMBase.cc:1377
Mdouble log(Mdouble Power)
LL< Log::ERROR > ERROR
Error log level.
Definition: Logger.cc:53
bool isEqual(Mdouble v1, Mdouble v2, Mdouble absError)
Compares the difference of two Mdouble with an absolute error, useful in UnitTests.
Mdouble cos(Mdouble x)
Definition: ExtendedMath.cc:64
const unsigned unsignedMax
Definition: GeneralDefine.h:46
file will not be created/read
Mdouble sin(Mdouble x)
Definition: ExtendedMath.cc:44
LL< Log::WARN > WARN
Warning log level.
Definition: Logger.cc:54
void getLineFromStringStream(std::istream &in, std::stringstream &out)
Reads a line from one stringstream into another, and prepares the latter for reading in...
Definition: Helpers.cc:424
const Mdouble pi
Definition: ExtendedMath.h:45
all data will be written into/ read from a single file called name_
void setMin(const Vec3D &min)
Sets the minimum coordinates of the problem domain.
Definition: DPMBase.cc:1079
void setSaveCount(unsigned int saveCount)
Sets File::saveCount_ for all files (ene, data, fstat, restart, stat)
Definition: DPMBase.cc:398
Mdouble getRealTime()
Definition: Helpers.cc:927
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...
Definition: Helpers.cc:374
void loadingTest(const ParticleSpecies *species, Mdouble displacement, Mdouble velocity, Mdouble radius, std::string name)
Used to test the loading/unloading properties of a contact law.
Definition: Helpers.cc:620
void gnuplot(std::string command)
Plots to a gnuplot window.
Definition: Helpers.cc:470
bool createDirectory(std::string)
Definition: Helpers.cc:971
virtual void actionsBeforeTimeStep()
A virtual function which allows to define operations to be executed before the new time step...
Definition: DPMBase.cc:1824
bool isNext(std::istream &is, const std::string name)
reads next value in stream as a string and compares it with name.
Definition: Helpers.cc:944
std::string to_string(const T &n)
Definition: Helpers.h:227
void setPosition(const Vec3D &position)
Sets the position of this BaseInteractable.
Mdouble round(const Mdouble value, unsigned precision)
Definition: Helpers.cc:600
bool readFromCommandLine(int argc, char *argv[], std::string varName)
Returns true if command line arguments contain varName, false else Usage example: if (readFromCommand...
Definition: Helpers.cc:876
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...
void setTimeStep(Mdouble newDt)
Sets a new value for the simulation time step.
Definition: DPMBase.cc:1195
bool compare(std::istream &is, std::string s)
Checks if the next argument in the input stream is a certain string.
Definition: Helpers.cc:860
This is a class defining walls.
Definition: InfiniteWall.h:47
#define GetCurrentDir
Definition: Helpers.cc:51
void writeCommandLineToFile(const std::string filename, const int argc, char *const argv[])
Writes a string to a file.
Definition: Helpers.cc:460
void setOrientation(const Quaternion &orientation)
Sets the orientation of this BaseInteractable.
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:503
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 ...
Definition: Helpers.cc:346
Implementation of a 3D matrix.
Definition: Matrix.h:37
virtual const Vec3D & getVelocity() const
Returns the velocity of this interactable.
Definition: Vector.h:49
T square(const T val)
squares a number
Definition: ExtendedMath.h:104
bool openFile(std::fstream &file, std::string filename, std::fstream::openmode mode)
Provides a simple interface for opening a file.
Definition: Helpers.cc:540
void checkTemplate(T real, T ideal, double error, std::string whatIsChecked)
Definition: Helpers.cc:891
const Quaternion & getOrientation() const
Returns the orientation of this BaseInteractable.
void 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...
Definition: Helpers.cc:774
bool addToFile(std::string filename, std::string filecontent)
Adds a string to an existing file.
Definition: Helpers.cc:485
void setAngularVelocity(const Vec3D &angularVelocity)
set the angular velocity of the BaseInteractble.
Mdouble getTime() const
Returns the current simulation time.
Definition: DPMBase.cc:797
Implementation of a 3D symmetric matrix.
void more(std::string filename, unsigned nLines=constants::unsignedMax)
Definition: Helpers.cc:582
Mdouble getTimeMax() const
Returns the maximum simulation duration.
Definition: DPMBase.cc:855
std::vector< double > readArrayFromFile(std::string filename, int &n, int &m)
Definition: Helpers.cc:561
virtual void setupInitialConditions()
This function allows to set the initial conditions for our problem to be solved, by default particle ...
Definition: DPMBase.cc:1961
void setSpecies(const ParticleSpecies *species)
Defines the species of the current wall.
Definition: BaseWall.cc:171
void setFileType(FileType fileType)
Sets File::fileType_ for all files (ene, data, fstat, restart, stat)
Definition: DPMBase.cc:449