一种基于蒙特卡洛模拟的边坡稳定可靠度敏感性分析方法

摘要

本发明提供了一种基于蒙特卡洛模拟的边坡稳定可靠度敏感性分析方法，包括：步骤1，构造不确定性参数的联合概率密度函数；步骤2，利用蒙特卡洛模拟方法获得边坡失效概率，并获得失效样本；步骤3，设计多种敏感性分析方案，并分别构造各敏感性分析方案下不确定性参数的联合概率密度函数；步骤4，获得各敏感性分析方案下的边坡失效概率；步骤5，根据各敏感性分析方案的边坡失效概率，获得边坡失效概率随不确定性参数的统计特征的变化趋势。本发明适用范围广、计算过程简单、计算效率高，可有效揭示边坡可靠度水平与不确定性参数的统计特征间的响应规律，对边坡风险控制、设计优化等具有一定的指导意义。

说明书

技术领域

本发明涉及一种边坡参数敏感性分析方法，尤其涉及一种基于蒙特卡洛模拟的边坡稳定 可靠度敏感性分析方法。

背景技术

基于可靠度理论的边坡敏感性分析能有效地表征输入参数的不确定性及其统计特征(均 值、方差等)对边坡稳定可靠度的影响，对边坡设计、加固等具有重要指导意义，近年来受 到广泛关注。蒙特卡洛模拟方法作为一种求解可靠度的随机模拟方法，因其概念简单、适用 性强在岩土工程可靠度领域得到广泛应用。然而蒙特卡洛模拟方法存在不确定性传播机理不 明确的缺点，无法直接获得失效概率随着输入参数不确定性变化的响应规律。因此，在利用 蒙特卡洛模拟进行不确定性参数的敏感性分析时，通常需要重复模拟，计算输入随机变量统 计特征的一系列离散点对应的失效概率，然后分析统计特征与失效概率之间的变化规律。在 工程实践中，一般边坡设计方案对应的失效概率较小(比如，低于0.001)，而对于小概率问 题，蒙特卡洛模拟方法计算效率较低，重复模拟需要耗费大量的时间和计算机资源，所需的 计算时间和资源随着边坡稳定性分析模型的复杂程度的增加而增加。

Wang^[1]等提出了一种基于蒙特卡洛模拟的敏感性分析方法，该方法能够有效地定量分析 不确定性参数对于失效概率的相对贡献以及计算不同不确定性参数值对应的条件失效概率， 但是未对失效概率与不确定性参数的统计特征间的响应规律进行研究。随后Wang^[2]提出了一 种结合贝叶斯理论的敏感性分析方法，该方法不需要重新执行蒙特卡洛模拟，就可以计算单 个不确定性参数的统计特征变化后边坡的失效概率，获得失效概率与不确定性参数的统计特 征之间的响应规律。但是，该方法需要制作频率直方图，计算过程较为复杂，且无法解决多 个变量统计特征同时变化的敏感性分析问题。

文中涉及如下参考文献：

[1]Wang Y,Cao Z J,Au S K.Efficient Monte Carlo Simulation of parameter sensitivity in  probabilistic slope stability analysis[J].Computers and Geotechnics,2010,37(7－8):1015－1022.

[2]Wang Y.Uncertain parameter sensitivity in Monte Carlo Simulation by sample  reassembling[J].Computers and Geotechnics,2012,46:39－47.

[3]Kiureghian A D,Liu P L.Structural reliability under incomplete probability information[J]. Journal of Engineering Mechanics,1986,112:86－114.

发明内容

针对现有技术存在的不足，本发明基于蒙特卡洛模拟提出了一种计算过程简便、计算效 率高的边坡稳定可靠度敏感性分析方法。

为解决上述技术问题，本发明采用如下的技术方案：

基于蒙特卡洛模拟的边坡稳定可靠度敏感性分析方法，包括步骤：

步骤1，确定不确定性参数的统计特征，并构造不确定性参数的联合概率密度函数f(x)；

步骤2，基于联合概率密度函数f(x)，利用蒙特卡洛模拟方法获得边坡失效概率P_f，并 获得失效样本；

步骤3，设计多种敏感性分析方案，根据敏感性分析方案对应的统计特征，分别构造各敏 感性分析方案下不确定性参数的联合概率密度函数f^k(x)；

所述的敏感性分析方案采用如下方法获得：基于步骤1确定的统计特征，考虑不确定参数 的实际变化范围，改变一个或多个不确定性参数的统计特征值，保持其它不确定性参数的统 计特征概率分布不变，即获得一种敏感性分析方案；

步骤4，获得各敏感性分析方案下的边坡失效概率，本步骤进一步包括：

