基于离散元的普洒崩塌过程分析研究

熊绍真，史文兵，彭雄武，王 勇

（1.贵州大学资源与环境工程学院，贵州贵阳 550025；
2.贵州大学教育部喀斯特地质资源与环境重点实验室，贵州贵阳 550025；
3.江西省地质局水文地质大队，江西南昌 330000）

地质灾害是世界上公认的仅次于地震的第二大破坏性自然灾害，且灾害类型多样［1-6］，有崩塌、滑坡、泥石流、地面塌陷等。同时，人类工程活动尤其是地下采矿活动加剧了地质灾害问题发生，往往造成重大人员伤亡和财产损失，一直以来受到国内外学者的重点关注。由地下采矿诱发的地质灾害往往具有规模大、突发性强、成因复杂、破坏力强等特点。1903年加拿大煤层开采引起的Frank滑坡规模达3.0×107m3，死亡76人［7］。1966年英国发生Aberfan滑坡由地下开采诱发，造成144人死亡［8］。1980年湖北盐池河磷矿区崩塌山体总体积100×104m3，造成284人罹难［9］。2017年贵州纳雍普洒采动崩塌体积49.1×104m3，造成35人死亡［10］。

采动作用诱发大规模灾害问题层出不穷，且此类灾害发生具有不可预见性，归根究底还是对该类地质灾害的变形破坏过程和形成机理认识有限。因此，很有必要研究采动斜坡变形破坏过程及其机理，以便为此类灾害的监测预警与防治提供可靠的科学依据。数值模拟方法是用来研究地下开采作用下斜坡变形破坏的机理和变形过程的重要手段［11］。贵州普洒崩塌发生后，许多学者用数值模拟方法对其进行研究，也有了一些研究成果，如肖锐铧等［12］认为降雨和采动作用是普洒斜坡变形破坏的重要因素；
钟祖良等［13］运用UDEC离散元方法探讨采动作用下岩溶山体崩滑失稳变形破坏特征；
郑光等［14］根据现场调查、无人机航拍及地质资料分析，初步分析了崩塌发生的成因机理；
Zhu等［15］利用地质信号及动态模型（DAN3D）进行数值模拟，分析了崩塌的形成机理和动力学特征；
余逍逍［16］运用底摩擦物理模拟和数值模拟研究多煤层开采对斜坡变形破坏的影响，认为地下采动对崩塌破坏具有显著作用，提出普洒崩塌的变形破坏机制为采空塌陷-拉裂-蠕滑-剪断。学者们都关注到地下采动作用是诱发普洒崩塌的主控因素［17］，其机制及破坏过程的研究对此类灾害认识具有重要意义。尽管学者们对普洒崩塌的变形破坏过程比较深入的研究，但是关注点多集中在机制研究上，而对于地下开采顺序、范围、历史等采矿状况影响高陡斜坡变形破坏过程的研究较少；
普洒崩塌规模大、运动距离远，已有的研究对于运动过程中复杂的下垫面地形的数值模拟与实际差别较大。因此，非常有必要借助三维模拟手段研究普洒崩塌全过程，为此类地质灾害的识别、预警和预报提供指导。

文中以贵州纳雍县普洒崩塌为研究对象，在野外工程地质条件调查基础上，查明了煤层开采历史、开采顺序、开采范围等开采信息，运用离散元方法（3DEC）模拟了崩塌的破坏过程，分析了采动作用下覆岩变形响应作用，并与崩塌实际过程进行了对比，由此验证了离散元方法研究崩塌过程的可靠性，为研究此类崩塌防治提供理论依据及参考。

1.1 地形地貌

研究区位于云贵高原向黔中山原过渡带的贵州高原第二阶梯黔西山原，属乌蒙山系溶蚀侵蚀地貌。区

图1 三维示意图Fig.1 3D schematic

