多年冻土边坡稳定性影响因素敏感性分析

褚志成，陈鹏辉，雷胜友，郭长志，陈虹宇，王汉文

(长安大学 公路学院，陕西 西安 710064)

0 引 言

多年冻土边坡易受极端气候周期性影响，从而产生局部冻融滑塌，尤其是在全球气候逐年变暖的大环境下，多年冻土地区边坡发生的热融滑塌失稳现象屡见不鲜。因此，针对多年冻土边坡稳定性影响因素进行研究，对该地区边坡防护治理和设计有很好的参考价值。

寒区土质边坡发生热融滑塌的表现形式根据滑塌部位可分为上部滑塌式、中部滑塌式和整体滑塌式[1]。暖季气温回升，边坡浅层土体和坡面积雪的融化使土体处于水分过饱和状态；冷季气温下降，水分向冻结锋面积聚，导致发生冻胀，边坡浅层土体承受反复的融化和冻结是造成寒区边坡发生失稳的主要原因。目前，已有大量关于寒区边坡稳定性分析方法的研究成果。刘红军等[2]分析了土质边坡冻融滑塌失稳的影响因素，表明边坡在春季融期冻结层的存在是冻土边坡发生失稳破坏的主要原因，并通过有效应力法推导了边坡的稳定性系数计算方法；牛富俊等[3]对热融滑塌边坡的地貌、地质及变形特征、影响因素进行了研究，利用静力平衡原理在考虑地下水沿坡面平行方向的渗流推导出热融滑塌型边坡的安全系数计算表达式；武鹤等[4]考虑滑坡面的坡端阻力作用，利用极限平衡理论推导出寒区土质边坡冻融滑塌稳定性的解析表达式；高樯等[5]利用有效应力极限平衡法在不考虑冻土融化过程中孔隙水压力增大的情况下，推导出斜坡稳定性系数的计算公式；孙国栋等[6]考虑冻融循环以及渗流对寒区高陡边坡的作用，提出基于冻融折线型滑移面的概念，并结合框架锚杆支护边坡，得到寒区高陡边坡稳定性评价方法。已有冻土边坡稳定性分析方法大多采用极限平衡法、有效应力等整体稳定性分析方法，但针对多年冻土边坡局部冻融的特殊性，上述方法并不完全适用。

在边坡影响因素敏感性分析方面，目前常用的边坡敏感性分析方法主要有正交试验分析、单因素分析、灰色关联分析等，而单因素分析在多因素交互影响下有一定的局限性。灰色关联理论自1982年提出以来，已在多个领域得到广泛应用。在寒区边坡稳定性分析中，现有研究成果相对较少，于晖等[7]通过组合权重灰色关联法评价多年冻土区路基稳定性，并以此建立了评价体系；范新宇等[8]结合均匀设计和灰色关联分析法研究寒区矿山边坡的影响因素，进一步说明灰色关联分析法对于边坡稳定性敏感度分析的可靠性。基于此，本文以东北大兴安岭地区某一级公路路堑土质边坡为研究对象，针对多年冻土边坡局部冻融的特殊性，提出基于有限元的温度场-应力场(热-力)耦合多年冻土边坡稳定性分析方法，并结合强度折减原理，计算边坡稳定性系数，进一步利用正交试验设计与改进的灰色关联模型，分析多年冻土边坡稳定性影响因素的敏感性大小，以期为该地区边坡防护治理提供一定的参考依据。

1 多年冻土边坡稳定性分析方法

1.1 基于有限元的热-力耦合边坡稳定性分析方法

冻融作用下多年冻土边坡稳定性是一个复杂的动态平衡问题，极端的气候条件使边坡土体性质具有明显的时效性，并且多年冻土边坡整体稳定性较好，坡体浅层易受极端气候反复冻融的影响，而永冻层和不冻层相对较稳定。因此，考虑土体冻融过程中环境温度和冻融循环对边坡土体的综合作用，利用多物理场仿真软件的函数编辑功能编写热-力耦合动态if函数，即利用“T=0 ℃”实时捕捉边坡潜在的冻融交界面，从而解决难以在数值计算中界定动态冻融交界面的问题。

强度折减法具有不需假设滑动面和能够适应各种复杂情况等优点，在边坡稳定性应用中(即在有限元计算中)，使岩土体的黏聚力和内摩擦角不断按比例减小，直至岩土体达到破坏的临界状态，将此时的减小比例即折减系数作为边坡的稳定性系数，即

cF=c/F，

φF=tan-1(tanφ/F)，

(1)

式中:cF，φF分别为折减后的黏聚力和内摩擦角；F为稳定性系数。

