Montel-Carlo模拟法在滑坡风险评价中的应用

(国核电力规划设计研究院有限公司，北京 100095)

针对山区边坡稳定性问题，常用的评价方法是基于详细的岩土体物理力学参数采用传统的力学公式进行计算，但在大范围的山区区域获得准确的岩土体力学计算参数所花费的代价往往非常昂贵。一般来说，在城市规划、环境评估、工程选址等工程建设前期阶段，需要对场地稳定性做出评价，分析不良地质作用的危害程度及其工程处理或避让方案，然而在该阶段所获得的岩土体的计算参数是非常有限的。因此，如何在有限的条件下选择合理的评价方法准确评价近场区的滑坡危害问题显得至关重要。

目前，工程中有很多定性或定量手段评价边坡的稳定性，但这些理论方法所得到正确结果主要基于正确的地质模型及准确的岩土体力学参数。然而，事实是我们通常采用有限数量的室内试验或原位测试方法所取得的岩土物理力学参数并不能完全反映大面积岩土体的工程特性，地质参数是呈非线性分布的随机变量。本文研究内容是将岩土体地质参数作为变量，利用蒙特－卡洛方法模拟参数的随机分布情况，基于GIS软件进行滑坡风险评价。

近年来，已开发的GIS软件功能日益强大，GIS技术在环境评价和地质灾害预测中发挥了越来越重要的作用，在边坡稳定性评价及滑坡预测中的一些适用方法也得到提高。另一方面，边坡稳定性分析结果因地质参数不准确造成的误差要远大于因地质模型不准确所产生的误差，这主要因为地质参数的随机性及缺乏满足符合实际的地质模型。因此，在现实中我们可以选择相对简单的模型利用GIS进行滑坡风险性评价。

本方法研究过程中，仍然采用在GIS滑坡危险性评价中使用最为广泛的无限斜坡模型，该模型具有较合理的GIS栅格特征。同时在研究过程中引入邻域影响权重矩阵，并进行邻域分析研究。本文研究采用了Arc gis 软件及Visual Basic 6.0软件编程进行计算。

1 Montel-Carlo模拟方法在边坡稳定性分析中的应用

边坡稳定性系数Fs等于滑面抗滑力与下滑力之比，由于岩土体地质参数具有随机性，所以Fs很难准确反映斜坡的安全性，另一方面，选取不同的地质模型也可能会得到不同的结果，同时针对不同的建筑物安全等级和地质条件复杂程度选取的安全性评价的标准也是不一致的。因此我们可以给出滑坡概率Pks来间接说明适用于反映岩土体随机地质参数获得的稳定性系数Fs的概率。

1.1 Montel－Carlo模拟方法

影响边坡安全的因子是具有一定随机分布的变量，显然边坡安全系数Fs也一个随机变量，其概率分布情况是由地质参数分布决定的。假定我们已知安全系数Fs的分布函数，那么滑坡概率Pks公式即可写成下式：

式中：Fs为边坡稳定性系数，Pks为边坡稳定性

系数小于设计安全系数的概率。

事实上，稳定性系数的分布函数f(Fs)是很难直接得到的，通常是利用Montel－Carlo模拟法计算失稳概率Pks。Montel－Carlo模拟法的原理是根据参数的随机分布情况得到参数的随机变量，然后代入求解稳定性系数Fs的计算公式中，经过多次计算得到f(Fs)的分布函数，并进一步求解得到失效概率Pks。

1.2 随机样本

根据在[0，1]区间均匀分布函数，我们能得到一定概率分布的随机变量样本。设Z为随机变量，其分布函数为F(u)：如果F(u)

为具有连续性且严格单调函数，则：

式(3)在[0，1]区间均匀分布。其中，X为在[0，1]区间均匀分布的随机变量。即：

式(4)是与分布函数F(u)相关的随机分布函数，式中F－1(X)为F(Z)的反函数，如果我们在[0，1]区间均匀分布的样本数量足够的大，即样本Xk={X1，X2，…}，根据公式：

即可求得样本分布函数F(Z)。

在计算过程中，主要考虑了四种随机变量分布：均匀分布、负指数分布、正态分布、对数正态分布。

在[0，1]区间均匀分布的随机变量采用LCG法：

1.3 安全系数计算

安全系数计算方法采用简单的无限斜坡计算公式，计算模型见图1。

图1 无限斜坡稳定性计算模型

无限斜坡稳定系系数计算公式为：

式中：c'为滑面有效粘聚力；φ'为滑面有效内摩擦角；β为坡面的倾角；γ为滑坡体重度；γW为水的重度；Z为滑块高度；m为地下水水头高度。

1.4 计算过程与相关问题讨论

Montel-Carlo模拟按照图2来完成。

图2 Montel-Carlo模拟方法流程图

Montel－Carlo模拟法循环次数由计算结果所要求的精度决定，根据从5～50000次的计算结果的发现，按照工程上对精度的要求，200次即可满足。