内总体地势南高北地，地势陡峭，高山耸立，山脉总体沿南西向延伸。研究区域（老鹰岩）三维示意图见图1，小老鹰岩山体最高标高为2 175 m，最低标高点位于斜坡下部普洒村的河沟处，标高为1 875 m，相对高差约300 m，坡面坡度为10°～25°，局部地段斜坡坡度可达55°（图2）。斜坡上部为陡崖，下部为缓坡，高差约200 m。下部煤系地层二叠系上统龙潭组含煤地层，露头部位地势平缓，一般标高在1 700～1 900 m之间，平均标高约1 850 m，大部分区域被第四系覆盖，主要为耕地及建房用地。研究区降水范围为1 270 mm左右，年均降水量为1 238.8 mm。

1.2 地层岩性

普洒崩塌位于张威背斜东南翼。区内发育F1、F2和F3三大断裂系。调查区域出露岩层由新至老分别为：第四系（Q），下三叠系夜郎组（T1y）、上二叠系长兴-大龙组（P3c+d）、龙潭组（P3l）。剖面图根据图2中1-1截面，见图3。

图2 研究区平面图Fig.2 The plane sketch about the study area

研究区内上部为三叠系下统夜郎组（T1y）青灰-灰褐色薄至中厚层状灰岩，发育三组节理面：（1）145°∠85°，（2）6°∠80°，（3）200°∠75°；
分布于斜坡顶部。中下部（T1y1）为三叠系紫红、灰色泥质粉砂岩，主要分布于斜坡中部。据调查，崩塌体后缘顶部有一层20～30 m的土黄色泥岩层出露，泥岩质地松散，手可挖动。岩性为二叠系上统长兴—大龙组（P3c+d）深灰、灰色泥质灰岩分布于洼地以上的陡坡。二叠系上统龙潭组（P3l）煤系地层出露于斜坡坡脚，斜坡坡度小于10°，岩性为泥质粉砂岩、泥岩等软岩以及煤层。第四系（Q）主要为粘土、砂质粘土，厚0～15 m。

1.3 采矿活动

矿区可采煤层分别为M6、M10、M14、M16、M18、M20层。崩塌前已开采3层煤层，分别为M16、M14以及M10，开采顺序为：M16-M14-M10。该矿开采历史主要分生产规模6万t/年和30万t/年2个阶段。采煤方式为走向长壁采煤法，炮采工艺，全部垮落法管理顶板。采空区位置如图2所示，红色表示M16煤层，蓝色表示M14煤层，紫色表示M10煤层。3层煤层有多个工作面，崩塌范围内M16有3个工作面，M14有3个工作面，M10有2个工作面。M16、M14、M10煤层倾角约为7°，顶板岩性为粉砂岩和粉砂质泥岩。M16平均层厚1.60 m，M14平均层厚1.23 m，M10平均层厚2.12 m。

2.1 崩塌破坏特征

崩塌发生后，源区岩体失稳，铲刮斜坡中部的裂隙化山体，而后进一步碎裂，最终以碎屑流的方式向300°～310°方向运动，堆积于斜坡坡脚及前缘平台。斜坡前缘桥边组及大树脚组部分居民房屋被碎屑流摧毁掩埋。堆积体形态受到斜坡前缘小山包阻挡呈不规则的扇形向两个方向堆积。崩塌陡壁后缘海拔约2 120 m，坡脚海拔约为1 922 m，相对高差约为200 m。崩塌方向的水平长度约为800 m左右，陡壁坡脚至堆积区前缘的水平长度约为660 m，堆积区平均厚度约为4m。崩塌源区方量约为49.1×104m3，铲刮方量约2.1×104m3，堆积区受铲刮的方量约为31×104m3，本次崩塌的总方量约为82.4×104m3。根据无人机航拍和现场调查，将普洒崩塌堆积区分为崩塌源区、下落铲刮区、流通停积区［13］，如图4（a）所示。

图4 分区特征图Fig.4 Partition feature picture

（1）崩塌源区

崩塌源区位于斜坡中上部，崩塌剪出口位于坡顶向下60～90 m处，剪出口附近的主要岩性是粉砂岩和泥灰岩。崩塌源区上部灰岩溶蚀发育，可见明显的凹岩腔。据无人机航拍及野外调查，崩塌源区崩落岩体平均高约85 m，宽约145 m，平均厚度约40 m，总方量约49.3×104m3。崩塌物质主要为三叠系下统夜郎组（T1y）灰岩、泥岩、粉砂岩夹泥灰岩，如图4（b）。

