基于MATLAB的基坑整體穩(wěn)定性分析

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

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

0 引言

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

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

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

1 基本原則

1.1 基坑穩(wěn)定性分析

對具有圍護結(jié)構(gòu)的粘性土邊坡或土質(zhì)邊坡分析，其目的是設計合理的埋設深度或檢查所提出的圍護結(jié)構(gòu)在給定條件下是否穩(wěn)定合理[2][6-7]?；觽?cè)壁的剪應力在自重和地面超載作用下超過土體抗剪強度時，土體將失去穩(wěn)定性，發(fā)生滑移。對于具有加固和支撐作用的粘性土，土會滑動形成滑動面。通常滑動面會被認為是圓弧形以方便計算，這種形式的破壞稱為整體剪切破壞?；诱w穩(wěn)定性分析是分析基坑側(cè)壁是否會發(fā)生整體剪切破壞。通常用圓弧滑動穩(wěn)定安全系數(shù)Ks表示，Ks由《國家行業(yè)標準—建筑基坑支護技術規(guī)程》(JGJ120-2012)規(guī)定。錨拉式圍護結(jié)構(gòu)、土釘墻、重力式水泥土墻的穩(wěn)定性均可采用圓弧滑片法分析，對于安全等級一級、二級和三級的基坑工程，安全系數(shù)Ks應分別不小于1.35、1.3、1.25。

1.2 圓弧滑動條分法

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

如圖1所示，ABC土體為分析對象，內(nèi)摩擦角為φ，黏聚力為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）計算出的滑動力矩更大，安全系數(shù)更高。

滑動土體的抗滑力矩為：

安全系數(shù)的計算公式如下：

根據(jù)多年的工程實踐經(jīng)驗，按照條分法計算安全系數(shù)一般比極限平衡法計算的安全系數(shù)低10%～20%。對于滑動面圓弧圓心角較大或孔隙水壓力較大時，安全系數(shù)將遠小于其他方法，對于實際工程，瑞典條分法的計算值是偏安全的。

2 編程

2.1 算法描述

在編寫程序之前描述算法，有助于在程序開發(fā)過程中保持邏輯清晰，使編程更高效。用流程圖描述了基坑整體穩(wěn)定性分析流程，先考慮最簡單的情況，即無荷載作用、無地下水的均質(zhì)土體一級放坡的土坡，且假設圓弧滑動面過坡腳。

如圖2所示，在輸入?yún)?shù)后，設置關于圓心坐標(ox，oy)的雙重循環(huán)，循環(huán)的范圍為 ox∈[X min，X max]，oy∈[Y min,Y max]。每次循環(huán)增量為step，step值設置的越小，搜索到圓心的位置就越精確。但相應的計算量會指數(shù)級上升，所以設置一個合理的大小即可。對于每一次循環(huán)中的圓心坐標(ox，oy)，計算所需與圓心坐標有關的參數(shù)，便可計算該圓心坐標對應的安全系數(shù)。安全系數(shù)所需的抗滑合力sum1和下滑合力sum2是通過土條從1到n循環(huán)，累計每個土條的抗滑分力clj+γbjhjcosθjtanφ和每個土條的下滑分力γbjhjsinθj得到。其中，n是土條劃分的條數(shù)，劃分的份數(shù)越多，抗滑合力sum1和下滑合力sum2計算的便越準確。

圖2 全局穩(wěn)定性算法流程圖

完成圓心坐標(ox，oy)的雙重循環(huán)結(jié)束后，得到一個存儲安全系數(shù)的二維數(shù)組。數(shù)組的索引值與圓心坐標有關，通過在數(shù)組中尋找最小值，并根據(jù)其索引值計算圓心坐標，便得到最危險滑動圓弧的位置和土坡的整體穩(wěn)定性安全系數(shù)。

2.2 坐標系的建立及參數(shù)計算

如圖3所示，為便于計算，選擇坡的頂點作為原點坐標，建立平面直角坐標系。

