地震对金沙江金坪子滑坡的稳定性影响分析

王 涛 周 葵 袁乾博 肖诗荣

（1.三峡大学 土木与建筑学院, 湖北 宜昌 443002；2.国家电网湖北省电力有限公司 兴山县供电公司, 湖北宜昌 443700）

地震地质灾害因其强大的致灾力而被各界学者广泛关注.地震不仅造成直接的灾害损失,由其引发的次生灾害如滑坡、崩塌等地质灾害也具有较强的破坏性,在地震活跃的西南一带尤为突出,故探讨地震作用如何影响滑坡体的稳定性这一问题极为重要.对于地震作用下的滑坡变形问题,我国很多学者对其进行大量的研究并取得了一定成果.孙崇绍等[1]通过对1 500年以来的国内地震和450年的历史资料记录整理,得出相应的地震滑坡的分布和规模与区域地质的联系.丁彦慧等[2]通过对中国194例地震触发滑坡、塌方事例的讨论,确定了地震滑坡与地震强度、震源距离彼此的关系.徐文杰等[3]利用ABAQUS软件,从三维空间出发考虑由地震引起的肖家桥滑坡失稳的相关性分析,该种方法克服了二维分析方法中的操作误差,且与现实问题较符合.许冲等[4-5]采用GIS结合确定性因子(CF)对汶川地震12个诱发因素进行了详细分析,提出了各因子有助于地震滑坡滑动的数值范围.范昭平等[6]在充分考虑岩体材料动力特性和地震特性的基础上,利用动力时程法对地震动载荷作用下的滑坡动力特性进行了分析,选取最小平均稳定系数作为评估结果.杨德旭[7]以GIS为技术手段,分析地震作用对滑坡影响因素的影响,并建立二元Logistic回归模型确定了地震对研究区的影响范围.樊晓一等[8]分析地震滑坡的优势坡度及其范围,指出滑坡平均坡度和地震抛掷作用受断层与坡向位置决定.

综上所述,本文以地理位置属于地震活跃的西南区域的金坪子滑坡Ⅱ区滑坡体为研究对象,依据乌东德水电站的抗震设防烈度,运用动力时程分析法,通过有限元岩土软件MIDAS-GTS建立金坪子滑坡Ⅱ区地质模型分别分析滑坡体在地震作用下的位移、速度、加速度曲线变化,得出地震对滑坡体稳定的影响.

1 滑坡概况

1.1 滑坡区自然地质条件

金坪子滑坡位于金沙江下游右岸,距乌东德水电站坝址约900 m,属典型峡谷地貌.金坪子滑坡所在地区地层物质组成为冲洪积物(Qal+pl)、崩塌堆积体(Qcol)、坡积物(Qdl)及滑坡堆积物(Qdel)、下伏中元古界会理群落雪组(Pt21)、黑山组(Pt2hs),下层震旦系上统观音崖组(Z2g)、灯影组(Z2d)等.区内断裂构造发育,东西、南北向断层交错切割.从历史记载来看,距离乌东德最近的强级震感为1985年禄劝县6.3级地震,震中距离20 km；小江断裂带为地震较频繁区域,历史上曾经发生6～7级共10次,7～8级一次,至近期仍有小震活动.根据我国地震烈度区划图以及乌东德水电站工程安全措施,金坪子所处区域一带地震基本烈度为Ⅶ度,滑坡区地震最大加速度为0.1g.

1.2 滑坡研究区基本特征

整个滑坡体在平面形态上呈“花瓶”状,坡体前缘瓶口处为滑坡剪出口,高程约880 m,坡体后缘界限相邻于Ⅰ区前缘,分布高程约1 480 m,宽度约450 m.金坪子滑坡Ⅱ区主滑方向126°,原始地形坡度约30°,滑坡前缘紧邻金沙江右岸,剪出口位于前缘基岩出露陡坎处.滑床倾向同坡向一致,倾角23°,滑带厚约2～5 m,滑带内微层理清晰,擦痕现象明显,运动擦痕方向与滑坡体主滑方向基本一致；滑体平均厚度约60～80 m,体积约2 700×104m3；滑体按物质组成由上至下可分为3层,第一层为白云岩块石碎石土层,厚度约10～15 m；第二层为碎石碎屑夹土层,碎石碎屑组成成分主要为千枚岩,含少量白云岩,厚约20～40 m；似层状碎石,碎屑岩成分为千枚岩,厚度约10 m,金坪子滑坡Ⅱ区代表性剖面如图1所示.