强度折减法仅针对边坡整体稳定性，对仅受到局部冻融的多年冻土边坡不能完全适用。本文在热-力耦合动态if函数的基础上，利用强度折减原理仅对发生冻融作用的浅层土体不断折减其黏聚力和内摩擦角，以坡顶位置处竖向位移发生突变作为边坡失稳破坏的判断标准，以此得到考虑各个综合影响因素下的边坡局部稳定性系数。

1.2 正交试验设计

正交试验利用正交表在考虑因素水平间交互作用的情况下，使用较少的试验次数得到最佳的试验搭配方案[9]。在多年冻土边坡稳定性分析中，因现场调研要求较高，正交试验的优点得到了充分运用。本文借助数理统计软件SPSS得到所需的正交试验表。

1.3 灰色关联分析法

1.3.1 灰色关联分析法原理

传统的数理统计分析常常因为需要大样本回归分析、样本数据须满足某个典型的概率分布的方差分析等特殊要求，在实际应用中受到很多限制，不能满足人们的实际需求。灰色关联分析法根据系统行为特征的发展态势将不同事物间的不确定性通过关联度表示出来，即将抽象事物量化。

1.3.2 传统灰色关联分析法[10]

(1)确定表征系统行为特征的映射量和影响系统行为的因素，即Xi=[xi(1)，xi(2)，…，xi(n)]因素序列(子序列X)，Yi=[yi(1)，yi(2)，…，yi(n)]为系统行为特征的映射量序列(母序列Y)。则X和Y可表示为

(2)

(3)

(2)当系统行为特征映射量和各因素间的量纲不一致时，为减小数值差异造成的误差，需要对各序列做进一步的标准化处理。常用的标准化处理措施有Min-max法、Z-score法等。已知各个序列的最小值和最大值，故采用Min-max标准化法，即

(4)

(5)

(3)根据标准化处理后的母序列和子系列矩阵，得到两序列的绝对差值序列，即

Δij=|x′ij-y′ij|。

(6)

(4)根据绝对差值序列，求得Δ中的Δmax和Δmin，即

Q=Δmax=maxΔij，

(7)

q=Δmin=minΔij。

(8)

(5)计算关联系数rij；根据式(9)求得各因素的关联系数，得到关联系数矩阵R，

rij=(q+ρQ)/(Δij+ρQ)，

(9)

式中，ρ为分辨率系数，ρ∈[0，1]。

(6)关联度ri为

(10)

1.3.3 分辨系数取值

分辨系数在灰色关联分析中可以增加各因素间关联系数的差异性、增大关联度之间分辨率。分辨系数的取值直接影响到关联度大小。而传统灰色关联分析中，分辨系数一般取0.5，显然弱化了分辨系数的影响。基于此，定义Δv(k)为因素绝对差值的均值，ε(k)为判断序列是否出现奇异值的指标，即

(11)

ε(k)=Δv(k)/maxΔij(k)，

(12)

以ρ=0.5为区分序列出现奇异值的临界分辨系数值，则分辨系数具体取值据如下确定[11]。

(1)1/ε(k)>3时，表示序列出现奇异值，取ρ∈[ε(k)，1.5ε(k)]，一般取1.5ε(k)。

(2)0≤1/ε(k)≤3时，序列较稳定，分以下3种情况增强maxΔij对关联度的作用，即，

(a)2≤1/ε(k)≤3时，ρ∈(1.5ε(k)，2ε(k)]。

(b)0<1/ε(k)<2时，ρ∈[0.8，1]。

(c)ε(k)=0时，表示ε(k)和ρ之间没有关系，则在[0，1]中取任意值。

确定ρ′后，式(9)变换为

r′ij=(q+ρ′Q)/(Δij+ρ′Q)。

(13)

1.3.4 权重确定

权重表示各影响因素对于边坡稳定性的相对重要程度。由于多年冻土边坡稳定性各影响因素之间差异性较小、交互影响，因此在评价过程中，如果以不考虑权重的关联度表征各因素的主次关系是不严谨的。指标权重确定方法(CRITIC)是一种以对比强度和指标间的冲突性为基础的客观权重赋权法[12]，通过相关系数表示各因素的相关性，依据式(14)可得到指标的客观权重。

(14)

其中lij为指标间的相关系数。

相关系数采用Pearson提出的方法得到，即

(15)

Gi表示各因素所包含的信息量，以因素类间标准差表示各个因素取值大小的差距，标准差越大表示各因素间的差距越大，Gi越大表示该因素的相对重要程度越大[13]。

Gi=σiti，

(16)

式中，σi为各因素的类间标准差。

最后将各因素的权重归一化处理，得到各因素的客观权重ωi，

(17)

