蒋水华 ，祁小辉 ，曹子君 ，李典庆

（1.武汉大学 水资源与水电工程科学国家重点实验室，湖北 武汉 430072；2.武汉大学 水工岩石力学教育部重点实验室，湖北 武汉 430072；3.南昌大学建筑工程学院，江西 南昌 330031）

1 引 言

土体是经过复杂地质作用而形成的天然材料。受复杂地质成因的影响，土体物理力学性质随着空间位置的变化而不同，表现出一定的层状分布特征。土体参数的这种层状分布特征对边坡稳定的安全性和可靠性有着显著的影响，使得边坡可能存在多种潜在失效模式[1－3]，采用单一失效模式（比如临界确定性滑动面或者临界概率滑动面[3]）进行边坡可靠度分析会低估边坡失效概率，以致于错误地评价边坡安全性。为了准确地评估边坡的安全性和可靠性，亟需发展可以同时考虑多失效模式的边坡系统可靠度分析方法。

目前，国内外许多学者在边坡系统可靠度问题上已进行了大量有益的研究。如张兴等[4]采用蒙特卡洛模拟（Monte Carlo simulation,MCS）研究了土体参数间相关性对边坡系统失效概率的影响。谭晓慧[5]将边坡视为含多条潜在滑动面的串联系统，采用Ditlevsen 上、下限法[6]分析了边坡系统可靠度。Oka 等[1]以及Chowdhury 等[2]提出，在边坡可靠度分析中应当合理地考虑边坡多失效模式，并采用了Ditlevsen 上、下限法进行边坡系统可靠度分析，得出边坡系统失效概率的可能变化范围。Low 等[7]采用一阶可靠度方法（first order reliability method,FORM）结合Ditlevsen 上、下限法计算边坡系统失效概率。以上这些研究都是基于部分代表性滑动面对边坡进行系统可靠度分析。这些代表性滑动面没有经过严格筛选，可能会忽略一些重要的失效模式，导致低估边坡系统失效概率。此外，Ching 等[8]在考虑大量潜在滑动面的情况下采用重要抽样方法进行边坡系统可靠度分析。Zhang 等[9]采用2 阶多项式对大量潜在滑动面构造其安全系数与输入参数间的响应面，再结合MCS 进行边坡系统可靠度分析。然而，实际工程边坡失效概率一般在10-6～10-3量级，基于MCS和大量潜在滑动面的边坡系统可靠度分析方法的计算量非常大，使得其在工程实际问题中的实用性受到一定的约束。

最近一些研究[9－16]提出，首先在大量潜在滑动面中筛选出代表性滑动面，然后基于代表性滑动面进行边坡系统可靠度分析以提高计算效率。Zhang 等[9,11]提出，从大量潜在滑动面中筛选出代表性滑动面，再采用2 阶多项式或者克里金方法对每条代表性滑动面构造响应面，并结合MCS 进行边坡系统可靠度分析。Ji 等[10]通过类似于Hassan和Wolff 方法[3]确定代表性滑动面，针对每条代表性滑动面采用2阶多项式构建响应面，再采用FORM 进行边坡系统可靠度分析。Cho[12]通过设障碍方法（barrier）确定边坡代表性失效模式，再采用多点一阶可靠度方法（multi-point FORM）计算边坡系统可靠度。Zhang等[13]基于Hassan和Wolff 方法[3]确定代表性滑动面，采用2 阶多项式对每条代表性滑动面构造响应面，再结合MCS 计算边坡系统失效概率。Li 等[14]根据大量潜在滑动面与临界确定性滑动面间的相关性筛选代表性滑动面，然后也采用MCS 计算边坡系统失效概率。目前代表性滑动面的确定方法通常需要计算大量潜在滑动面之间的相关性，计算过程较为繁琐。此外，摩擦/黏性边坡安全系数通常是土体参数的非线性隐式函数，当采用相关非正态分布表征土体参数分布特征时，常用的2 阶多项式展开可能难以准确地近似潜在滑动面安全系数与土体参数间的高阶非线性隐式函数关系。

本文提出了一套基于随机响应面法的边坡系统可靠度分析方法。首先从大量潜在滑动面中有效地筛选出代表性滑动面。然后针对每条代表性滑动面，采用Hermite 多项式展开建立其安全系数与土体参数间的非线性显式函数关系(即随机响应面)，以提高安全系数的计算效率。在此基础上通过直接MCS计算边坡系统失效概率。本文首先介绍了提出的基于随机响应面法的边坡系统可靠度分析方法及其计算流程，最后通过两个算例验证了该方法的有效性，并探讨了不同阶响应面以及土体参数分布类型和相关性对边坡系统可靠度的影响。

