Simulation of Landslide Movement Process by Discontinuous Deformation Analysis
-
摘要: 坡是一个动态过程, 滑坡体的运动是一个集滑动、转动、拉张等运动方式的复杂运动过程.传统的极限平衡计算和有限元分析均无法描述滑坡的运动学特点和运动过程.非连续变形分析(DDA) 是最近发展起来的一种新的离散数值分析方法, 该方法基于块体的运动学理论及数值分析, 可以开展块体的静力和动力学计算.应用非连续变形分析方法对长江三峡库区新滩滑坡的运动全过程进行了数值模拟研究, 模拟方案充分依据该滑坡的地质、地形特征, 按不同岩土体和地质结构面类型进行块体单元的划分, 共划分成5 0 4个块体单元.模拟结果表明, 新滩滑坡是以斜坡中部姜家坡一带的局部破坏为其运动的开始阶段, 并进一步牵引上部滑体和推动下部滑体.代表性块体单元的位移变化曲线和滑动速度变化曲线反映了滑动过程中滑坡体块体系统的变形是非连续的, 各处块体的运动形态各异.从而很好地再现了新滩滑坡的整个动态过程, 揭示了滑坡的运动机制.Abstract: Landslide is a dynamic process. The movement of landslide body is a complex dynamic process composed of sliding, rotational, and tensile movements. Either the traditional limit equilibrium calculation or finite element analysis cannot be used to describe the characteristics and process of landslide movements. The discontinuous deformation analysis (DDA), a new discrete numerical simulation method, incorporates both the massive dynamic theory and the numerical analysis, so as to make a static and dynamic calculation of the massive body. In this paper, the discontinuous deformation analysis is applied to the numerical simulation of the dynamic movement process of Xintan landslide in the Three Gorges reservoir. This simulation classifies the landslide profile as 504 units in terms of geology, topography, soil and rock body types and geological structural surface types. The simulation results show that the Xintan landslide started with the partial failure of the Jiangjiapo belt in the middle of the landslide, and then the upper part of the landslide was further trailed upward and the lower part of the landslide was further pushed downward. Both the displacement variation curves and the sliding speed variation curves of the representative massive units show that the deformation of the massive body system of the sliding slope is discontinuous in the sliding process, and that the movements of different masses are quite different in shape, reflecting well the whole dynamic process of the landslide and revealing the dynamic mechanism of the landslide.
-
Key words:
- landslide /
- simulation /
- discontinuous deformation analysis
-
滑坡问题是典型的不连续变形块体的动力学问题, 目前的研究大多集中在边坡的稳定性分析方面, 而对滑坡的动力学过程研究还很少.传统方法中, 极限平衡法是进行边坡稳定性评价分析的主要方法, 但该方法通常需要预先假定破坏模式和滑动面, 而且不考虑岩土体的变形, 因此很难计算滑坡渐进破坏产生的真实位移场; 基于连续性假设的有限元方法应用于滑坡稳定性分析时, 虽然不必预先假定滑动面, 但它不能求解滑坡滑动实际可能产生的大变形和大位移; 基于牛顿运动定律的离散单元法是一种非连续分析方法, 该方法可以充分考虑离散块体之间的滑动、转动和张开等运动形态, 能得到大变形、大位移解, 但离散单元法在确定阻尼和能量耗散方面还存在不易解决的问题.
非连续变形分析(discontinuous deformation analysis, 简称DDA) 是最近发展起来的一种新的离散数值分析方法[1].它可作静力、动力计算, 包括正分析和反分析.DDA方法有完备的运动学理论及数值实现、全一阶位移近似、严格的平衡假定、正确的能量消耗和较高的计算效率.由于其刚体位移和块体变形采用统一的有限元格式求解, 不仅允许块体本身有位移和变形, 而且还允许块体间有滑动、转动、张开等运动形式, 能得到大变形、大位移解, 因此, DDA法兼有有限元和离散元的优点, 又克服了它们的缺点, 具有广泛的应用前景.非连续变形方法自提出以来, 便受到国内外广泛的关注和深入的研究[2~7].本文将这一方法应用到新滩滑坡动力学过程的数值仿真模拟研究, 有利于进一步揭示滑坡产生的运动机制.
1. DDA方法的基本原理
1.1 块体位移及变形
DDA分时步进行计算, 大位移和大变形是由小位移、小变形累加而成.假定每个时间步满足小位移、小变形条件.设每个块体处具有常应力、常应变, 块体的运动及变形由6个独立的变形参数确定:
(1) 块体中任意点(x, y) 的位移可由变形变量[Di]表示:
(2) 式中:
(3) 此处, (u0, v0) 指块体质心(x0, y0) 的刚体平移; γ0指块体绕质心(x0, y0) 的转角; (εx, εy, γxy) 指块体的正应变和剪应变; [Ti]为块体位移转换矩阵.
由于线性位移函数(2) 的使用, 当块体发生大的刚体转动时容易产生体积膨胀[6, 7], 为了消除这个误差, 每个时步计算完成后所有块体位移利用式(4) 重新修正:
(4) 1.2 总体平衡方程
通过块体之间的约束和作用在各块体上的位移约束条件, 把若干个单独的块体连接起来构成一个块体系统.假设所定义的块体系统有n个块体, 联立平衡方程具有如下形式:
(5) 因为每个块体有6个自由度, 所以式(5) 的系数矩阵中, 每个元素kij是6×6阶子矩阵, Di和Fi是6×1阶子矩阵, Di代表块体i的变形变量, Fi是块体i分配到6个变形变量的荷载.子矩阵[kii]与块体i的材料特性有关, [kij] (i≠j) 则由块体i与块体j之间的接触和约束条件决定.
由给定的常应变位移模式, 根据变分原理, 总体平衡方程可以通过对应力及外力(包括惯性力和阻尼力) 作用下的总势能Πp求导得出, 即总势能最小:
(6) (7) 1.3 块体系统运动学
DDA方法以天然存在的不连续面, 如节理、裂隙等切割岩体形成块体系统, 块体单元可凸可凹, 系统中块体之间的接触类型有4种: 凸角与凸角的接触、凸角与凹角的接触、凸角与棱边的接触及棱边与棱边的接触.块体之间的接触界面是位移间断面, 呈理想的刚塑性关系.块体系统运动变形时, 接触界面间需满足不嵌入和无张拉条件.该条件的满足是通过在各接触位置加上或去掉相应的刚硬弹簧, 即开-闭迭代来实现的[1].
2. 滑坡动态过程的仿真分析
滑坡的失稳破坏, 是一个动态过程, 存在岩土体的断裂、滑动、转动及剪胀等复杂过程, 具有宏观上的不连续性和单个块体运动的随机性.为了研究滑坡失稳破坏的动力学过程, 我们选取较为典型的新滩滑坡作为研究对象.
2.1 地质建模和网格剖分
新滩滑坡位于三峡大坝坝址上游26 km处的长江北岸.姜家坡—新滩镇堆积斜坡的形成主要是岸坡巨型岩崩发生后, 大量崩塌碎石在重力作用下, 不断沿斜坡向下堆积的结果, 故而斜坡堆积物具有倾斜层理的构造.由于岩崩作用的不断发展和古滑动面的起伏, 在姜家坡堆积斜坡前缘形成坡角达50°~60°的斜坡, 该段斜坡稳定性相对最差(图 1).
选择滑坡主滑剖面作为DDA计算剖面, 其长度为1 945 m, 高度为834 m.由于斜坡堆积物具有倾斜层理的构造, 根据地质资料, 在计算中考虑4种不同类型的层理结构面, 并作为划分块体单元的主要依据.块体尺寸大小和疏密程度则根据斜坡岩土体的结构, 引入竖直的假想不连续面剖分确定.在坡体内部, 基岩的完整性较好, 属于滑床部分, 块体划分得较大.DDA计算网格共划分块体单元504个.计算域约束条件为: 地面自由, 左边界为水平方向不动, 而底部边界为垂直方向不动.
2.2 计算参数的选取
(1) 材料参数和荷载条件.滑坡块体系统主要包括:
古滑坡体、古滑动面、滑坡体内部的次级滑动面和滑动面以下的志留系泥岩, 各材料的主要性质概括见表 1.新滩滑坡发生时, 地下水的作用并不明显, 且滑坡体中的地下水较少, 故作用荷载为自重体积力.计算中由位移速度变化引起的惯性力, 在计算过程中自动形成.
表 1 计算参数的取值Table Supplementary Table Parameters for calculation(2) 计算参量.DDA计算中需要输入以下参量:
计算时步、时间步长和最大位移比.计算时步规定计算过程分几个时间增量步进行, 本文计算取1 000时步, 时间步长为0.001 s.最大位移比是一个无量纲数, 用来判别接触的可能性, 本文计算时的取值为0.001.
2.3 DDA计算成果分析
不连续变形分析方法的特点之一是变形的不连续, 另一特点是在计算中引入了时间因数, 考虑变形有一时间过程, 因此它是动态的.只有在小位移情况且岩体很稳定时才可忽略时间因数, 用静力方法计算.根据以上所建立的力学模型, 采用DDA方法对新滩滑坡发展演化的全过程进行了模拟研究, 图 2给出了滑坡体滑动的全过程, 其中图 2a, 2b, 2c分别为200时步、500时步和1 000时步滑坡体的瞬时状态, 即块体系统的变形情况.从图 2a中可以看出, 在整体性滑坡发生的前期阶段, 主要是新滩姜家坡至广家崖危险体变形急剧, 后缘大幅度下陷拉裂, 前缘出现崩滑, 并且在姜家坡坡脚处发生明显的剪出鼓胀现象, 这与新滩滑坡的位移观测资料和实地观察结果是非常吻合的.当姜家坡发生滑动后, 整个滑坡体的后缘出现下滑, 顶部块体与基岩面逐渐拉开, 块体之间有较大的相互错动, 表现为张开、滑动和倾倒破坏等运动形式(图 2b), 新滩滑坡出现整体性滑移.滑坡最后阶段的快速活动过程中块体之间发生相互碰撞, 能量和速度不断交换, 运动形式极其复杂, 新滩滑坡前缘最终滑入长江(图 2c).
图 3, 4分别表示了滑坡体滑动时, 后缘、前缘及中间姜家坡地带的3个代表性块体单元的滑动速度变化曲线和位移变化曲线.从图中可以看出, 滑坡体块体系统在滑动过程中变形是非连续的, 因此各处块体的运动形态各异.由于姜家坡堆积斜坡前缘坡角较大, 块体崩滑后具有高势能特征, 故而此处块体位移和速度随时间变化较为急剧.在整体性滑坡的初期, 滑坡体后缘块体的速度和位移都大于前缘块体, 说明姜家坡崩滑后, 后缘块体变形加剧, 再加上崩塌块石的冲击和加载堆积, 促使斜坡发生整体性滑动, 故滑坡机理上具有推动式的趋势.
3. 结论
本文应用非连续变形分析方法对新滩滑坡的运动全过程进行了数值模拟研究, 模拟方案充分依据该滑坡的地质、地形特征.模拟结果表明, 新滩滑坡是以斜坡中部姜家坡一带的局部破坏为其运动的开始阶段, 并进一步牵引上部滑体和推动下部滑体.代表性块体单元的位移变化曲线及滑动速度变化曲线反映了滑动过程中滑动坡体块体系统的变形是非连续的, 各处块体的运动形态各异.对新滩滑坡运动全过程所进行的数值仿真模拟研究, 不仅再现了该滑坡大变形的整个动态过程, 还有助于对该滑坡机理进行深入的描述.
-
表 1 计算参数的取值
Table 1. Parameters for calculation
-
[1] Shi G H. Discontinuous deformation analysis — a new numerical model for the statics and dynamics of block system[D]. Berkeley: University of Califo rnia, 1988. [2] 邬爱清, 任放, 董学晟. DDA数值模型及其在岩体工程上的初步应用研究[J]. 岩石力学与工程学报, 1997, 16 (5): 411-417. doi: 10.3321/j.issn:1000-6915.1997.05.003WU A Q, REN F, DONG X S. Numerical model of DDA and its application in rock engineering [J]. Chinese Journal of Rock Mechanics and Engineering, 1997, 16(5) : 411-417. doi: 10.3321/j.issn:1000-6915.1997.05.003 [3] 裴觉民, 石根华. 岩石滑坡体的块体动态稳定和非连续变形分析[J]. 水利学报, 1993, (3): 28-34. doi: 10.3321/j.issn:0559-9350.1993.03.004PEI J M, SHI G H. Dynamic stability of rock slide and discontinuous deformation analysis [J]. Shuili Xuebao, 1993, (3) : 28-34. doi: 10.3321/j.issn:0559-9350.1993.03.004 [4] 殷坤龙, 姜清辉, 汪洋. 新滩滑坡运动全过程的非连续变形分析与仿真模拟[J]. 岩石力学与工程学报, 2002, 21 (7): 959-962. doi: 10.3321/j.issn:1000-6915.2002.07.005YIN K L, JIANG Q H, WANG Y. Numerical simulation on the movement process of Xintan landslide by DDA method [J]. Chinese Journal of Rock Mechanics and Engineering, 2002, 21(7) : 959-962. doi: 10.3321/j.issn:1000-6915.2002.07.005 [5] 姜清辉, 丰定祥. 三维非连续变形分析方法中角-面接触模型的研究[J]. 岩石力学与工程学报, 2000, 19 (增刊): 930-935. https://www.cnki.com.cn/Article/CJFDTOTAL-YSLX2000S1023.htmJIANG Q H, FENG D X. A study on angle-plane model in 3-dimensional discontinuous deformation analysis [J]. Chinese Journal of Rock Mechanics and Engineering, 2000, 19 (Suppl) : 930-935. https://www.cnki.com.cn/Article/CJFDTOTAL-YSLX2000S1023.htm [6] Koo C Y, Chern J C. Modification of the DDA method for rigid block problems [J]. Int J Rock Mech Min Sci, 1998, 35(6) : 683-693. doi: 10.1016/S0148-9062(97)00319-7 [7] Te-Chih Ke. The issue of rigid-body rotation in DDA[A]. Proceedings of the first international forum on discontinuous deformation analysis (DDA) and simulations of discontinuous media [C]. Berkeley, CA: [s. n. ], 1996.318-325. -