4.1根据联合概率密度函数f^k(x)和失效样本，获得失效样本的权重指标f(x_j)和f^k(x_j)分别为第j个失效样本的权重指标、不确定性参数联合概率密度函数和第 k个敏感性分析方案下的联合概率密度函数；

4.2根据权重指标获得敏感性分析方案的边坡失效概率为第k个敏感 性分析方案的边坡失效概率，n_s为失效样本数，N为步骤2蒙特卡洛模拟方法中产生的随机 向量样本数；

步骤5，根据各敏感性分析方案的边坡失效概率，获得边坡失效概率随不确定性参数统计 特征的变化趋势，从而识别影响边坡稳定可靠度的关键不确定性参数。

步骤1进一步包括子步骤：

1.1根据不确定性参数的试验数据，或现有文献中不确定性参数的统计特征取值，确定 不确定性参数的统计特征以及各不确定性参数间的相关系数；

1.2根据不确定性参数的统计特征及各不确定性参数间的相关系数，构造不确定性参数 的联合概率密度函数f(x)。

子步骤1.2具体为：

(1)若不确定性参数全部相互独立，则不确定性参数的联合概率密度函数为各不确定性 参数概率密度函数的乘积；

(2)若不确定性参数全部相关，则基于Copula理论或Nataf变换法构造不确定性参数的联 合概率密度函数；

(3)若不确定性参数部分独立、部分相关，首先，采用步骤(1)中方法获得相互独立 的不确定性参数的联合概率密度函数，记为其中，q为相互独立的不确定性参数个 数；然后，采用步骤(2)中方法获得相关的不确定参数的联合概率密度函数f(x_q+1,…,x_n)； 则不确定性参数的联合概率密度函数

步骤2进一步包括子步骤：

2.1产生服从于联合概率密度函数f(x)的随机向量样本；

2.2对随机向量样本进行边坡稳定性分析，判断各随机向量样本对应的边坡是否失效；

2.3统计失效样本数，获得边坡失效概率其中，n_s为失效样本数，N为随机向 量样本数；所述的失效样本即失效边坡对应的随机向量样本；

2.4输出并保存失效样本。

与现有技术相比，本发明具有以下优点和有益效果：

1、适用范围广：

本发明基于蒙特卡洛模拟，概念清晰、简单易懂。此外，本发明能够有效地解决单一不 确定性参数和多种不确定性参数的边坡稳定可靠度敏感性分析问题，同时也能较好地适用于 功能函数表达式为隐式的复杂边坡稳定可靠度敏感性分析问题，适用范围广。

2、计算过程简单，计算效率高：

传统蒙特卡洛模拟方法在计算敏感性分析方案所对应的边坡失效概率时，需重复执行蒙 特卡洛模拟，计算效率较低，尤其对于复杂边坡的稳定可靠度敏感性分析而言，蒙特卡洛模 拟中边坡确定性分析所需时间长，重复执行蒙特卡洛模拟进行敏感性分析计算时间长，对计 算资源要求高。本发明充分利用不确定性参数原概率分布蒙特卡洛模拟的计算信息，在求解 敏感性分析方案失效概率过程中，只需计算失效样本的权重指标，避免了蒙特卡洛模拟的重 复执行，简化了计算过程，提高了计算效率，是一种简单高效的边坡稳定可靠度敏感性分析 方法。尤其对于上述复杂边坡的稳定可靠度敏感性分析而言，本发明对计算效率的提高显得 尤为突出，有利于推动蒙特卡洛模拟方法在边坡可靠度工程的应用。

3、揭示不确定性传播机理，为边坡工程设计提供重要的参考依据：

本发明能有效地揭示边坡可靠度水平与不确定性参数的统计特征间的响应关系，从而明 确了蒙特卡洛模拟过程中不确定性的传播机理，对边坡风险控制、设计优化等具有一定的指 导意义。

综上，本发明可简化计算过程，提高计算效率，揭示蒙特卡洛模拟过程中不确定性的传 播机理，对推动岩土工程可靠度与风险理论在工程实践中应用有着重要的实用价值。

附图说明

图1为本发明的具体流程示意图；

图2为基于蒙特卡洛模拟的边坡可靠度分析流程图；

图3为实施例中边坡的剖面图；

图4为实施例中不同敏感性方案的不确定性参数的概率分布图，以硬质粘土厚度T_cr为例；

图5～7为实施例中敏感性分析结果。

具体实施方式

下面将结合附图和具体实施方式进一步说明本发明。

本发明的具体流程见图1，具体步骤如下：

步骤1，确定不确定性参数的统计特征，并构造不确定性参数的联合概率密度函数f(x)， 即不确定性参数原分布的联合概率密度函数f(x)，X为表征不确定性参数的随机向量。

