一种节理岩质边坡可靠度分析下限法的制作方法

本发明涉及一种节理岩质边坡可靠度分析下限法，属于岩质边坡工程可靠度分析技术领域。

背景技术：

节理岩质边坡的边坡稳定分析一直是岩土工程界的研究热点，一般通过求解边坡的安全系数或极限荷载来预测边坡的安全度。边坡稳定性分析方法一般可分为确定性分析方法和非确定性分析方法，前者是将所有的参数取确定值进行计算，而后者是将某些计算参数当成随机变量进行计算分析。总所周知，岩质边坡的很多计算参数都是不确定性的，即实际工程边坡的任何性质都不具有唯一的确定性，如岩体容重、节理的分布、岩体材料的抗剪强度、地下水水位等，这些参数都具有很大的不确定性和变异性；我们常规计算所采用的确定值一般很难代表参数的真实值。因此，随机数值计算方法的发展和计算机性能的提升，对边坡进行可靠度分析越来越被工程师们所重视。

在当前节理岩质边坡的可靠度分析方法的研究成果中，主要分成两类：一类是基于极限平衡法的可靠度分析，另一类是基于有限单元法的可靠度分析。极限平衡法具有概念简单、计算效率高、工程应用便捷的特点，特别是一些极限平衡法有显式安全系数的计算表达式，其不足是必须事先假定滑面的位置。有限单元法比极限平衡法通用性更强、能适用于复杂的边坡稳定性分析，并且能够获得边坡的应力应变分布，采用有限元法的边坡可靠度分析可以得到较为合理的结果；但其不足是材料本构方程不完善、计算较大计算效率低。

在进行边坡的可靠度分析时，与边坡稳定性分析方法：极限平衡法、有限单元法配套的随机模拟方法有蒙特卡洛模拟、重要性抽样模拟等，这些随机模拟方法可以适用于复杂的功能函数，可以求解高纬度的复杂问题，因此这些方法被广泛应用于边坡的可靠度分析中。

在岩土工程极限承载力分析中，极限分析下限法是一种高效的求解工具，根据下限定理可建立岩土体的静力许可应力场并采用数学规划方法求解极限荷载。岩土体的下限分析法主要包括解析法、基于数值方法的下限法等。特别是近三十年来，下限法与有限元离散思想、块体单元离散思想的结合，使得其在构造大规模岩土结构的静力许可应力场成为可能；而且下限法求解得到的下限解小于真解(偏于安全)，因此下限解对工程也最有用。下限法的优点是：下限法在理论上比刚体极限平衡法更严密；同时，由于不考虑材料的本构关系下限法在数值计算时比有限元法更高效；下限法可以同时获得安全系数和应力场。使用下限法对岩质边坡进行确定性分析的研究成果已经非常多，但采用下限法对岩质边坡的可靠度分析的研究成果还非常少。

鉴于此，本发明基于国家自然科学基金项目(51564026)的研究工作提出了一种新的节理岩质边坡可靠度分析的下限法。

技术实现要素：

本发明提供了一种节理岩质边坡可靠度分析下限法，为节理岩质边坡可靠度计算提供一种新的方法。

本发明的技术方案是：一种节理岩质边坡可靠度分析下限法，所述方法具体步骤如下：

步骤1、根据节理岩质边坡的实际情况拟定节理岩质边坡的基本信息；

步骤2、生成结构面的抗剪参数的蒙特卡洛随机数；

步骤3、建立求解节理岩质边坡可靠度的下限法随机规划模型；

步骤4、节理岩质边坡可靠度计算的下限法随机规划模型的数值求解，得到节理岩质边坡的极限状态函数；

步骤5、根据步骤4的结果计算节理岩质边坡的可靠度。

所述节理岩质边坡的基本信息包括：节理岩质边坡的几何参数、完整岩体的物理力学参数、结构面地质参数、边坡外荷载信息。

所述步骤2具体为：分别确定结构面的凝聚力ck的均值μc、标准差σc和变异系数δc，摩擦角的均值标准差和变异系数同时设凝聚力ck、摩擦角的概率分布类型是正态分布；使用蒙特卡洛随机数生成函数生成结构面k的凝聚力和摩擦角的随机数

