Update Hertz-Mindlin model and docs
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
andIP2_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}
andcs=2 \cdot \beta_s \cdot \sqrt{mbar \cdot k_s}
, wherek_n
,k_s
the tangent normal and shear stiffness values, as in [Thornton2015]. I should note that forc_n
specifically, [PFC3D], [Tsuji1992] and [Antypov2011] consider the secant stiffnessk_n=4/3 \cdot E \cdot \sqrt{R \cdot u_n}
and not the tangent onek_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.