基于上限分析的渣土边坡抗滑桩桩位优化研究

肖林超，陈文强2，姚文敏2，熊 爽2，曾江波

(1.深圳市勘察测绘院有限公司，广东 深圳 518028； 2.中国地质大学(武汉) 工程学院，武汉 430074)

1 研究背景

近20 a来我国城市化进程发展迅速，将城市建设活动中产生的大量以弃土弃渣为主的余泥渣土进行堆填形成新的城市用地，已成为一种使城市可持续发展的新思路，如天津市历时3 a建成的以人工渣土堆填山体堆积而成的占地40万m2的南翠屏公园[1]。然而余泥渣土的物理力学性质差，加之多种外部因素的影响，渣土堆填边坡若处理不当极易诱发滑坡灾害，给邻近的城市建筑和人民生命财产带来巨大的损失。如2015年12月20日，深圳光明新区红坳淤泥渣土收纳场发生的特大滑坡事故，造成了73人死亡、4人失踪，33栋建筑物被损毁或掩埋，直接经济损失人民币8.8亿余元[2]。因此，对余泥渣土堆填边坡进行稳定性评价和治理工程的相关研究具有重大意义。

作为一种常用的加固余泥渣土堆填边坡的抗滑措施，在实际工程中，抗滑桩的桩位设计通常借鉴规范中的经验参数，这样往往使得设计偏于保守或者不合理，从而造成实际治理工程造价的浪费。国内外许多学者对边坡抗滑桩的桩位优化问题进行了相关研究，研究方法目前主要包括理论计算、数值模拟和物理模型试验3种方法[3]。理论计算方法主要分为极限平衡法、极限分析法。Ito等[4]、Poulos[5]借助传统极限平衡方法得出抗滑桩布设在边坡中部附近对提高边坡稳定性最有利。基于极限分析方法，国内外许多学者对均质边坡进行了桩位优化研究，结果表明当抗滑桩布设在边坡坡脚的一定范围内时，对于给定设计安全系数条件下的抗滑桩的受力情况最有利，即抗滑桩此时受力最小[6-9]。在数值模拟方面，雷文杰等[10]采用有限元强度折减法研究得出单排桩布置在边坡中部时，边坡稳定性最好，但是抗滑桩的受力最大，并且边坡破坏模式随抗滑桩桩位变化而变化。Wei和Cheng[11]、夏永成[12]分别采用FLAC3D、ABAQUS和Geo-slope软件进行研究也得出与雷文杰等相似的结论。而在物理模型试验方面，许多学者也开展了一些颇有意义的相关研究[13-15]。

从上述研究可以看出，目前关于抗滑桩加固边坡的合理桩位尚无一致的结论，另外对于加固渣土边坡的抗滑桩最优桩位的研究相对较少。将现有的极限分析上限法应用于分析抗滑桩加固边坡问题时，常需要较复杂的程序求解或者借助非线性规划方法[6-8]，但目前采用的大多数非线性规划算法为序列二次规划(Sequence Quadratic Program，SQP)算法[9]，该算法对计算初值具有很强的依赖性，要求出真实可靠的解比较费时费力。本文采用遗传算法与非线性规划SQP算法相结合的方法，可以较好地解决该问题。

基于强度折减原理与极限分析上限法，本文对单排抗滑桩加固渣土边坡的最优桩位问题进行了分析研究。采用遗传算法与非线性规划SQP算法相结合的方法，方便快捷地求解出抗滑桩加固渣土边坡的极限分析上限解。选取深圳市某渣土堆填边坡为研究对象，对某一设计安全系数下的单排抗滑桩需提供的抗滑力与布设桩位的关系进行了求解分析，结合不同桩位处抗滑桩滑裂面以上的桩长、滑裂面变化规律以及局部稳定性变化情况，综合确定案例渣土边坡的抗滑桩最优布设桩位。

2 抗滑桩加固渣土边坡稳定性上限分析方法

2.1 强度折减原理