岩土工程可靠度分析中，土体参数可分为两类参数：确定性参数和不确定性参数，不确 定性参数包括土体粘聚力、内摩擦角、土体容重、土层厚度等。对于不确定性参数，其不确 定性的量化可根据有无试验数据分为两种情况：一，若边坡设计中有现场或室内试验结果， 可根据不确定性参数的试验数据，采取统计方法量化其不确定性。比如对于不确定性参数土 体重度而言，在现场取多个测点试样，根据各测点试样的土体重度，利用统计方法确定其统 计特征，例如均值、方差等，然后根据获得的统计特征量化土体重度的不确定性。二，若没 有关于不确定性参数的试验资料，可根据现有文献不确定性参数的统计特征取值，对所研究 边坡中相应不确定性参数的统计特征进行合理的设定。

本步骤进一步包括以下子步骤：

1.1根据不确定性参数的试验数据，或现有文献不确定性参数的统计特征取值，确定不 确定性参数的统计特征以及各不确定性参数间的相关系数ρ。

1.2根据不确定性参数的统计特征以及各不确定性参数间的相关系数，构造不确定性参 数的联合概率密度函数f(x)。

联合概率密度函数f(x)的构造可分为三种情况：(1)不确定性参数全部相互独立情况下 联合概率密度函数的构造；(2)不确定性参数全部相关情况下联合概率密度函数的构造；(3) 不确定性参数部分独立、部分相关的情况下联合概率密度函数的构造。

对于相互独立的不确定性参数而言，n维不确定性参数的联合概率密度函数f(x)为各不 确定性参数的概率密度函数的乘积，因此f(x)可以表示为：

$f\left(x\right)=f(x_1,x_2,...,x_n)=\underset{i=1}{\overset{n}{\Pi }}f\left(x_i\right)---\left(1\right)$

式(1)中：f(x_i)为第i个不确定性参数X_i的概率密度函数，其中x_i为X_i的取值，n为不 确定性参数个数，Π表示连乘运算。

对于相关的不确定性参数而言，联合概率密度函数的构造可采用Copula理论、Nataf变 换法等。本具体实施中采用较常用的Nataf变换法构造相关的不确定性参数的联合概率密度 函数。现以二维相关不确定性参数X₁和X₂为例，介绍Nataf变换法的基本原理。

将相关的不确定性参数X₁和X₂的概率密度函数分别记为f(x₁)和f(x₂)，其中，x₁和x₂分 别表示X₁和X₂的取值，其累计分布函数分别记为F₁(x₁)和F₂(x₂)，不确定性参数X₁和X₂间 的皮尔逊相关系数为ρ。对不确定性参数X₁和X₂进行如下等概率变换：

$\left(\begin{array}{c}\Phi \left(y_i\right)=F_i\left(x_i\right)\\ y_i={\Phi }^{-1}\left(F_i\right(x_i))\end{array}\right),i=1,2---\left(2\right)$

式(2)中：Φ(y_i)为标准正态变量Y_i的累计分布函数，Φ^-1(·)为标准正态变量累计分布 函数Φ(·)的反函数，Y_i为X_i等概率变换后对应的标准正态分布变量，x_i和y_i表示X_i和Y_i的取 值。

根据Nataf变换法，可推导不确定性参数X₁和X₂的联合概率密度函数f(x)为：

$f\left(x\right)=f(x_1,x_2)=f\left(x_1\right)f\left(x_2\right)\frac{{\phi }_2(y_1,y_2;{\rho }_0)}{\phi \left(y_1\right)\phi \left(y_2\right)}---\left(3\right)$

式(3)中：ρ₀为变量Y₁和Y₂的皮尔逊相关系数；φ(·)为标准正态分布的概率密度函数； φ₂(·)为二维标准正态分布联合概率密度函数，其表达式为：

${\phi }_2(y_1,y_2;{\rho }_0)=\frac{1}{2\pi \sqrt{1-{{\rho }_0}^2}}\exp \left\{\frac{-1}{2(1-{{\rho }_0}^2)}\right[{y_1}^2-2{\rho }_0{y}_1{y}_2+{y_2}^2\left]\right\}---\left(4\right)$

根据相关系数定义，ρ和ρ₀具有关系：

$\rho ={\int }_{-\infty }^{+\infty }{\int }_{-\infty }^{+\infty }\left(\frac{x_1-{\mu }_1}{{\sigma }_1}\right)\left(\frac{x_2-{\mu }_2}{{\sigma }_2}\right)f\left(x_1\right)f\left(x_2\right)\frac{{\phi }_2(y_1,y_2;{\rho }_0)}{\phi \left(y_1\right)\phi \left(y_2\right)}{\mathrm{dx}}_1{\mathrm{dx}}_2---\left(5\right)$

