基于MATLAB的基坑整体稳定性分析

周建平，李银，任勇，王明阳，马海春

（1.中国水利水电第七工程局有限公司，四川 成都 231400；
2.合肥工业大学资源与环境工程学院，安徽 合肥 230009）

随着社会的发展和进步，高层建筑和地下工程的数量与日俱增。紧随而来基坑开挖及支护成为的问题。在空旷区域放坡开挖是最经济、最简单的施工方法[13-14]。在基坑支护设计中，基坑土体的整体稳定性尤为重要。分析粘性土整体稳定性的方法有很多，极限平衡法已被用于评估斜坡的稳定性[1][11-12][15]。

目前，边坡工程中常用的安全计算方法是费伦纽斯（W.Fellenius）在1927～1937年提出的瑞典条分法[8][10]。瑞典条分法具有原理简单、公式简洁等优点，但瑞典条分法的计算过程繁琐，计算量巨大[3-4]，需要反复试算圆心位置与半径，直至逼近最危险滑动圆弧。

本文利用MATLAB数学软件对基坑整体稳定性进行编程，结合瑞典条分法理论综合考虑土层、地下水、地面荷载等多种影响基坑边坡开挖整体稳定性的因素，大大减少了基坑设计人员的计算工作，并结合具体实例验证了结果的可行性，提高了基坑稳定性计算效率。

1.1 基坑稳定性分析

对具有围护结构的粘性土边坡或土质边坡分析，其目的是设计合理的埋设深度或检查所提出的围护结构在给定条件下是否稳定合理[2][6-7]。基坑侧壁的剪应力在自重和地面超载作用下超过土体抗剪强度时，土体将失去稳定性，发生滑移。对于具有加固和支撑作用的粘性土，土会滑动形成滑动面。通常滑动面会被认为是圆弧形以方便计算，这种形式的破坏称为整体剪切破坏。基坑整体稳定性分析是分析基坑侧壁是否会发生整体剪切破坏。通常用圆弧滑动稳定安全系数Ks表示，Ks由《国家行业标准—建筑基坑支护技术规程》(JGJ120-2012)规定。锚拉式围护结构、土钉墙、重力式水泥土墙的稳定性均可采用圆弧滑片法分析，对于安全等级一级、二级和三级的基坑工程，安全系数Ks应分别不小于1.35、1.3、1.25。

1.2 圆弧滑动条分法

圆弧滑动条分法有两个基本假设（Navarrov，2014）：土体的剪切面是贯穿土体内圆弧面，安全系数是圆弧包含的剪切面上提供的抗滑力矩和滑动力矩之比；
土条间作用力对边坡整体稳定性影响不大[5]。以最典型的均质土坡为例，对瑞典条分法的计算过程进行说明。

如图1所示，ABC土体为分析对象，内摩擦角为φ，黏聚力为c。滑动面为圆弧AC，圆心为O，半径为R。将滑动土ABC分成n个土条，随意取出一条土条j。土条j的宽度为bj，土条中心高度为hj，土条的滑动弧长为lj，对土条j的受力分析，土条受重力Gj作用，下部土提供的支撑力Nj，下部土提供的摩擦阻力Tj。重力向正方向分解导致土条的滑移力Gjsinθj，对土体下部的压力 Gjcosθj等于支撑力，与土条的压力相反。

图1 瑞典条分法

滑动土体下滑力矩为：

在公式（1a）中，Lj为土条j重心到滑动圆心距离。由于每个土条中心位置都不一样，给计算带来困难。为了计算简便，采用半径R作为下滑力到圆心的距离。这样，运用公式（1b）计算出的滑动力矩更大，安全系数更高。

滑动土体的抗滑力矩为：

安全系数的计算公式如下：

根据多年的工程实践经验，按照条分法计算安全系数一般比极限平衡法计算的安全系数低10%～20%。对于滑动面圆弧圆心角较大或孔隙水压力较大时，安全系数将远小于其他方法，对于实际工程，瑞典条分法的计算值是偏安全的。

2.1 算法描述

在编写程序之前描述算法，有助于在程序开发过程中保持逻辑清晰，使编程更高效。用流程图描述了基坑整体稳定性分析流程，先考虑最简单的情况，即无荷载作用、无地下水的均质土体一级放坡的土坡，且假设圆弧滑动面过坡脚。