（2）下落铲刮区

崩塌源区失稳岩体突然崩落后，在巨大的重力势能及动能作用下沿着坡面不断向下运动，沿途铲刮坡表碎裂岩体及坡面堆积物，铲刮厚度约1.5 m，铲刮区宽约180 m，高80 m，则铲刮方量约2.1×104m3［13］，如图4（c）。

（3）流通停积区

崩塌堆积区位于斜坡前缘缓坡区域。主要为村寨和荒地，坡脚前缘约200 m处有一小山包。崩塌体失稳倾倒，经碰撞、解体、碎裂，最终以碎屑流的方式高速运动，摧毁了部分居民房屋，造成重大的人员伤亡。堆积区形态呈现“鱼尾形”，中部较宽，前后缘较窄，延伸方向为SW-NE方向，如图4（d）。

2.2 崩塌过程分析

崩塌体主要发生于2017年8月28日上午10点30分许，且崩塌发生前坡表不时有块石崩落。据斜坡失稳过程的视频资料显示，崩塌发生整个过程可以清晰的反映出来。崩塌发生后，作者收集了监测员用无人机拍摄了崩塌的全过程的视频。根据现场视频影像可将斜坡的变形破坏过程描述为5个阶段（图5）。

图5 崩塌过程图Fig.5 Collapse process picture

视频中清楚地呈现了崩塌前大约7 min的山体变形破坏现象。首先是崩滑灾害初始阶段，见图5（a）。在0′22″时，斜坡左侧后缘有石块挤出，并发生了小型坍塌。之后1 min内斜坡顶部及右侧有多次小崩落，并在2 min时斜坡中部有破碎块石坠落。而后是坡面局部崩塌阶段。2′24″，斜坡顶部有一小部分岩体塌落，并伴有声响。经过崩塌源区右侧坡面多次小崩塌后，崩塌源区中部岩体进一步的碎裂，进而发生小规模整体垮塌，见图5（b）、（c）。坡面由右侧小规模崩塌进一步加剧到向坡面中部挤压变形，说明老鹰岩山体在不断碎裂，应力在不断调整，开始进入大变形阶段，见图5（d）、（e）。进入整体失稳阶段，见图5（f），在4′26″，斜坡顶部再次发生大方量崩塌，而后老鹰岩山体变形迅速向斜面右侧、左侧、中部推进，并有进一步加剧的趋势，山体极不稳定。由于斜坡坡面变形加剧，崩塌源区两侧不断有崩落体滚落，撞击坡脚激起白烟。崩塌源区从南侧到斜坡中部区域发育有裂缝，这些裂缝相互连接，向下逐步扩展。同时，斜坡中部发育有许多垂直裂缝，逐渐在崩塌源区下部形成一条挤压破碎带。最终进入斜坡整体溃散阶段，见图5（g），坡体位移快速增加，应力发生系统性的转变，斜坡在重力作用下滑移-剪出，一泻而下，带有高势能的崩滑体迅速冲向坡脚。崩滑体并没有仅仅停留在坡脚，而是转化为碎屑流高速推进到坡脚居民区，摧毁了大量的农田和房屋。在6′21″的时候，该事件结束，仅残留直立的陡崖壁，见图5（h）。

3.1 数值模型的建立

离散元方法可以用来分析不连续反倾斜坡变形破坏过程，能够揭示斜坡演化过程和机制。文中结合已有地质资料及现场调查，建立符合实际的三维离散元数值模型（图6）。模型范围长550 m，宽350 m，高335 m。考虑到模型实际范围太大及计算的可操作性，文中计算范围为老鹰岩及坡脚采空区一带。岩层赋本构模型为M-C弹塑性模型，结构面本构模型为摩尔库伦模型。只考虑自重应力场，不考虑构造应力，取重力加速度为10 m·s-2，边界条件为固定前后（X方向）、两侧（Y方向）和底面（Z方向），而模型老鹰岩顶部边界和坡面为自由边界。模型中Y轴为正表示斜坡东侧，X轴为正表示斜坡内部，Z轴为正表示斜坡顶部。

