基于Newmark法的设定地震滑坡危险性评估

林高聪, 潘书华, 叶振南

(中国地质调查局水文地质环境地质调查中心, 河北 保定 071051)

0 引 言

地震是触发滑坡发生的重要因素, 地震滑坡滑移速度快、滑移远的特点使之难以做到快速预警或预报[1-2], 国内外重大地震事件都引发了大量的滑坡灾害[3-10], 造成巨大的损失和隐患。因此, 开展地震滑坡危险性研究已成为国内外滑坡风险评价的热点之一, 尤其是近年来在全球地震活动进入相对活跃的背景下[11], 关注地震诱发区域性滑坡地质灾害的研究更具有显著的社会效应和科学意义。

关于地震滑坡危险性评价的方法可分为基于数学统计模型和确定力学模型的评价方法, 以震后滑坡发育规律与发震构造和地质环境指标影响因子间的关系研究最为广泛, 常采用基于数学方法的统计模型进行分析评价, 有: 二元统计法, 如频率比法[12-13]; 多元统计方法, 如逻辑回归模型[14-17]; 自学习识别模型, 如神经网络模型[12, 18-19]、支持向量机模型[18, 20-22]等; 还有与主观决策相结合的评价方法, 如层次分析法[24]、模糊模型等。此外, 在滑坡易发性或危险性评价[24]中常用的评价方法如信息量法、证据权法、确定系数法等, 也常应用于震后同震滑坡危险性评价中。

基于数理和统计模型评价的前提是指标要素彼此间是相互独立条件, 其侧重对现有滑坡与其影响因素间的宏观调查与统计, 获得其统计规律, 并不侧重滑坡形成机理的严格数理意义表达。评价的结果高度依赖于对已有滑坡的调查编录程度, 同时评价指标体系的结果也局限应用于地理地质环境条件相似区域。

采用力学模型基于地震强度或地面加速度对斜坡的作用而进行滑坡危险性评价, 其中以Newmark模型应用最为广泛, 如王涛等[25]、杨志华等[26]采用基于Newmark模型的震后同震滑坡危险性快速评估研究等。这些方法适用发震构造清楚或地震动强度和参数记录比较完整的情况, 注重的是震后对地震滑坡影响要素分析和危险性反演评估。

对于潜在地震区的滑坡危险性预测评估, 有研究者采用概率地震危险性分析方法和地震滑坡危险性评价模型结合的方法对潜在地震区的地震滑坡危险性进行预测评估。如陈晓利等[27]采用概率地震与综合指标相结合的方法, 对汶川地区潜在地震滑坡进行了评估; 刘甲美等[28]提出采用概率地震危险性分析及Newmark累积位移评估模型对天水地区进行潜在地震滑坡危险性区划。但是概率地震危险性评估方法对衰减关系的不确定性校正会混淆对低概率地震危险性估计的正确性[29], 对于局地尺度、有大震发生的情况, 需要加强关于设定地震的考虑。

以往滑坡危险性评价中, 鲜有注重对已有滑坡体单元的评价, 同时针对局地尺度、有潜在大震发生的情况, 基于设定地震的滑坡危险性评估较少。本文提出采用历史地震重现的方法, 组合多次历史地震对场地影响的结果，评价其地震危险性, 结合Newmark模型, 在缺少地震动参数记录条件下对地震滑坡危险性进行评价。以往区域滑坡危险性评价多采用栅格数据作为评价单元, 本文同时采用滑坡单元及栅格单元, 对已有滑坡复活和斜坡失稳的问题进行了评价, 解决了调查滑坡分布与危险性评估脱节的问题。同时，充分考虑工程地质因素, 区分不同类型滑坡的滑体厚度、力学强度参数等, 精细了评价结果的精度，其成果可为潜在震区地震诱发滑坡的危险性评估研究提供参考, 及为地质灾害防灾减灾管理提供支持和依据。

1 基于Newmark位移模型的地震滑坡危险性评估原理

Newmark模型的理论基础是极限平衡理论, 认为滑块的永久变形是在地震荷载作用下块体沿着滑动面发生破坏并不断累积位移所致, 当施加于滑动面处的加速度超过其极限平衡状态下的临界加速度时, 块体即沿破坏面发生滑动, 对荷载加速度a(t)与临界加速度ac的差值部分进行二次积分即得到累积位移值Dn。

在实际情况中由于缺乏完整的地震动加速度记录曲线, 无法有效获得加速度的二次积分。有研究者通过搜集海量地震记录[30-35], 针对任一记录可以计算得到地震强度IA、相应的临界加速度ac和累计位移量Dn, 最后拟合得到位移经验式, 用于区域评估。Dn通过ac、Arias强度IA以及地面峰值加速度PGA(即amax)之间的函数关系进行计算[31-32]