强度折减法是一种用来求解边坡稳定安全系数的常用方法，该方法最早由Bishop[16]提出。强度折减技术的原理是将土体抗剪强度参数(c，φ)同时除以一个比较小的初始折减系数(Fs)，这样就得到了一组新的抗剪强度参数(cf，φf)；然后采用岩土体新的参数进行试算边坡的稳定性，不断增大折减系数(Fs)直至达到某一值时，若折减系数再大一些则边坡就将发生失稳破坏，那么发生破坏之前的那个折减系数就为边坡的稳定性系数。强度折减原理用数学表达式可表示为

(1)

式中：cf，φf分别为经强度折减后的土体黏聚力和内摩擦角；c，φ分别为土体黏聚力和内摩擦角。

2.2 渣土边坡抗滑桩稳定性上限分析

本文在基于强度折减原理与极限分析上限法对抗滑桩加固的渣土边坡桩位优化分析时，引入如下假定：①渣土边坡足够长，将本问题简化为平面应变问题来进行分析；②渣土边坡填土均质，各向同性；③渣土边坡滑动土体为理想刚塑性体，破坏时满足摩尔-库伦屈服条件和相关联流动法则；④抗滑桩长度足够长，假定其为刚性桩。

当考虑抗滑桩的加固效应时，单排抗滑桩加固渣土边坡的破坏模式采用对数螺旋线旋转破坏模式[17]，抗滑桩桩侧有效土压力分布模式采用水平方向分布模式[6，8]，对数螺旋线旋转破坏模式示意图见图1。图1中：H为边坡坡高；Lx为边坡坡面AC在水平方向上的投影长度；xF为抗滑桩与边坡坡趾间的水平距离；L为边坡坡顶与坡面的交点A和边坡潜在破坏面与坡顶的交点B间的距离；D为边坡坡趾C点和潜在破坏面与坡脚的交点C′间的距离。

图1 单排桩加固的渣土边坡旋转破坏机构示意图[8]Fig.1 Rigid rotation failure mechanism for construction solid waste landfill slope reinforced with one row of piles[8]

图1中在自重条件下，刚性块体ABC′C绕旋转中心O点相对对数螺旋破坏面BC′以下的静止土体作刚性转动，旋转角速度为Ω。采用对数螺旋破坏模式时，速度间断面BC′上任意点的旋转半径r都满足对数螺旋线方程，用数学表达式表示为[17]

r=r0e(θ-θ0)tanφf。

(2)

式中：r，θ为速度间断面上任一点关于旋转中心O点的旋转半径和旋转角；r0，θ0为破坏面与坡顶的交点B点关于O点的初始旋转半径和旋转角。

由图1中几何关系，不难得出H/r0，L/r0，D可通过变量θ0，θh，β′表示为[17]：

sinθ0+α] ；

(3)

[sinθh+αeθh-θ0tanφf-sinθ0+α] ；

(4)

(5)

式中：α，β，β′分别为边坡坡顶AB、边坡坡面AC、边坡A点与C′点的连线和水平线的夹角；θh为旋转半径rh关于O点的旋转角。

根据极限分析上限定理，如图1所示的抗滑桩加固渣土边坡旋转破坏机构中，整个体系能耗计算主要包括外力做功功率(包括刚性块体ABC′C重力做功功率WG和抗滑桩有效抗力做功功率WF)和发生在速度间断面上的内能耗散功率Wint。刚性块体ABC′C重力做功功率WG的数学表达式为[17]

(6)

式中：γ为土体重度；f1，f2，f3，f4均为关于θ0，θh，β′，Fs，φ，α，β的函数，具体表达式如下[17]。

(7)

(8)

(9)

(10)

如图1所示，当抗滑桩布置在与渣土边坡坡趾水平距离为xF的位置时，其与坡面的交点E与其与破坏面的交点F间的距离h即为抗滑桩在潜在破坏面内的长度。由几何关系，可得h与抗滑桩桩位xF间的关系为[6]：

(11)

xF=rFcosθF-rhcosθh-D。

(12)

式中rF，θF分别为抗滑桩与破坏面的交点F在旋转破坏机构中所对应的旋转半径和旋转角。

对于抗滑桩桩侧有效土压力分布模式采用水平方向分布模式且考虑抗滑桩有效抗力的弯矩效应时，可计算得抗滑桩有效抗力做功功率WF的表达式为[6]

WF=Fp(sinθFrF-mh)Ω。

(13)