图6 数值模型图Fig.6 Numerical model picture

3.2 参数确定

模型中岩体及结构面力学参数主要是根据现场采样进行室内试验所获取，同时也参考了已有的成果及地区工程经验综合获得［18-21］。各岩层力学参数都是室内试验及参考前人成果获得，结构面参数也是室内试验获得相应的力学参数，再做相应的折减获得，其中刚度参数是先通过结构面直剪试验得到的力学参数再参考经验公式拟合得到［22］，如表1～表2所示。

表1 岩体力学参数取值表Table 1 Rock mechanics parameter value

表2 结构面力学参数取值表Table 2 Mechanical parameter value of structural plane

3.3 模拟结果分析

天然情况下，斜坡上部灰岩节理裂隙发育，其中一组陡直节理与坡面倾向差不多，把临空破碎岩体和后缘岩体分隔开，在降雨、风化、温差、自重等作用下易形成贯通裂缝，进而斜坡失稳破坏。灰岩岩溶发育，有落水洞、溶蚀裂隙等，利于雨水渗入弱化节理强度。斜坡中部厚层泥质粉砂岩节理也较发育，可见明显多组节理。同时，斜坡陡峭，加上岩体被节理切割破碎，易成为崩塌堆积体的物质来源。节理裂隙附近容易形成应力集中区，加剧斜坡失稳破坏。天然情况下，斜坡各岩体、节理、层面赋予参数，再进行迭代计算，斜坡能够达到平衡，没有出现大的沉降，位移没有大的变化，证明斜坡天然情况下是稳定的，应力分布总体平衡，不会发生斜坡失稳现象。

根据3DEC数值模拟的结果，崩塌的失稳过程可以很好的反映出来，可以分为以下4个阶段：坡顶裂缝形成阶段、山体局部变形失稳阶段、斜坡整体失稳破坏阶段和崩滑堆积阶段。

（1）坡顶裂缝形成阶段

斜坡节理裂隙发育，多组结构面相互交错切割岩体，使得岩体呈碎块状，且斜坡三面临空，这些都是斜坡变形破坏的基本地质条件。但是，多煤层的开采是斜坡变形破坏的主控因素。由图7和图8可知，随着煤层的开采，上覆岩体不断被扰动，岩体内部应力进行重分布，坡体的应力状态也进行重新调整，采空区顶板岩层产生卸荷作用，在重力作用下岩层会出现轻微变形。随后多煤层开挖，采空区的扩大，斜坡内部应力再次重新分布，导致采空区顶板变形加剧，引起坡体表面以及坡顶出现多条裂缝，岩体趋向碎裂化，这与图7中现场斜坡变形情况比较吻合，坡顶形成了多条裂缝。坡顶出现明显的变形，最大位移达到了5 m，且坡顶位移大小不一，说明坡顶岩体产生了不同程度的错动，坡体开始出现不均匀沉降现象。

（2）山体局部变形失稳阶段

煤层开采完后，斜坡岩体应力会重新调整，煤层顶板因失去支撑向下下沉，逐渐形成采空塌陷区。这个阶段重力起着重要的作用，其对岩体的作用力可以理解为垂直方向和水平方向的合力。垂直方向分力垂直于岩层层面，使得上覆岩体产生不均匀沉陷变形，上覆岩体强度降低到小于重力垂直方向分力，导致裂缝产生，岩体趋向于破碎，从而整个坡体向下发生不均匀沉降，坡顶也发生强烈变形，出现错动现象。同时，重力的水平分力推动岩体向临空面方向运动，坡体裂缝进一步变宽边长，向下扩展，坡表岩体更加破碎。整个坡体不断的发生不均匀沉降，且临空面附近沉降量值大，岩块有向临空面倾倒趋势，往坡内沉降值相对减小。

图7 斜坡变形图Fig.7 Slope deformation picture

图8 斜坡变形整体位移图Fig.8 The overall displacement of slope deformation

图9 坡面变形图Fig.9 Slope deformation picture

图10 山体局部变形阶段整体位移图Fig.10 The overall displacement of the hollowing collapse phase