图1 金坪子滑坡Ⅱ区代表性剖面详图

2 基于刚体极限平衡分析法的滑坡稳定性计算

运用极限平衡法计算天然状态下和地震状态下滑坡体的安全稳定性系数,初步判断滑坡体的稳定性.极限平衡法运用假定虚拟状态分析滑坡体稳定性,该方法虽忽略了土体自身应力应变关系,不能精确描述滑坡体稳定性,但此方法经多年实践,积累较多成功经验,结果偏向于安全.目前常用的极限平衡分析方法有Morgenstein-Price法、Janbu法、Bishop法以及Spencer法等[9].

根据金坪子滑坡Ⅱ区滑坡的地质情况,以图1所示的典型代表性工程地质剖面建立滑坡计算模型.模型长约2733 m,高约1423 m,滑坡物质主要由白云岩碎块石土、千枚岩碎屑夹土等构成.山体模型如图2所示,岩土参数见表1.

图2 滑坡剖面模型示意

表1 金坪子滑坡Ⅱ区滑坡体物理力学参数

运用4种极限平衡分析理论分别计算滑坡体在天然状态和地震状态下的安全系数.模拟采用乌东德水电站的Ⅶ度抗震设防烈度,加速度峰值为0.1g,故可在模型边界施加0.1g的震动加速度来模拟地震.计算结果见表2.

表2 安全系数计算结果

由表2可知,自然状态下经极限平衡分析得出的安全系数皆大于1.2,即处于稳定状态.当施加地震后,4种极限平衡分析法计算得到的安全系数均小于0.8,可见在加速度峰值0.1g的作用下滑坡体不稳定可能产生滑坡现象.

3 基于动力时程分析法的滑坡稳定性分析

为进一步研究地震对滑坡体的影响,运用有限元岩土软件MIDAS-GTS建立金坪子滑坡Ⅱ区地质模型,采用动力时程分析法分别分析滑坡体在加速度峰值为0.1g地震作用下的位移、速度、加速度曲线变化,研究地震对滑坡稳定性的影响.

3.1 计算原理

时程分析法的基本原理为根据材料的弹性、非弹性特性对其动力运动方程逐步迭代进行积分求解.在求解计算时,不仅需要考虑地震的频率、振幅和持时,同时也需要考虑滑坡体的地形地貌、地质构造条件,因此时程法可以进行边坡非线性分析,从而得到坡体内各点在任意时刻的位移、速度、加速度等变化,根据其变化影响,分析边坡的稳定性.

时程分析法的动力方程如下：

式中：[M]代表质量矩阵；[C]代表阻尼矩阵；[K]代表刚度 矩阵；F(t)代 表 动 力 荷 载.其 中 ：¨u(t)、˙u(t)、u(t)分别代表加速度、速度、相对位移.

在已知地震阻尼的前提下,运用时程分析法中的振型叠加法进行求解.该方法的使用条件是阻尼矩阵[C]可以由刚度矩阵[K]、质量矩阵[M]线性组合,表达式如下：

式中：α,β为瑞利(Rayleigh)系数.

阻尼是边坡动力分析的条件特性之一,在时程法分析中常用Rayleigh阻尼.瑞利阻尼的理论是利用两个振型不相同的阻尼α、β产生的等效阻尼比和等于另一个对应振型的阻尼比,从而求出阻尼系数α、β.求解阻尼系数具体公式如下：

式(3)～(6)中：ξ1、ξ2分别为振型相同的边坡阻尼比,ω1、ω2代表边坡的自振圆频率.通过上述方程的求解,就可得出瑞利阻尼系数α、β.