圖3 斜率計算模型

與圓心坐標(ox，oy）有關的參數(shù)計算如下。

坡腳C坐標：

其中kb為放坡坡度，寬度與高度之比；h為基坑深度。

滑動弧半徑：根據(jù)勾股定理計算，得

滑動圓弧與地面交點A橫坐標：

聯(lián)立滑動圓弧方程與地面方程

解得：

土條寬度：為滑動土體寬度1/n

與條數(shù)j有關的參數(shù)計算如下。

土條滑動面中心坐標：

土條滑動面弧長：近似取土條端點連線長度

土條重量：從滑弧中心至地面的完整土條的重量減去被開挖部分挖去的楔形土條的重量。

2.3 程序編寫與參數(shù)調(diào)整

2.3.1 編程

根據(jù)流程圖所示的算法，可以編寫此程序。對于計算結(jié)果輸出，我們可以繪制圓弧滑動安全系數(shù)與圓心坐標的關系圖像。Ks=?(ox,oy)由于Ks值是若干散點值，因此可以使用內(nèi)置函數(shù)mesh?gird()創(chuàng)建繪圖網(wǎng)格，并使用surf()函數(shù)繪制三維曲面。我們還可以對圖像進行調(diào)整，增加坐標軸名稱，增加網(wǎng)格，色彩插值處理，標記最小圓心位置等。

需輸入的參數(shù)包括放坡坡度與基坑幾何深度，土體重度值，內(nèi)聚力及內(nèi)摩擦角。需要調(diào)整的參數(shù)為：搜索范圍ox∈[X min,X max]，oy∈[Y min,Y max]、搜索步長step和土條數(shù)n。這些參數(shù)可以根據(jù)實際需要進行更改。

搜索范圍過大會導致無意義的計算，搜索范圍過小的可能會遺漏Ks最小值。根據(jù)理論分析，滑動面的中心只會出現(xiàn)在坐標系的第一象限。經(jīng)過調(diào)試，將搜索范圍設置為 ox∈[0,10]，oy∈[0,10]，即可保證最小值Ks在此范圍內(nèi)，并且能較好的表現(xiàn)出三維表面變化的趨勢。

搜索步長step越小，計算精度越高，搜索到的圓心位置越精確，對Ks影響不大，但同時計算所需時間越長。step值縮小為原來的N倍，計算量約為原來的N2倍。經(jīng)調(diào)試，搜索步長step設置為0.1～0.01即可滿足計算要求。

土條劃分數(shù)量n越小，對土條各項參數(shù)的計算結(jié)果越準確。土條劃分數(shù)量n值對Ks有影響。土條劃分數(shù)量n越大，算得Ks越大。這是由于土條滑動弧長近似取土條底面端點連線長度，比實際弧長短，計算的抗滑力小于實際的抗滑力，計算結(jié)果偏安全。土條劃分數(shù)量n值越大，計算的抗滑力越接近真值。此外，土條劃分數(shù)量n的值越大，三維曲面越光滑，圓心坐標也會發(fā)生變化。土條劃分數(shù)量n值的增加對計算時間的增加沒有step值的敏感，因為n只涉及一重循環(huán)，且循環(huán)內(nèi)的代碼較少。經(jīng)調(diào)試，土條劃分數(shù)量n設為1000片是較為合適的值。

2.3.2 功能擴展

在完成基坑放坡整體穩(wěn)定性分析最簡單的程序之后，可以考慮更復雜的情況，使該程序能夠進行更多種類的基坑放坡整體穩(wěn)定性分析。

①土體分層

絕大多數(shù)基坑場地土體是分層的，每層土的厚度、重度、內(nèi)聚力和內(nèi)摩擦角都不同，需要設置并輸入這些參數(shù)。土體分層影響土條受到的重力，應該按照土層分層計算，計算公式為 Gj=Σhjiγibj，hji是土條 j在 i層土內(nèi)的高度，γi是 i層土的重度。土體分層還影響土條受到的抗滑力，根據(jù)滑動面中心縱坐標判斷滑動面在哪一土層內(nèi)，取該土層的內(nèi)聚力和內(nèi)摩擦角。

②坡頂有超載

大多數(shù)基坑周圍存在建筑荷載、車輛荷載和建筑材料堆載。為了計算簡便，設地面超載為均布荷載。輸入均布荷載參數(shù)qj，可以發(fā)現(xiàn)均布荷載作用于地面，即y軸左側(cè)，可以根據(jù)土條滑動面中心橫坐標是否小于0判斷。若xc(j)<0，土條受到地面荷載是qj=qj，若xc(j)>0，土條受到地面荷載是qj=0。地面超載作用點在土條頂面，與圓心的距離和重力與圓心的距離相同，且作用方向也與重力相同。為計算方便，可將土條受到的地面超載加到重力上計算。

考慮坡頂超載，安全系數(shù)為：

③存在地下水

當土體中存在有地下水時，地下水位以下的粉土、砂土、碎石土會受到地下水的浮力作用，粘土和地下水位以上土體無浮力作用。因此，需要設置進入每一層土壤的類別cli，如果土層i是粘土，cli=1。如果土層i不是粘土，cli=0。由于基坑降水措施，基坑內(nèi)外水位不同，需要設置并輸入基坑內(nèi)外水位參數(shù)值?；觾?nèi)水位為hwn，基坑外水位為hww。浮力作用即孔隙壓力減小了土條與下部土體之間的內(nèi)聚力，造成抗滑力降低，而土條壓力不變，安全系數(shù)變小，基坑偏向于危險。根據(jù)土條滑動面中心的橫坐標xc(j)判斷是基坑內(nèi)還是基坑外，根據(jù)土條滑動面中心的縱坐標yc(j)判斷滑動面所處土層，若該土層為粘土或高于水位，則孔隙壓力為0，否則計算孔隙壓力，在基坑內(nèi) uj=[hwn－yc(j)]γw，基坑外uj=[hwn－yc(j)]γw，其中，γw為地下水重度。