图11 斜坡失稳图Fig.11 Slope instability picture

图12 斜坡失稳整体位移图Fig.12 Slope instability overall displacement picture

据调查，2016年老鹰岩山体受到地下开采扰动，山体变形剧烈，从而山体右侧发生小规模崩塌，坡面也不时出现落石。由图9和图10所示，当迭代时步增加到27 000步时，坡顶变形更加剧烈，坡顶位移在2.5～7 m之间，不均匀沉降明显，出现了多条裂缝，岩体更加碎裂化。斜坡整体也发生了不均匀沉降，采空塌陷区影响范围不断向坡顶延伸，导致坡顶拉裂缝不断增多变宽，坡体整体开始向下发生不均匀的沉降，坡表岩体沿着拉裂缝不断向临空方向移动，易产生倾倒破坏，这与图9中斜坡实际变形范围比较一致。图10也显示坡表岩体产生位移，岩体趋向碎裂化，开始脱离坡表，沿着坡面掉落，形成小崩塌，说明坡表变形剧烈，这与崩塌前发生过小崩塌现象相吻合。

（3）斜坡整体失稳破坏阶段

采空区顶板塌陷后，其影响范围逐渐向坡顶方向发展，采空区被充满，坡体的不均匀沉降也在不断扩大直至结束。不均匀沉降影响范围逐渐扩大，山体局部变形进一步加剧，使得向外倾倒的趋势愈来愈明显，坡顶及坡面间岩体多组节理裂隙与层面相互作用发生错动，导致节理裂隙张开度变宽并趋向碎裂化。当运行至52 000步，图11和图12显示破碎岩体与后缘稳定岩体被裂缝分割开，岩体向临空面发生倾倒破坏，坡顶裂缝不断被拉开，拉开间距不断变大，渐渐与坡体内部节理贯通，使得岩体变形剧烈。同时，坡体岩体破碎，并开始产生大位移，最大可达94 m，岩块大规模溃散、松动，由表及里，向下滚落或坠落，这与崩塌视频中影像相一致，见图5（g）。当山体向临空面倾倒变形不断加剧时，在重力的作用下，山体失稳掉落，斜坡内部岩体也逐渐失稳，坡顶岩体也失去支撑掉落，斜坡发生变形破坏，形成大规模崩塌。崩落体在向下滚动的过程中相互碰撞，逐渐解体，堆积于坡表前缘。形成陡崖面。

（4）崩滑堆积阶段

山体变形剧烈，斜坡失稳溃散，上部岩体崩落，堆积于坡脚。由于斜坡高陡，崩塌体重力势能转化为巨大的动能，斜坡快速失稳滑动，倾滑而下，铲刮斜坡中下部的松散坡积物，向斜坡前缘高速运动。当迭代到350 000步时，图13和图14显示，崩塌块体间受到冲击和碰撞作用，崩塌体进一步解体、碎裂，而后以碎屑流方式高速运动，位移值也发生激增，最终呈现“马鞍形”堆积于斜坡前缘，中部较宽，前后缘较窄。由于坡脚缓坡处有一个小山包，崩塌体运动到此处出现分流情况，大部分崩塌体向两个方向呈不规则扇形堆积，少部分崩塌体出现爬坡现象堆积于小山包上。由图13所示，数值模拟结果与现场调查图比较吻合。

图13 崩滑堆积图Fig.13 Collapse stacking picture

图14 堆积阶段位移图Fig.14 Stacking stage displacement picture

运用3DEC数值模拟分析普洒崩塌，可以很好的揭示崩塌的形成，再现崩塌发生整个过程。结合现场崩塌视频和数值模拟结果，可以发现3DEC模拟的崩塌过程与现场视频影像分析结果具有较好的一致性。同时，该类灾害具有一定的特征：

（1）高陡的地形为斜坡提供良好的临空条件。老鹰岩山体高陡，坡度为上陡下缓，三面临空，易产生向临空面倾倒破坏。