式(5)中：μ₁和μ₂分别为不确定性参数X₁和X₂的均值；σ₁和σ₂分别为不确定性参数X₁和X₂的标准差。

变量Y₁和Y₂的相关系数ρ₀可通过式(5)迭代求解，也可以根据Der Kiureghian和Liu^[3]提供的针对不同概率统计特征的经验公式求解，还可以根据高斯-埃尔米特积分求解。但是一 般情况下，ρ₀和ρ相差不大，在近似求解中可以直接取ρ₀＝ρ。

需要指出的是，岩土体不确定性参数通常服从正态分布。若不确定性参数X₁和X₂服从 正态分布，此时根据式(3)和(4)，不确定性参数X₁和X₂的联合概率密度函数f(x)可以表 达为：

$\begin{array}{c}f\left(x\right)=f(x_1,x_2;\rho )=\frac{1}{2{\mathrm{\pi \sigma }}_1{\sigma }_2\sqrt{1-{\rho }^2}}\times \\ \exp \left\{\frac{-1}{2(1-{\rho }^2)}\right[{\left(\frac{x_1-{\mu }_1}{{\sigma }_1}\right)}^2-2\rho \left(\frac{x_1-{\mu }_1}{{\sigma }_1}\right)\left(\frac{x_2-{\mu }_2}{{\sigma }_2}\right)+{\left(\frac{x_2-{\mu }_2}{{\sigma }_2}\right)}^2\left]\right\}\end{array}---\left(6\right)$

对于部分独立、部分相关的不确定性参数而言，根据公式(1)获得独立的不确定性参数 的联合概率密度函数，记为根据上述Nataf变换方法或Copula理论可获得相关的 不确定参数的联合概率密度函数f(x_q+1,…,x_n)，故部分独立部分相关情况下不确定性参数的 联合概率密度函数其中，f(x_i)表示独立不确定性参数 X_i,(i＝1,...q)的概率密度函数，X_q+1、X_q+2、…X_n表示相关的不确定性参数。

步骤2，根据步骤1构建的联合概率密度函数f(x)，利用蒙特卡洛模拟方法计算边坡失效 概率P_f，并输出和保存失效样本x_j，j＝1,2,...,n_s。

本步骤采用传统蒙特卡洛模拟方法实现，具体流程见图2，主要包括如下子步骤：

2.1根据步骤1构建的联合概率密度函数f(x)，产生服从于联合概率密度函数f(x)的随 机向量样本x_i，i＝1,2,...N，i表示随机向量样本编号。

2.2构建边坡稳定性分析模型，计算随机向量样本x_i对应的边坡功能函数值 g(x_i)＝FS(x_i)-1，并根据功能函数值g(x_i)判断边坡是否失效。

目前常用的边坡稳定性分析法有极限平衡法和有限元法，其中极限平衡法包括瑞典条分 法、简化毕肖普法、摩根斯坦-普赖斯法等。用于边坡稳定性的商业软件主要有GeoStudio(极 限平衡法)、ANSYS(有限元法)，ABAQUS(有限元法)等。本发明适用任何一种边坡稳定性 分析法。

对于具体的边坡工程，根据要采用的边坡稳定性分析方法，运用相应的商业软件建立边 坡稳定性分析模型。随后，代入随机向量x_i，计算其对应的边坡功能函数值g(x_i)＝FS(x_i)-1_，当边坡功能函数值g(x_i)>0，即边坡安全系数FS(x_i)小于1时，表示边坡处于稳定状态；当 g(x_i)＝0时，表示边坡处于临界状态；当g(x_i)<0时，表示边坡失稳，即边坡失效。

循环执行子步骤2.1和子步骤2.2，直至i＝N。

2.3统计失效样本数目n_s，计算边坡失效概率

满足g(x_i)<0的随机向量样本x_i为失效样本。统计失效样本数量计为n_s，失效样本记为x_j， j＝1,2,...,n_s。根据蒙特卡洛模拟方法原理，边坡失效概率P_f可以表示为：

$\widehat{P_f}=\frac{1}{N}\underset{i=1}{\overset{N}{\Sigma }}I\left[g\left(x_i\right)\right]=\frac{n_s}{N}---\left(7\right)$

式(7)中，表示边坡失效概率P_f的估值，I[g(x)]为指示函数，当g(x)<0时，I[g(x)] 取1；当g(x)≥0时，I[g(x)]取为0。

2.4输出并保存失效样本x_j，j＝1,2,...,n_s。

步骤3，基于步骤1确定的不确定性参数统计特征，通过增加和减小不确定性参数统计特 征获得一系列敏感性分析方案，根据各敏感性分析方案对应的不确定性参数的统计特征构造 联合概率密度函数。

