基于MATLAB的基坑整體穩定性分析

周建平，李銀，任勇，王明陽，馬海春

（1.中國水利水電第七工程局有限公司，四川 成都 231400；2.合肥工業大學資源與環境工程學院，安徽 合肥 230009）

0 引言

隨著社會的發展和進步，高層建筑和地下工程的數量與日俱增。緊隨而來基坑開挖及支護成為的問題。在空曠區域放坡開挖是最經濟、最簡單的施工方法[13-14]。在基坑支護設計中，基坑土體的整體穩定性尤為重要。分析粘性土整體穩定性的方法有很多，極限平衡法已被用于評估斜坡的穩定性[1][11-12][15]。

目前，邊坡工程中常用的安全計算方法是費倫紐斯（W.Fellenius）在1927～1937年提出的瑞典條分法[8][10]。瑞典條分法具有原理簡單、公式簡潔等優點，但瑞典條分法的計算過程繁瑣，計算量巨大[3-4]，需要反復試算圓心位置與半徑，直至逼近最危險滑動圓弧。

本文利用MATLAB數學軟件對基坑整體穩定性進行編程，結合瑞典條分法理論綜合考慮土層、地下水、地面荷載等多種影響基坑邊坡開挖整體穩定性的因素，大大減少了基坑設計人員的計算工作，并結合具體實例驗證了結果的可行性，提高了基坑穩定性計算效率。

1 基本原則

1.1 基坑穩定性分析

對具有圍護結構的粘性土邊坡或土質邊坡分析，其目的是設計合理的埋設深度或檢查所提出的圍護結構在給定條件下是否穩定合理[2][6-7]?；觽缺诘募魬υ谧灾睾偷孛娉d作用下超過土體抗剪強度時，土體將失去穩定性，發生滑移。對于具有加固和支撐作用的粘性土，土會滑動形成滑動面。通?；瑒用鏁徽J為是圓弧形以方便計算，這種形式的破壞稱為整體剪切破壞?；诱w穩定性分析是分析基坑側壁是否會發生整體剪切破壞。通常用圓弧滑動穩定安全系數Ks表示，Ks由《國家行業標準—建筑基坑支護技術規程》(JGJ120-2012)規定。錨拉式圍護結構、土釘墻、重力式水泥土墻的穩定性均可采用圓弧滑片法分析，對于安全等級一級、二級和三級的基坑工程，安全系數Ks應分別不小于1.35、1.3、1.25。

1.2 圓弧滑動條分法

圓弧滑動條分法有兩個基本假設（Navarrov，2014）：土體的剪切面是貫穿土體內圓弧面，安全系數是圓弧包含的剪切面上提供的抗滑力矩和滑動力矩之比；土條間作用力對邊坡整體穩定性影響不大[5]。以最典型的均質土坡為例，對瑞典條分法的計算過程進行說明。

如圖1所示，ABC土體為分析對象，內摩擦角為φ，黏聚力為c?；瑒用鏋閳A弧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 編程

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。由于基坑降水措施，基坑內外水位不同，需要設置并輸入基坑內外水位參數值?；觾人粸閔wn，基坑外水位為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 工程實例

3.1 工程和地質概況

一個住宅項目位于道路交叉口，本工程主要包括一棟21層高層建筑及周邊多棟低層建筑，以及場地西北側的地下車庫。

場地北側設置廠房，西側、南側為道路，東側為空地，場地周邊無重要管線，周邊環境較為平坦開闊。

場地土層第四系地貌類型屬江淮丘陵區地貌單元，地下第一層為雜填土，d=1m,γ=18.8kN/m3,c=0kPa,φ=15o。第二層為棕黃色，黃灰色，可塑-硬粘土，d=3m，γ =19.4kN/m3，c=57.5kPa，φ =13.7o。第三層為棕黃色或灰黃色，硬塑至硬態粘土，d=10m，γ=19.6kN/m3，c=71.9kPa，φ=14.5o。僅場地上層為積水，且分布極不均勻，無連續潛水面。

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，整體穩定性良好。

4 結論

①利用MATLAB編程分析邊坡開挖工程基坑整體穩定性，應用程序基于瑞典切片法，根據規定的臺階尺寸逐點搜索滑動弧點中心在規定范圍內的位置，確定基坑邊坡最危險滑移面的位置，可以求解出基坑邊坡安全系數。

②編程從最簡單的邊坡模型分析，再考慮土層分層、坡頂超載、地下水的作用、土釘加固、施工條件等，并逐步完善程序，使程序能夠檢查一般基坑邊坡開挖方案的整體穩定性，并能考慮施工條件，施工過程中的整體穩定性也更高。

③工程師可根據工程需要和當地工程經驗，設置搜索范圍，搜索步驟和土條切片數量，達到所需計算精度，從而提高效率，減少誤差，具有一定的實用價值。