lgDn=1.521lgIA-1.993lgac-1.546;

(1)

lgDn=0.067lgIA-3.719ac/amax+0.852±0.365。

(2)

根据累计位移量与斜坡失稳之间的概率函数关系[32], 对地震条件下的滑坡危险性进行评价, 主要步骤包括(图1): ① 依据岩土体c、φ值和地形参数，计算斜坡的静态稳定系数Fs; ② 求解极限平衡状态下临界加速度ac; ③IA强度与PGA地震动参数的场地分布计算; ④ 斜坡累积位移量Dn与滑坡失稳概率计算, 分析评价滑坡危险性。

图1 基于Newmark位移模型的地震滑坡危险性评价流程

2 研究区概况及滑坡灾害发育特征

研究区位于甘肃省天水市秦州区境内, 温带半湿润气候区, 地理坐标105°19′34″E—105°30′05″E、34°30′N—34°40′N。

2.1 地质环境

研究区地处陇中黄土高原与西秦岭山地过渡地带, 水系发达, 呈现山间河谷盆地地貌和黄土梁峁地貌。地层岩性以第四系上更新世黄土、砂黏土和新近系红色泥岩及黏土岩为主, 占全区面积77%以上, 其余地段出露古近系粗粒碎屑岩、含砾泥岩和泥盆系上统粗粒碎屑岩及下古生界牛头河群变质岩。第四系黄土、砂黏土与新近系泥岩、黏土岩为该区域易滑岩组, 该类岩组发育的滑坡占全部滑坡总数的90%以上。

2.2 地质构造与新构造运动

研究区位于祁吕贺兰山造山带前弧, 西秦岭北缘构造带中段和陇西旋卷构造带的交汇地带, 地质构造发育强烈, 西秦岭北缘断裂从研究区域穿过(图2), 并发育多条走向NW、NWW的断裂和褶皱, 且多为逆冲和走滑断裂, 倾角大于70°。新构造运动强烈, 具有较强的、持续性活动的变形特点。

图2 研究区地震与活动断裂分布图

研究区内有两条地震带穿过, 即天水-兰州-河西走廊地震带和武都-天水地震带, 曾发生143年甘谷西7级地震、1604年礼县6级地震、1654年天水72级地震、1885年天水南6级地震等历史强震[36]。

根据中国地震台网中心记录的数据, 自1970—2015年, 研究区中心周围300 km内发生4级以上地震592次, 6～7级地震10次, 7级以上地震1次; 研究区周围150 km内4级以上地震57次, 6级以上地震1次。根据南北地震带中段5级以上地震活动记录显示, 目前该区域处于自1879年开始的活跃期, 未来有发生地震的可能。

2.3 滑坡灾害发育特征

结合地质资料、DOM影像(0.25 m分辨率)和DEM数据(2 m分辨率), 开展野外调查及滑坡地质灾害填图工作, 研究区230 km2范围内共计发育滑坡447处(图3), 据此建立了区域滑坡数据库。

图3 研究区滑坡编目图

依据滑坡滑体物质成分及滑面发育位置的不同, 将研究区内滑坡分为黄土型、泥岩型、砂岩及砂砾岩型、板岩型、变质岩类和侵入岩类6类, 依据调查数据, 对上述6类滑坡滑体厚度统计获得不同物质成分类型滑体厚度的平均值, 详见表1。

表1 不同物质成分类型的滑体平均厚度

滑坡在空间分布上沿河流水系方向呈带状分布, 遍及整个区域的河谷及支沟两岸, 其中基岩型滑坡主要发育在新近系泥岩中, 滑面位于深部的强弱风化界面或基岩内部中, 发育多个滑面, 体积规模上可达108m3, 其形成可能与地震有关, 这些已发生滑坡的部位受到地震的作用会更加容易破坏, 为危险性评价的重点地段。

3 研究区地震危险性分析

以往研究多采用概率地震危险性评估地震危险性, 但是该方法对衰减关系的不确定性校正会混淆对低概率地震危险性估计的正确性, 对于局地尺度、有大震发生的情况, 需要加强关于设定地震的考虑, 本文采用历史地震法评价潜在震区的地震危险性, 假定发生过强震的位置在未来时间可能再次发生相同强度的地震。选取以研究区中心300 km范围以内震级大于4级的地震数据作为评价地震危险性的初始筛选数据, 搜集到公元前780年—2015年4月期间的663次地震数据。

将不同震级单位间的地震数据源转化为同一震级, 采用有关对不同震级间的相互转化经验公式[37-38]，统一将地震标为面波震级MS

(4)

其中：ML为体波震级;MW为矩震级。