2 边坡系统可靠度分析

对于多层边坡稳定性问题[9－13]，可能存在大量潜在滑动面，只要任何一条潜在滑动面安全系数小于1.0，边坡就会沿该滑动面滑动进而导致整个边坡失稳。因此，边坡安全性可定义为一串联系统可靠度问题[5,14]。根据串联系统失效概率的定义[9,14]，边坡系统失效概率pfs计算表达式为

式中：P(·)为某一事件或系统的失效概率；E[Sj]为边坡沿第j 条滑动面失稳事件，j=1,2,…,Ns；Ns为所有潜在滑动面数目。式（1）可进一步表示为

由于岩土工程中试验数据通常非常有限，难以获得较为完整的土体参数概率信息等，同时式（2）的被积分项是一个复杂的函数。此外，一般情况下土体参数服从相关非正态分布，每条潜在滑动面安全系数FS 是土体参数X 的非线性隐式函数关系，因此，对式（2）直接积分计算的难度较大。然而，式（2）可以通过MCS 方便地计算，根据MCS 式（2）进一步简化为

式中：N为MCS 中产生的样本数目。为了准确地计算pfs，在MCS 中通常需要产生大量的样本进行边坡稳定性分析，并基于大量潜在滑动面计算每组样本所对应的最小安全系数，计算量非常大。

3 随机响应面

为了提高计算效率，本文将在大量潜在滑动面中筛选出Nr条代表性滑动面，并对每条代表性滑动面采用Hermite 随机多项式展开（Hermite polynomial chaos expansion，HPCE）建立其安全系数FSj与随机土体参数X 间的显式函数关系（即第j 重随机响应面）[15－19]。然后，将构建的代表性滑动面的随机响应面用在MCS 中计算每组样本所对应的最小安全系数。对于第j 条代表性滑动面，其安全系数表达式为[15,17]

式中：j=1,2,…,Nr，Nr为代表性滑动面数目；n为随机变量数目；为待定系数；ξ=（ξi1,ξi2,…,ξin）为独立标准正态随机向量；是自由度为n 的Hermite 多项式展开[15,17]。然后根据独立标准正态随机向量ξ与土体参数X 间的等概率变换函数关系[17]为

式中：FXi(x)和 fXi(x)分别为输入随机变量Xi的累积分布函数和概率密度函数；Φ-1(·)为一维标准正态变量的累积分布函数的逆函数。

因此，根据式（5）、（6）便可得到第j 条代表性滑动面安全系数与土体参数X 间的显式函数关系。其中关键问题是确定式（5）中多项式展开的待定系数a：通过选择Np组独立标准正态随机样本ξ，基于样本ξ 进行Np次边坡稳定性分析，然后对于第j 条代表性滑动面就有Np个安全系数，再根据样本ξ和第j 条代表性滑动面所对应的Np个安全系数，建立线性代数方程组求解多项式展开待定系数a，进而建立第j 条代表性滑动面安全系数所对应的随机响应面。类似地可以构建Nr重随机响应面，最后给定一组土体参数X 的随机样本，根据式（5）、（6）就可以得到Nr条代表性滑动面的安全系数，将其中的最小值作为边坡稳定性分析的系统输出响应，进而通过式（4）MCS 计算边坡系统失效概率。

4 代表性滑动面确定方法

图1 边坡系统可靠度计算流程图Fig.1 Flow chart for system reliability analysis of slope

边坡系统可靠度分析中重要的一步是需要确定代表性滑动面。本文边坡代表性滑动面是在随机响应面构建过程中确定的，首先在独立标准正态空间进行概率配点或者拉丁超立方抽样，选取Np组独立标准正态样本点ξ，将这Np组样本点通过式（6）分别转换到原始空间中得到Np组X 的随机样本值。然后，将其中每一组随机样本值作为输入参数，通过边坡稳定性分析（如简化毕肖普法）计算得到Ns条潜在滑动面的安全系数，并确定其中最小安全系数对应的滑动面，即临界确定性滑动面。对于Np组X 的随机样本值，理论上可从中获得Np条临界确定性滑动面。然而由于基于某两组或者多组X 的随机样本值所搜索到的临界确定性滑动面可能相同，所以只可筛选出Nr条不同的临界确定性滑动面，可见Nr≤Np。本文将这Nr条不同的临界确定性滑动面作为代表性滑动面进行边坡系统可靠度分析。此外，其他方法[9－14]也可用于确定边坡代表性滑动面，然而当考虑参数空间变异性所需离散随机变量数目较多时，它们的计算量较大。与之相比，本文方法不需要另外专门计算潜在滑动面安全系数之间的相关系数，简化了计算过程，提高了计算效率。