式中：Fp为抗滑桩有效抗力合力；m为表征弯矩效应的经验参数，其值与桩侧有效土压力分布模式有关。一般来说，桩侧有效土压力分布模式为三角形、矩形、梯形时，m分别取1/3，1/2，7/20[18]，特殊地，当m取0时则表示只考虑抗滑桩与潜在破坏面相交处的剪力作用。

发生在速度间断面BC′上的内能耗散功率Wint可由黏聚力c与速度间断面上的微分面积dθ/cosφ和切向速度Vcosφ之乘积沿该面积分得到[6]，具体表达式为

式中f5为关于θ0,θh,Fs，φ的函数。

根据上限定理，使体系外力做功功率等于内能耗散功率，便可得渣土边坡达到临界失稳状态时抗滑桩桩侧有效抗力Fp的上限表达式[6]，即

(15)

式(15)即为使渣土边坡稳定性提高到某一设计安全系数下，抗滑桩需提供的有效抗滑力的表达式。可以看出当安全系数给定情况下桩侧有效抗力Fp为关于变量θ0，θh，θF，β′,Fs的函数；不难看出Fp与Fs两者之间可以相互转化，即当Fp已知时，相应的Fs也可以求解出来。由上限定理可知，当渣土边坡稳定性达到临界状态时，有效抗滑力取得最大值，也即上限值。为求得单排抗滑桩所需提供的有效抗力上限解，本文通过选取有效抗力Fp为目标函数，通过遗传算法与非线性规划SQP算法相结合的方法，编写了相应的优化求解程序，采用的数学规划表达式为：

maxFp=Fpθ0,θh,θF,β′ ；

(16)

式中：xF为依据式(12)计算得出的抗滑桩布设桩位理论值，而xFact为实际给定的抗滑桩布设桩位。通过式(16)、式(17)求得变量θ0，θh，θF，β′的上限解后，代入式(3)—式(5)和式(12)便可求出单排抗滑桩加固后的渣土边坡潜在破坏面与抗滑桩在滑面以内的桩长。

3 案例分析

3.1 工程概况

深圳市某余泥渣土受纳场地原始地貌为丘陵和壑沟，因场地大规模的人工堆填余泥渣土，现状地形已发生巨大改变，场地中间现地形较为平缓，大致呈自西南向东北倾伏，在东侧、西侧及北侧边缘顺接原来冲沟的位置形成了边坡，地势相对较低。

根据钻探成果揭露，场地内地层自上而下依次为以余泥渣土为主的人工填土层(Qml)、以黏土为主的第四系全新统冲洪积层(Q4al+pl)、以含砂粉质黏土为主的第四系上更新统坡积层(Q3dl)、以砂质粉质黏土为主的第四系中更新统残积层(Q2el)，以及不同风化程度的花岗片麻岩组成的燕山期基岩层。

场地内分布有雨源型沟流及面流，场地内径流主要为临时性的地表径流。受地形地貌及气候影响，具有源短流急、流水暴涨暴落特性。场地内上层滞水大多分布在第四系土层中，场地内其储量较少，接受降雨补给，受季节性影响；基岩裂隙水主要赋存在风化岩中，水量一般，主要间接补给来源为大气降水。场地内地下水位埋深较深。

本文选取渣土受纳场东侧一高约10 m，坡度为25°的边坡为研究对象，该边坡沿走向长达200余米，且边坡所受外力皆垂直于走向，沿边坡走向改变较小，可认为案例边坡满足平面应变条件假设，边坡地质剖面图如图2所示。坡体主要由深厚淤泥渣土组成，碾压较差，结构较松散，整体上目前稳定性较差。根据《建筑边坡工程技术规范》(GB 50330—2013)标准划分，该边坡安全等级可划分为一级。

图2 渣土受纳场案例边坡地质剖面图Fig.2 Geological profile of the case slope of a construction solid waste landfill

图3 渣土受纳场案例边坡概化模型Fig.3 Schematic diagram of the simplified calculation model of the case slope

3.2 抗滑桩加固的渣土边坡概化模型