3.2 计算模型与边界条件

以金坪子滑坡Ⅱ区滑坡地质资料中的岩土体物理力学性质参数为依据,充分模拟该滑坡的实际特征,以图1所示的典型代表性工程地质剖面建立滑坡计算模型,如图3所示.

计算模型长约2733 m,高约1423 m,网格划分时采用三角形+四边形格式,整个模型剖面共有3 938个单元,8 756个节点.边坡模型底部为固定边界,同时施加X、Y方向的约束；模型左右两边均施加竖直方向的约束,固定Y方向.虽然底部边界由于波的反射作用在一定程度上会放大滑坡上部动力响应,但鉴于底部边界为基岩,放大作用影响较小,基本符合约束条件.

为了更好地研究滑坡体在地震作用下的变形情况,分别在滑坡体的后缘上部、下部,坡体前缘的上部、下部设置4个观察点A、B、C、D,以供后文进行相应的数据分析.相关参数见表1.

3.3 地震波的选取

由于地震波在边坡破坏动力响应中占主导地位的是横波,因此本文仅考虑横波,在MIDAS-GTS软件输入相应水平加速度,因此在上述模型的左右边界仅固定Y方向.选取的地震波为EL-Centro波水平方向分量,持续时间20 s,烈度Ⅶ度,最大加速度峰值为0.1g.加速度曲线如图4所示.

图4 加速度曲线图

3.4 位移分析

在地震发生后,坡体内各点不同时刻均开始发生了不同程度的震动位移,滑带处水平方向实时动态位移曲线如图5所示.图中各线代表滑带上不同高程的监测点,相对位移代表不同时刻加速度所产生的位移值.

图5 滑带X向不同时刻相对位移

由图5可以看出,滑带处的水平位移以坐标原点为竖直方向的对称轴左右对称分布,且水平位移振幅由滑坡体的前缘至后缘呈逐渐增大的趋势,坡体后缘处位移最大,最大值达0.20 m,坡体前缘较小.同时图形中部分时刻位移曲线成较大弧形,坡体前缘与后缘的位移运动方向相差较大,不在同一时步上.说明滑坡体受到地震作用时,坡体对地震荷载作用的响应存在一个高程放大效应作用.随着高程的增大,坡体受到地震波效应越强烈,震动幅度越大.虽然坡体前缘最早受地震波震动,但其幅度较小.同时由于地震波在坡体内部传播时,坡体前缘与后缘所受波的作用大小以及传播时间均存在时间差,从而导致坡体的变形位移前缘与后缘出现了不一致的情况,坡体变形曲线呈弧形.

由计算可得4个观察点(A、B、C、D)的水平位移变化,如图6所示.

图6 观察点X向相对位移

分析图6可以发现,滑坡体前缘上下部C、D点的变形位移小于滑坡体后缘上部记录点A、B的位移.前缘变形最大水平相对位移值0.15 m,后缘最大累积变形相对位移值达到0.26 m,后缘变形是前缘的2.3倍.进一步分析发现,滑坡体表部观察点A、C变形水平相对位移值远大于滑坡体下部观察点B、D的水平相对位移值,滑坡体表部与下部水平位移最大变形皆位于坡体后缘的A、B点处,上部变形是下部的2.5倍.滑坡体上4个观察点的变形位移曲线时程运动基本一致,都随地震持续时间的增加,滑坡体变形位移沿着滑动方向呈递增的趋势.说明由于滑坡体表面为自由边界,坡体表部受到的两侧间接约束较小,而坡体内部所受到的间接约束较大,同时由于坡体表部处于临空状态,更利于地震力的释放,所以导致当滑坡体受到水平方向的地震加速度时,坡体表部在动力作用下其变形相对水平位移远大于坡体内部.滑坡体后缘观察点A、B的变形水平相对位移均大于坡体前缘的观察点C、D,综合判断是因为高程的放大效应作用.同时由于发生地震时,坡体内部的地震力会产生相应惯性力,其惯性力在坡体内沿着坡体前缘与后缘之间做循环运动,从而不断打破了岩土体的原始平衡状态,导致岩土体之间的黏聚力呈逐渐减小的状态.所以随着地震持续时间的增加,滑坡体的变形累积位移逐渐增大.