5 计算流程

本文提出的边坡系统可靠度分析方法的计算流程图如图1 所示，主要步骤如下：

（1）取土体参数随机变量均值建立边坡稳定性分析模型。

（2）产生可覆盖整个边坡可能失稳区域的大量潜在滑动面，本文首先采用SLOPE/W 模块[20]剪入、剪出方法生成大量圆弧型滑动面，为保证计算精度要求，建议潜在滑动面总数Ns＞104较为合适。

（3）将边坡稳定性分析模型存为slope.xml 源文件，并对其进行边坡稳定性分析自动搜索临界确定性滑动面。

（4）选取独立标准正态空间概率配点或者拉丁超立方样本点ξNp×n（维度为Np×n，Np为随机样本数目），根据随机变量统计特性（均值、变异系数和分布类型等），通过式（6）等概率变换方法得到Np组原始空间中土体参数X 随机样本值，再用X 的样本值分别代替边坡源模型文件slope.xml 中对应随机变量的均值，从而生成Np个新的slope.xml 计算文件。

（5）借助批处理软件分别对Np个slope.xml 文件进行边坡稳定性分析批次计算，如Winbatch[21]软件。Winbatch 是一微软Windows 脚本语言，其运行环境由解释器、代码编辑器、对话框设计和一个用来创建可执行文件的编译器组成，其主要思想是通过Winbatch 执行文件*.wbt 控制SLOPE/W 软件自动启动、打开与求解计算文件以及自动关闭，给出了相应的Winbatch 软件与SLOPE/W 结合程序如下：

（6）Np次边坡稳定性分析可自动搜索出Nr条不同的临界确定性滑动面，将这Nr条不同的临界确定性滑动面作为代表性滑动面进行边坡系统可靠度分析，并分别从生成的Np个计算结果文件slope.fac中提取每条代表性滑动面对应的Np个安全系数。

（7）根据独立标准正态样本点ξNp×n和Np个安全系数基于式（5）建立线性方程组，求解多项式展开系数，进而建立边坡每条代表性滑动面安全系数与土体参数X 间的显式函数关系（即随机响应面）。

（8）最后根据式（4）～（6）采用直接MCS 计算边坡系统失效概率。

本文提出的边坡系统可靠度分析方法巧妙地实现了边坡确定性分析与概率分析的有机结合，而且充分地利用了现有商业软件。本文以GEOSTUDIO软件的SLOPE/W 模块为例，予以实现边坡系统可靠度分析，所提出的方法同样可以拓展到其他边坡稳定性分析软件。本文提出方法只需要对Nr条代表性滑动面进行边坡系统可靠度分析，计算效率较高。当边坡每条代表性滑动面安全系数均表示为土体参数的显式函数关系后，最终只需要将一些数学表达式视为功能函数进行直接MCS，不再需要进行原始的边坡稳定性分析，从而可以快速地获得边坡系统失效概率。

6 算 例

6.1 两层不排水饱和黏土边坡

首先以一个两层不排水饱和黏土边坡系统可靠度问题为例验证提出方法的有效性，Low 等[7]、Ching 等[8]、Ji 等[10]、Zhang 等[11]和Cho[12]均对该边坡稳定进行了可靠度分析，其相应计算结果可用来验证本文方法的有效性。边坡计算模型如图2 所示，坡高为24 m，坡度为0.75:1。根据以上文献，将两黏土层的不排水抗剪强度cu1和cu2视为随机变量，均值分别为120、160 kPa，变异系数均为0.3，并且都服从对数正态分布。土体重度γsat＝19 kN/m3。