所述步骤3具体为：(1)采用块体单元离散节理岩质边坡，将节理岩质边坡离散成块体单元+结构面的几何系统；(2)建立目标函数：将容重超载系数的随机变量设为目标函数；(3)建立块体单元的下限法约束条件；(4)建立节理岩质边坡可靠度计算的极限状态函数；(5)根据目标函数、下限法的约束条件，并结合结构面的抗剪参数的蒙特卡洛随机数以及节理岩质边坡可靠度计算的极限状态函数，得到求解节理岩质边坡可靠度的下限法随机规划模型；其中块体单元的下限法约束条件包括块体单元的平衡条件、结构面的屈服条件、块体单元的静力边界条件。

所述步骤4具体为：将凝聚力随机数摩擦角随机数从1至N逐次代入下限法随机规划模型；对每一个确定的随机规划模型变成一个系数矩阵是定值的线性规划问题，采用单纯形法求解对应的目标函数λ(t)；将求解得到的N组λ(t),(t＝1,…,N)代入极限状态函数，得到节理岩质边坡的极限状态函数Z＝λ(t),(t＝1,…,N)。

所述节理岩质边坡的可靠度包括节理岩质边坡容重超载系数的均值、标准差、可靠度指标和节理岩质边坡的失效概率。

所述下限法随机规划模型为：

式中，t＝1,…,N，N是生成的随机数的数量；

节理岩质边坡可靠度计算的极限状态函数Z＝λ(t)，目标函数λ(t)表示容重超载系数的随机变量；

表示块体单元的平衡条件，mi是块体单元i中结构面的数量，αk是结构面k的nk轴与x轴的夹角且逆时针为正，Nk是结构面k沿外法线nk方向的法向力且拉正压负，Vk是结构面k沿切向sk方向的剪力且使块体产生逆时针转动为正，λ(t)是第t个超载系数随机量，γ是块体单元材料的容重，Ai是块体单元i的面积，fxi是第i个块体单元的形心处沿x方向的等效外力，fyi是第i个块体单元的形心处沿y方向的等效外力，nb是节理岩质边坡中块体单元的数量；

表示结构面的屈服条件，是结构面k凝聚力的第t个随机数，是结构面k摩擦角的第t个随机数，lk是结构面k的长度，ns是整个节理岩质边坡中结构面的数量；

表示块体单元的静力边界条件，是已知边界条件的结构面的内力向量，b表示节理岩质边坡中已知的边界外力等于0的结构面，nc是整个岩质边坡中边界上外力等于0的结构面的数量；

表示凝聚力和摩擦角蒙特卡洛随机数的生成函数；分别是结构面凝聚力和摩擦角的均值，分别是结构面凝聚力、摩擦角的标准差。

本发明的有益效果是：

1、本发明为节理岩质边坡的极限承载力可靠度分析提供一种新方法，首次将下限法理论、块体单元离散思想、随机规划方法和蒙特卡洛法结合起来分析节理岩质边坡的可靠度，可求解岩质边坡的可靠度(容重超载系数的均值、标准差、可靠度指标，以及节理岩质边坡的失效概率)。

2、本发明方法可获得与结构面凝聚力、摩擦角的蒙特卡洛随机数(正态分布)对应的容重超载系数的正态分布情况，可以由容重超载系数＜1来计算边坡的失效概率；

3、本发明方法中采用蒙特卡洛方法求解随机规划问题，具有精度高、计算速度快等特点；

4、本发明方法采用单纯形法求解线性规划问题，编制计算程序简单。

附图说明

图1为本发明的流程图；

图2为节理岩质边坡块体单元受力示意图；

图3为节理岩质边坡相邻块体单元之间的结构面示意图；

图4为实施例1节理岩质边坡的几何形状示意图；

图5为实施例1结构面AC凝聚力的随机数正态分布示意图；

图6为实施例1结构面AC摩擦角的随机数正态分布示意图；

图7为实施例1节理岩质边坡块体单元离散示意图；

图8为实施例1节理岩质边坡容重超载系数正态分布示意图；

图9为实施例2节理岩质边坡的几何形状示意图；

图10为实施例2结构面凝聚力的随机数正态分布示意图；

图11为实施例2结构面摩擦角的随机数正态分布示意图；