图7为滑坡体整体变形水平方向相对位移云图,由云图可知坡体后缘的变形最为剧烈,最大累积水平变形位移达26 cm,前缘变形较小.这与上述观察点的变化情况一致,同时也说明滑坡体后缘已经开始出现较大变形裂缝,坡体内部塑性变形累积增大,此时坡体可能处于临界破坏状态.

图7 滑坡整体X向相对位移

3.5 速度时程分析

滑坡体前后缘速度变化反映了滑坡体内能量的变化.滑坡体速度变化时程曲线如图8所示.由图8可知,滑坡体后缘速度变化值最大,坡体前缘速度变化小于后缘,说明由于地震中高程的放大效应,导致随着高程的增大,坡体的速度放大效应也逐步增大.由图可知滑坡体后缘上部A点的速度变化最大值为0.7 m/s,下部B点最大速度0.4 m/s,坡体前缘上部C点最大变化速度0.45 m/s,下部D点最大变化速度0.2 m/s,对比其相应数值可知,滑坡体表部速度变化大于坡体内部,表部坡体速度最大值是下部的1.75,后缘是前缘的3.5倍.

图8 滑坡体水平速度变化时程曲线

由上述分析可知,由于坡体各区域速度变化不一致,导致其岩土体之间会产生彼此的不协调运动,滑坡体上岩土体各层物质在互相摩擦期间,会不断减小岩土体物质的粗糙度,使其逐渐碾碎趋向于圆滑,从而降低了岩土体之间的粘结力.另一方面,随着地震持续时间的增加,岩土体速度变化的增大,使其原先处于静止状态的滑动面慢慢会变成动摩擦状态.在物理界中,动摩擦系数远小于静摩擦系数.从而导致岩土体在运动中可能会出现小型崩塌,一旦滑坡体的滑动面逐步趋向动摩擦状态,此时可能原先稳定的岩土体已变成不稳定状态.

图9为滑坡体整体变形水平方向速度变化云图,由云图可知变形运动速度最大的位置集中于滑坡体后缘岩土体,最大水平运动速度达0.9 m/s,前缘岩土体运动速度较小.这与上述观察点的变化情况一致,同时也说明滑坡体的岩土体在地震作用下,速度变化随高程放大效应而增大.

图9 滑坡整体X向速度云图

3.6 加速度分析

滑坡体各观察点水平加速度变化时程曲线如图10所示.

对比各点加速度曲线发现,随着高程的增加,滑坡体整体岩土体的加速度成放大效应.加速度最大值坡体后缘均大于坡体前缘,滑坡体后缘上部A点最大水平加速度为2.5g,后缘下部B点最大加速度达1.5g,滑坡体前缘上部C点最大水平加速度为0.9g,下部最大为0.7g.分析比较数值发现,滑坡体后缘水平加速度是前缘的2.8倍,坡体表部水平加速度是坡体内部的3.6倍.说明滑坡体表部受地震波水平加速度的影响远大于坡体内部.分析认为：由于滑坡体的岩土体受原始地震水平加速度影响,导致岩土体内加速度随着地震持续时间的增长,其加速度变化也在增大.由牛顿第二定律可知,当坡体内的加速度逐渐变大时,岩土体内的惯性力也会逐渐增大；同时由于加速度变化呈现循环往复的规律,导致坡体内惯性力也呈现循环变化状态,致使岩土体在惯性力作用下,不断在坡体后缘与前缘之间循环变化,致使岩土体出现弧形破坏模式.持续的震动,会使弧形破坏的深度与长度加大,间接为地下水的渗入创造了条件.当地下水在震动的过程中,渗入滑坡体的岩土体内时会产生相应的超孔隙水压力,导致岩土体的有效应力降低,土体的抗剪强度逐渐减小.

图10 滑坡体各观察点水平加速度变化时程曲线

