Skip to content

Update Hertz-Mindlin model and docs

Vasileios Angelidakis requested to merge updateHertzMindlin into master

In this MR I submit the following changes to our Hertz-Mindlin (no-slip solution) implementation:

  • Effective shear modulus G* calculation
  • Shear stiffness coefficient Kno calculation
  • Viscous damping model harmonisation.
  • Doc changes to reflect these
  • I did these also for the Ip2_FrictMat_FrictMat_MindlinCapillaryPhys and IP2_PartialSatMat_PartialSatMat_MindlinPhys which use the same model. @robcaulk do let me know if you see any issue with that.

For the viscous damping model specifically, I made the following modifications:

  • I corrected the relationship between betan-en and betas-es, by applying Eqs. (B6)-(B7) from [Thornton2013], which work for the Hertzian model where the end of a contact is considered when the normal force is zero and not when the normal displacement is zero (i.e. when attractive normal forces are not allowed).

  • When the coefficients of restitution (en, es) are defined, I calculated c_n=2 \cdot \beta_n \cdot \sqrt{mbar \cdot k_n} and cs=2 \cdot \beta_s \cdot \sqrt{mbar \cdot k_s}, where k_n,k_s the tangent normal and shear stiffness values, as in [Thornton2015]. I should note that for c_n specifically, [PFC3D], [Tsuji1992] and [Antypov2011] consider the secant stiffness k_n=4/3 \cdot E \cdot \sqrt{R \cdot u_n} and not the tangent one k_n=2 \cdot E \cdot \sqrt{R \cdot u_n}, in contrast to [Thornton2015]. I used the tangent one as in [Thornton2015] but I am happy to change that if you think otherwise.

So far, when the user defined the coefficients of restitution (en, es) as input, only en was being used and it was assumed that es=en and also that cs=cn, as in Tsuji et al (1992). I believe that both assumptions are not a good idea. Especially the latter, since the elastic shear forces are calculated with the shear stiffness, and the shear damping coefficient cs should be also proportional to that stiffness, and not to the normal stiffness, to avoid overdamping collisions (since kn>ks).

  • I should note that [Thornton2015] assumed betan=betas (=gamma in his book/papers) but I used his Eqs. (B6)-(B7) for different en/es, as we already allow different viscous damping ratios betan/betas (also PFC3D does).

  • With my current analysis, the user can define either en/es or betan/betas and get to the same results, which was not achieved so far.

Edited by Vasileios Angelidakis

Merge request reports