考虑多年冻土边坡稳定性各影响因素不同分辨系数和各因素客观权重后，得到改进的灰色关联分析模型，从而评价关联度明显优于式(10),即传统灰色关联模型得到的关联度。改进的灰色关联模型关联度计算公式为

(18)

2 试验因素与范围的选取

2.1 影响边坡稳定性的因素

某一级公路位于东北大兴安岭(内蒙古境内)，据气象资料数据，该区年平均气温-3 ℃，最大冻结深度3.1 m，根据中国季节性冻土标准冻深线，该区属于多年冻土区。根据课题组在该线某段钻孔实测数据表明，地表大多覆盖0.5～2 m的草炭土，以下多为粉质黏土和砂砾。边坡土质为黏性土，整条路线人工边坡数量较大，边坡坡率大多集中在1∶1～1∶3。

与普通地区边坡相比，多年冻土边坡滑坡失稳是长期、极其复杂的过程，不仅受到土质和开挖削坡等影响，而且还受到极端气候引起的反复冻融作用。本文从引起多年冻土边坡发生失稳破坏的角度，总结出影响因素有以下几个。

2.1.1 自身因素

对于土质边坡而言，边坡坡度、密度、土体黏聚力和内摩擦角在工程实际应用中起着极其重要的作用，是影响土质边坡稳定性的基本因素。

2.1.2 环境因素

(1)季节活动层冻融深度。冻融过程中，土体内的水分不断发生迁移：随着气温下降，季节活动层土体逐渐冻结，土层水分开始出现相变、迁移等现象；到了暖季，气温升高使水分在冻融锋面聚集[14]。对于多年冻土边坡，气温升高直接影响坡体温度场的分布，此时边坡处于融化期，季节活动层水分增加，造成边坡浅层土体承载力急剧下降，边坡出现滑塌。所以季节活动层冻融深度在多年冻土边坡稳定性分析中是不可忽略的影响因素。将边坡浅层最大融化深度定义为季节活动层冻融深度，由于冻融深度数据不易采集，本文借助有限元仿真模拟得到边坡的冻融深度。

(2)冻融损伤系数。冻土边坡在温度周期循环作用下坡体浅层经历着反复冻融。大量研究表明，冻融循环对土体力学性质的改变较为突出，长期冻融作用下，土体强度出现不同程度的损伤[15-16]。葛琪[16]给出了冻融作用下，土体抗剪强度指标损伤系数与冻融循环次数的关系:

kc=0.368×e-0.28n+0.64，

kφ=0.379×e-0.39n+0.63，

(19)

式中:kc为黏聚力冻融损伤系数；kφ为内摩擦角冻融损伤系数；n为冻融循环次数。

2.2 基于有限元的边坡稳定性计算

2.2.1 热-力耦合模型的建立

假设边坡土体为均质、各向同性的弹塑性材料，土体屈服准则选取Drucker-Prager准则匹配Mohr-Coulomb条件即DP4准则，温度场采用带相变的非稳态各向同性传热模型。

(1)边界条件。环境温度的冷暖交替直接影响到地层温度场分布，主要取决于地表与大气的对流换热过程。在不考虑大气逐渐升温的情况下，地面温度边界条件为[17]

(20)

式中:t为时间，s；Ta为气温；Ts为地表温度；TM为年平均气温；A为年平均气温较差。

年平均气温与年较差根据当地气象数据取为-3 ℃和50 ℃。模型底部考虑地热流的影响，施加第二类边界条件-λ(∂T/∂n)=0.02；模型两侧设为绝热边界条件:(∂T/∂n)=0；应力分析模型侧面为法向约束，底部设置为全位移固定约束。

(2)初始条件。初始温度场对整个坡体状态的影响较为重要，本文将实测年平均地温作为整个坡体初始温度场计算的初始值，计算30 a温度场，稳定后提取此时温度场作为热-力耦合阶段的初始温度场。

(3)热-力耦合机制。土体冻结时，其力学参数抗剪强度、泊松比和弹性模量与土体温度的关系为(21)[18]，

(21)

式中:ai，bi均为试验常数i=1,2,3,4，土体温度T>0 ℃时，bi=0；m=0.6；ai的取值对应下文的不同工况，其中b1=27，b2=-0.008，b3=0.09，b4=0.6。

2.2.2 天然温度场验证

(1)参数选取。根据实测钻孔数据，该段天然地层从上依次为0.3 草炭土、2.2 粉质黏土、2.2 碎石土，以下为砂砾。按照实际地层情况建立80 m×40 m的矩形几何模型。限于篇幅，几何模型不再展示。土体热学参数密度ρ、比热C和导热系数的取值见表1。

