1.一种基于牵引式斜坡变形破坏机理的稳定性分析和预测预警方法，其特征在于：

1)牵引式斜坡的变形机理分析，具体地讲，斜坡前缘为破坏后区、临界状态面及斜坡后缘稳定状态区的分析，牵引式斜坡高度与变形的关系呈现出“S”曲线特征，存在三种可能的破坏形式分析，提出了牵引式斜坡破坏的决定方法，各种物理量相对应的关系曲线图；

2)本构方程及参数决定方法；

3)稳定分析中的不平衡拉力法条块划分及计算和拉破坏的决定，牵引式斜坡的综合下滑力-抗滑力计算法、牵引式斜坡的主拉力法和牵引式斜坡变形稳定系数法的计算步骤及计算方法；

4)牵引式斜坡有限元滑面边界计算法的计算模型、计算步骤、收敛标准及位移、应力边界的赋值方法，所述应力边界即边界赋值满足本发明提出的本构关系；

5)主滑方向定义：主滑方向所指的对象为已滑动的坡体，简称已滑体，走向为已滑体的最大剖面方向，倾向为已滑体剖面所具有的倾向，倾角为已滑体下滑力的矢量和方向。

2.一种基于牵引式斜坡变形破坏机理的稳定性分析和预测预警方法，其特征在于包含如下步骤：(1)分析牵引式斜坡的变形机理，将斜坡前缘以临界滑面为基准，定义为破坏后区、临界状态面及稳定状态区，牵引式斜坡高度与变形的关系呈现出S曲线特征；分析可能发生破坏的三种情况，情况I为沿最弱的滑面发生破坏，情况II为斜坡发生新的剪切破坏，情况III为坡体或坡体上的节理或裂隙发生拉破坏；

(2)计算获得滑面各点的临界应变值，并计算出现行滑面临界应力状态点及位移值；

(3)针对步骤(1)的三种情况分别计算稳定系数，进行稳定性分析；

(4)利用边坡稳定性计算的滑面边界法获取斜坡坡面位移，进行预测预警。

3.如权利要求2所述的基于牵引式斜坡变形破坏机理的稳定性分析和预测预警方法，其特征在于：所述情况II的剪切破坏由滑面的剪应力与剪应变的全过程曲线决定，情况III的拉破坏由坡体的拉应力与拉应变全过程曲线决定。

4.如权利要求2所述的基于牵引式斜坡变形破坏机理的稳定性分析和预测预警方法，其特征在于所述步骤(2)的具体计算过程如下：步骤(2.1)分析斜坡的基本形态、特征，试验获得滑面的基本物理力学参数G、S、m、ρ、C、φ、a1，a2，a3，ξN和坡体的弹性模量力学参数，计算相对应的位移场和应力场，通过应力场、位移场决定相对应的稳定系数；

m ρ

步骤(2.2)将步骤(2.1)得到的参数代入公式τ＝Gγ[1+γ/S] ，式中τ、γ分别表示材料的剪应力和剪应变，τ、G的单位为MPa或kPa或Pa，S、m和ρ为无单位参数，且-1<ρ≤0和1+mρ≠0；

peak peak

临界应力空间τ 采用摩尔库伦准则τ ＝C+σntanφ，式中C为凝聚力，σn为法向应力，C和σn的单位为MPa或kPa或Pa，φ为滑面摩擦角；

临界应变空间γpeak采用公式 式中σn为法向应力，单位为MPa或kPa或Pa；

peak ρ

临界应力空间与临界应变空间的关系为τ /γpeak＝G[1-1/(1+mρ)] ，临界应变空m间符合公式S+(1+mρ)γpeak＝0；

c ζ

参数ρ＝ρ0/(1+(ρ0/ρc-1)(σn/σn) )，式中ρ0为法向应力σn为零值的ρ值，cρc为σn等于σn 时的ρ值，ζ为常系数。

5.如权利要求2所述的基于牵引式斜坡变形破坏机理的稳定性分析和预测预警方法，其特征在于所述步骤(3)利用不平衡拉力法计算稳定系数：(3.1)先对选定坡体的滑面，结合坡面特征点决定滑面特征点，对滑面特征点的切线作垂线，交滑体坡面形成不同特征点；