图12为实施例2节理岩质边坡块体单元离散示意图；

图13为实施例2节理岩质边坡容重超载系数正态分布示意图。

具体实施方式

实施例1：如图1-13所示，一种节理岩质边坡可靠度分析下限法，具体方法步骤如下：

(一)、拟定节理岩质边坡的基本信息

实施例1是含有一组软弱结构面的岩质边坡，其结构如图4所示，其主要计算参数如下：节理岩质边坡的几何参数：高度55.94m、顶宽36.6m、临空面角度45度；完整岩体的容重设为定值取25kN/m³；结构面(节理面)AC的倾角是30度；边坡只受岩体自重荷载作用。

(二)、生成结构面的抗剪参数的蒙特卡洛随机数

结构面AC的凝聚力的均值取μc＝100kPa、变异系数δc分别取0.10、0.15、0.20、0.25，摩擦角的均值取变异系数分别取0.10、0.15、0.20、0.25，凝聚力和摩擦角的概率分布类型均为正态分布。采用式(1)生成了结构面AC抗剪参数的随机数。本实施例中取N＝3000，即分别生成了凝聚力和摩擦角3000个随机数在变异系数δc＝0.2的情况下结构面AC凝聚力的随机数的正态分布如图5所示，在摩擦角变异系数的情况下结构面AC摩擦角的随机数的正态分布如图6所示。

(三)、建立求解节理岩质边坡可靠度的下限法随机规划模型

将实施例1边坡采用块体单元离散如图7所示，共离散为2个块体单元。根据式(9)建立实施例1求解节理岩质边坡可靠度的下限法随机规划模型。

(四)、节理岩质边坡可靠度计算的下限法随机规划模型的数值求解

采用如图1所示的节理岩质边坡可靠度下限法数值计算流程计算3000组容重超载系λ(t),(t＝1,…,3000)。在结构面凝聚力的均值μc＝100kPa、变异系数δc＝0.20，摩擦角的均值取变异系数的条件下，计算得到的实施例1边坡容重超载系数λ(t),(t＝1,…,3000)的正态分布情况如图8所示。

(五)、计算节理岩质边坡的可靠度。

按式(10)、(11)、(12)、(14)求解得到实施例1节理岩质边坡的可靠度如表1所示。

表1实施例1节理岩质边坡的可靠度统计表

在所有参数都是定值的情况下，求解此边坡的容重超载系数是一个静定问题，对应的定值λ的极限值可以由解析法求出，在结构面AC的凝聚力ck＝100000kPa、摩擦角的情况下，得到的容重超载系数的解析是1.259，由本发明方法计算得到的各种不同变异系数下均值与解析解的误差均很小(在1.1％以内)，且两者的误差随着变异系数的增加而增大。

边坡的可靠度指β标随着变异系数的增大而减小，说明抗剪强度参数离散程度越高则边坡的可靠度越低。

边坡的失效概率Pf随着变异系数的增大而增大，说明抗剪强度参数离散程度越高则边坡的失效概率越高。

实施例2：一种节理岩质边坡可靠度分析下限法，所述方法具体步骤如下：

(一)、拟定节理岩质边坡的基本信息

实施例2是含有两组成组节理的岩质边坡，其几何信息如图9所示，其主要计算参数如下：节理岩质边坡的几何参数：高度30m、顶宽45m、临空面角度80度；完整岩体的容重设为定值取27kN/m³，两组结构面(节理面)的间距为6m、两组结构面的倾角分别是25度和80度；边坡只受岩体自重荷载作用。

(二)、生成结构面的抗剪参数的蒙特卡洛随机数

实施例2：所有结构面的凝聚力的均值取μc＝50kPa、变异系数δc分别取0.10、0.15、0.20、0.25，所有结构面的摩擦角的均值取变异系数分别取0.10、0.15、0.20、0.25。采用式(1)生成了实施例2边坡所有结构面的抗剪参数的随机数。本实施例中取N＝3000，即分别生成了节理凝聚力和摩擦角的3000个随机数

在变异系数δc＝0.15的情况下所有结构面凝聚力的随机数的正态分布如图10所示，在变异系数的情况下所有结构面摩擦角的随机数的正态分布如图11所示。

