From a420d04418b967523b1cf3f311d76f41a4dbcc18 Mon Sep 17 00:00:00 2001 From: jtclemm Date: Tue, 26 Nov 2024 11:28:40 -0700 Subject: [PATCH] Unifying desc/code for tsuji and coeff restitution --- doc/src/pair_granular.rst | 20 ++++++++++---------- src/GRANULAR/gran_sub_mod_damping.cpp | 17 +---------------- src/GRANULAR/gran_sub_mod_damping.h | 3 +-- 3 files changed, 12 insertions(+), 28 deletions(-) diff --git a/doc/src/pair_granular.rst b/doc/src/pair_granular.rst index b4341961cc..8558832464 100644 --- a/doc/src/pair_granular.rst +++ b/doc/src/pair_granular.rst @@ -253,26 +253,26 @@ of the normal contact model parameters should be between 0 and 1, but no error check is performed on this. The *coeff_restitution* model is useful when a specific normal coefficient of -restitution :math:`e` is required. In these models, the normal coefficient of -restitution :math:`e` is specified as an input in place of the usual :math:`\eta_{n0}` -value in the normal model. Following the approach of -:ref:`(Brilliantov et al) `, when using the *hooke* normal model, -*coeff_restitution* then calculates the damping coefficient as: +restitution :math:`e` is required. It operates much like the *Tsuji* model +but, the normal coefficient of restitution :math:`e` is specified as an input +in place of the usual :math:`\eta_{n0}` value in the normal model. Following +the approach of :ref:`(Brilliantov et al) `, when using the *hooke* +normal model, *coeff_restitution* then calculates the damping coefficient as: .. math:: - \eta_n = \sqrt{\frac{4m_{eff}k_n}{1+\left( \frac{\pi}{\log(e)}\right)^2}} , + \eta_n = \sqrt{\frac{4m_{eff}k_{nd}}{1+\left( \frac{\pi}{\log(e)}\right)^2}} , +where :math:`k_{nd}` is the same stiffness defined in the above *Tsuji* model. For any other normal model, e.g. the *hertz* and *hertz/material* models, the damping coefficient is: .. math:: - \eta_n = -2\sqrt{\frac{5}{6}}\frac{\log(e)}{\sqrt{\pi^2+(\log(e))^2}}(R_{eff} \delta_{ij})^{\frac{1}{4}}\sqrt{\frac{3}{2}k_n m_{eff}} , + \eta_n = -2\sqrt{\frac{5}{6}}\frac{\log(e)}{\sqrt{\pi^2+(\log(e))^2}}(R_{eff} \delta_{ij})^{\frac{1}{4}}\sqrt{\frac{3}{2}k_{nd} m_{eff}} , -where :math:`k_n = \frac{4}{3} E_{eff}` for the *hertz/material* model. Since -*coeff_restitution* accounts for the effective mass, effective radius, and -pairwise overlaps (except when used with the *hooke* normal model) when calculating +Since *coeff_restitution* accounts for the effective mass, effective radius, +and pairwise overlaps (except when used with the *hooke* normal model) when calculating the damping coefficient, it accurately reproduces the specified coefficient of restitution for both monodisperse and polydisperse particle pairs. This damping model is not compatible with cohesive normal models such as *JKR* or *DMT*. diff --git a/src/GRANULAR/gran_sub_mod_damping.cpp b/src/GRANULAR/gran_sub_mod_damping.cpp index 9dd68060ac..3aa888372c 100644 --- a/src/GRANULAR/gran_sub_mod_damping.cpp +++ b/src/GRANULAR/gran_sub_mod_damping.cpp @@ -152,9 +152,8 @@ double GranSubModDampingTsuji::calculate_forces() ------------------------------------------------------------------------- */ GranSubModDampingCoeffRestitution::GranSubModDampingCoeffRestitution(GranularModel *gm, LAMMPS *lmp) : - GranSubModDamping(gm, lmp) + GranSubModDampingTsuji(gm, lmp) { - allow_cohesion = 0; } /* ---------------------------------------------------------------------- */ @@ -171,17 +170,3 @@ void GranSubModDampingCoeffRestitution::init() damp /= sqrt(MY_PI * MY_PI + logcor * logcor); } } - -/* ---------------------------------------------------------------------- */ - -double GranSubModDampingCoeffRestitution::calculate_forces() -{ - // in case argument <= 0 due to precision issues - double sqrt1; - if (gm->delta > 0.0) - sqrt1 = MAX(0.0, gm->meff * gm->Fnormal / gm->delta); - else - sqrt1 = 0.0; - damp_prefactor = damp * sqrt(sqrt1); - return -damp_prefactor * gm->vnnr; -} diff --git a/src/GRANULAR/gran_sub_mod_damping.h b/src/GRANULAR/gran_sub_mod_damping.h index c931e385cc..98c31d680a 100644 --- a/src/GRANULAR/gran_sub_mod_damping.h +++ b/src/GRANULAR/gran_sub_mod_damping.h @@ -87,11 +87,10 @@ namespace Granular_NS { /* ---------------------------------------------------------------------- */ - class GranSubModDampingCoeffRestitution : public GranSubModDamping { + class GranSubModDampingCoeffRestitution : public GranSubModDampingTsuji { public: GranSubModDampingCoeffRestitution(class GranularModel *, class LAMMPS *); void init() override; - double calculate_forces() override; }; /* ---------------------------------------------------------------------- */