diff --git a/doc/src/pair_granular.rst b/doc/src/pair_granular.rst index 696c87b0ab..7998f11cc5 100644 --- a/doc/src/pair_granular.rst +++ b/doc/src/pair_granular.rst @@ -44,7 +44,7 @@ Examples pair_coeff * * hertz 1000.0 50.0 tangential mindlin 1000.0 1.0 0.4 heat area 0.1 pair_style granular - pair_coeff * * mdr 5e6 0.4 1.9e5 2.0 0.5 0.5 tangential linear_history 940.0 0.0 0.7 rolling sds 2.7e5 0.0 0.6 damping mdr 1 + pair_coeff * * mdr 5e6 0.4 1.9e5 2.0 0.5 0.5 tangential linear_history 940.0 1.0 0.7 rolling sds 2.7e5 0.0 0.6 damping mdr 1 Description """"""""""" @@ -89,7 +89,7 @@ and their required arguments are: 4. *dmt* : E, :math:`\eta_{n0}` (or :math:`e`), :math:`\nu`, :math:`\gamma` 5. *jkr* : E, :math:`\eta_{n0}` (or :math:`e`), :math:`\nu`, :math:`\gamma` 6. *mdr* : :math:`E`, :math:`\nu`, :math:`Y`, :math:`\Delta\gamma`, - :math:`\psi_b`, :math:`\eta_{n0}`, :math:`d_{type}` + :math:`\psi_b`, :math:`\eta_{n0}` Here, :math:`k_n` is spring stiffness (with units that depend on model choice, see below); :math:`\eta_{n0}` is a damping prefactor (or, in its @@ -206,9 +206,6 @@ The model requires the following inputs: 6. *Damping coefficent* :math:`\eta_{n0} \ge 0` : The damping coefficient is a tunable parameter that controls damping in the normal direction. - 7. *Damping type* :math:`d_{type} =` 1 or 2 : The damping type specfies the - damping model defined for the *mdr* damping class described below. - .. note:: The values for :math:`E`, :math:`\nu`, :math:`Y`, and :math:`\Delta\gamma` (i.e., @@ -263,7 +260,7 @@ The *mdr* model currently only supports *fix wall/gran/region*, not any *fix wall/gran/region* commands must also use the *mdr* model. Additionally, the following *mdr* inputs must match between the *pair_style* and *fix wall/gran/region* definitions: :math:`E`, -:math:`\nu`, :math:`Y`, :math:`\psi_b`, and :math:`e`. The exception +:math:`\nu`, :math:`Y`, :math:`\psi_b`, and :math:`\eta_{n0}`. The exception is :math:`\Delta\gamma`, which may vary, permitting different adhesive behaviors between particle-particle and particle-wall interactions. @@ -329,7 +326,7 @@ for the damping model currently supported are: 3. *viscoelastic* 4. *tsuji* 5. *coeff_restitution* -6. *mdr* +6. *mdr* (class) : :math:`d_{type}` If the *damping* keyword is not specified, the *viscoelastic* model is used by default. @@ -420,9 +417,8 @@ restitution for both monodisperse and polydisperse particle pairs. This damping model is not compatible with cohesive normal models such as *JKR* or *DMT*. The *mdr* damping class contains multiple damping models that can be toggled between -by specifying different integer values for the :math:`d_{type}` parameter in the *mdr* -normal model. This damping option is therefore only compatible with the normal *mdr* -contact model. +by specifying different integer values for the :math:`d_{type}` input parameter. This +damping option is only compatible with the normal *mdr* contact model. Setting :math:`d_{type} = 1` is the suggested damping option. This specifies a damping model that takes into account the contact stiffness :math:`k_{mdr}` calulated diff --git a/examples/granular/in.tableting.200 b/examples/granular/in.tableting.200 index 02eb6743de..58fc8dde3e 100644 --- a/examples/granular/in.tableting.200 +++ b/examples/granular/in.tableting.200 @@ -50,8 +50,8 @@ variable mu_roll equal 0.6 variable k_roll equal 2.25*${mu_roll}*${mu_roll}*${YoungsModulus}*${atomRadius} variable gamma_roll equal 0.0 -pair_coeff * * mdr ${YoungsModulus} ${PoissonsRatio} ${YieldStress} ${SurfaceEnergy} ${psi_b} ${damp} ${damp_type} & - damping mdr & +pair_coeff * * mdr ${YoungsModulus} ${PoissonsRatio} ${YieldStress} ${SurfaceEnergy} ${psi_b} ${damp} & + damping mdr ${damp_type} & tangential linear_history ${kt} ${xgammat} ${mu_s} & rolling sds ${k_roll} ${gamma_roll} ${mu_roll} @@ -60,7 +60,7 @@ pair_coeff * * mdr ${YoungsModulus} ${PoissonsRatio} ${YieldStress} ${SurfaceEne variable disp_upper equal 0.0 variable disp_lower equal 0.0 -variable wall_contact_string string "granular mdr ${YoungsModulus} ${PoissonsRatio} ${YieldStress} ${SurfaceEnergyWall} ${psi_b} ${damp} ${damp_type} damping mdr tangential linear_history ${kt_wall} ${xgammat} ${mu_s_wall} rolling sds ${k_roll} ${gamma_roll} ${mu_roll}" +variable wall_contact_string string "granular mdr ${YoungsModulus} ${PoissonsRatio} ${YieldStress} ${SurfaceEnergyWall} ${psi_b} ${damp} damping mdr ${damp_type} tangential linear_history ${kt_wall} ${xgammat} ${mu_s_wall} rolling sds ${k_roll} ${gamma_roll} ${mu_roll}" variable dieHeight2 equal 2*${dieHeight} diff --git a/examples/granular/in.triaxial.compaction.12 b/examples/granular/in.triaxial.compaction.12 index 6a948157d3..0a35ab3c0a 100644 --- a/examples/granular/in.triaxial.compaction.12 +++ b/examples/granular/in.triaxial.compaction.12 @@ -38,8 +38,8 @@ variable kt equal 2/7*${YoungsModulus}*${atomRadius} variable xgammat equal 0.0 variable mu_s equal 0.5 -pair_coeff * * mdr ${YoungsModulus} ${PoissonsRatio} ${YieldStress} ${SurfaceEnergy} ${psi_b} ${damp} ${damp_type} & - damping mdr & +pair_coeff * * mdr ${YoungsModulus} ${PoissonsRatio} ${YieldStress} ${SurfaceEnergy} ${psi_b} ${damp} & + damping mdr ${damp_type} & tangential linear_history ${kt} ${xgammat} ${mu_s} ######################################### ADD IN PLANES ################################################ @@ -57,7 +57,7 @@ region plane_xz_neg plane 0 -${halfBoxWidth} 0 0 1 0 side in move NULL v_plane_d region plane_xy_pos plane 0 0 ${halfBoxWidth} 0 0 -1 side in move NULL NULL v_plane_disp_neg units box region plane_xy_neg plane 0 0 -${halfBoxWidth} 0 0 1 side in move NULL NULL v_plane_disp units box -variable wall_contact_string string "granular mdr ${YoungsModulus} ${PoissonsRatio} ${YieldStress} ${SurfaceEnergy} ${psi_b} ${damp} ${damp_type} damping mdr tangential linear_history ${kt} ${xgammat} ${mu_s} " +variable wall_contact_string string "granular mdr ${YoungsModulus} ${PoissonsRatio} ${YieldStress} ${SurfaceEnergy} ${psi_b} ${damp} damping mdr ${damp_type} tangential linear_history ${kt} ${xgammat} ${mu_s} " fix plane_yz_pos all wall/gran/region ${wall_contact_string} region plane_yz_pos contacts fix plane_yz_neg all wall/gran/region ${wall_contact_string} region plane_yz_neg contacts diff --git a/src/GRANULAR/gran_sub_mod_damping.cpp b/src/GRANULAR/gran_sub_mod_damping.cpp index c35d63d668..b68a720413 100644 --- a/src/GRANULAR/gran_sub_mod_damping.cpp +++ b/src/GRANULAR/gran_sub_mod_damping.cpp @@ -190,6 +190,13 @@ void GranSubModDampingCoeffRestitution::init() GranSubModDampingMDR::GranSubModDampingMDR(GranularModel *gm, LAMMPS *lmp) : GranSubModDamping(gm, lmp) { + num_coeffs = 1; +} + +void GranSubModDampingMDR::coeffs_to_local() +{ + damp_type = coeffs[0]; // damping type 1 = mdr stiffness or 2 = velocity + if (damp_type != 1 && damp_type != 2) error->all(FLERR, "Illegal MDR normal model, damping type must an integer equal to 1 or 2"); } /* ---------------------------------------------------------------------- */ @@ -200,7 +207,6 @@ void GranSubModDampingMDR::init() error->all(FLERR, "Damping mdr can only be used with mdr normal model"); damp = gm->normal_model->get_damp(); - damp_type = gm->normal_model->get_damp_type(); } /* ---------------------------------------------------------------------- */ diff --git a/src/GRANULAR/gran_sub_mod_damping.h b/src/GRANULAR/gran_sub_mod_damping.h index c73a4a8e5e..274c3d9e29 100644 --- a/src/GRANULAR/gran_sub_mod_damping.h +++ b/src/GRANULAR/gran_sub_mod_damping.h @@ -42,7 +42,7 @@ namespace Granular_NS { protected: double damp_prefactor; double damp; - int damp_type; + int damp_type; // damping type is only used by normal mdr model }; /* ---------------------------------------------------------------------- */ @@ -99,6 +99,7 @@ namespace Granular_NS { class GranSubModDampingMDR : public GranSubModDamping { public: GranSubModDampingMDR(class GranularModel *, class LAMMPS *); + void coeffs_to_local() override; void init() override; double calculate_forces() override; }; diff --git a/src/GRANULAR/gran_sub_mod_normal.cpp b/src/GRANULAR/gran_sub_mod_normal.cpp index aa53f5fcf6..68cb853d9d 100644 --- a/src/GRANULAR/gran_sub_mod_normal.cpp +++ b/src/GRANULAR/gran_sub_mod_normal.cpp @@ -441,7 +441,7 @@ GranSubModNormalMDR::GranSubModNormalMDR(GranularModel *gm, LAMMPS *lmp) : { if (lmp->citeme) lmp->citeme->add(cite_mdr); - num_coeffs = 7; + num_coeffs = 6; contact_radius_flag = 1; size_history = 27; nsvector = 1; @@ -474,7 +474,6 @@ void GranSubModNormalMDR::coeffs_to_local() gamma = coeffs[3]; // effective surface energy psi_b = coeffs[4]; // bulk response trigger based on ratio of remaining free area: A_{free}/A_{total} damp = coeffs[5]; // coefficent of restitution - damp_type = coeffs[6]; // damping type 1 = mdr or 2 = velocity if (E <= 0.0) error->all(FLERR, "Illegal MDR normal model, Young's modulus must be greater than 0"); if (nu < 0.0 || nu > 0.5) error->all(FLERR, "Illegal MDR normal model, Poisson's ratio must be between 0 and 0.5"); @@ -482,7 +481,6 @@ void GranSubModNormalMDR::coeffs_to_local() if (gamma < 0.0) error->all(FLERR, "Illegal MDR normal model, effective surface energy must be greater than or equal to 0"); if (psi_b < 0.0 || psi_b > 1.0) error->all(FLERR, "Illegal MDR normal model, psi_b must be between 0 and 1.0"); if (damp < 0.0) error->all(FLERR, "Illegal MDR normal model, damping coefficent must be greater than or equal to 0"); - if (damp_type != 1 && damp_type != 2) error->all(FLERR, "Illegal MDR normal model, damping type must an integer equal to 1 or 2"); G = E / (2.0 * (1.0 + nu)); // shear modulus kappa = E / (3.0 * (1.0 - 2.0 * nu)); // bulk modulus diff --git a/src/GRANULAR/gran_sub_mod_normal.h b/src/GRANULAR/gran_sub_mod_normal.h index 0b7bb664d5..8142d88cea 100644 --- a/src/GRANULAR/gran_sub_mod_normal.h +++ b/src/GRANULAR/gran_sub_mod_normal.h @@ -41,7 +41,6 @@ namespace Granular_NS { int get_cohesive_flag() const { return cohesive_flag; } double get_damp() const { return damp; } - double get_damp_type() const { return damp_type; } double get_emod() const { return Emod; } double get_fncrit() const { return Fncrit; } int get_material_properties() const { return material_properties; } @@ -52,7 +51,6 @@ namespace Granular_NS { protected: double damp; // argument historically needed by damping // typically (but not always) equals eta_n0 - int damp_type; // damping type is only used by normal mdr model double Emod, poiss; double Fncrit; int material_properties, cohesive_flag;