由于案例边坡的稳定性主要受深厚的人工填土层控制，因此本节仅考虑人工填土层的稳定性，基于如图3所示的边坡概化模型，对案例边坡进行抗滑桩桩位优化分析。考虑到目前国内外关于暴雨条件下的抗滑桩相关分析的研究相对较少，因此本节将对天然和暴雨工况下渣土边坡抗滑桩桩位进行优化分析。需要说明的是案例边坡坡体厚度较小，地下水埋藏较深，可认为暴雨工况下当该边坡岩土体由于降雨入渗达到整体饱和状态时即为案例边坡最危险状态；此外，由于降雨入渗机制较为复杂，研究理论尚不完善，因此本文对暴雨工况下的案例边坡进行抗滑桩优化分析研究时，只通过采用岩土体饱和状态参数(包括重度和抗剪强度指标)来代表暴雨工况的影响，作为一种简化分析方法，此种方法也被许多研究者所采用[19-23]。计算中采用的素填土物理力学参数为：天然重度γ=19.0 kN/m3，天然土体黏聚力和内摩擦角分别为c=10 kPa，φ=15°；饱和重度γsat=20.0 kN/m3，饱和土体黏聚力和内摩擦角分别为csat=8 kPa,φsat=12°。案例边坡安全等级可划分为一级，故边坡设计安全系数取1.35。

运用本文编写的优化程序，可求得该案例渣土边坡天然和暴雨工况下的稳定性系数分别为1.203和0.939，即该边坡天然和暴雨工况下的稳定性状况分别为欠稳定与不稳定，因此有必要采用抗滑桩的方式来对该边坡进行加固。本文引入无量纲抗滑桩桩位系数ε=xF/Lx(如图3所示，Lx为边坡坡趾距坡肩的水平距离)来代表抗滑桩的布设桩位。当抗滑桩布设在太靠坡脚或坡顶时，边坡易发生局部失稳，造成单排抗滑桩加固效果不明显，因此本文重点研究当ε介于0～0.8区间内的抗滑桩桩位对边坡加固效果的影响，并且局部失稳情况将在讨论环节中进行讨论。

图4 抗滑桩桩位系数ε与桩侧有效抗力Fp和 滑裂面以上桩长h的关系曲线Fig.4 Relation of pile location coefficient(ε) versus effective side resistance(Fp)and pile length(h) above slip surface

3.3 渣土边坡抗滑桩桩位优化上限结果分析

对不同工况下布设在不同桩位处的抗滑桩加固的案例边坡，运用本文所提出的渣土边坡抗滑桩桩位优化分析模型进行分析，可获得抗滑桩布设桩位与给定安全系数下抗滑桩所需提供的抗滑力(桩侧有效抗力)、最可能滑裂面以上抗滑桩桩长以及最可能滑裂面之间的变化关系，如图4和图5所示。

图5 抗滑桩桩位对最可能滑裂面参数(L和D)的影响Fig.5 Impact of pile location on parameters (L and D) of the critical slip surface

从图4(a)可看出抗滑桩桩位系数ε与桩侧有效抗力Fp间的变化关系很好地符合指数函数关系，2种工况下的拟合优度R2都高达0.999。随着抗滑桩桩位系数ε由0增大到0.4左右，抗滑桩需提供的抗滑力变化较小；当ε>0.4后，抗滑桩需提供的抗滑力急剧增大。因此，从抗滑桩所需提供的抗滑力来看，当抗滑桩桩位系数ε介于0～0.4之间时，抗滑桩取得较合适的桩位。另外，我们也可以看出，暴雨工况下抗滑桩需提供的抗滑力要明显大于天然工况下的抗滑桩需提供的抗滑力，而抗滑桩需提供的抗滑力与抗滑桩的截面设计密切相关，因此在实际抗滑桩设计工程中需要注意该问题。

由图4(b)可知，抗滑桩滑裂面以上桩长随抗滑桩桩位系数ε的增大呈先增大后减小的变化趋势。这种变化关系不仅与对数螺旋线型滑裂面的自身特点有关，而且还与抗滑桩桩位系数ε的改变所引起的案例边坡最可能滑动面的改变有关。抗滑桩滑裂面以上桩长与抗滑桩桩位系数ε的变化关系较好地符合抛物线变化规律，2种工况下的拟合优度R2都高达0.983，抗滑桩滑裂面以上桩长在ε=0.5～0.6取得最大值。因此，从抗滑桩滑裂面以上桩长来看，较合适的抗滑桩布设桩位应在案例边坡中部附近，越靠近坡脚，抗滑桩滑裂面以上桩长越小，抗滑桩造价也越低。另外，相同抗滑桩布设桩位下，暴雨工况下抗滑桩滑裂面以上桩长比自然工况下的桩长大一半左右。

