MercuryDPM  Trunk
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
SinterNormalSpecies.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 
27 #include "SinterNormalSpecies.h"
29 #include "Logger.h"
30 #include "Species/BaseSpecies.h"
31 
32 class BaseParticle;
33 
34 class BaseInteractable;
35 
37  : BaseNormalForce()
38 {
40  loadingStiffness_ = 0.0;
42  cohesionStiffness_ = 0.0;
44  dissipation_ = 0.0;
45  sinterAdhesion_ = 0.0;
46  sinterRate_ = 0.0;
48  complianceZero_ = 0.0;
49  surfTension_ = 0.0;
50  constantC1_ = 0.0;
51  separationDis_ = 0.0;
52 #ifdef DEBUG_CONSTRUCTOR
53  std::cout<<"SinterNormalSpecies::SinterNormalSpecies() finished"<<std::endl;
54 #endif
55 }
56 
61  : BaseNormalForce(p)
62 {
76 #ifdef DEBUG_CONSTRUCTOR
77  std::cout<<"SinterNormalSpecies::SinterNormalSpecies(const SinterNormalSpecies &p) finished"<<std::endl;
78 #endif
79 }
80 
82 {
83 #ifdef DEBUG_DESTRUCTOR
84  std::cout<<"SinterNormalSpecies::~SinterNormalSpecies() finished"<<std::endl;
85 #endif
86 }
87 
91 void SinterNormalSpecies::write(std::ostream& os) const
92 {
93  os << " loadingStiffness " << loadingStiffness_;
94  os << " maxUnloadingStiffness " << unloadingStiffnessMax_;
95  os << " cohesionStiffness " << cohesionStiffness_;
96  os << " maxPenetration " << penetrationDepthMax_;
97  os << " dissipation " << dissipation_;
98  os << " sinterAdhesion " << sinterAdhesion_;
99  os << " inverseSinterViscosity " << inverseSinterViscosity_;
100  os << " sinterRate " << sinterRate_;
101  os << " complianceZero " << complianceZero_;
102  os << " surfTension " << surfTension_;
103  os << " constantC1 " << constantC1_;
104  os << " separationDis " << separationDis_;
105  os << " sinterType " << (unsigned) sinterType_;
106 }
107 
111 void SinterNormalSpecies::read(std::istream& is)
112 {
113  std::string dummy;
114  is >> dummy >> loadingStiffness_;
115  is >> dummy >> unloadingStiffnessMax_;
116  is >> dummy >> cohesionStiffness_;
117  is >> dummy >> penetrationDepthMax_;
118  is >> dummy >> dissipation_;
119  is >> dummy >> sinterAdhesion_;
120  is >> dummy >> inverseSinterViscosity_;
121  is >> dummy >> sinterRate_;
122  is >> dummy >> complianceZero_;
123  is >> dummy >> surfTension_;
124  is >> dummy >> constantC1_;
125  is >> dummy >> separationDis_;
126  unsigned type;
127  is >> dummy >> type;
128  sinterType_ = (SINTERTYPE) type;
129 }
130 
135 {
136  return "Sinter";
137 }
138 
145 {
157 
158 }
159 
166 void SinterNormalSpecies::setPlasticParameters(Mdouble loadingStiffness, Mdouble unloadingStiffnessMax,
167  Mdouble cohesionStiffness, Mdouble penetrationDepthMax)
168 {
169  if (loadingStiffness <= 0 || unloadingStiffnessMax < loadingStiffness || cohesionStiffness < 0 ||
170  penetrationDepthMax < 0)
171  {
172  logger(ERROR, "arguments of setPlasticParameters do not make sense");
173  }
174  setLoadingStiffness(loadingStiffness);
175  setUnloadingStiffnessMax(unloadingStiffnessMax);
176  setCohesionStiffness(cohesionStiffness);
177  setPenetrationDepthMax(penetrationDepthMax);
178 }
179 
184 {
185  return loadingStiffness_;
186 }
187 
192 {
193  return unloadingStiffnessMax_;
194 }
195 
200 {
201  return cohesionStiffness_;
202 }
203 
208 {
209  return penetrationDepthMax_;
210 }
211 
216 {
217  loadingStiffness_ = loadingStiffness;
218 }
219 
224 {
225  unloadingStiffnessMax_ = unloadingStiffnessMax;
226 }
227 
232 {
233  cohesionStiffness_ = cohesionStiffness;
234 }
235 
240 {
241  penetrationDepthMax_ = penetrationDepthMax;
242 }
243 
249 {
250  if (unloadingStiffnessMax_ / (.5 * mass) < mathsFunc::square(dissipation_ / mass))
251  {
252  logger(ERROR, "Dissipation too high; max. allowed %", sqrt(2.0 * unloadingStiffnessMax_ * mass));
253  }
254  return 0.02 * constants::pi /
255  std::sqrt(unloadingStiffnessMax_ / (.5 * mass) - mathsFunc::square(dissipation_ / mass));
256 }
257 
263 {
264  if (dissipation < 0)
265  {
266  logger(ERROR, "setDissipation(%) argument should be non-negative!", dissipation);
267  }
268  dissipation_ = dissipation;
269 }
270 
276 {
277  if (sinterAdhesion < 0)
278  {
279  logger(ERROR, "setSinterAdhesion(%) argument should be non-negative!", sinterAdhesion);
280  }
281  sinterAdhesion_ = sinterAdhesion;
282 }
283 
289 {
290  if (sinterRate < 0)
291  {
292  logger(ERROR, "setSinterRate(%) argument should be non-negative!", sinterRate);
293  }
294  sinterRate_ = sinterRate;
295 }
296 
298 {
299  if (complianceZero < 0)
300  {
301  logger(ERROR, "complianceZero(%) argument should be non-negative!", complianceZero);
302  }
303  complianceZero_ = complianceZero;
304 }
305 
307 {
308  if (surfTension < 0)
309  {
310  logger(ERROR, "SurfTension(%) argument should be non-negative!", surfTension);
311  }
312  surfTension_ = surfTension;
313 }
314 
316 {
317  if (complianceOne < 0)
318  {
319  logger(ERROR, "ComplianceOne(%) argument should be non-negative!", complianceOne);
320  }
321  constantC1_ = complianceOne;
322 }
323 
325 {
326  if (separationDis < 0)
327  {
328  logger(ERROR, "SeparationDis(%) argument should be non-negative!", separationDis);
329  }
330  separationDis_ = separationDis;
331 }
332 
333 
334 
335 
341 {
342  sinterType_ = sinterType;
344  logger(INFO, "Sintertype set to CONSTANT_RATE");
346  logger(INFO, "Sintertype set to PARHAMI_MCKEEPING");
348  logger(INFO, "Sintertype set to TEMPERATURE_DEPENDENT_FRENKEL");
350  logger(INFO, "Sintertype set to REGIME_SINTERING");
351  else
352  logger(ERROR, "Sintertype not understood");
353 }
354 
360 {
361  //assertOrDie(inverseSinterViscosity < 0, "inverseSinterViscosity should be non-negative!");
362  if (inverseSinterViscosity < 0)
363  {
364  logger(ERROR, "setInverseSinterViscosity(%) argument should be non-negative!", inverseSinterViscosity);
365  }
367  inverseSinterViscosity_ = inverseSinterViscosity;
368 }
369 
370 void SinterNormalSpecies::setSinterForceAndTime(Mdouble adhesionForce, Mdouble sinterTime, Mdouble radius)
371 {
374  setSinterAdhesion(adhesionForce / radius);
375  setInverseSinterViscosity(1.0 / (93.75 * adhesionForce * sinterTime / std::pow(radius, 5)));
376  logger(INFO, "Set sintering parameters: adhesion force f_a=%*r, sinter viscosity is nu = contactRadius^4/%",
378 }
379 
385  (Mdouble alpha, Mdouble beta, Mdouble atomicVolume /*Omega*/, Mdouble surfaceEnergy /*gamma_s*/,
386  Mdouble thicknessDiffusion /*deltaB*D0B*/, Mdouble activationEnergy /*QB*/, Mdouble temperature /*T*/)
387 {
388  setSinterType(SINTERTYPE::PARHAMI_MCKEEPING);
389  const Mdouble boltzmannConstant /*k_B*/ = 1.38064852e-23;
390  const Mdouble gasConstant /*R_g*/ = 8.314459848;
391  const Mdouble thicknessDiffusionVacancy /*DB*/ =
392  thicknessDiffusion * exp(-activationEnergy / gasConstant / temperature);
393  const Mdouble diffusionParameter /*DeltaB*/ =
394  atomicVolume / boltzmannConstant / temperature * thicknessDiffusionVacancy;
395  setInverseSinterViscosity(constants::pi / (2.0 * beta * diffusionParameter));
396  setSinterAdhesion(alpha / beta * constants::pi * surfaceEnergy);
397  logger(INFO,
398  "Set sintering parameters: the adhesion force f_a=%*r, the rate of change of the plastic overlap is d(delta0)/dt=f_ep/(%*a^4)",
399  getSinterAdhesion(),
400  getInverseSinterViscosity());
401 }
402 
407 {
408  setLoadingStiffness(new_.k);
409  setDissipation(new_.disp);
410 }
411 
416 {
417  return dissipation_;
418 }
419 
424 {
425  return sinterAdhesion_;
426 }
427 
432 {
434 }
435 
440 {
441  return sinterRate_;
442 }
443 
445 {
446  return complianceZero_;
447 }
448 
450 {
451  return surfTension_;
452 }
453 
455 {
456  return constantC1_;
457 }
458 
460 {
461  return separationDis_;
462 }
463 
465 {
466  return sinterType_;
467 }
468 
470 {
471  return temperatureDependentSinterRate_(temperature);
472 }
473 
474 std::function<double(double temperature)> SinterNormalSpecies::getTemperatureDependentSinterRate() const
475 {
477 }
478 
480  std::function<double(double temperature)> temperatureDependentSinterRate)
481 {
482  temperatureDependentSinterRate_ = temperatureDependentSinterRate;
483 }
484 
485 
494 {
495  if (eps == 0.0)
496  {
498  dissipation_ = std::sqrt(2.0 * mass * loadingStiffness_);
499  }
500  else
501  {
502  dissipation_ = -mass / tc * std::log(eps);
504  .5 * mass * (mathsFunc::square(constants::pi / tc) + mathsFunc::square(dissipation_ / mass));
505  }
507 }
508 
516 {
517  loadingStiffness_ = stiffness;
518  if (eps == 0.0)
519  {
520  dissipation_ = std::sqrt(2.0 * mass * stiffness);
521  }
522  else
523  {
524  dissipation_ =
525  -std::sqrt(2.0 * mass * stiffness / (constants::sqr_pi + mathsFunc::square(log(eps)))) * log(eps);
526  }
528 }
529 
532 {
533  if (mass <= 0)
534  {
535  logger(ERROR, "Error in getCollisionTime(%) mass is not set or has an unexpected value, (getCollisionTime(%))",
536  mass, mass);
537  }
538  if (loadingStiffness_ <= 0)
539  {
540  logger(ERROR, "Error in getCollisionTime(%) stiffness=% is not set or has an unexpected value "
541  "(getCollisionTime(%), with stiffness=%)",
542  mass, loadingStiffness_, mass, loadingStiffness_);
543  }
544  if (dissipation_ < 0)
545  {
546  logger(ERROR, "Error in getCollisionTime(%) dissipation=%"
547  " is not set or has an unexpected value, (getCollisionTime(%), with dissipation=%)",
548  mass, dissipation_, mass, dissipation_);
549  }
550  Mdouble tosqrt = loadingStiffness_ / (.5 * mass) - mathsFunc::square(dissipation_ / mass);
551  if (tosqrt <= -1e-8 * loadingStiffness_ / (.5 * mass))
552  {
553  logger(ERROR, "Error in getCollisionTime(%) values for mass, stiffness and dissipation would lead to an "
554  "overdamped system, (getCollisionTime(%), with stiffness=% and dissipation=%)",
555  mass, mass, loadingStiffness_, dissipation_);
556  }
557  return constants::pi / std::sqrt(tosqrt);
558 }
void setParhamiMcKeeping(Mdouble alpha, Mdouble beta, Mdouble atomicVolume, Mdouble surfaceEnergy, Mdouble thicknessDiffusion, Mdouble activationEnergy, Mdouble temperature)
Sets the sinterAdhesion_ and inverseSinterViscosity_ based on the Parhami-McKeeping parameters...
Mdouble unloadingStiffnessMax_
the maximum elastic constant (k_2^max) for plastic deformations
SINTERTYPE sinterType_
sinterType options determin how the rate of sintering, d(equilibriumOverlap)/dt is computed PARHAMI_M...
SinterNormalSpecies()
The default constructor.
void setSeparationDis(Mdouble separationDis)
void mix(SinterNormalSpecies *S, SinterNormalSpecies *T)
creates default values for mixed species
void setDissipation(Mdouble dissipation)
Sets the linear dissipation coefficient of the linear plastic-viscoelastic normal force...
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
return type specifically for fuctions returning k and disp at once
Definition: Helpers.h:45
LL< Log::INFO > INFO
Info log level.
Definition: Logger.cc:55
void setCohesionStiffness(Mdouble cohesionStiffness)
Sets the cohesive stiffness of the linear plastic-viscoelastic normal force.
void setTemperatureDependentSinterRate(std::function< double(double temperature)> temperatureDependentSinterRate)
Mdouble exp(Mdouble Exponent)
Definition: ExtendedMath.cc:84
Mdouble loadingStiffness_
(normal) spring constant (k_1)
Mdouble cohesionStiffness_
the adhesive spring constant (k^c) for plastic deformations
Mdouble beta(Mdouble z, Mdouble w)
This is the beta function, returns the approximation based on cmath's implementation of ln(gamma) ...
void read(std::istream &is)
Reads the species properties from an input stream.
void setPenetrationDepthMax(Mdouble penetrationDepthMax)
Sets the maximum penetration depth of the linear plastic-viscoelastic normal force.
void setPlasticParameters(Mdouble loadingStiffness, Mdouble unloadingStiffnessMax, Mdouble cohesionStiffness, Mdouble penetrationDepthMax)
Sets all parameters of the linear plastic-viscoelastic normal force at once.
Mdouble getSinterAdhesion() const
Accesses sinterAdhesion_.
Mdouble computeTimeStep(Mdouble mass)
Returns the optimal time step to resolve a collision of two particles of a given mass.
void setStiffnessAndRestitutionCoefficient(Mdouble k_, Mdouble eps, Mdouble mass)
Sets k, disp such that it matches a given tc and eps for a collision of two copies of P...
void setSinterAdhesion(Mdouble sinterAdhesion)
Sets sinterAdhesion_.
Mdouble penetrationDepthMax_
the depth (relative to the normalized radius) at which k_2^max is used (phi_f)
Mdouble getCollisionTime(Mdouble mass)
Calculates collision time for two copies of a particle of given disp, k, mass.
Mdouble log(Mdouble Power)
void setCollisionTimeAndRestitutionCoefficient(Mdouble tc, Mdouble eps, Mdouble mass)
Set k, disp such that is matches a given tc and eps for a collision of two different masses...
SINTERTYPE getSinterType() const
void setInverseSinterViscosity(Mdouble inverseSinterViscosity)
Sets inverseSinterViscosity_.
void write(std::ostream &os) const
Writes the species properties to an output stream.
void setSinterRate(Mdouble sinterRate)
Sets sinterRate_.
void setSinterForceAndTime(Mdouble adhesionForce, Mdouble sinterTime, Mdouble radius)
void setSurfTension(Mdouble complianceZero)
MERCURY_DEPRECATED void setLoadingStiffnessAndDissipation(helpers::KAndDisp new_)
Allows the spring and dissipation constants to be changed simultaneously.
Mdouble disp
Definition: Helpers.h:49
const Mdouble pi
Definition: ExtendedMath.h:45
Mdouble dissipation_
linear dissipation coefficient
std::function< double(double temperature)> temperatureDependentSinterRate_
static Mdouble average(Mdouble a, Mdouble b)
Returns the harmonic mean of two variables.
Definition: BaseSpecies.cc:110
Mdouble getComplianceZero() const
Mdouble getLoadingStiffness() const
Returns the loading stiffness of the linear plastic-viscoelastic normal force.
Mdouble getCohesionStiffness() const
Returns the cohesive stiffness of the linear plastic-viscoelastic normal force.
void setConstantC1(Mdouble constantC1_)
Mdouble getSeparationDis() const
~SinterNormalSpecies()
The default destructor.
void setLoadingStiffness(Mdouble loadingStiffness)
Sets the loading stiffness of the linear plastic-viscoelastic normal force.
std::string getBaseName() const
Used in Species::getName to obtain a unique name for each Species.
void setComplianceZero(Mdouble complianceZero)
Mdouble getPenetrationDepthMax() const
Returns the maximum penetration depth of the linear plastic-viscoelastic normal force.
SinterNormalSpecies contains the parameters used to describe a plastic-cohesive normal force (Stefan ...
std::function< double(double temperature)> getTemperatureDependentSinterRate() const
void setUnloadingStiffnessMax(Mdouble unloadingStiffnessMax)
Sets the maximum unloading stiffness of the linear plastic-viscoelastic normal force.
Defines the basic properties that a interactable object can have.
Mdouble getUnloadingStiffnessMax() const
Returns the maximum unloading stiffness of the linear plastic-viscoelastic normal force...
Mdouble getDissipation() const
Allows the normal dissipation to be accessed.
Mdouble getSurfTension() const
void setSinterType(SINTERTYPE sinterType)
Sets sinterRate_.
T square(const T val)
squares a number
Definition: ExtendedMath.h:106
Mdouble getInverseSinterViscosity() const
Accesses inverseSinterViscosity_.
Mdouble getConstantC1() const
const Mdouble sqr_pi
Definition: ExtendedMath.h:47
Mdouble getSinterRate() const
Accesses sinterRate_.