考慮地下水，安全系數(shù)為：

④用土釘和錨桿對邊坡進行加固

基坑工程中，常常使用錨桿或土釘對邊坡進行加固。加固機制是利用錨桿或土釘與土體之間的錨固力和錨桿或土釘自身的抗拉強度，將滑動土體與下部土體連接在一起。錨桿或土釘對滑動土體的拉力提供額外的抗滑力，從而提高基坑的安全系數(shù)。要想在程序中增加土釘提供的抗滑力，則需要另外對土釘進行循環(huán)計算，算法流程圖如圖3所示。

⑤土釘共設置有m道

第k道土釘?shù)目估瓨O限承載力：

其中，fyk為第k道土釘或錨桿的屈服強度標準值；As為第k道土釘或錨桿抗拉材料的截面積。

第k道土釘錨固抗拔極限承載力：

圖4 計算土釘抗滑力的算法

其中，fptk為第k道土釘或土錨與土體的摩阻力；Ap為第k道土釘或錨桿與土體的接觸面積。

值得注意的是，滑動土體中的土釘錨固抗拔極限承載力和底層土體中的土釘錨固抗拔極限承載力大小不同。應取兩者最小值作為土釘錨固抗拔極限承載力Rk1。

第k道土釘?shù)臉O限拉力值Rk取Rk0和Rk1當中最小值。

計算系數(shù)：

其中，θk為第k道土釘或錨桿與滑動圓弧交點處的法線與垂線的夾角；αk為第k道土釘或錨桿的傾角；φ為第k道土釘或錨桿與滑動圓弧交點處的土體內(nèi)摩擦角。

第k道土釘或錨桿提供的抗滑力：

其中，sxk為第k道土釘或錨桿的水平間距。

土釘或錨桿提供的抗滑合力：

考慮土釘或錨桿，安全系數(shù)為：

⑥考慮施工進度