(3.2)连接各特征点形成直线，形成不平衡拉力法的计算条块；

(3.3)不平衡拉力法的假设同不平衡推力法假设；

(3.4)计算时，按条块编号从下往上进行；

(3.5)计算方法为使斜坡体最上面的一个条块推力为零，从而获得不平衡拉力法的各条块下滑力、阻力、剩余拉力和稳定系数；

上述不平衡拉力法的条块划分是为了决定坡体各特征点对应的最小断面，比较各断面所对应的剩余拉力Pi和滑体相应断面对应的最大抗拉力 大小，从而决定各断面是否发生拉破坏，当 时，滑体发生拉破坏；

t

式中σc，Sd分别为滑体或节理、裂隙的抗拉强度和对应的断面面积乘积。

6.如权利要求2所述的基于牵引式斜坡变形破坏机理的稳定性分析和预测预警方法，其特征在于所述步骤(3)利用综合下滑力-抗滑力法计算稳定系数：i

对于情况I或情况II选定的滑体的每个条块，对沿滑面的单元的下滑力Ps 和临界摩阻力Ti分别求其矢量和，形成综合下滑力Ps和综合摩阻力T，且分别与水平轴形成的最小夹角为αs和αf，定义稳定性系数为：式(1)的物理意义为单位下滑力所产生的抗力；

对于情况III，假定计算所得e点为临界状态点，选取可能破坏的滑体为研究对象，对i于选定滑体的每个条块，对沿滑面的单元的下滑力Ps 和可能破坏面的剩余拉力与临界摩阻力Ti和可能破坏面的最大拉力P 分别求其矢量和，形成综合下滑力Ps和综合抗滑力T，且分别与水平轴形成的最小夹角为αs和αf，定义稳定性系数为式(1)；

对于情况III，以最小稳定系数作为可能发生拉破坏的稳定系数；在此基础上分析下一步可能的破坏形式。

7.如权利要求2所述的基于牵引式斜坡变形破坏机理的稳定性分析和预测预警方法，其特征在于所述步骤(3)利用牵引式斜坡的主拉力法计算稳定系数：对于情况I或情况II选定的滑体的每个条块，m-1条块为临界状态条块，求1～m-1条i m m块滑面的剩余拉力Ps 矢量和，形成综合剩余拉力P，且与X轴的夹角为αs ；求m至n每个i m-n条块的下滑力Ps 与临界摩阻力Ti差值的矢量和，形成综合剩余抗滑力T ，其与水平轴的m-n夹角为αf ，定义富余稳定性系数为：

式(2)的物理意义为沿主滑方向单位下滑力的剩余抗力；

对于情况III，假定计算所得e点为临界状态点，选取可能破坏的滑体为研究对象，m-1i m条块为临界状态条块，求1～m-1条块滑面的剩余拉力Ps 矢量和，形成综合剩余拉力P，且m i与X轴的夹角为αs ；求m条块至bf断面对应的每个条块的下滑力Ps 与临界摩阻力Ti差m-n值，以及bf断面最大抗拉力与此时具有的拉力差值矢量和，形成综合剩余抗滑力T ，其与m-n水平轴的夹角为αf ，定义富余稳定性系数为式(2)。

8.如权利要求2所述的基于牵引式斜坡变形破坏机理的稳定性分析和预测预警方法，其特征在于所述步骤(3)利用牵引式斜坡变形稳定系数法计算稳定系数：定义牵引式边坡现状稳定系数，X轴方向的现状稳定系数为边坡后缘关键块发生破坏时沿滑面在X轴方向的位移矢量和Sp-t除以现状边坡前缘至滑面现状临界状态单元或条块在X轴方向的位移矢量和Sc-t，采用相同的方法定义Y、Z轴方向的现状稳定系数；综合稳定系数定义为边坡后缘关键块发生破坏时沿滑面在X、Y、Z轴方向的位移矢量和除以现状边坡前缘至滑面现状临界状态单元或条块在X、Y、Z轴方向的位移矢量和；存在三个不同方x x y y向的稳定系数，X、Y、Z轴三方向的稳定系数分别为Fs-x＝Sc-t/Sp-t、Fs-y＝Sc-t/Sp-t、Fs-z＝z zSc-t/Sp-t。

