FEM/DEM法在典型柱状危岩体破坏过程数值分析中的应用

陈小婷，黄波林(防灾减灾湖北省重点实验室/三峡大学，湖北 宜昌 443002)

目前开展斜坡变形破坏数值模拟的方法有基于有限单元的方法(FEM)、基于离散单元的方法(DEM)和连续离散单元方法(FEM/DEM耦合方法)。在处理复杂材料小变形机制方面，FEM有着明显的优势[1]；而DEM在处理离散介质的大变形问题上，具有非常重大的优势[2]。FEM/DEM耦合分析方法结合了连续和非连续数值方法的优点，既可计算连续介质的小变形，又可在裂隙间开展大变形，同时在采用断裂力学的相关准则后，可以对裂纹的产生和扩展进行模拟[3]。本文选用该方法对典型柱状危岩体变形破坏进行数值模拟分析。

1 FEM/DEM方法研究现状

自1989年Munjiza[4]提出FEM/DEM耦合方法以来，这一方法被应用在了广泛的领域，其有效性和灵活性也得了极大的提高。Munjiza等、Owen等[5]、Xiang J等[6]和Latham等[7]推动了在FEM/DEM模型框架中的裂隙扩展算法。在这个算法中，裂隙以有限单元网格的形式被植入。在有限元法中引入损伤力学和断裂力学理论，当单元网格受力超过强度准则后，连续单元变成离散单元。离散单元的开裂发生在有限单元的边界上，并采用离散元法处理开裂后块体间的接触。基于FEM/DEM方法，Morris等开发了Livermore distinct element code (LDEC)来进行地质体的数值模拟[8]，Rockfield公司则推出了ELFEN。国内LI S H等[9]、魏怀鹏等[10]、冯春等[11]发展了连续—非连续离散单元法(CDEM)。Eberhardt等[12]利用EFLFEN对1991年Swiss的Randa岩质崩滑体进行了数值分析，研究了岩桥的破坏、崩滑体的运动及堆积过程。Vyazmensky等[13]使用基于FEM/DEM方法的ELFEN软件研究了岩体开挖引起地面沉降的机理，指出该法可以加深对开挖过程中岩体的力学行为理解，并改进沉降预测的方法，使用FEM/DEM方法的ELFEN软件对露天采矿引发的边坡破坏及寻找合适开挖方式进行了研究。Mahabadi等[14]使用FEM/DEM方法构建了Y-GUI程序，并对动力条件下的巴西圆盘试验进行数值模拟，模拟结果与试验相吻合。Chang Xiaolin等[15]利用FEM/DEM方法和物理试验同时开展了巴西圆盘试验、单轴抗压试验、单边切口梁试验、平面滑移和撞击试验，模拟结果与试验结果吻合程度较高。

本次研究选用了基于MunjizaY程序的Y-GUI程序开展箭穿洞危岩体的变形破坏过程研究。Y-GUI程序包括Mahabadi等建立的FEM/DEM前处理模块和FEM/DEM方法的Y程序计算模块。该模型包括摩尔库伦模型、岩石断裂模型和摩擦模型，可以用于模拟分析岩质斜坡的变形破坏过程。

2 典型柱状危岩体概况

箭穿洞危岩体位于三峡库区重庆市巫山县长江左岸，上距巫山县城约15 km，位于巫峡著名的神女峰西侧坡脚处，距下游神女溪口和下游青石居民聚集区1.3 km。在构造上，箭穿洞斜坡位于神女峰背斜的核部。神女峰背斜核部出露的最老地层为三叠系大冶组三段(T1d3)，然后依次出露大冶组四段(T1d4)和三叠系嘉陵江组1-4段(T1j1-4)。

箭穿洞斜坡在地貌上呈现陡崖与中缓坡交替的台阶状地貌，共有三级台阶。第一级台阶下发育箭穿洞危岩体，其台阶面为T1d4。由于岩层位于背斜核部，地层较平缓，产状为260°∠6°。因此，箭穿洞危岩体斜坡总体为近水平层状岸坡，岩性组成为高程155 m以下的中层状含泥质条带灰岩(T1d3)和高程155 m以上的灰色厚层砾屑、砂屑灰岩(T1d4)。

箭穿洞危岩体的上游侧边界为冲沟，下游侧边界为一上陡下缓的大型裂隙，其上部产状为324°∠65°，中段局部近直立转缓倾，下部产状324°∠45°。危岩体从临空面陡崖至后缘共发现有卸荷裂缝，产状为276°～260°∠75°～85°。后缘高程为278～305 m，基座高程为155 m，平均高差为135 m。平均横宽约55 m，平均厚度约50 m，体积约36×104m3，主崩方向为260°(图1)，其典型剖面图见图2。危岩体基座较破碎，但主柱体完整。

图1 箭穿洞危岩体平面图Fig.1 Geological engineering map of the Jianchuandong rockmass1—裂隙；2—岩体边界；3—地层界线；4—岩层和结构面产状

泥质条带灰岩位于箭穿洞危岩体的基座部位，相对厚层灰岩较为软弱。基座岩体长期受压，较破碎。2008年三峡库区175 m蓄水后，泥质条带灰岩处于库水波动带。水位周期性波动造成该基座岩体处于浸泡-风干周期性作用下，造成基座岩体加速劣化[10]。基座的加速变形破坏，导致危岩体稳定性降低，发生基座压裂式崩塌。

图2 箭穿洞危岩体工程地质剖面图Fig.2 Geo-engineering section map of the Jianchuandong dangerous rockmass