在基坑開挖過程中，可能也會出現(xiàn)整體剪切破壞。既要保證基坑開挖支護后基坑的安全穩(wěn)定，又要維持基坑開挖支護施工安全穩(wěn)定。施工進度模擬可以選擇基坑施工過程中的幾個關鍵節(jié)點，并根據(jù)節(jié)點狀態(tài)進行驗算。程序設計以開挖支護為主，改變基坑深度的值，深度由淺至深，模擬開挖過程。模擬支護過程可設置一個變量sg(k)，代表第k道土釘錨桿的施工情況。sg(k)=1代表第k道土釘錨已施工，按正常計算；sg(k)=0代表第k道土釘錨桿尚未施工，Rk值設置為0。

⑦幾種特殊情況

在調(diào)試過程中，部分區(qū)域存在異常值，即負值和極大數(shù)值，影響安全系數(shù)最小值的取值和變化趨勢的觀察。經(jīng)分析，出現(xiàn)異常的原因為：圓心在該區(qū)域內(nèi)的圓弧與土坡坡頂?shù)孛鏌o交點，即滑土體極小或滑動圓弧不存在，為了避免該區(qū)域異常值對其他區(qū)域變化趨勢的觀察，可以在循環(huán)的開始時增加判斷。如果圓弧與土坡坡頂?shù)孛鏌o交點，可跳過本次循環(huán)，并且將該圓心位置的安全系數(shù)設為4。

滑動面可能與有些土釘沒有交點，此時這道土釘極限拉力值應直接設置為零，并且跳過極限拉力值計算過程。

3 工程實例

3.1 工程和地質(zhì)概況

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

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

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

3.2 基坑支護設計方案

本工程一部分開挖深度為6.1m，另一部分開挖深度為3m，支護結(jié)構(gòu)采用放坡結(jié)合土釘墻加固?；訙\部采用一級自然放坡，放坡系數(shù)取0.8，放坡高度取3m，放坡寬度取2.4m。放坡剖面圖如圖5所示。

圖5 放坡剖面圖一

淺基坑和深基坑交接部分也采用一級自然放坡，坡高為3.1m，放坡系數(shù)為0.81，坡寬2.5m，邊坡剖面圖如圖6所示。

圖6 放坡剖面圖二

深基坑采用放坡加土釘墻，坡高為6.1m，放坡系數(shù)為0.85，坡寬為5.2m，邊坡剖面圖如圖7所示。

圖7 放坡土釘墻剖面圖

設置兩道土釘，土釘參數(shù)如下表所示。

土釘設計參數(shù)表

3.3 基坑整體穩(wěn)定性驗算

3.3.1 3m自然放坡段

圖8 開挖至2m計算結(jié)果

圖9 開挖至3m計算結(jié)果

3.3.2 3.1m自然放坡段

3m的放坡與3.1m的放坡的基坑深度和坡度基本相同，但后者安全系數(shù)遠高于前者，這是由于后者放坡所在土層的粘聚力和內(nèi)摩擦角高于前者，各基坑開挖深度越大，安全深度越低，驗算結(jié)果符合實際情況。

圖10 開挖至2m計算結(jié)果

圖11 開挖至3.1m計算結(jié)果

圖12 開挖至2m，土釘均未施工計算結(jié)果

圖13 開挖至4m，土釘均未施工計算結(jié)果

3.3.3 6.1m放坡土釘墻段

第一道土釘施工后，土釘提供抗滑力，安全系數(shù)增大。

圖15 開挖至-6.1m，僅設置第一道土釘計算結(jié)果

繼續(xù)開挖至-6.1m，安全系數(shù)變小，此時為整個基坑開挖最危險的時刻，進行土釘施工，并且加強監(jiān)測。

圖16 開挖至-6.1m，兩道土釘均未施工計算結(jié)果

通過計算，在放坡加土釘墻開挖支護條件下，基坑邊坡安全系數(shù)為2.96，整體穩(wěn)定性良好。

4 結(jié)論

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

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

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