综合分析发现,加速度变化值在6 s处振幅最大,而当滑坡体受到如此大的地震波幅度时,其相应速度、位移也在此时刻附近取得最大速度和位移,虽然有较小的时差,但变化基本与加速度一致.随着后期地震震动的持续,滑坡体的加速度、速度、位移也在不断变化着,但变化值均趋向于逐渐减小.说明6 s时刻是地震的转折点,此时刻滑坡体已经形成较大位移,但还未形成永久变形,随着震动的持续,直至最后振幅趋向于平稳时,滑坡体的塑性变形已经形成,地震最大永久位移达26 cm,滑坡体已经基本处于不稳定状态.

通过数值计算,滑坡体整体水平方向加速度变化如图11所示.由云图可知滑坡体惯性力最大的位置集中于滑坡体后缘岩土体,最大水平运动加速度达2.8 m/s2,滑坡体前缘岩土体运动加速度较小,因此其惯性力也较小.这与上述监测点的变化情况一致,同时也说明滑坡体的岩土体在地震作用下,加速度变化也随高程放大效应而增大,坡体内的惯性力由坡体后缘至前缘做循环变化运动.

图11 滑坡整体X向加速度云图

4 治理建议及讨论

4.1 治理建议

经分析表明金坪子滑坡Ⅱ区坡体在地震作用下极不稳定,一旦坡体发生失稳,会直接对乌东德坝区产生危害,故有必要对其采取相应措施.综合金坪子滑坡Ⅱ区滑坡体地形地貌、物质组成以及坝区安全、施工操作方便、经济合理等方面考虑,建议采取抗滑桩+削坡反压+地表排水措施治理方案,提高滑坡体的稳定性.

4.2 讨论

1)两种分析方法比较

时程分析法可进行边坡非线性分析,得到坡体内各点在任意时刻的位移、速度、加速度等变化,从微观角度对滑坡的变形进行精确分析,但从宏观角度看,该方法在实际工程中运用较少,对滑坡稳定性的计算没有极限平衡法直观,且极限平衡法在实际工程中运用娴熟.但极限平衡法多用于计算滑坡的安全系数,宏观判定滑坡是否稳定,无法对滑坡进行微观角度的精准分析.如何结合上述两种方法分析滑坡稳定性是值得探讨的理论与实践问题.

2)研究结论的可靠性及其检验

通过结合时程分析法与极限平衡法对滑坡的对比分析可知,地震烈度Ⅶ度地震作用下,滑坡最大变形达0.26 m,坡体安全系数在0.7～0.8之间,似乎由此可以得出结论,地震可导致金坪子Ⅱ区滑坡体变形失稳.由于分析方法的多种条件限制、概化和假定,此研究结论的可靠性需要在滑坡运行实践中检验和修正.

5 结 语

本文探讨了地震作用对金坪子滑坡Ⅱ区稳定性的影响,利用MIDAS-GTS进行数值模拟探讨在加速度峰值为0.1g的地震作用下,滑坡体的位移、速度、加速度变化情况,从而得到滑坡体在地震作用下的永久变形.根据分析,得到如下结论：

1)计算表明,天然状态下,4种极限平衡分析方法计算得到的安全系数均大于1.2,在加速度峰值为0.1g的地震作用下安全系数急速下降,均小于0.8,证明坡体在遭遇地震时可能发生失稳.

2)坡体受到地震作用时,坡体对地震荷载的响应存在一个高程放大效应,坡体变形的位移、速度、加速度随高程的增加而增大,滑体后缘变形大于前缘,根据分析计算可知,坡体变形加速度振幅值在6 s时刻最大,随之开始衰减；滑坡体最终地震永久变形位移达26 cm.

3)滑坡体在地震波水平加速度影响下,因其加速度、速度、位移的振幅呈循环往复变化,坡体内的惯性力也会由坡体后缘至前缘呈不协调往返运动,不仅促使岩土体之间的运动呈现不协调规律、降低岩土体物质的粗糙度,导致产生弧形破坏,同时也会在运动过程中为地下水渗入提供路径.从而产生孔隙水压力,土体强度减小,导致滑坡的产生.