9.如权利要求2所述的基于牵引式斜坡变形破坏机理的稳定性分析和预测预警方法，其特征在于步骤(4)中所述边坡稳定性计算的滑面边界法，对于牵引式斜坡，取滑体作为计算对象，潜在滑动面上的力与位移以边界条件加以处理，滑面以梁单元或其它单元加以处理，基本假设为梁单元沿横向没有变形或变形遵循不同本构关系，坡体采用现行的各种单元；具体计算步骤如下：第一步，选取坡体作为有限元计算对象，除滑面的正应力与剪应力外，坡体表面边界条件计算方法与常规一致，以滑面底边梁单元为第1单元，沿着滑面底边直至滑面顶面的最后一个梁单元依次排序定义为第N单元；

第二步，在滑面垂直方向赋应变为零，并在此基础上进行斜坡的第一次计算，从而获得滑面梁单元底面第一次计算的法向应力、切向应力和应变，如式(3，4)，i i

式中Tn，s，Sn，s 分别表示滑面第i次计算所得的法向应力、切向应力和法向应变、切向m，i m，i m，i m，i应变矩阵，σn ，τs ，εn 和εs 分别表示第m单元第i次计算所得的法向应力、切向应力、法向应变和切向应变；n表示法向，s表示切向，i表示计算次数；

比较每一单元第一次计算所得的切向应力与计算所得的法向应力相对应的临界摩阻力的大小及方向，当切向应力与临界摩阻力方向相反时，且切向应力绝对值大于临界摩阻力时，相应单元切向应力作为第二次计算时的边界条件取为，临界摩阻力加上该单元计算所得切向力作为该单元下次计算的初始切向应力边界条件，假定1～k单元切向应力与摩阻力方向相反；对于某单元计算所得切向应力方向与摩阻力方向一致或相反，但绝对值小于临界摩阻力，下次计算时不作为初始条件，假定为k+1～N单元，如式(5)，k，1 k，1 k，1△τs ＝τcrit +τs ，L∈(1，k) (5)第三步，进行第二次计算，相对应1～k及1～N个单元力的初始应力及应变边界条件分别如下式，按式(6)和式(7)，进行第二次计算所得的滑面单元的法向应力及切向应力矩阵见式(8)，首先对1～N单元进行判定，按式(8)的法向应力以下式计算对应的临界摩阻力：i，2 i，2

τcrit ＝Ci+σn tanφi i∈(1，N) (9)i，j i，j

式中，τcrit 为第i单元第j次计算所得临界摩阻力，σn 为第i单元第j次计算所得法向应力；

对于已经施加了临界摩阻力的滑面对应单元，比较法向应力大小，当此时计算的法向应力之差大于某定值D1时，重新对该单元切向应力按式(11)赋值，L，2 L，1

|σn -σn |≥D1，L∈(1，k) (10)L，2 L，1 L，1 L，2 L，1

△τs ＝τs +τcrit +(τcrit -τcrit )/2，L∈(1，k) (11)

对于未赋临界摩阻力的单元，当切向应力与临界摩阻力方向相反时，且绝对值大于临界摩阻力时，对这样的单元施加临界摩阻力如下式(12)，i，2 i，2

如第i单元第2次计算，施加的切向应力((τcrit +τs )，假定i∈(k，m))；

i，2 i，2 i，2

△τs ＝τcrit +τs ，i∈(k，m) (12)按照上述计算，在进行第三次计算时，1～N个滑面应力及应变边界条件分别为：第四步：重复第三步，计算直至多次，对于已经施加了临界摩阻力的滑面对应单元，比较法向应力大小，当此时与上一次计算的法向应力之差最大值小于某定值D1时，则可以结束计算，L，MM+1 L，MM

|σn -σn |max

i，MM+1 i，MM+1 i，

此时记下法向应力σn 和切向应力τs ，其中i∈(1，N)，以及切向应变εsMM+1

，其中i∈(1，N)。