首先取参数均值采用SLOPE/W 模块剪入、剪出方法[20]生成可覆盖整个边坡可能失稳区域的大量潜在滑动面，根据文献[9,13－14]可知，滑动面总数Ns＞104可以满足精度要求，边坡系统失效概率基本上不再随着Ns的增加而发生变化。本文为保证计算精度，共生成了10 164 条圆弧型滑动面，其位置如图2 所示。然后与文献[2,8,11－12]一样，取参数均值采用简化毕肖普法计算得到边坡安全系数为1.995，与Ji 等[10]采用普通条分法得到的1.997基本一致。显然自动搜索到的临界确定性滑动面包含在所有潜在滑动面中，如图2 所示。

图2 边坡计算模型和10 164 条潜在滑动面Fig.2 Slope geometry and 10 164 potential slip surfaces

图3 边坡3 条代表性滑动面Fig.3 Three representative slip surfaces of slope

获得边坡所有潜在滑动面之后，边坡系统可靠度分析的重要一步是从中筛选出代表性滑动面。以4 阶HPCE 方法为例，由于4 阶Hermite 多项式展开待定系数数目为15，如要获得多项式展开系数至少需进行15 次边坡稳定性分析。通过这15 次边坡稳定性分析可自动搜索到3 条不同的临界确定性滑动面作为代表性滑动面如图3 所示，这3 条代表性滑动面与文献[7,10－12]所确定的代表性滑动面位置非常吻合，其中包含了临界确定性滑动面。然后利用Hermite 随机多项式展开分别建立每条代表性滑动面的随机响应面（即安全系数与随机变量cu1和cu2之间的显式函数关系）。以图3 所示的临界确定性滑动面为例，首先采用4 阶HPCE 可拟合得到该滑动面安全系数FS 与独立标准正态变量ξ1和ξ2间的显式函数关系如下：

式中：ξ1=(lncu1-4.744)/0.294；ξ2=(lncu2-5.032)/0.294。据此便可得到临界确定性滑动面FS 与随机变量cu1和cu2间的显式函数关系。最后对式（7）采用20×104次直接蒙特卡洛模拟计算得到临界确定性滑动面的失效概率为1.43×10-3。基于以上3 条代表性滑动面的随机响应面，通过直接MCS 计算得到边坡系统失效概率为4.11×10-3，如表1 所示。可见多层边坡系统失效概率明显大于临界确定性滑动面的失效概率，仅采用单一临界确定性滑动面计算边坡失效概率会低估边坡的失稳风险。

表1 边坡可靠度结果的比较Table 1 Comparison of reliability results of the slope

为了进一步验证本文基于代表性滑动面进行边坡系统可靠度分析的有效性。采用4 阶Hermite多项式展开构建所有（10 164 条）潜在滑动面的随机响应面，在此基础上通过MCS 计算得到边坡系统失效概率为4.11×10-3。与基于代表性滑动面的计算结果完全吻合，说明了本文方法可有效地识别边坡的代表性滑动面，而且其计算效率明显高于基于所有潜在滑动面的计算方法。此外，根据基于所有潜在滑动面的边坡系统可靠度分析结果，可以确定最大失效概率对应的滑动面（即临界概率滑动面），如图3 所示。可见该临界概率滑动面也包含在本例所确定的3 条代表性滑动面中。与临界确定性滑动面类似，临界概率滑动面失效概率（2.65×10-3，见表1）也小于边坡系统失效概率（4.11×10-3）。

此外，表1 给出了已有文献不同方法计算得到的边坡系统失效概率及其采用的代表性滑动面数目。本文的计算结果（4.11×10-3）与Low 等[7]采用Ditlevsen 上、下限法（[4.32×10-3,4.41×10-3]）、Ching等[8]采用重要性抽样方法（4.1×10-3）、Ji 等[10]采用分层响应面法结合一阶可靠度方法（[4.02×10-3,4.11×10-3]）、Zhang 等[11]采用克里金响应面结合MCS 方法（4.58×10-3）以及Cho[12]采用多点一阶可靠度方法（4.36×10-3）得到的结果基本一致，这进一步说明本文所提方法的有效性。