地震对场地造成的破坏与震中距呈衰减关系, 沿发震断裂的平行方向破坏程度大, 垂直于断裂方向地震动衰减快, 采用适合于甘肃省境内的长轴方向地震烈度衰减的经验公式[39], 计算和对比了663次地震在研究区中心造成的场地烈度值(图4), 其中有13处历史地震能够造成研究区场地烈度大于Ⅵ, 满足引发区域群发地震滑坡的条件。

图4 300 km范围内历史地震在研究区中心的烈度值图

I=4.846+1.464M-1.783ln(Ra+22)，

(5)

式中:I为场地地震烈度;Ra为场地与震源的距离。

筛选出3处对研究区影响程度最大的历史地震, 作为评价研究区地震危险性的重现地震事件, 分别是来自NNW方向的143年甘谷西7级地震、E方向的734年天水7级地震和SE方向的1654年天水72级地震, W方向地震引起研究区场地烈度均小于Ⅵ, 利用地震烈度与Arias强度之间的经验公式, 分别计算上述3处地震在研究区内Arias强度IA。

lgIA=0.75M-2lgRa-2.35。

(6)

前人在用一个设定地震综合表达对研究区或场地周围多个地震产生的影响时, 曾选取对研究区或场址产生大于或等于给定参数值的全部可能地震的震级、震中距的期望值作为设定地震的震级和震中距[40]。本文在对研究场地引起地震滑坡的有效地震中, 考虑了地震作用方向及其沿距离衰减特征, 期望值并不能表达出场地未来遭遇地震影响的严重程度, 因此借助ArcGIS平台对比筛选出3次地震作用下场内部各点遭遇到最强地震动参数值, 采取全部可能地震引起的最大值作为研究场址未来地震作用下的最危险性情况, 将其组合综合表达为设定地震危险性Arias强度评价结果(图5)。

图5 研究区Arias强度筛选流程图

4 研究区滑坡危险性分析

4.1 参数赋值

前人鲜有开展采用滑坡为单元的滑坡危险性评估, 致使在评估环节无法充分利用调查成果。为此, 本文考虑地震作用下已有滑坡的复活情况, 使二者有效衔接, 可以兼顾评估滑坡复活和斜坡失稳的可能性。根据研究区地层岩性和岩土工程地质特性, 在调查和滑坡填图的基础上对研究区滑坡与斜坡分别进行了岩组区划, 将已有滑坡视为一种特殊的工程地质岩组, 按照滑体的物质成分和结构特征分为黄土类滑坡、泥岩滑坡、砂岩及砂砾岩滑坡、板岩滑坡、变质岩滑坡和侵入岩滑坡岩组, 对于斜坡体按照地层岩性和结构等各划分为单一岩组(图6)。

图6 研究区岩组分布图

考虑斜坡失稳前后岩土体力学强度的变化, 将岩组的天然强度与残余强度[41-42](表2)分别赋予斜坡体岩组和滑坡体岩组, 并对其进行区划。

表2 研究区岩土体力学强度参数

对于滑坡岩组, 依据调查结果, 将每个滑坡滑体的深度作为其在地震作用下复活的深度; 对于研究区内不同岩性的斜坡体, 依据滑坡发育特征, 取同种岩性下滑体的平均厚度作为斜坡在地震作用下失稳时的深度。

4.2 稳定性计算与易发性评价

采用Newmark滑块理论对整个区域所有滑坡和斜坡体的稳定性系数(Fs)进行计算。

(7)

式中:h代表滑块的竖向深度(m);h′代表滑块被水浸入的厚度与滑块厚度比值;γ代表岩土体的容重；γw为水的容重(kN/m3);c为黏聚力(kPa);φ为内摩擦角(°);α表示潜在滑动面与水平方向的夹角(°)。

根据稳定性系数Fs和临界加速度ac之间的关系函数

ac=(Fs-1)gsinα，

(8)

计算临界加速度值, 其值反映的物理意义是坡体抵御外部荷载阻止斜坡失稳的能力, 临界加速度值越大说明滑(斜)坡失稳需要的能量和地震动越大, 临界加速度值越小则越容易被地震触发。因此，临界加速度值可作为评价其滑坡易发性的指标。

在野外调查过程中发现，研究区斜坡体浅层及已有滑坡体岩土体内部含水量很低, 未见坡体中有地下水出露, 与前人“黄土梁峁地区的缓斜坡地段不含地下水, 土质湿度不大[1]”的结论相符, 因此在稳定性计算中暂不考虑地下水对斜坡的影响, 仅考虑岩土体的黏聚力和内摩擦角对安全系数的影响。