如图2所示，在输入参数后，设置关于圆心坐标(ox，oy)的双重循环，循环的范围为 ox∈[X min，X max]，oy∈[Y min,Y max]。每次循环增量为step，step值设置的越小，搜索到圆心的位置就越精确。但相应的计算量会指数级上升，所以设置一个合理的大小即可。对于每一次循环中的圆心坐标(ox，oy)，计算所需与圆心坐标有关的参数，便可计算该圆心坐标对应的安全系数。安全系数所需的抗滑合力sum1和下滑合力sum2是通过土条从1到n循环，累计每个土条的抗滑分力clj+γbjhjcosθjtanφ和每个土条的下滑分力γbjhjsinθj得到。其中，n是土条划分的条数，划分的份数越多，抗滑合力sum1和下滑合力sum2计算的便越准确。

图2 全局稳定性算法流程图

完成圆心坐标(ox，oy)的双重循环结束后，得到一个存储安全系数的二维数组。数组的索引值与圆心坐标有关，通过在数组中寻找最小值，并根据其索引值计算圆心坐标，便得到最危险滑动圆弧的位置和土坡的整体稳定性安全系数。

2.2 坐标系的建立及参数计算

如图3所示，为便于计算，选择坡的顶点作为原点坐标，建立平面直角坐标系。

图3 斜率计算模型