（2）地层岩性为上硬下软型地层组合。普洒斜坡上部为灰岩，属于硬岩；
下部为泥岩、泥质粉砂岩、煤层，属于软岩，具有典型的上硬下软“二元结构”。下部软岩的物理力学性质较低，而上部硬岩力学强度高，容易导致上覆岩体挤压下部岩体产生变形破坏。

（3）岩体发育有三组陡倾结构面，切割山体。岩体内发育3组节理面：第1组结构面将后缘稳定岩体与源区失稳岩体切割开；
第2组结构面将南侧的小老鹰岩陡壁与崩塌源区失稳岩体切割开；
第3组结构面构成崩塌源区的原始临空面，3组结构面的切割，使斜坡内部岩体被切割成破碎的不规则的“三棱柱体”，构成块体的来源。且结构面陡倾，不利于岩体的稳定性。

（4）地下采空加速上覆岩层弯曲变形，岩层节理裂隙扩张，上覆岩层出现“悬臂效应”。地下采动作用下，斜坡顶部岩体逐渐向临空方向变形，使得裂缝沿结构面追踪发展，逐渐拉开，向下发展，形成贯通性的潜在滑面。

文中结合现场调查和3DEC数值模拟发现，地下开采是诱发普洒崩塌的主要因素。地下开采引起坡体内部的应力重分布，上覆岩体产生塑性变形，逐渐形成塌陷带，随着时间推移，坡顶出现明显变形，拉裂缝变大变宽，岩体愈加破碎，加剧岩体向临空方向倾倒破坏，斜坡失稳破坏。而斜坡高陡，上部岩体失稳溃散后倾滑而下快速堆积于下部缓坡处，崩塌体经过冲击作用、铲刮作用和重力荷载突然增大产生巨大的动能，以高速碎屑流的形式失稳滑坡堆积于斜坡前缘。因此，普洒崩塌的变形破坏过程可以分为4个阶段，即坡顶裂缝形成阶段、山体局部变形失稳阶段、斜坡整体失稳破坏阶段和崩滑堆积阶段。

此外，普洒斜坡变形迹象比较复杂，崩塌发生前发生过多次小型崩塌。地下开采影响斜坡变形和后缘裂缝发展过程都还是需要解决的问题，需要进一步重点研究。

（1）现场视频中清楚地呈现了崩塌发生前约7 min的山体变形破坏现象。根据现场视频影像可将崩塌过程分为5个阶段；
崩滑灾害初始阶段、坡面局部崩塌阶段、斜坡大变形阶段、斜坡失稳阶段和斜坡整体溃散阶段。

（2）地下开采是诱发普洒崩塌的主要因素。采用3DEC模拟斜坡在地下开采作用下的变形破坏后运动情况，清楚地展现了斜坡失稳溃散的整个过程。且后崩塌过程与现场视频影像以及现场调查结果具有较好的一致性。地下开采引起采空区顶板塌陷，导致斜坡出现不均匀沉降，坡顶变形剧烈，产生拉裂缝；
坡表岩体也愈发破碎，山体局部变形失稳，岩体开始脱离坡表，沿着坡面掉落；
而后山体变形剧烈，使得斜坡发生倾倒式破坏，斜坡失稳溃散，形成大规模崩塌。

（3）普洒崩塌的主要过程为：煤矿的开采过程中，斜坡受到扰动，加剧岩体产生变形，上覆岩体出现开裂塌陷现象，坡顶出现拉裂缝；
采空区形成后，整个坡体出现不均匀沉降，裂缝进一步向下扩展，岩体破碎趋向临空面倾倒，开始脱离坡表，形成小崩塌；
坡顶岩体变形剧烈，斜坡整体失稳，发生大规模崩塌；
崩塌块体间冲击-铲刮作用，崩塌体以高速碎屑流的方式向下滑动，并堆积于坡脚前缘。

（4）普洒崩塌堆积体呈“马鞍形”堆积于斜坡前缘。由于崩塌块体间冲击作用和铲刮作用，崩塌体产生巨大的动能，以碎屑流方式沿坡面高速向下运动；
在斜坡前缘遇到一个小山包，崩塌体就沿两侧方向滑动，堆积于两侧，具分流现象，堆积于山包两侧，只有少数崩塌体爬坡，堆积于小山包上。