针对所研究的不确定性参数，通过改变其统计特征值，在确保其不超出相应不确定性参 数的实际变化范围前提下，设置M种敏感性分析方案；根据各敏感性分析方案对应的统计特 征，分别构造相应的不确定性参数联合概率密度函数f^k(x)，k＝1,2,...M，M为敏感性分析 方案数量，k表示敏感性分析方案编号，f^k(x)表示第k个敏感性分析方案相应的不确定性参 数的联合概率密度函数。

利用蒙特卡洛模拟方法分析不确定性参数敏感性时，需要通过计算输入不确定性参数统 计特征(如均值，方差)的一系列离散点的边坡失效概率，以探究边坡失效概率和不确定性 参数的统计特征间的响应规律。

如对土体粘聚力参数进行敏感性分析时，可通过增加和减小粘聚力的均值μ或标准差σ， 设置M种粘聚力的不同分布作为敏感性分析方案。然后，根据各敏感性分析方案下的统计特 征，构造各敏感性分析方案下的不确定性参数联合概率密度函数f^k(x)。联合概率密度函数 f^k(x)的具体构造方法同子步骤1.2。

这里需要指出的是，本发明能够对多个不确定性参数同时进行敏感性分析，如可以设置 一系列同时改变多个不确定性参数的均值μ和标准差σ的敏感性分析方案。

步骤4，利用步骤3获得的联合概率密度函数f^k(x)和步骤2获得的失效样本x_j， j＝1,2,...,n_s，求解各敏感性分析方案下的边坡失效概率P_f^k，根据边坡失效概率P_f^k随不确定 性参数统计特征的变化规律对不确定性参数的统计特征和变异系数进行敏感性分析。

根据边坡失效概率通用公式，第k个敏感性分析方案下的边坡失效概率P_f^k可表示为：

$\begin{array}{c}{P_f}^k=\underset{g\left(x\right)<0}{\int }{f}^k\left(x\right)dx=\int I\left[g\left(x\right)\right]f^k\left(x\right)dx\\ =E_{f^k\left(x\right)}\left\{I\right[g\left(x\right)\left]\right\}\end{array}---\left(8\right)$

式(8)中：为在联合概率密度函数f^k(x)下的数学期望算子。

结合不确定性参数原分布的联合概率密度函数f(x)，即步骤1构造的联合概率密度函数 f(x)，第k个敏感性分析方案下的边坡失效概率P_f^k可进一步表示为：

$P_j^k=\int I\left[g\left(x\right)\right]\frac{f^k\left(x\right)}{f\left(x\right)}f\left(x\right)dx=E_{f\left(x\right)}\left\{I\right[g\left(x\right)\left]{\omega }^k\right\}---\left(9\right)$

式(9)中，ω^k为本发明自定义的权重指标，其计算表达式如下：

${\omega }^k=\frac{f^k\left(x\right)}{f\left(x\right)}---\left(10\right)$

通过权重指标ω^k可以把求解敏感性分析方案条件下的边坡失效概率P_f^k转化为求解 I[g(x)]ω^k在f(x)分布下的数学期望。因此，可以直接利用步骤2.1中蒙特卡洛模拟的随机样 本x_i,i＝1,2,···,N，计算边坡失效概率P_f^k。故基于联合概率密度函数f(x)产生的N组随机 样本求解第k个敏感性分析方案的边坡失效概率P_f^k可表示为：

${{\hat{P}}_f}^k=\frac{1}{N}\underset{i=1}{\overset{N}{\Sigma }}\left(I\right[g\left(x_i\right)\left]{\omega }_i^k\right)---\left(11\right)$

考虑到指示函数I[g(x)]的取值特征，仅失效样本x_j,j＝1,2,···,n_s会对失效概率P_f^k产 生贡献，公式(11)可进一步表达为：

${{\hat{P}}_f}^k=\frac{1}{N}\underset{j=1}{\overset{n_s}{\Sigma }}\left({\omega }_j^k\right)---\left(12\right)$

P_f^k估计值的方差可表示为：

$\begin{array}{c}Var\left[{\hat{P}}_f^k\right]\approx \frac{1}{N}\left\{\frac{1}{N-1}\right\{\underset{i=1}{\overset{N}{\Sigma }}{\left(I\right[g\left(x_i\right)\left]\frac{f^k\left(x_i\right)}{f\left(x_i\right)}\right)}^2-N\left[\frac{1}{N}{\underset{i=1}{\overset{N}{\Sigma }}\left(I\right[g\left(x_i\right)\left]\frac{f^k\left(x_i\right)}{f\left(x_i\right)}\right)}^2\right]\left\}\right\}\\ \approx \frac{1}{N-1}\left\{\frac{1}{N}\underset{i=1}{\overset{N}{\Sigma }}{\left(I\right[g\left(x_i\right)\left]{\omega }_i^k\right)}^2-{\left({\hat{P}}_f^k\right)}^2\right\}\\ =\frac{1}{N-1}\left\{\frac{1}{N}\underset{i=1}{\overset{n_s}{\Sigma }}{\left({\omega }_j^k\right)}^2-{\left({\hat{P}}_f^k\right)}^2\right\}\end{array}---\left(13\right)$