与圆心坐标(ox，oy）有关的参数计算如下。

坡脚C坐标：

其中kb为放坡坡度，宽度与高度之比；
h为基坑深度。

滑动弧半径：根据勾股定理计算，得

滑动圆弧与地面交点A横坐标：

联立滑动圆弧方程与地面方程

解得：

土条宽度：为滑动土体宽度1/n

与条数j有关的参数计算如下。

土条滑动面中心坐标：

土条滑动面弧长：近似取土条端点连线长度

土条重量：从滑弧中心至地面的完整土条的重量减去被开挖部分挖去的楔形土条的重量。

2.3 程序编写与参数调整

2.3.1 编程

根据流程图所示的算法，可以编写此程序。对于计算结果输出，我们可以绘制圆弧滑动安全系数与圆心坐标的关系图像。Ks=ƒ(ox,oy)由于Ks值是若干散点值，因此可以使用内置函数mesh⁃gird()创建绘图网格，并使用surf()函数绘制三维曲面。我们还可以对图像进行调整，增加坐标轴名称，增加网格，色彩插值处理，标记最小圆心位置等。

需输入的参数包括放坡坡度与基坑几何深度，土体重度值，内聚力及内摩擦角。需要调整的参数为：搜索范围ox∈[X min,X max]，oy∈[Y min,Y max]、搜索步长step和土条数n。这些参数可以根据实际需要进行更改。

搜索范围过大会导致无意义的计算，搜索范围过小的可能会遗漏Ks最小值。根据理论分析，滑动面的中心只会出现在坐标系的第一象限。经过调试，将搜索范围设置为 ox∈[0,10]，oy∈[0,10]，即可保证最小值Ks在此范围内，并且能较好的表现出三维表面变化的趋势。

搜索步长step越小，计算精度越高，搜索到的圆心位置越精确，对Ks影响不大，但同时计算所需时间越长。step值缩小为原来的N倍，计算量约为原来的N2倍。经调试，搜索步长step设置为0.1～0.01即可满足计算要求。

土条划分数量n越小，对土条各项参数的计算结果越准确。土条划分数量n值对Ks有影响。土条划分数量n越大，算得Ks越大。这是由于土条滑动弧长近似取土条底面端点连线长度，比实际弧长短，计算的抗滑力小于实际的抗滑力，计算结果偏安全。土条划分数量n值越大，计算的抗滑力越接近真值。此外，土条划分数量n的值越大，三维曲面越光滑，圆心坐标也会发生变化。土条划分数量n值的增加对计算时间的增加没有step值的敏感，因为n只涉及一重循环，且循环内的代码较少。经调试，土条划分数量n设为1000片是较为合适的值。

2.3.2 功能扩展

在完成基坑放坡整体稳定性分析最简单的程序之后，可以考虑更复杂的情况，使该程序能够进行更多种类的基坑放坡整体稳定性分析。

①土体分层

绝大多数基坑场地土体是分层的，每层土的厚度、重度、内聚力和内摩擦角都不同，需要设置并输入这些参数。土体分层影响土条受到的重力，应该按照土层分层计算，计算公式为 Gj=Σhjiγibj，hji是土条 j在 i层土内的高度，γi是 i层土的重度。土体分层还影响土条受到的抗滑力，根据滑动面中心纵坐标判断滑动面在哪一土层内，取该土层的内聚力和内摩擦角。

②坡顶有超载

大多数基坑周围存在建筑荷载、车辆荷载和建筑材料堆载。为了计算简便，设地面超载为均布荷载。输入均布荷载参数qj，可以发现均布荷载作用于地面，即y轴左侧，可以根据土条滑动面中心横坐标是否小于0判断。若xc(j)<0，土条受到地面荷载是qj=qj，若xc(j)>0，土条受到地面荷载是qj=0。地面超载作用点在土条顶面，与圆心的距离和重力与圆心的距离相同，且作用方向也与重力相同。为计算方便，可将土条受到的地面超载加到重力上计算。

考虑坡顶超载，安全系数为：

③存在地下水

当土体中存在有地下水时，地下水位以下的粉土、砂土、碎石土会受到地下水的浮力作用，粘土和地下水位以上土体无浮力作用。因此，需要设置进入每一层土壤的类别cli，如果土层i是粘土，cli=1。如果土层i不是粘土，cli=0。由于基坑降水措施，基坑内外水位不同，需要设置并输入基坑内外水位参数值。基坑内水位为hwn，基坑外水位为hww。浮力作用即孔隙压力减小了土条与下部土体之间的内聚力，造成抗滑力降低，而土条压力不变，安全系数变小，基坑偏向于危险。根据土条滑动面中心的横坐标xc(j)判断是基坑内还是基坑外，根据土条滑动面中心的纵坐标yc(j)判断滑动面所处土层，若该土层为粘土或高于水位，则孔隙压力为0，否则计算孔隙压力，在基坑内 uj=[hwn－yc(j)]γw，基坑外uj=[hwn－yc(j)]γw，其中，γw为地下水重度。

考虑地下水，安全系数为：

④用土钉和锚杆对边坡进行加固

基坑工程中，常常使用锚杆或土钉对边坡进行加固。加固机制是利用锚杆或土钉与土体之间的锚固力和锚杆或土钉自身的抗拉强度，将滑动土体与下部土体连接在一起。锚杆或土钉对滑动土体的拉力提供额外的抗滑力，从而提高基坑的安全系数。要想在程序中增加土钉提供的抗滑力，则需要另外对土钉进行循环计算，算法流程图如图3所示。

⑤土钉共设置有m道

第k道土钉的抗拉极限承载力：

其中，fyk为第k道土钉或锚杆的屈服强度标准值；
As为第k道土钉或锚杆抗拉材料的截面积。

第k道土钉锚固抗拔极限承载力：

图4 计算土钉抗滑力的算法

其中，fptk为第k道土钉或土锚与土体的摩阻力；
Ap为第k道土钉或锚杆与土体的接触面积。

值得注意的是，滑动土体中的土钉锚固抗拔极限承载力和底层土体中的土钉锚固抗拔极限承载力大小不同。应取两者最小值作为土钉锚固抗拔极限承载力Rk1。

第k道土钉的极限拉力值Rk取Rk0和Rk1当中最小值。

计算系数：

其中，θk为第k道土钉或锚杆与滑动圆弧交点处的法线与垂线的夹角；
αk为第k道土钉或锚杆的倾角；
φ为第k道土钉或锚杆与滑动圆弧交点处的土体内摩擦角。

第k道土钉或锚杆提供的抗滑力：

其中，sxk为第k道土钉或锚杆的水平间距。

土钉或锚杆提供的抗滑合力：

考虑土钉或锚杆，安全系数为：

⑥考虑施工进度

在基坑开挖过程中，可能也会出现整体剪切破坏。既要保证基坑开挖支护后基坑的安全稳定，又要维持基坑开挖支护施工安全稳定。施工进度模拟可以选择基坑施工过程中的几个关键节点，并根据节点状态进行验算。程序设计以开挖支护为主，改变基坑深度的值，深度由浅至深，模拟开挖过程。模拟支护过程可设置一个变量sg(k)，代表第k道土钉锚杆的施工情况。sg(k)=1代表第k道土钉锚已施工，按正常计算；
sg(k)=0代表第k道土钉锚杆尚未施工，Rk值设置为0。

⑦几种特殊情况

在调试过程中，部分区域存在异常值，即负值和极大数值，影响安全系数最小值的取值和变化趋势的观察。经分析，出现异常的原因为：圆心在该区域内的圆弧与土坡坡顶地面无交点，即滑土体极小或滑动圆弧不存在，为了避免该区域异常值对其他区域变化趋势的观察，可以在循环的开始时增加判断。如果圆弧与土坡坡顶地面无交点，可跳过本次循环，并且将该圆心位置的安全系数设为4。

滑动面可能与有些土钉没有交点，此时这道土钉极限拉力值应直接设置为零，并且跳过极限拉力值计算过程。

3.1 工程和地质概况

一个住宅项目位于道路交叉口，本工程主要包括一栋21层高层建筑及周边多栋低层建筑，以及场地西北侧的地下车库。

场地北侧设置厂房，西侧、南侧为道路，东侧为空地，场地周边无重要管线，周边环境较为平坦开阔。

场地土层第四系地貌类型属江淮丘陵区地貌单元，地下第一层为杂填土，d=1m,γ=18.8kN/m3,c=0kPa,φ=15º。第二层为棕黄色，黄灰色，可塑-硬粘土，d=3m，γ =19.4kN/m3，c=57.5kPa，φ =13.7º。第三层为棕黄色或灰黄色，硬塑至硬态粘土，d=10m，γ=19.6kN/m3，c=71.9kPa，φ=14.5º。仅场地上层为积水，且分布极不均匀，无连续潜水面。

3.2 基坑支护设计方案

本工程一部分开挖深度为6.1m，另一部分开挖深度为3m，支护结构采用放坡结合土钉墙加固。基坑浅部采用一级自然放坡，放坡系数取0.8，放坡高度取3m，放坡宽度取2.4m。放坡剖面图如图5所示。

图5 放坡剖面图一

浅基坑和深基坑交接部分也采用一级自然放坡，坡高为3.1m，放坡系数为0.81，坡宽2.5m，边坡剖面图如图6所示。

图6 放坡剖面图二

深基坑采用放坡加土钉墙，坡高为6.1m，放坡系数为0.85，坡宽为5.2m，边坡剖面图如图7所示。

图7 放坡土钉墙剖面图

设置两道土钉，土钉参数如下表所示。

土钉设计参数表

3.3 基坑整体稳定性验算

3.3.1 3m自然放坡段

图8 开挖至2m计算结果

图9 开挖至3m计算结果

3.3.2 3.1m自然放坡段

3m的放坡与3.1m的放坡的基坑深度和坡度基本相同，但后者安全系数远高于前者，这是由于后者放坡所在土层的粘聚力和内摩擦角高于前者，各基坑开挖深度越大，安全深度越低，验算结果符合实际情况。

图10 开挖至2m计算结果

图11 开挖至3.1m计算结果

图12 开挖至2m，土钉均未施工计算结果

图13 开挖至4m，土钉均未施工计算结果

3.3.3 6.1m放坡土钉墙段

第一道土钉施工后，土钉提供抗滑力，安全系数增大。

图15 开挖至-6.1m，仅设置第一道土钉计算结果

继续开挖至-6.1m，安全系数变小，此时为整个基坑开挖最危险的时刻，进行土钉施工，并且加强监测。

图16 开挖至-6.1m，两道土钉均未施工计算结果

通过计算，在放坡加土钉墙开挖支护条件下，基坑边坡安全系数为2.96，整体稳定性良好。

①利用MATLAB编程分析边坡开挖工程基坑整体稳定性，应用程序基于瑞典切片法，根据规定的台阶尺寸逐点搜索滑动弧点中心在规定范围内的位置，确定基坑边坡最危险滑移面的位置，可以求解出基坑边坡安全系数。

②编程从最简单的边坡模型分析，再考虑土层分层、坡顶超载、地下水的作用、土钉加固、施工条件等，并逐步完善程序，使程序能够检查一般基坑边坡开挖方案的整体稳定性，并能考虑施工条件，施工过程中的整体稳定性也更高。

③工程师可根据工程需要和当地工程经验，设置搜索范围，搜索步骤和土条切片数量，达到所需计算精度，从而提高效率，减少误差，具有一定的实用价值。

推荐访问:基坑 稳定性 分析