柱状危岩体的变形破坏过程取决于危岩体的应力应变场。同时，由于危岩体入江会产生涌浪，入江方式不同会造成涌浪的大小差异。因此，需要对该柱状危岩体的变形破坏过程进行分析。

3 危岩体FEM/DEM模型构建

野外调查表明，箭穿洞危岩体存在缓倾层面和较少的近垂向裂隙，除此以外没有其他裂隙发育。危岩体整体可划分为完整厚-块状岩体结构，要模拟这样一种岩体结构的变形破坏过程，需要能实现大变形和形成新裂隙的数值工具。本文采用FEM/DEM模型的Y-GUI程序，对箭穿洞危岩体进行数值模型构建、参数设置和破坏全过程计算。

图3为箭穿洞危岩体的FEM/DEM模型，模型长384 m，高222 m。模型总单元数为700个，其中危岩体及基座岩体单元数为250个，它们是可离散单元；斜坡岩体单元数为450个，它们是不可离散单元。模型中岩体和结构面所用参数见表1，模型采用弹性Mohr-Coulomb准则和拉张破坏准则。

4 危岩体破坏过程分析

从初期的危岩体受力情况来看，基座和基座附近岩体明显有应力集中现象(图4)。危岩体基座的拉应力主要来源于压力，即压致拉裂效应。从此时基座岩体的运动矢量图来看(图5)，基座岩体受压向两侧和向下运动明显。

图3 箭穿洞危岩体的FEM/DEM模型Fig.3 FEM/DEM model of the Jianchuandong dangerous rockmass

表1 模型中岩体和结构面参数Table 1 Rockmass and structure face parameters for the numerical model

图4 基座岩体的应力集中现象Fig.4 Stress concentration in the basic rockmass注：红色为应力阴影

图5 基座岩体的运动方向Fig.5 Moving direction of the basic rockmass注：红色为向x-方向，蓝色向x+方向

裂隙最开始新生在基座泥质条带灰岩以及与泥质条带灰岩接触的灰岩处。基座泥质条带灰岩区由于压力作用裂隙虽然出现破坏但没有分离，因此显示为局部出现空洞，而不是裂隙(红色区域)，基座裂隙最先出现在两端(图6)。上部危岩体的裂隙大多是从底部开始往上延伸，其形成机制是典型的压致拉裂，且裂隙最开始出现的位置靠近内侧(图7)。

图6 基座岩体最初破坏位置图Fig.6 Original failure location map of the basic rockmass注：黑色圈为破坏位置，蓝色为压破坏区域，红色为拉张破坏区域

图7 危岩体初始破坏图Fig.7 Original failure of the main rockmass注：红色线为岩体出现张开的裂隙

当基座岩体和危岩体岩体开始出现新裂隙后，裂隙网络迅速发展，最终将危岩体分成3个相互关联的岩柱。裂隙网络的发展和基座的破坏导致了单个岩柱出现破坏。破坏从临空面开始，基座开始垮塌，然后上部的岩柱或滑动或垮落或倾倒。当临空侧的岩柱倒塌开始后，剩余的岩柱也加速运动，最后整个危岩体解体失稳。整个破坏过程，有滑移、坠落、倾倒等状态，与坐落过程类似。图8展示了裂隙网络的发展、危岩体的失稳过程。

FEM/DEM方法预演了箭穿洞危岩体可能的变形破坏过程。从上述分析来看，箭穿洞危岩体今后裂隙的主要形成机制为压致拉裂，这与定性的变形失稳模式分析一致。压致拉裂主要是基座岩体被压缩导致的拉张破坏。破坏的基座岩体开始向临空面运动，直接导致了危岩体开始垮落。值得注意的是，基座岩体虽然压裂了，并未全部挤出，只有部分基座岩体离开了原位。图8显示这部分运动离开的岩体约占全部基座岩体的50%。

5 FEM/DEM讨论

箭穿洞危岩体案例分析充分展示了FEM/DEM方法的独特之处，它通过明确地考虑断裂和碎裂过程来模拟从连续到不连续的转变过程。小变形和应力应变由有限元进行分析；显式时间推进方案被用来整合离散单元牛顿运动方程，因此可以完全动态地模拟破坏过程。FEM/DEM方法对帮助分析大小变形混合的斜坡破坏与运动十分有利，适用于岩质危岩体或滑坡的破坏分析。

图8 箭穿洞危岩体变形破坏全过程数值模拟图(S为时步)Fig.8 Failure process map of the Jianchuandong dangerous rockmass

同时，需要注意的是，当前的FEM/DEM分析框架仅考虑了自重作用，并没有考虑水及振动等因素。对水因素的叠加，在FEM/DEM框架中可采用增加外部压力和增减重度的方式进行；对振动因素的增加，可采用脉冲式加速度来进行。增加考虑水和振动工况后，将使得FEM/DEM的适用范围更为广泛，其分析结果将更接近工程实际。

6 结论与建议

(1)柱状危岩体的基座处于应力集中状态，既有拉应力集中区也有压应力集中区，应力集中造成柱状危岩体的基座附近边界岩体易破坏。

(2)柱状危岩体裂隙的延展大部分是从下而上发展的，也有一部分是从上向下发展的，裂隙的延展导致了危岩体的解体和最终的破坏。

(3)FEM/DEM方法预演了岩体压致拉裂的裂隙网络形成和破碎岩体运动过程。数值模拟验证了箭穿洞危岩体会发生基座压裂型崩塌这一定性分析。基于该危岩体未来发展的认识，建议对该柱状危岩体基座进行补强防护。

(4)FEM/DEM耦合了有限元法和离散元法的优点，有利于分析岩体破裂及运动过程，是未来数值模拟的发展方向之一。