的变异系数可表示为：

$COV\left[{\hat{P}}_f^k\right]=\frac{\sqrt{Var\left[{\hat{P}}_f^k\right]}}{E\left[{\hat{P}}_f^k\right]}\approx \frac{\sqrt{\frac{1}{N-1}\left\{\frac{1}{N}\underset{j=1}{\overset{n_s}{\Sigma }}{\left({\omega }_j^k\right)}^2-{\left({\hat{P}}_f^k\right)}^2\right\}}}{{\hat{P}}_f^k}---\left(14\right)$

的变异系数能够有效地表征计算结果的准确性，为本发明方法计算结果的可 靠性提供依据。

由式(12)可以看出，在求解敏感性分析方案的边坡失效概率时，本发明只需要计算原 分布产生的失效样本处的权重指标并对其求和，无需重新进行蒙特卡洛模拟，计算过程较为 简单、计算效率高。

需要说明的是，本步骤中计算失效概率所用联合概率密度函数为失效样本范围内的联合 概率密度函数，有别于步骤1到3中全体样本范围内的联合概率密度函数。

本步骤的具体子步骤如下：

本步骤的具体子步骤如下：

4.1根据失效样本的联合概率密度函数f(x_j)和f^k(x_j)，获得第j个失效样本x_j的权重指 标f(x_j)和f^k(x_j)分别为第j个失效样本的权重指标、不确定性参数联合概 率密度函数和第k个敏感性分析方案下的联合概率密度函数。

4.2根据权重指标计算敏感性分析方案的边坡失效概率

步骤5，根据敏感性分析方案的边坡失效概率，绘制边坡失效概率随不确定性参数统计特 征变化的趋势曲线，识别影响边坡稳定可靠度的关键不确定性参数。

实施例1

一、工程简介

詹姆斯湾堤坝位于加拿大魁北克地区的一项水利工程，堤坝长50km，其中一种设计方案 剖面如图3所示，堤坝高12m，中间有一个56m宽平台，边坡坡角为18.4°，坡比为3:1。坝 体以下土层从上至下依次为硬质粘土、海成粘土、湖积粘土和冰碛层。

参考文献[1]，本实施例考虑6个随机变量，即待研究的不确定性参数，分别为硬质粘土 厚度T_cr、海成粘土不排水抗剪强度S_uM、湖积粘土的不排水抗剪强度S_uL、坝基到冰碛层的距 离D_Till、筑坝填土的摩擦角φ_Fill和重度γ_Fill。不确定性参数均服从正态分布且相互独立，相应 的统计特征见表1。其它确定性参数如下：海成粘土的厚度为8m，硬质粘土不排水抗剪强度 为41kPa，硬质粘土、海成粘土、湖积粘土的重度分别为19kN/m³、19kN/m³、20.5kN/m³。 此外，由几何关系可知湖积粘土厚度T_L＝D_Till-T_cr-8。

为简化计算，文献[1]对最危险滑动面增加两个限制：(1)最危险滑动面与冰碛层相切 并过点(x,y)，其中，x＝4.9m，y＝36m；(2)最危险滑动面圆心横坐标x＝85.9m。这种 处理方法使得最危险滑动面的位置仅与坝基到冰碛层的距离D_Till有关。

表1不确定性参数的统计特征

二、具体实施流程

步骤1，根据表1中各不确定性参数的统计特征，构造不确定性参数原分布的联合概率密 度函数f(φ_Fill,γ_Fill,T_cr,S_uM,S_uL,D_Till)。

由于本实施例中各不确定性参数服从正态分布且相互独立，因此联合概率密度函数 f(φ_Fill,γ_Fill,T_cr,S_uM,S_uL,D_Till)可以表示为：

f(φ_Fill,γ_Fill,T_cr,S_uM,S_uL,D_Till)＝f(φ_Fill)f(γ_Fill)f(T_cr)f(S_uM)f(S_uL)f(D_Till)   (15)

其中：