计算结果显示(图7), 处于欠稳定和高易发的区域有两部分: 第1种是部分现有滑坡地段处于非稳定状态, 静态稳定性系数小于1的有200余处; 第2种是藉河北岸分水岭以北渭河一级支流和金河南岸一级支流内沟谷两岸坡体处于非稳定状态和高易发区, 尤其是坡体中下部斜坡的稳定性最差, 这与在野外调查中观察到的中小型滑坡集中发育在河流支沟岸坡下部的情况相符。在研究区的北部地区, 藉河北岸分水岭以北的渭河一级支流地区, 有呈纵带状集中不稳定的区域, 这些区域主要是支沟的两岸, 属于坡体的中下部, 被切割深度较大, 并堆积着大量的黄土, 滑坡灾害发育较密集。

图7 研究区斜坡临界加速度值分布图

较稳定和低易发区主要分布支沟的分水岭地区, 属于坡体的上部, 坡度相对小一些, 坡体的表层黄土已被剥蚀, 出露的为古近纪砂岩和含砾砂岩, 胶结好、致密坚硬, 不易发生失稳, 同时部分已有地貌平而宽缓的滑坡, 评价结果显示其处于低易发区。

4.3 潜在地震诱发滑坡危险性评价

使用上述地震危险性评价的结果和累积位移量间拟合的经验关系(式(2))，计算获得潜在地震作用下非稳定坡体的位移量。对研究区累积位移量的分析比较发现，地震作用下大部分坡体累积位移量小于2 cm；同时，计算结果也显示部分已有滑坡在地震作用下会出现较明显位移的现象，表明在地震作用下，研究区内会出现老滑坡复活的可能(图8)；并且在山脊线、梯田陡坎、冲沟边界附近出现较高的累积位移量，且呈条带状分布，这些部位在野外调查中也是常见的斜坡表层变形部位，尽管在空间上分布不连续，但在地震作用下也可能会演化为斜坡陡坎的局部变形破坏、小规模滑塌。

图8 设定地震作用下典型滑坡单元位移及地形转折处的条带状位移示例

Newmark模型计算得到的潜在滑块累积位移量与斜坡失稳并不存在直接的关系,累积位移量较大的地区说明在地震作用下斜坡发生失稳的可能较大, Jibson提出采用Weibull概率密度函数，拟合得到滑坡永久位移量与发生概率之间的经验公式[33]

式中：P为斜坡累积位移下失稳的概率；k、a、b为拟合函数中的常数。由式(4)和累积位移量Dn的分布图可以计算获得研究区内潜在地震诱发滑坡发生的概率分布图，即滑坡危险性区划图(图9)。分区结果显示：河岸地带为滑坡危险性最高的区域，并且支沟上游危险性明显高于支沟的沟口位置，表明该区域的滑坡地质灾害除受地震诱发外还受地形地貌控制；研究区北部渭河一级支沟黄土丘陵区同为研究区的滑坡高风险区，其原因为黄土类岩组为该区域的滑坡易发岩组，地震作用加剧了其向灾害发展的程度。

图9 研究区滑坡危险性分布图

5 结 论

通过采用Newmark模型结合历史地震法对天水潜在震区地震诱发滑坡的危险性进行评价，采用滑坡单元，建立调查和评估的继承关系，强化了评估模型中对灾害地质因素的考虑，对天水藉河上游地区地震作用下的滑(斜)坡位移进行了计算，研究了滑坡复活和斜坡失稳等问题，结果表明基于设定地震的Newmark位移法能够反映研究区的滑坡危险性。

计算结果和滑坡危险性区划显示，在藉河及其支流河岸地带为滑坡的高危险区，支沟上游危险性高于支沟的沟口位置，北部渭河一级支沟黄土区为研究区的连片高危险区，体现了该区域河谷侵蚀、溯源侵蚀造成的梁峁地貌与黄土作为滑坡灾害易发岩组对该区域地质灾害空间分布的控制作用。

对滑坡单元稳定性和危险性计算显示，有半数滑坡处于欠稳定状态，在地震作用下，现有稳定状态的滑坡出现了明显的位移现象，表明该区域已有滑坡存在地震诱发复活的可能性；同时部分已有地貌平而宽缓的老滑坡区，评价结果显示其处于较低易发区，在设定地震作用下也未出现较大的位移现象，显示出了较好的稳定状态。该结果可为天水地区地震诱发滑坡的危险性评估研究和地质灾害防灾减灾管理提供参考。

与此同时，笔者也注意到计算评价过程中采用的经验公式源于不同地区、不同时期，具有一定的局限性，尤其是基于累积位移量失稳的经验公式，缺乏不同岩性组合、岩土体强度等条件下建立的失稳准则。此外，对于Newmark模型平面滑动对弧形滑面滑坡模拟的局限性，也是今后需要加强研究的地方。