(2)结果验证。图1为天然地表钻孔实测温度与数值模拟地温值对比曲线。由图1可知，模拟值与实测值随深度的变化曲线吻合度较高，除了两个不同时间点在接近地表位置差异较小(可能是接近地表位置测温钢管内的温度受到外界空气的影响较大所致)外，在冻土上限以下，基本一致，证明了本文模型的适用性和合理性。也可得该段年平均地温-1 ℃，年变化深度在地表以下9 m处。

表1 地基土体热学参数

图1 温度实测值与模拟值对比曲线

2.2.3 季节活动层冻融深度的确定

通过对该段某处边坡的简化处理，建立边坡分析模型，边坡高度15 m，坡度45°，密度1.816×103kg/m3，黏聚力34.1 kPa，内摩擦角20.9°，弹性模量20 MPa，泊松比0.3，几何模型见图2。

图2 边坡几何模型

下文以坡度45°、密度1.816×103kg/m3的边坡阐述冻融深度的计算过程。图3为边坡达到最大融化深度时的温度场云图(“黑线”为“0 ℃等温线”，即融化深度线)。图4为边坡融化深度和环境温度随时间变化曲线，可看出坡体融化深度与环境温度有很强的相关性，但在时间上与环境温度存在一定的滞后性。9月份时坡体达到最大融化深度2.86 m，将此时的融化深度作为季节活动层冻融深度。基于此，通过有限元仿真得到研究区段不同形态边坡冻融深度的取值范围,为2.70～3.15 m。

图3 边坡达到最大融化深度时温度场云图

图4 融化深度与环境温度变化曲线

综上，本文主要研究边坡坡度X1、密度X2、季节活动层冻融深度X3、土体黏聚力X4、内摩擦角X5和冻融损伤系数X6对多年冻土质边坡稳定性的影响。参数具体取值的范围如表2所示。

表2 边坡参数取值范围

3 正交试验结果的灰色关联分析

为有效获得上述6种影响因素的正交试验设计表，每个因素选取5个水平，通过SPSS软件得到L25(56)正交试验表。表3为各影响边坡稳定性因素及水平,正交试验结果如表4所示。

表3 影响因素水平

表4 正交试验设计结果

续表4

在正交试验设计数据的基础上，根据式(2)～(3)得到子序列X和母序列Y，

按照式(4)～(5)将X序列和Y序列标准化处理后，通过式(6)得到差值序列Δ矩阵为

将式(9)代入处理后的分辨系数，得到关联系数矩阵

最后，得到各因素分辨系数与权重，见表5，关联度结果见表6。

表5 各因素分辨系数与权重

表6 关联度计算结果

根据关联度计算结果可知，多年冻土边坡影响因素敏感程度大小依次为黏聚力>冻融损伤系数>季节活动层冻融深度>内摩擦角>坡度>密度。从该排序可以看出，土体黏聚力对多年冻土边坡稳定性最为重要，其次为季节活动层冻融深度和冻融损伤系数，而边坡坡度、内摩擦角和土体重度对边坡稳定性的影响相对较小。此结论也表明在多年冻土边坡防护治理中，除了加强防排水措施外，应通过有效的保温措施减少气候环境对边坡的影响，例如植物防护、土工格室植被防护等措施既能有效避免坡面长期暴露，又可利用植物根系的加筋作用提高土体抗剪强度，增加多年冻土边坡的稳定性。

4 结 论

(1)考虑冻融作用对边坡土体的影响，在影响边坡稳定性自身因素的基础上，引入季节活动层冻融深度和冻融损伤系数，为冻土边坡影响因素敏感性分析提供重要参考依据。

(2)针对多年冻土边坡局部冻融的特殊性，根据实测数据验证模型的正确性。基于有限元提出多年冻土边坡稳定性分析方法，并结合强度折减原理得到多年冻土边坡局部稳定性系数，为多年冻土边坡稳定性分析提供一种新方法。

(3)考虑分辨系数和权重对各因素间差异性和系统整体性的影响，将边坡稳定性系数作为母序列，边坡稳定性影响因素作为子序列，改进传统的灰色关联模型，以正交设计结果为源数据，较好地运用到多年冻土边坡稳定性影响因素敏感度分析中，并得到较合理的分析结果。

(4)依据灰色关联分析结果，得到影响多年冻土边坡稳定性的6个因素，其敏感程度大小排序依次为黏聚力，冻融损伤系数，冻融深度，内摩擦角，坡度，密度；在多年冻土边坡的防护治理中，应从热学角度考虑，减少外界极端气候对边坡的影响。本文方法和结论可为多年冻土边坡稳定性分析和治理提供参考。