$\left(\begin{array}{c}f\left({\phi }_{Fill}\right)=\frac{1}{\sqrt{2\pi }{\sigma }_1}\exp [-\frac{1}{2{{\sigma }_1}^2}{({\phi }_{Fill}-{\mu }_1)}^2]\\ f\left({\gamma }_{Fill}\right)=\frac{1}{\sqrt{2\pi }{\sigma }_2}\exp [-\frac{1}{2{{\sigma }_2}^2}{({\gamma }_{Fill}-{\mu }_2)}^2]\\ f\left(T_{cr}\right)=\frac{1}{\sqrt{2\pi }{\sigma }_3}\exp [-\frac{1}{2{{\sigma }_3}^2}{(T_{cr}-{\mu }_3)}^2]\\ f\left(S_{uM}\right)=\frac{1}{\sqrt{2\pi }{\sigma }_4}\exp [-\frac{1}{2{{\sigma }_4}^2}{(S_{uM}-{\mu }_4)}^2]\\ f\left(S_{uL}\right)=\frac{1}{\sqrt{2\pi }{\sigma }_5}\exp [-\frac{1}{2{{\sigma }_5}^2}{(S_{uL}-{\mu }_5)}^2]\\ f\left(D_{Till}\right)=\frac{1}{\sqrt{2\pi }{\sigma }_6}\exp [-\frac{1}{2{{\sigma }_6}^2}{(S_{Till}-{\mu }_5)}^2]\end{array}\right)---\left(16\right)$

式中，μ₁、μ₂、μ₃、μ₄、μ₅、μ₆和σ₁、σ₂、σ₃、σ₄、σ₅、σ₆分别为参数φ_Fill、γ_Fill、T_cr、 S_uM、S_uL、D_Till的均值和标准差，其具体数值见表1。

步骤2，根据不确定性参数原分布的联合概率密度函数f(φ_Fill,γ_Fill,T_cr,S_uM,S_uL,D_Till)，利用 蒙特卡洛模拟法获得詹姆斯湾堤坝的边坡失效概率P_f，并输出和保存失效样本x_j， j＝1,2,...,n_s。

2.1利用VBA程序产生服从于联合概率密度函数f(φ_Fill,γ_Fill,T_cr,S_uM,S_uL,D_Till)的随机向量 样本x_i，i＝1,2,...N，i表示随机向量样本编号，N＝2×10⁵。

2.2构建边坡稳定性分析模型，计算随机向量样本x_i对应的边坡功能函数值 g(x_i)＝FS(x_i)-1，并根据边坡功能函数值g(x_i)判断詹姆斯湾堤坝是否失效。

首先，采用简化毕肖普法，利用EXCEL表单构建詹姆斯湾堤坝稳定性分析模型。然后， 将随机向量样本作为输入参数带入詹姆斯堤湾坝的稳定性分析模型中。最后，根据稳定性分 析模型计算随机向量变量x_i对应的边坡安全系数FS(x_i)和边坡功能函数值g(x_i)。

循环执行子步骤2.1和子步骤2.2，直至i＝N。

2.3统计失效样本数n_s，共有失效样本n_s＝452个。根据式(7)计算詹姆斯湾堤坝的边坡 失效概率P_f＝452/200000＝2.26×10^-3。

2.4将失效样本x_j储存在EXCEL表单内，供敏感性分析使用。

以上基于不确定性参数原分布的蒙特卡洛模拟法(20万次)计算时间约10小时，所采用 计算机配置为内存4GB、CPU Intel Core i3和主频3.3GHz。因此，如果采用传统的蒙特卡洛 模拟法求解詹姆斯湾堤坝的敏感性分析问题，对于每种敏感性分析方案，在相同的模拟次数 下所需要的计算时间都约为10个小时。由此可见，采用传统蒙特卡洛模拟法求解边坡稳定可 靠度敏感性分析问题，计算效率非常低。

步骤3，通过改变所研究不确定性参数的统计特征值，设置一系列敏感性分析方案。根据 各敏感性分析方案对应的统计特征，构造各敏感性分析方案下不确定性参数的联合概率密度 函数f^k(φ_Fill,γ_Fill,T_cr,S_uM,S_uL,D_Till)，f^k(φ_Fill,γ_Fill,T_cr,S_uM,S_uL,D_Till)表示第k种敏感性分析方案下 的不确定性参数联合概率密度函数。

本步骤通过改变某不确定性参数的均值、方差或其他统计特征，来获得多种敏感性分析 方案，并通过分析各敏感性分析方案的边坡失效概率，获得该不确定性参数的统计特征对边 坡稳定可靠度的影响程度，从而得到影响边坡稳定可靠度的关键不确定性参数。

