[论文解读] Numerical Improvement of the Discrete Element Method applied to Shear of Granular Media

Andrés A. Peña, Pedro G. Lind+2 more|arXiv (Cornell University)|Sep 24, 2007

Geotechnical and Geomechanical Engineering参考文献 3被引用 5

一句话总结

本文指出，在颗粒松弛过程中，标准离散元法（DEM）的时间积分步长通常过大，导致数值不收敛，尤其是在剪切模拟中。为此，本文提出一种基于几何重叠和杆件矢量更新的改进摩擦力计算方法，使模拟在远大于以往认为可能的时间步长下仍能保持稳定与精确，恢复了传统方法失效时的收敛性。

ABSTRACT

We present a detailed analysis of the bounds on the integration step in Discrete Element Method (DEM) for simulating collisions and shearing of granular assemblies. We show that, in the numerical scheme, the upper limit for the integration step, usually taken from the average time $t_c$ of one contact, is in fact not sufficiently small to guarantee numerical convergence of the system during relaxation. In particular, we study in detail how the kinetic energy decays during the relaxation stage and compute the correct upper limits for the integration step, which are significantly smaller than the ones commonly used. In addition, we introduce an alternative approach, based on simple relations to compute the frictional forces, that converges even for integration steps above the upper limit.

研究动机与目标

• 确定DEM颗粒剪切模拟中时间步长的真实上限，该上限远小于通常使用的值。
• 解决在使用大时间步长时，标准DEM在松弛阶段出现的数值不稳定与不收敛问题。
• 开发一种替代的摩擦力计算方法，即使在较大时间步长下也能保持精度与收敛性。
• 将改进方法的有效性从球形颗粒扩展到多边形颗粒，后者接触点动力学更为复杂。
• 为计算效率与数值精度均至关重要的慢速剪切场景，提供可靠的DEM模拟框架。

提出的方法

• 基于松弛过程中的动能衰减，推导出新的时间步长上限，表明其远小于通常使用的 $ t_c $（接触时间）估计值。
• 提出一种几何方法，通过追踪连续时间步长间重叠区域的中心与杆件矢量，计算切向摩擦力。
• 利用预测-校正步骤，基于重叠区域中心与杆件矢量的变化，计算切向力增量的平动与转动贡献。
• 将改进的力模型应用于圆盘与多边形颗粒，考虑旋转与平动引起的接触点位置变化。
• 采用预测-校正格式，并在每一步更新几何量，即使在较大时间步长下也能保持二阶精度。
• 通过比较不同时间步长下的松弛时间，验证该方法，结果表明使用新方法后恢复了收敛性。

实验结果

研究问题

• RQ1 DEM模拟中确保颗粒松弛过程数值收敛的时间步长真实上限是什么？
• RQ2 为何标准Cundall弹簧模型在使用大时间步长时无法维持收敛性？
• RQ3 如何重新定义摩擦力，以在时间步长超过传统极限时仍保持数值精度与稳定性？
• RQ4 该改进方法在多边形颗粒中是否依然有效，其中接触点因旋转而发生位移且重叠几何形状复杂？
• RQ5 新方法是否能在不牺牲精度或收敛性的前提下，使慢速剪切模拟的时间步长显著增大？

主要发现

• DEM的时间步长上限远小于通常使用的 $ t_c $，实际上限取决于摩擦力计算的精度与自由度数量。
• 标准Cundall弹簧模型用于切向力时，误差随时间步长线性增长，导致整体方案精度仅能保持 $ O(\Delta t^2) $，即使采用预测-校正格式亦然。
• 所提出的几何方法计算摩擦力，可确保误差与预测-校正格式一致，即使在较大 $ \Delta t $ 下也能恢复二阶精度。
• 改进方法稳定了模拟，消除了松弛时间对时间步长的依赖性，该结果在圆盘与多边形颗粒系统中均得到验证。
• 该方法在多边形颗粒中依然有效，通过动态更新的杆件矢量与重叠区域中心，显式考虑了因旋转导致的接触点位移。
• 初步结果表明，增大 $ k_t/k_n $ 的比值可提高时间步长的有效上限，表明数值稳定性具有可调性。