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: %, "
386  "timestep: %\n Arguments need to be positive",
387  numberOfSaves, timeMax, timeStep);
388  }
389 }
390 
391 //seems to be unused, consider taking out \author weinhartt
392 //unsigned int helpers::getSaveCountFromNumberOfSavesPerTimeUnitAndTimeStep(unsigned int numberOfSaves, Mdouble time step)
393 //{
394 // if (numberOfSaves > 1 && time step > 0)
395 // {
396 // return static_cast<unsigned int>(ceil(1.0 / time step / static_cast<double>(numberOfSaves - 1)));
397 // }
398 // else
399 // {
400 // std::cerr << "Error in getSaveCountFromNumberOfSavesPerTimeUnitAndTimeStep (" << numberOfSaves << "," << time step << ")" << std::endl;
401 // std::cerr << "Arguments need to be positive" << std::endl;
402 // exit(-1);
403 // }
404 //}
405 
423 void helpers::getLineFromStringStream(std::istream& 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 }
430 
445 bool helpers::writeToFile(std::string filename, std::string filecontent)
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 }
458 
459 void helpers::writeCommandLineToFile(const std::string filename, const int argc, char * const argv[])
460 {
461  std::stringstream ss;
462  for (int i=0; i<argc; ++i) {
463  ss << argv[i] << ' ';
464  }
465  writeToFile(filename,ss.str());
466 }
467 
468 
469 void helpers::gnuplot(std::string command)
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 }
483 
484 bool helpers::addToFile(std::string filename, std::string filecontent)
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 }
497 
502 bool helpers::fileExists(std::string strFilename)
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 }
530 
539 bool helpers::openFile(std::fstream& file, std::string filename, std::fstream::openmode mode)
540 {
541  file.open(filename.c_str(), mode);
542  if (file.fail())
543  return false;
544  else
545  return true;
546 }
547 
556 {
557  return mass0 * mass1 / (mass0 + mass1);
558 }
559 
560 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 }
579 
580 void helpers::more(std::string filename, unsigned nLines)
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 }
597 
598 Mdouble helpers::round(const Mdouble value, unsigned precision)
599 {
600  const Mdouble logValue = log10(value);
601  const int factor = std::pow(10, precision - std::ceil(logValue));
602  return std::round(value * factor) / factor;
603 }
604 
605 std::string helpers::to_string(const Mdouble value, unsigned precision)
606 {
607  std::ostringstream stm;
608  stm << round(value,precision);
609  return stm.str();
610 }
611 
618 void helpers::loadingTest(const ParticleSpecies* species, Mdouble displacement, Mdouble velocity, Mdouble radius,
619  std::string name)
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);
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 }
685 
693  Mdouble tangentialDisplacement, Mdouble velocity, Mdouble radius,
694  std::string name)
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);
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 }
765 
772 void helpers::objectivenessTest(const ParticleSpecies* species, Mdouble displacement, Mdouble tangentialDisplacement,
773  Mdouble velocity, Mdouble radius, std::string name)
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);
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 }
857 
858 bool helpers::compare(std::istream& is, std::string s)
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 }
873 
874 bool helpers::readFromCommandLine(int argc, char *argv[], std::string varName)
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 }
886 
887 std::vector<Mdouble> helpers::linspace(Mdouble Min, Mdouble Max, int numberOfBins)
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 }
902 
903 template<class T>
904 void checkTemplate(T real, T ideal, double error, std::string whatIsChecked)
905 {
906  logger.assert_always(mathsFunc::isEqual(real, ideal, error),
907  whatIsChecked + ": % (should be %) ", real, ideal);
908  logger(INFO, whatIsChecked + ": % (correct)", real);
909 }
910 
911 void helpers::check(double real, double ideal, double error, std::string whatIsChecked)
912 {
913  checkTemplate(real,ideal,error,whatIsChecked);
914 }
915 
916 void helpers::check(Vec3D real, Vec3D ideal, double error, std::string whatIsChecked)
917 {
918  checkTemplate(real,ideal,error,whatIsChecked);
919 }
920 
921 void helpers::check(MatrixSymmetric3D real, MatrixSymmetric3D ideal, double error, std::string whatIsChecked)
922 {
923  checkTemplate(real,ideal,error,whatIsChecked);
924 }
925 
926 void helpers::check(Matrix3D real, Matrix3D ideal, double error, std::string whatIsChecked)
927 {
928  checkTemplate(real,ideal,error,whatIsChecked);
929 }
930 
931 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 }
939 
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 }
948 
957 bool helpers::isNext(std::istream& is, const std::string name) {
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 }
968 
969 template<>
970 std::string helpers::readFromCommandLine<std::string>(int argc, char *argv[], std::string varName, std::string value)
971 {
972  for (unsigned i=0; i<argc-1; ++i) {
973  if (varName == argv[i]) {
974  value = argv[i+1];
975  logger(INFO, "readFromCommandLine: % set to % ", varName.substr(1), value);
976  return value;
977  }
978  }
979  //if the variable is not found
980  logger(INFO, "readFromCommandLine: % set to default value % ", varName.substr(1), value);
981  return value;
982 }
983 
984 bool helpers::createDirectory(std::string path) {
985  //see https://stackoverflow.com/questions/20358455/cross-platform-way-to-make-a-directory
986  mode_t nMode = 0733; // UNIX style permissions
987  int nError = 0;
988 #if defined(_WIN32)
989  nError = _mkdir(path.c_str()); // can be used on Windows
990 #else
991  nError = mkdir(path.c_str(),nMode); // can be used on non-Windows
992 #endif
993  if (nError != 0) {
994  // handle your error here
995  }
996  return false;
997 }
998 
999 Mdouble helpers::getRayleighTime(Mdouble radius, Mdouble shearModulus, Mdouble poisson, Mdouble density) {
1000  return constants::pi*radius*sqrt(density/shearModulus)/(0.1631*poisson+0.8766);
1001 }
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:445
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:75
std::string lower(std::string s)
Definition: Helpers.cc:55
A spherical particle is the most simple particle used in MercuryDPM.
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:870
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:555
Logger< MERCURY_LOGLEVEL > logger("MercuryKernel")
Definition of different loggers with certain modules. A user can define its own custom logger here...
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:931
void setMax(const Vec3D &max)
Sets the maximum coordinates of the problem domain.
Definition: DPMBase.cc:1073
LL< Log::INFO > INFO
Info log level.
Definition: Logger.cc:55
Mdouble getRayleighTime(Mdouble radius, Mdouble shearModulus, Mdouble poisson, Mdouble density)
Definition: Helpers.cc:999
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:911
const std::complex< Mdouble > i
Definition: ExtendedMath.h:51
void setParticleDimensions(unsigned int particleDimensions)
Sets the particle dimensionality.
Definition: DPMBase.cc:1439
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:692
void setSystemDimensions(unsigned int newDim)
Sets the system dimensionality.
Definition: DPMBase.cc:1408
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:423
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:1109
void setSaveCount(unsigned int saveCount)
Sets File::saveCount_ for all files (ene, data, fstat, restart, stat)
Definition: DPMBase.cc:406
Mdouble getRealTime()
Definition: Helpers.cc:940
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:618
void gnuplot(std::string command)
Plots to a gnuplot window.
Definition: Helpers.cc:469
bool createDirectory(std::string)
Definition: Helpers.cc:984
virtual void actionsBeforeTimeStep()
A virtual function which allows to define operations to be executed before the new time step...
Definition: DPMBase.cc:1855
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:957
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:598
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:874
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...
std::vector< Mdouble > linspace(Mdouble a, Mdouble b, int N)
Definition: Helpers.cc:887
void setTimeStep(Mdouble newDt)
Sets a new value for the simulation time step.
Definition: DPMBase.cc:1225
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
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:459
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:502
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:106
bool openFile(std::fstream &file, std::string filename, std::fstream::openmode mode)
Provides a simple interface for opening a file.
Definition: Helpers.cc:539
void checkTemplate(T real, T ideal, double error, std::string whatIsChecked)
Definition: Helpers.cc:904
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:772
bool addToFile(std::string filename, std::string filecontent)
Adds a string to an existing file.
Definition: Helpers.cc:484
void setAngularVelocity(const Vec3D &angularVelocity)
set the angular velocity of the BaseInteractble.
Mdouble getTime() const
Returns the current simulation time.
Definition: DPMBase.cc:805
Implementation of a 3D symmetric matrix.
void more(std::string filename, unsigned nLines=constants::unsignedMax)
Definition: Helpers.cc:580
Mdouble getTimeMax() const
Returns the maximum simulation duration.
Definition: DPMBase.cc:885
std::vector< double > readArrayFromFile(std::string filename, int &n, int &m)
Definition: Helpers.cc:560
virtual void setupInitialConditions()
This function allows to set the initial conditions for our problem to be solved, by default particle ...
Definition: DPMBase.cc:1989
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:457