为了探讨不同响应面方法对边坡系统可靠度的影响，本文也采用了不含交叉项的普通2 阶多项式展开[9]、2和3 阶HPCE 构建代表性滑动面的响应面，并计算边坡系统失效概率。如表1 所示，基于不含交叉项普通2 阶多项式的计算结果（4.11×10-3）与4 阶HPCE 方法的计算结果完全一致；3 阶HPCE方法的计算结果（5.73×10-3）与4 阶HPCE 方法的计算结果略有差别；2 阶HPCE 方法的计算结果（3.90×10-4）与4 阶HPCE 方法的计算结果有显著差别。此外，Zhang 等[11]采用2 阶经典响应面法得到边坡系统失效概率为1.7×10-3与本文4 阶HPCE方法的计算结果（4.11×10-3）也有明显差别。可见采用不含交叉项2 阶多项式和4 阶HPCE 方法可以得到相对较为准确的计算结果，然而2 阶经典响应面法、2 阶和3 阶HPCE 方法的计算精度不够。其原因在于对于不排水饱和黏土边坡，基于简化毕肖普法得到的安全系数FS 与不排水抗剪强度参数之间是线性函数关系[2]，Zhang 等[9]采用不含交叉项2阶多项式展开可以准确拟合FS 与不排水抗剪强度参数间的函数关系。然而，由于不排水抗剪强度参数服从对数正态分布，根据式（6）可知，不排水抗剪强度与独立标准正态变量间呈非线性指数关系（如cu1=exp(4.744+0.294ξ1)），故低阶（如2 阶）Hermite 多项式展开难以准确地拟合潜在滑动面安全系数与独立标准正态随机变量间的高阶非线性函数关系。如式（7）所示，除2 阶多项式展开之外，3 阶和4 阶多项式展开对准确地建立安全系数FS 与cu1和cu2间的函数关系起到一定的作用。本文提出方法可采用高阶多项式展开准确地近似边坡安全系数与土体参数间的显式函数关系，具有较好的计算精度。同时基于代表性滑动面进行边坡系统可靠度分析极大地提高了计算效率。确定代表性滑动面时无需单独计算滑动面间的相关性，简化了计算过程。

6.2 芝加哥国会街切坡

下面以芝加哥国会街切坡系统可靠度问题为例，进一步验证本文提出方法的有效性。Oka 等[1]、Chowdhury 等[2]、Ching 等[8]和Zhang 等[11]也对该边坡稳定性进行了可靠度分析。芝加哥国会街切坡计算模型如图4 所示，包含4个土层（顶部砂土层及其下部3个黏土层），坡高为13.8 m，上、下两层坡角分别为36.3°和36°。根据文献[2]，将3个黏土层土体黏聚力c和内摩擦角φ 视为随机变量，均服从独立正态分布，它们的统计特征如表2 所示；顶层砂土黏聚力c为0 kPa和内摩擦角φ为30°，均视为确定量。此外，各层土体重度γ 也视为确定量，取为18.5 kN/m3。同样基于SLOPE/W 剪入、剪出方法生成10 164条可覆盖整个边坡可能失稳区域的潜在滑动面，其位置如图4 所示。各层土体参数取均值采用简化毕肖普法得到边坡安全系数FS为1.484和相应的临界确定性滑动面如图4 所示。

图4 边坡计算模型和10 164 条潜在滑动面Fig.4 Slope geometry and 10 164 potential slip surfaces

表2 土体参数的统计特征Table 2 Statistics of soil parameters

本例采用2 阶HPCE 方法进行250 次边坡稳定性分析，得到35 条不同的临界确定性滑动面作为代表性滑动面，如图5 所示。针对每条代表性滑动面，采用2 阶HPCE 构建随机响应面，再通过100×104次MCS 计算得到边坡系统失效概率为3.41×10-2，如表3 所示。同样采用2 阶Hermite 多项式展开构建了所有（10 164 条）潜在滑动面的随机响应面，并通过MCS 计算得到边坡系统失效概率也为3.41×10-2，与基于代表性滑动面的计算结果完全吻合。再次说明了本文提出方法不仅能够有效地识别边坡的代表性滑动面，而且其计算效率明显高于基于所有潜在滑动面的计算方法。

如表3 所示，本文采用2 阶HPCE 方法得到的失效概率（3.41×10-2）高于Chowdhury 等[2]采用Ditlevsen 上、下限法得到的边坡系统失效概率([1.92×10-2,2.78×10-2])。其主要原因为Chowdhury等[2]在计算中只考虑了4 条代表性滑动面，却忽略了经过第1个黏土层的一些重要失效模式，从而低估了边坡失效概率。基于本文识别的35 条代表性滑动面，可用于准确地估计边坡系统失效概率。

图5 抗剪强度参数服从独立正态分布时边坡35 条代表性滑动面Fig.5 35 representative slip surfaces of slope with independent normal shear strength parameters

表3 抗剪强度参数服从独立正态分布时边坡可靠度结果的比较Table 3 Comparison of reliability results of the slope with independent normal shear strength parameters