(三)、求解节理岩质边坡可靠度的下限法随机规划模型

将实施例2边坡采用块体单元离散如图12所示，共离散为37个块体单元。根据式(9)建立实施例2节理岩质边坡可靠度的下限法随机规划模型。

(四)、节理岩质边坡可靠度计算的下限法随机规划模型的数值求解

采用如图1所示的节理岩质边坡可靠度下限法数值计算流程计算实施例2的3000组容重超载系λ(t),(t＝1,…,3000)。

在结构面凝聚力的均值μc＝50kPa、变异系数δc＝0.15，摩擦角的均值变异系数的条件下，计算得到的实施例2边坡容重超载系数λ(t),(t＝1,…,3000)的正态分布情况如图13所示。

(五)、计算节理岩质边坡的可靠度。

按式(10)、(11)、(12)、(14)求解得到实施例2节理岩质边坡的可靠度如表2所示。

表2实施例2节理岩质边坡的可靠度统计表

节理岩质边坡的容重超载系数的均值和标准差都随变异系数的增大而增大，可靠度指标β随变异系数的增大而减小。

节理岩质边坡的失效概率Pf随着变异系数的增大而减增大，当变异系数取0.1时失效概率为17.867％，当变异系数取0.25时失效概率增加至35.7％。

本发明的工作原理是：

以节理岩质边坡为研究对象，将蒙特卡洛方法、下限法理论、块体单元离散思想和随机规划方法结合起来，求解节理岩质边坡的可靠度；首先，采用蒙特卡洛方法生成结构面抗剪参数的随机数；然后，将节理岩质边坡离散成块体单元+结构面的几何体；根据下限法理论构建节理岩质边坡同时满足块体单元平衡条件、结构面屈服条件和静力边界条件的静力许可应力场，并以容重超载系数为目标函数；建立了求解节理岩质边坡可靠度的下限法随机规划模型，同时提出了随机规划模型的求解方法；通过计算可获得节理岩质边坡的可靠度(包括容重超载系数的均值、标准差、可靠度指标以及边坡的失效概率)。

本发明的技术路线图如图1所示。

本发明的求解节理岩质边坡可靠度的下限法的技术方案依次按以下步骤进行：

一、拟定节理岩质边坡的基本信息

根据节理岩质边坡的实际情况，拟定进行可靠度分析的主要计算参数，包括：节理岩质边坡的几何参数、完整岩体的物理力学参数(容重等)，结构面地质参数、边坡外荷载信息等。

二、生成结构面的抗剪参数的蒙特卡洛随机数