从图5可看出，抗滑桩桩位对案例边坡最可能滑裂面有较显著影响。抗滑桩桩位系数ε增大，最可能滑裂面参数L和D都会减小，说明最可能滑裂面的范围逐渐减小。此外，抗滑桩桩位对暴雨工况下的边坡最可能滑裂面的影响大于自然工况。

4 讨 论

由上述结果分析可知，当抗滑桩布设在案例渣土边坡越靠近坡趾的位置时，抗滑桩为使案例渣土边坡稳定性提高到给定安全系数水平下所需提供的抗滑力越小，滑裂面以上抗滑桩桩长也越小。但是，当抗滑桩布设太靠近坡脚时，边坡很可能会发生抗滑桩的越顶破坏[10]或桩前土体的局部失稳破坏。Li 等[24]研究发现当抗滑桩布设在与边坡坡脚距离达到0.7倍坡面水平距离左右时，桩前土体的局部失稳破坏才先于越顶破坏发生。因此，本节重点研究抗滑桩桩位系数ε=0～0.6时抗滑桩加固渣土边坡的越顶破坏局部稳定性和考虑局部失稳的整体稳定性变化。

从图6可看出，抗滑桩加固渣土边坡的越顶破坏局部稳定性与抗滑桩未加固时2种不同工况下稳定性差异相对应，都表现为暴雨工况下小于天然工况下。2种工况下的局部稳定性系数Fslocal都随抗滑桩桩位系数ε增大呈下凸型曲线增大，而整体稳定性系数也随抗滑桩桩位系数ε增大由低于设计安全系数水平逐渐转变为达到设计安全系数。当ε=0.3和ε=0.6时(更确切地说，由局部稳定性系数变化曲线来看，ε=0.25和ε=0.55左右时)，天然和暴雨工况下的渣土边坡整体稳定性系数达到设计安全系数条件，满足安全要求。因此，为使边坡整体稳定性达到安全要求，并且抗滑桩需提供的抗滑力、滑裂面以上桩长、最可能滑裂面也取得较优值，对于案例渣土边坡而言，天然和暴雨工况下抗滑桩最优桩位分别位于ε=0.3和ε=0.6附近。

图6 考虑局部失稳的抗滑桩桩位对案例渣土边坡局部 和整体稳定性系数的影响Fig.6 Impact of pile location on the local and holistic safety of factor of the case slope in consideration of local instability

5 结 论

本文基于强度折减法和极限分析上限法，通过采用遗传算法结合非线性规划SQP算法，编写了求解极限分析上限解的优化程序，提出抗滑桩加固的渣土边坡桩位优化方法，所得结论如下：

(1)给定安全系数下，抗滑桩桩位对抗滑桩需提供的抗滑力、滑裂面以上桩长、最可能滑裂面有显著影响。抗滑桩桩位系数ε增大，抗滑桩需提供的抗滑力呈指数式增大，滑裂面以上桩长呈抛物线形先增大后减小的变化趋势，渣土边坡最可能滑裂面逐渐减小。

(2)为使案例渣土边坡稳定性提高到给定安全系数，暴雨工况下抗滑桩需提供的抗滑力和滑裂面以上桩长都大于自然工况下的结果，暴雨工况下渣土边坡最可能滑裂面范围也比自然工况下的最可能滑裂面范围大。

(3)综合考虑抗滑桩需提供的抗滑力、滑裂面以上桩长、最可能滑裂面范围以及局部稳定性问题，天然和暴雨工况下加固案例渣土边坡的抗滑桩最优桩位分别位于桩位系数ε=0.3和ε=0.6附近。

(4)基于案例研究，计算结果表明本文所提出的遗传算法结合非线性规划SQP算法的求解方法，能有效而快速地计算出基于极限分析上限法的渣土边坡稳定性系数和最可能滑动面，为渣土边坡防治设计优化奠定基础。