为了进一步探讨不同响应面方法对边坡系统可靠度的影响，本文采用不含交叉项的普通2 阶多项式[9]、3 阶和4 阶HPCE 多项式构建代表性滑动面的响应面，并计算其系统失效概率，计算结果如表3 所示。采用不含交叉项的普通2 阶多项式、3 阶和4 阶HPCE 方法所得的边坡系统失效概率分别为3.59×10-2，3.47×10-2和3.68×10-2，与2 阶HPCE 方法的计算结果基本一致。采用不含交叉项的普通2阶多项式和2～4 阶HPCE 方法得到的边坡功能函数的累积分布函数（cumulative distribution function,CDF）非常接近，只在较小失效概率区域（如小于10-3）略有不同，如图6 所示。这说明2 阶多项式展开可以准确地拟合边坡安全系数和土体参数间的非线性隐式函数关系。尽管基于简化毕肖普法的摩擦/黏性土坡安全系数是抗剪强度参数c和φ 的非线性隐式函数关系[2]，但由于随机变量均服从正态分布，边坡安全系数与独立标准随机变量间的非线性程度较低，故2 阶多项式展开可以满足精度要求。然而，由于土体抗剪强度参数不能取为负值，因此，岩土工程可靠度分析时常用对数正态分布、Gamma 分布、Beta 分布或者截尾指数等非正态分布来表征抗剪强度参数的分布特征[22]。为了探讨土体参数分布类型对边坡系统失效概率的影响，进一步将3个黏土层土体黏聚力c和内摩擦角φ 视为相互独立的对数正态随机变量，并且保持均值和变异系数不变（见表2），再次采用本文提出方法进行系统可靠度分析。如采用3 阶HPCE方法进行300 次边坡稳定性分析，可自动搜索得到34 条不同的临界确定性滑动面作为代表性滑动面，如图7 所示。采用3 阶Hermite 多项式展开构建每条代表性滑动面的随机响应面，进而通过MCS 得到边坡系统失效概率为4.79×10-3，明显低于土体参数服从独立正态分布时的边坡失效概率（3.47×10-2，见表3）。

图6 抗剪强度参数服从独立正态分布时边坡功能函数概率分布函数曲线Fig.6 Cumulative distribution functions of performance function of slope with independent normal shear strength parameters

图7 抗剪强度参数服从独立对数正态分布时边坡34 条代表性滑动面Fig.7 34 representative slip surfaces of slope with independent lognormal shear strength parameters

图8 抗剪强度参数服从独立对数正态分布时边坡功能函数概率分布函数曲线Fig.8 Cumulative distribution functions(CDF)of performance function of slope with independent lognormal shear strength parameters

图8给出了采用不含交叉项的普通2阶多项式[9]、2～4 阶HPCE 方法得到的边坡功能函数的累积分布函数。在相对较大的失效概率区域（如大于10-2），4 种多项式展开方法的CDF 曲线十分接近。然而在相对较小的失效概率区域（如小于10-3），不同多项式展开方法得到的CDF 有明显差异。随着多项式展开阶数增加，2 阶和3 阶HPCE 方法得到的CDF 曲线逐渐接近于4 阶HPCE 方法的CDF 曲线。根据Hermite 随机多项式的自身收敛性判别方法[15－19]，可将基于1 000 次拉丁超立方抽样的4 阶HPCE 方法的CDF 曲线可视为精确解。可见对于该含非正态抗剪强度参数的摩擦/黏性边坡系统可靠度问题，2阶多项式展开的计算精度显然不够。

此外，土体抗剪强度参数c和φ 之间通常存在着负相关关系[23]。为了探讨抗剪强度参数间相关性对边坡可靠度的影响，计算过程中仍采用对数正态分布表征土体参数分布特征，并且均值和变异系数保持不变（见表2），另外考虑各黏土层的黏聚力和内摩擦角间的相关系数为-0.7。如采用3 阶HPCE方法进行300 次边坡稳定性分析，可自动搜索得到38 条不同的临界确定性滑动面作为代表性滑动面。采用3 阶Hermite 多项式展开构建每条代表性滑动面的随机响应面，进而通过MCS 得到边坡系统失效概率为8.16×10-4。这明显低于土体参数服从独立对数正态分布时的失效概率（4.79×10-3），可见抗剪强度参数间的相关性对边坡系统可靠度有重要影响，忽略c和φ 间的负相关性会高估边坡系统失效概率。