节理岩质边坡的承载能力主要受控于结构面的强度，同时结构面的两个抗剪强度参数(凝聚力ck、摩擦角具有较大的不确定性和变异性，因此本发明设结构面的凝聚力ck、摩擦角为相互独立的随机量，并符合正态分布；同时设岩体的容重等其他参数为定值。本发明采用蒙特卡洛方法生成凝聚力和摩擦角随机量的随机数，具体过程如下：

1.分别确定结构面的凝聚力ck的均值μc、标准差σc和变异系数δc，摩擦角的均值标准差和变异系数同时设凝聚力ck、摩擦角的概率分布类型是正态分布；

2.使用蒙特卡洛随机数生成函数生成结构面k的凝聚力和摩擦角的随机数，凝聚力和摩擦角蒙特卡洛随机数生成函数如下：

上式中:t＝(1,…,N)，N是生成的随机数的数量，一般取2000～5000；是结构面凝聚力的第t个随机数(上标r表示变量为随机数变量)；是结构面摩擦角的第t个随机数(上标r表示变量为随机数变量)；分别是结构面凝聚力和摩擦角的均值，分别是结构面凝聚力、摩擦角的标准差。

三、建立求解节理岩质边坡可靠度的下限法随机规划模型

塑性极限分析下限法可以叙述为：与所有静力许可应力场对应的荷载中，极限荷载最大。同时满足平衡条件、屈服条件以及静力边界条件的应力场称为静力许可应力场。在一个节理岩质边坡结构体内，可以构造出无数个符合上述条件的静力场；一一寻求其对应的荷载，其中最大的荷载必然是最接近真实的极限荷载。下限解必然小于真实解。

为了构建节理岩质边坡的静力许可应力场，本发明采用块体单元离散节理岩质边坡，将边坡外荷载的超载系数设定为目标函数，构建块体单元同时满足平衡方程、屈服条件和边界条件的静力许可应力场，建立求解节理岩质边坡极限荷载的下限法随机规划模型。具体步骤如下：

1.采用块体单元离散节理岩质边坡

本发明采用块体单元离散节理岩质边坡，块体单元受力如图2所示，相邻块体单元的结构面k受力如图3所示。其中(x,y)是总体坐标系；(nk,sk)是相邻块体单元i、j之间结构面k的局部坐标系，如图3所示nk为结构面k的外法线方向，sk为结构面的切向；块体单元i形心Ci上作用有等效荷载力向量fxi是第i个块体单元的形心处沿x方向的等效外力，fyi是第i个块体单元的形心处沿y方向的等效外力。；相邻块体单元i、j之间的结构面k形心Pk处作用有内力向量结构面另一侧作用有其反作用力向量Nk是沿结构面k沿外法线nk方向的法向力(正负号规定：拉正压负)，Vk是沿结构面k沿切向sk方向的剪力(正负号规定：使块体产生逆时针转动为正)，N'k是Nk的反作用力，Vk'是Vk的反作用力。

采用块体单元离散节理岩质边坡以后，节理边坡变成块体单元+结构面的几何系统，为了简化计算和便于构建下限法的静力许可应力场的约束条件，本发明作如下假设：(1)假设块体单元只发生平动，不会发生转动；(2)假设块体单元为刚体，块体单元不会发生任何变形和破坏，边坡的破坏只会发生在相邻块体单元之间的结构面上；(3)假设结构面材料为理想刚塑性材料，(4)假设岩体中的结构面仅发生剪切破坏。

2.建立目标函数

本发明采用容重超载的方式使节理岩质边坡达到极限状态，设边坡岩块的实际容重为γ，即当容重γ逐步增大至γ'时边坡处于破坏的极限状态。本发明定义容重超载系数如下：

上式中：λ是节理岩质边坡的容重超载系数；γ'是节理岩质边坡超载时块体单元材料的极限容重；γ是节理岩质边坡的块体单元材料的当前实际容重。

根据下限法理论，节理岩质边坡在达到极限状态时，需求容重γ的最大值γ'；根据式(1)结构面的凝聚力随机变量和摩擦角随机变量分别被离散为N个随机数，则容重超载系数λ也对应变为与和相关的随机变量如下：

上式中：t＝(1,…,N)，λ(t)是节理岩质边坡容重超载系数的随机变量。

根据下限法理论，本发明将容重超载系数的随机变量λ(t)设为目标函数如下：

Maximize:λ(t)(4)

3.建立块体单元的下限法约束条件

构建静力许可应力场需要块体单元同时满足平衡方程、屈服条件和静力边界条件。

(1)块体单元的平衡条件

对于图2所示的任意一个块体单元i受自重、外荷载、相邻岩块的作用而保持平衡，其力的平衡方程如下：

上式中：nb是节理岩质边坡中块体单元的数量；mi是块体单元i中结构面的数量；λ(t)是第t个超载系数随机量，Ai是块体单元i的面积，γ是块体单元材料的容重，αk是结构面k的nk轴与x轴的夹角(正负号规定：逆时针为正)，Nk是沿结构面k沿外法线nk方向的法向力(正负号规定：拉正压负)，Vk是沿结构面k沿切向sk方向的剪力(正负号规定：使块体产生逆时针转动为正)，fxi是第i个块体单元的形心处沿x方向的等效外力，fyi是第i个块体单元的形心处沿y方向的等效外力。

(2)结构面的屈服条件

对于如图3所示的任意一条结构面k，由于结构面的凝聚力和摩擦角的是N组随机数因此结构面k的剪切滑移破坏的Mohr-Coulomb屈服条件可用凝聚力和摩擦角的随机数写成如下形式：

上式中：ns是整个节理岩质边坡中结构面的数量；lk是结构面k的长度；是结构面k的凝聚力的第t个随机数(上标r表示变量为随机数变量)；是结构面k摩擦角的第t个随机数(上标r表示变量为随机数变量)；Nk是结构面k的法向力(符号规定：拉正压负)，Vk是结构面k的剪力；N是凝聚力随机数和摩擦角随机数的数量。

(3)块体单元的静力边界条件

根据下限法理论，静力许可应力场必须满足已知的静力边界条件。考虑节理岩质边坡中已知的边界外力等于0的结构面b，其边界条件表达式为：

上式中：nc是整个岩质边坡中边界上外力等于0的结构面的数量；是已知边界条件的结构面的内力向量。

4.建立节理岩质边坡可靠度计算的极限状态函数

为了求解节理岩质边坡容重超载系数的数学分布规律，根据式(3)本发明设节理岩质边坡可靠度计算的极限状态函数为：

上式中：Z是节理岩质边坡的极限状态函数。

5.建立节理岩质边坡可靠度计算的下限法随机规划模型

根据静力许可应力场的目标函数式(4)以及下限法的约束条件式(5)、(6)、(7)，并结合结构面的凝聚力随机数和摩擦角随机数生成函数式(1)以及节理岩质边坡可靠度计算的极限状态函数式(8)，可到求解节理岩质边坡可靠度的下限法随机规划模型如下：

上式中：(t＝1,…,N)。

四、节理岩质边坡可靠度计算的下限法随机规划模型的数值求解

式(9)求解目标是：由结构面的凝聚力随机数和摩擦角随机数的数学特征求解容重超载系数随机变量λ(t)的数学特征。式(9)的下限法模型的矩阵的系数中包含结构面的凝聚力和摩擦角的随机数其是一个典型的随机规划问题，因此企图直接通过的数学特征精确地求解λ(t)的数学特征是非常困难的，只能采用近似的数值计算方法。

本发明采用如下的数值方法求解节理岩质边坡可靠度的下限法随机规划模型：

(1)将t从t＝1到t＝N循环，逐次将凝聚力随机数摩擦角随机数代入式(9)；

(2)对于每一个确定的式(9)变成一个系数矩阵是定值的线性规划问题，此时采用单纯形法求解对应的目标函数λ(t)；

(3)将求解得到的N组λ(t),(t＝1,…,N)代入极限状态函数，可计算得到节理岩质边坡的极限状态函数Z＝λ(t),(t＝1,…,N)。

本发明的数值求解流程如图1所示。

五、计算节理岩质边坡的可靠度。

根据计算得到的N组λ(t)计算极限状态函数Z的基本统计量可得到节理岩质边坡的可靠度，包括：节理岩质边坡容重超载系数的均值、标准差、可靠度指标和节理岩质边坡的失效概率。

节理岩质边坡容重超载系数的均值计算公式如下：

上式中：μλ是节理岩质边坡容重超载系数的均值。

节理岩质边坡容重超载系数的标准差计算公式如下：

上式中：σλ是节理岩质边坡容重超载系数的标准差。

节理岩质边坡可靠度指标计算公式如下：

上式中：μλ是节理岩质边坡容重超载系数的均值、σλ是节理岩质边坡容重超载系数的标准差、β是节理岩质边坡的可靠度指标。

当边坡的容重超载系数≥1时边坡是安全的、当边坡的容重超载系数＜1时边坡发生失稳，因此节理岩质边坡的失效功能函数如下：

上式中：I(t)是节理岩质边坡的失效功能函数。

则边坡的失效概率P可按下式计算：

上式中：Pf是节理岩质边坡的失效概率。

本发明的特点是：本发明将蒙特卡洛方法、塑性极限分析下限法、块体单元离散思想以及随机规划结合起来，提出了一种节理岩质边坡可靠度计算的下限法，可以计算得到边坡容重超载系数的均值、标准差、可靠度指标，同时可以得到节理岩质边坡的失效概率。相对传统的极限平衡可靠度分析法其具有更高的计算精度；相比有限元可靠度分析法，由于本发明方法忽略了材料的本构关系，因此计算速度和效率更高。

上面结合附图对本发明的具体实施方式作了详细说明，但是本发明并不限于上述实施方式，在本领域普通技术人员所具备的知识范围内，还可以在不脱离本发明宗旨的前提下做出各种变化。