在本研究中我们选取岩土参数c，φ，γ和地下水深度hw作为随机变量进行Montel－Carlo方法模拟。一般来讲，地质参数的分布符合正态分布或对数正态分布，如果选取的c，φ，γ和地下水深度hw为正态分布或对数正态分布，则安全系数Fs计算结果的分布与参数一致。

计算发现在其他地质参数保持不变的情况下，边坡的失稳概率随坡角的增加呈指数增长，见图3。

图3 边坡坡角与斜坡滑坡概率的关系

通过对参数的敏感性研究，主要针对参数变化对概率值的影响。由图4中可知，对概率值影响最大的参数为地下水埋深hw，参数粘聚力c和土体重度γ影响较小。对于不同的工程，由于每个参数相对变化值，得到的敏感变化程度是不同的。

图4 岩土参数与斜坡滑坡概率的关系

参数偏差对滑坡概率的影响见图5，随着参数标准偏差的增大，在参数平均值条件下，安全系数Fs大于1和Fs小于1两种状况下失稳概率均趋向于50%。

图5 标准偏差与斜坡滑坡概率的关系

2 GIS在西南山区特高压输电线路选线滑坡灾害评价中的应用

雅安—武汉1000 kV特高压直流输电线路工程路径穿过我国西南山区，该区域地形起伏大，常年降雨量大，地质灾害频发，尤其是滑坡灾害严重危害了当地生产生活。该线路自西向东穿越整个山区，为避开滑坡灾害对沿线塔基的影响，工程选线阶段工作至关重要。由于线路很长，大范围的滑坡评价工作量非常大。因此，工程中采用了GIS方法对线路途经区域进行大面积滑坡调查评价。

2.1 基本参数设置和计算程序流程

本次计算时选取面积约为1720×1540(m)作为实例，选用的栅格模型大小为10×10(m)，利用GIS软件通过Arc/View功能可较容易取得区域内某一点的地质参数，根据已有点的工程地质条件通过差值取得滑坡体的大致厚度。

同时，导致自然边坡失稳的外界因素主要包括诸如降雨、地震、人类活动等因素，在本研究过程中，降雨因素按地下水埋深条件考虑，然后对边坡模型施加外部不利条件进行稳定性计算。

滑坡稳定性按无限边坡公式利用VB程序计算，基本计算流程见图6。

图6 编写程序流程图

2.2 计算结果

通过计算，我们获得了每个划分网格区域的安全系数及失稳概率。需要说明的是，在计算中将滑坡概率设置成稳定性系数小于1时的概率。图7所示为该区域经计算的滑坡概率分布示意图，表1是滑坡概率对应的安全等级。

表1 滑坡概率对应安全等级

由图7所示，我们得到了评价区域范围内危险地段及其失稳概率，如图中某些地段的滑坡概率超过了90%，说明该区域已经出现了大面积滑坡迹象。根据GIS软件中显示危险斜坡的分布情况，然后针对危险区域再进一步进行详细勘察评价，使现场工作更加具体，节省了大量的外业工作。

图7 蒙特－卡洛模拟法斜坡失稳概率

2.3 邻域关系分析

在使用GIS过程中，我们是将每一网格单元假定为一个无限斜坡体，但实际上单元间是相互联系影响的，通过邻域分析，利用权重矩阵分析网格单元间影响关系大小。计算见式(8)：

式中：wij为单元间影响关系的权重值。

每个计算单元的滑坡概率计算公式可整理成式(9)：

式中：Pks为每个网格单元的滑坡概率；n为网格数，其中n值得大小仅与网格划分大小及可能滑坡面积有关。

计算时，利用邻近网格单元间影响的衰减值a(a＜1)推断权值wij，例如假定n=5，相邻网格单元影响值得权重值矩阵见式(10)：

通过计算得知，改变a和n值的大小使滑坡的范围大小随之改变，但是滑坡位置及安全性等级并未改变，可见采用GIS软件计算，邻域间关系影响仅对滑坡范围有影响，滑坡的地点及危险程度基本没有改变。

3 结论

通过对基于GIS的Montel－Carl模拟方法计算区域内斜坡危险性概率的研究，得到如下结论：

(1)在工程前期阶段进行滑坡评价，由于计算参数的不确定性，比较好的评价方法是选用简单滑坡计算模型，基于参数的随机分布，可较容易得到滑坡失稳Fs＜[Fs]的概率。

(2)在大型山区滑坡灾害评价研究中，GIS发挥了越来越重要作用，尤其在工程前期阶段。利用其栅格单元数据，可很容易的计算出每个栅格单元的滑坡概率，并可以图示形式显示出来。利用无限斜坡模型得出了每个计算单元的安全系数及滑坡的概率分布。

(3)利用权重矩阵模型，分析了网格间关系对滑坡概率的影响，通过计算，利用邻域间影响权重矩阵计算公式可减小计算结果的随机性和不协调性。

综上所述，基于GIS的强大功能，可较容易的评价山区线路工程塔基附近的滑坡灾害危险性并能展示其分布情况，通过初步评价结果，便于有针对性的对危险地段进行详细勘察评价。