在考虑c和φ 负相关性的情况下，图9 也给出了不同方法得到的边坡功能函数CDF 曲线。在失效概率较小区域（如小于10-3）基于不含交叉项普通2阶多项式和2 阶HPCE 方法得到的功能函数的CDF与3 阶和4 阶HPCE 方法的CDF 有显著差别。相应的边坡系统失效概率如表4 所示。可见此时2 阶多项式展开方法的计算精度非常差，基于不含交叉项普通2 阶多项式和2 阶HPCE 方法的计算结果（9.1×10-5和1.2×10-5）与4 阶HPCE 方法的计算结果（8.1×10-4）相差一个数量级。3 阶HPCE 方法的计算结果（8.16×10-4）与4 阶HPCE 方法的计算结果非常吻合。这主要是由于当随机变量服从相关对数正态分布时，边坡安全系数与独立标准正态随机变量之间的非线性程度较高，2 阶多项式展开难以准确地近似这种高阶非线性函数关系。因此，当c和φ服从相关非正态分布时，2 阶多项式展开方法计算精度可能不够，应当采用相对高阶（如3 阶和4 阶）多项式展开构建代表性滑动面的响应面，从而更为准确地估计边坡系统失效概率。此外，本文提出方法可以有效地评价该低失效概率水平（10-5～10-3）的边坡系统可靠度。

图9 抗剪强度参数服从相关对数正态分布时边坡功能函数概率分布函数曲线Fig.9 Cumulative distribution functions of performance function of slope with correlated lognormal shear strength parameters

表4 抗剪强度参数服从相关对数正态分布时边坡可靠度结果Table 4 Comparison of reliability results of the slope with correlated lognormal shear strength parameters

7 结 论

（1）提出基于随机响应面法的边坡系统可靠度分析方法，能够有效地识别边坡代表性滑动面，并采用Hermite 多项式展开构建每条代表性滑动面的随机响应面，然后通过MCS 计算边坡系统失效概率，具有较高的计算效率。此外，本文提出的代表性滑动面确定方法，无需单独计算潜在滑动面间的相关性，简化了计算过程。该方法为实际工程边坡系统可靠度问题提供了一条有效的分析途径。

（2）基于单一临界确定性滑动面或者部分代表性滑动面进行边坡系统可靠度分析均会低估边坡系统失效概率，进而导致低估滑坡风险。本文提出方法能够有效地、合理地识别代表性滑动面，准确地计算低失效概率水平边坡系统可靠度，为制定合理的边坡加固措施提供了参考依据。

（3）本文提出方法能够有效地分析含相关非正态土体参数的边坡系统可靠度问题。土体抗剪强度参数分布类型及其相关性对边坡系统失效概率具有显著的影响。土体参数服从相关对数正态分布时的边坡系统失效概率远小于服从独立正态分布时的边坡系统失效概率。对于含相关非正态土体参数的边坡系统可靠度问题，2 阶多项式展开精度可能不够，建议采用相对高阶（如3 阶或4 阶）多项式展开，建立代表性滑动面的响应面，以提高计算精度。此外，本文假定边坡所有潜在滑动面为圆弧型，如何有效地分析任意形状滑动面的边坡系统可靠度问题需要进一步研究。

[1]OKA Y,WU T H.System reliability of slope stability[J].Journal of Geotechnical Engineering,ASCE,1990,116(8):1185－1189.

[2]CHOWDHURY R N,XU D W.Geotechnical system reliability of slopes[J].Reliability Engineering and System Safety,1995,47(3):141－151.

[3]HASSAN A M,WOLFF T F.Search algorithm for minimum reliability index of earth slopes[J].Journal of Geotechnical and Geoenvironmental Engineering,ASCE,1999,125(4):301－308.

[4]张兴,廖国华.多滑面边坡的破坏概率[J].岩土工程学报,1990,12(6):55－62.ZHANG Xing,LIAO Guo-hua.Failure probability of slope with multiple slip surfaces[J].Chinese Journal of Geotechnical Engineering,1990,12(6):55－62.

[5]谭晓慧.多滑面边坡的可靠性分析[J].岩石力学与工程学报,2001,20(6):822－825.TANG Xiao-hui.Reliability analysis on a slope with several surfaces[J].Chinese Journal of Rock Mechanics and Engineering,2001,20(6):822－825.

