自然界中的岩体常年受复杂水文地质条件等环境因素影响，因此岩体内部结构通常存在着裂缝等不连续结构面，表现出较强的非连续性。岩体裂缝的力学效应作为岩体力学的中心内容，且岩体裂缝扩展和贯通的原理是预防和控制重大工程地质灾害的基础，所以研究岩体裂缝的力学效应具有十分重要的现实意义和理论价值。然而，基于连续介质的方法只能够模拟裂缝形成之前的断裂过程，无法模拟裂缝形成后岩石块体间的碰撞、运动过程[1]。相反，基于不连续介质的数值方法，例如非连续变形分析方法(discontinuous deformation analysis，简称DDA)[2]、数值流行方法(numerical manifold method，简称NMM)[3]、离散单元法(discrete element method，简称DEM)[4-5]和光滑粒子流体动力学方法(smoothed particle hydro- dynamics，简称SPH)[6]可以作为有效的替代方法弥补此类不足，已广泛应用于分析岩体之间的相互作用[7]。其中，离散单元法(DEM)因其理论的成熟性，已广泛应用于研究离散颗粒的力学行为[8-9]。
球体单元因其简单有效的接触模式常被用于离散单元法的数值模拟。但是，它不能准确模拟不规则形状的块体。多面体单元因其形状和岩体具有相似性，在离散元中常被用于模拟不规则块体。但是，对于不同的接触模式，多面体单元缺少一个简单统一的数学模型。在计算块体单元之间的接触力之前，需要确定块体之间的接触模式，例如点-面接触、点-边接触、点-点接触、面-面接触、边-边接触和边-面接触[10]。此外，由于顶点处法向接触力的不确定性，它不能解决点-点接触模式的问题。这也导致了法向接触力的不一致，进而带来能量失衡和数值不正确的问题。
Munjiza等[11]提出了有限-离散单元耦合分析方法(combined finite-discrete element method，简称FDEM)，与传统离散单元法相比，能够在不考虑多变的接触模式的情况下实现接触力的计算。在该方法中，接触力是由分布力代替中心接触力进行计算，通过对嵌入区域的势函数进行积分，实现接触力的计算，既方便又统一，并保证了能量平衡和动量守恒。大量的应用实践已经证实了FDEM的优越性[12-13]。但FDEM的计算过程是基于复杂的布尔运算，这会大幅度提高计算成本。此外，它对单元类型有严格的限制，复杂模型只能离散成四面体单元而不能离散成其他多面体单元，这会导致单元数量巨增[14]，大幅度降低计算效率。
相比于FDEM，圆化多面体离散单元法具有计算效率高、接触力计算形式简单等优点。圆化多面体是由多面体和球体的闵可夫斯基和生成[15]，是一个和多面体具有同等尺寸但角点已经圆化的多面体[16]，因此两个多面体的接触类似于两球体之间的接触，简化接触模式的判断，提高了计算效率。虽然圆化多面体离散单元法是可行的，但由于该方法认为离散单元服从刚体假设，因此它并不能模拟软颗粒或弹性颗粒的变形。此外，该方法还缺少一个合适的模型来模拟断裂的全过程，包括岩体的变形、断裂和接触相互作用。颗粒的变形会影响岩体的应力分布，进而影响断裂破坏过程。准确模拟岩体的断裂全过程需要既能反映颗粒的相互作用，又能反映颗粒自身的变形。
为解决上述问题，提出一种可用于模拟断裂全过程的三维可变形圆化多面体离散单元法，并建立了一个三维断裂模型用于该方法预测断裂的演化过程。该方法将岩体离散成有限单元，并沿着有限单元的边界嵌入无厚度节理单元。采用虚拟裂缝模型[17]和Mohr-Coulomb准则来判断节理单元的破坏状态。圆化多面体离散单元法解决了离散岩体之间的接触力问题，单元的应力应变问题则由有限单元法求解。与现有的圆化多面体离散单元法相比，该方法能够反映颗粒的应力应变，并能够模拟岩体的断裂全过程，包括变形、断裂和接触相互作用。与FDEM相比，该方法在计算接触力时避免了复杂的布尔运算，大幅度提高了计算效率。此外，它还突破四面体单元的限制，能够应用于更多类别的单元类型。 结 论
(1) 本文提出了一种用于模拟断裂全过程的三维可变形圆化多面体离散单元法。该方法将圆化多面体离散单元法与有限单元法相结合，有限单元法计算颗粒的应力应变，圆化多面体离散单元法计算颗粒之间的接触相互作用，在有限单元的边界嵌入无厚度节理单元，通过判断节理单元的破坏状态来捕捉裂缝萌生与扩展，可以有效模拟岩体断裂的全过程。
(2) 通过5个算例来验证该方法的有效性和可行性，其中颗粒堆积试验模拟结果表明，该方法比FDEM更高效；巴西圆盘劈裂试验验证了该方法在处理拉伸问题的准确性；单缺口梁试验证明了该方法能够捕捉到拉剪复合应力下的断裂裂缝扩展；最后，四点弯曲梁试验和球体冲击试验表明，该方法能够准确捕捉岩体受冲击载荷后裂缝的萌生和扩展。因此，该方法为模拟准脆性材料的断裂全过程提供了一个简单高效的工具。