由于本实施例中不确定性参数较多，此处仅选取硬质粘土厚度T_cr、冰碛层厚度D_Till和湖 积粘土的不排水抗剪强度S_uL，研究T_cr、D_Till、S_uL的均值变化对詹姆斯湾堤坝边坡失效概率 的影响。如表2所示，本实施例中通过改变不确定性参数T_cr、D_Till、S_uL的均值，共设置18种 敏感性分析方案：方案1～6中，硬质粘土厚度T_cr的均值分别设为u₃-3σ₃、u₃-2σ₃、u₃-σ₃、 u₃+σ₃、u₃+2σ₃、u₃+3σ₃，u₃、σ₃分别为T_cr原分布的均值和标准差，其概率分布曲线见图4； 其他不确定性参数的统计特征值均保持原分布的统计特征值。类似的，分别对参数S_uL和D_Till做同样的处理，可以得到方案7～12和方案13～18。

各敏感性分析方案的联合概率密度函数f^k(φ_Fill,γ_Fill,T_cr,S_uM,S_uL,D_Till)的构造与案例分析步 骤1类似，在此不再赘述。

表2敏感性分析方案

步骤4，根据求得的联合概率密度函数f^k(φ_Fill,γ_Fill,T_cr,S_uM,S_uL,D_Till)和步骤2中输出的失效 样本x_j，j＝1,2,···,n_s，求解各敏感性分析方案的边坡失效概率P_f^k。

根据式(10)，第k个敏感性分析方案对应的权重指标ω^k可表示为：

${\omega }^k=\frac{f^k({\phi }_{Fill},{\gamma }_{Fill},T_{cr},S_{uM},S_{uL},D_{Till})}{f({\phi }_{Fill},{\gamma }_{Fill},T_{cr},S_{uM},S_{uL},D_{Till})}---\left(17\right)$

如当k取1时，权重指标ω¹可以表示为

${\omega }^1=\frac{f^1({\phi }_{Fill},{\gamma }_{Fill},T_{cr},S_{uM},S_{uL},D_{Till})}{f({\phi }_{Fill},{\gamma }_{Fill},T_{cr},S_{uM},S_{uL},D_{Till})}=\frac{f^1\left(T_{cr}\right)}{f\left(T_{cr}\right)}---\left(18\right)$

失效样本的权重指标值ω_j^k带入式(12)，则可以计算出各敏感性分析方案所对应的边坡 失效概率。

在此需要指出的是，权重指标的计算可以结合EXCEL内置的概率密度函数(如正态分布 概率密度函数NORM.DIST)，计算十分简便。同时，由步骤3～4可以看出，在求解敏感性分析 方案的过程中，由于充分利用了原概率分布蒙特卡洛模拟的失效样本，直接利用失效样本的 权重指标求解失效概率，无需重新执行蒙特卡洛模拟，计算过程非常简单高效，计算一种敏 感性分析方案对应的边坡失效概率所需时间仅为5.4×10^-3s。

步骤5，根据敏感性分析方案的边坡失效概率，绘制失效概率随不确定性参数的统计特征 变化的趋势曲线，识别影响詹姆斯湾堤坝稳定可靠度的关键不确定性参数。

图5～7分别为詹姆斯湾堤坝边坡失效概率随不确定性参数T_cr、S_uL和D_Till均值变化的趋势 曲线。如图5所示，随着硬质粘土厚度T_cr均值从2.56m逐步增加到5.44m、失效概率从8×10^-3缓慢下降到1×10^-3左右。如图6所示，当湖积粘土粘聚力S_uL均值的从12.27kN/m²增加50.13 kN/m²时，詹姆斯湾堤坝失效概率则从10^-1左右急剧下降到10^-7，说明湖积粘土粘聚力S_uL的均 值的增加会导致詹姆斯湾堤坝可靠度水平急剧降低。如图7所示，随着冰碛层厚度D_Till均值从 15.5m逐步增加21.5m，詹姆斯湾堤坝的失效概率从10^-5增加到10^-2左右，可见詹姆斯湾堤坝的 失效概率会随着不确定性参数冰碛层厚度D_Till均值的增加而增加。

通过上述分析可以得出，在本实施例所研究三个参数之中，詹姆斯湾堤坝的失效概率对 于冰碛层厚度D_Till、湖积粘土粘聚力S_uL较敏感，其中，湖积粘土粘聚力S_uL均值的变化对詹 姆斯湾堤坝稳定可靠度影响最大。故在工程上，应该采取相应的工程技术获得更多关于湖积 粘土粘聚力S_uL的信息，从而更加准确的评估边坡的安全性。

上述实施例表明，本发明所提出的蒙特卡洛模拟的边坡稳定可靠度敏感性分析方法，能 够高效地求解敏感性分析方案的失效概率，有效地揭示边坡可靠度水平与不确定参数的统计 特征之间的响应关系，并识别出影响边坡稳定可靠度的关键不确定性参数，为边坡风险控制、 稳定加固，优化设计等提供了参考依据。