[6]DITLEVSEN O.Narrow reliability bounds for structural systems[J].Journal of Structural Mechanics,1979;7(4):453－472.

[7]LOW B K,ZHANG J,TANG W H.Efficient system reliability analysis illustrated for a retaining wall and a soil slope[J].Computers and Geotechnics,2011,38(2):196－204.

[8]CHING J,PHOON K K,HU Y G.Efficient evaluation of reliability for slopes with circular slip surfaces using importance sampling[J].Journal of Geotechnical and Geoenvironmental Engineering,ASCE,2009,135(6):768－777.

[9]ZHANG J,ZHANG L M,TANG W H.New methods for system reliability analysis of soil slopes[J].Canadian Geotechnical Journal,2011,48(7):1138－1148.

[10]JI J,LOW B K.Stratified response surfaces for system probabilistic evaluation of slopes[J].Journal of Geotechnical and Geoenvironmental Engineering,ASCE,2012,138(11):1398－1406.

[11]ZHANG J,HUANG H W,PHOON K K.Application of the Kriging-based response surface method to the system reliability of soil slopes[J].Journal of Geotechnical and Geoenvironmental Engineering,ASCE,2013,139(4):651－655.

[12]CHO S E.First-order reliability analysis of slope considering multiple failure modes[J].Engineering Geology,2013,154:98－105.

[13]ZHANG J,HUANG H W,JUANG C H,et al.Extension of Hassan and Wolff method for system reliability analysis of soil slopes[J].Engineering Geology,2013,160:81－88.

[14]LI L,WANG Y,CAO Z J,CHU X S.Risk deaggregation and system reliability analysis of slope stability using representative slip surfaces[J].Computers and Geotechnics,2013,53:95－105.

[15]JIANG S H,LI D Q,CAO Z J,et al.Efficient system reliability analysis of slope stability in spatially variable soils using Monte Carlo simulation[J].Journal of Geotechnical and Geoenvironmental Engineering,ASCE,2015,141(2):04014096.

[16]蒋水华,李典庆,曹子君,等.考虑参数空间变异性的边坡系统可靠度分析[J].应用基础与工程科学学报,2014,22(5):841－855.JIANG Shui-hua,LI Dian-qing,CAO Zi-jun,et al.System reliability analysis of slopes considering spatial variability of soil properties[J].Chinese Journal of Basic Science and Engineering,2014,22(5):841－855.

[17]LI D Q,CHEN Y F,LU W B,et al.Stochastic response surface method for reliability analysis of rock slopes involving correlated non-normal variables[J].Computers and Geotechnics,2011,38(1):58－68.

[18]李典庆,周创兵,陈益峰,等.边坡可靠度分析的随机响应面法及程序实现[J].岩石力学与工程学报,2010,29(8):1513－1523.LI Dian-qing,ZHOU Chuang-bing,CHEN Yi-feng,et al.Reliability analysis of slope stability using stochastic response surface method and code validation[J].Chinese Journal of Rock Mechanics and Engineering,2010,29(8):1513－1523.

[19]蒋水华,冯晓波,李典庆,等.边坡可靠度分析的非侵入式随机有限元法[J].岩土力学,2013,34(8):2347－2354.JIANG Shui-hua,FENG Xiao-bo,LI Dian-qing,et al.Reliability analysis of slope using non-intrusive stochastic finite element method[J].Rock and Soil Mechanics,2013,34(8):2347－2354.

[20]GEO-SLOPE International Ltd.Stability modeling with SLOPE/W 2007 Version:an engineering methodology[M].Calgary:GEO-SLOPE International Ltd.,2010.

[21]Wilson Window Ware Inc.WIL Reference Manual[M].Seattle:Wilson Window Ware Inc.,2012.http://techsupt.winbatch.com/techsupt/PDFs/wil-lang20120411.pdf.

[22]陈立宏,陈祖煜,刘金梅.土体抗剪强度指标的概率分布类型研究[J].岩土力学,2005,26(1):38－45.CHEN Li-hong,CHEN Zu-yu,LIU Jin-mei.Probability distribution of soil strength[J].Rock and Soil Mechanics,2005,26(1):38－45.

[23]范明桥,盛金保.土强度指标φ,c 的互相关性[J].岩土工程学报,1997,19(4):100－104.FAN Ming-qiao,SHENG Jin-bao.Cross correlation of soil strength indexes φ and c[J].Chinese Journal of Geotechnical Engineering,1997,19(4):100－104.