基于數(shù)值流形法的反傾層狀巖質(zhì)邊坡傾倒破壞分析1

王歡歡 郭明珠

（北京工業(yè)大學(xué)，城市建設(shè)學(xué)部，北京 100124）

引言

反傾層狀巖質(zhì)邊坡指巖層走向與邊坡走向接近一致，而傾向與坡面傾向相反的邊坡，傾倒變形破壞是巖質(zhì)邊坡典型破壞失穩(wěn)形式（李明霞，2015）。傾倒破壞模式通常分為塊體傾倒、塊體彎曲傾倒和彎曲傾倒（Goodman 等，1976）。青藏高原地區(qū)是板塊運(yùn)動(dòng)最活躍的地區(qū)，該地區(qū)地殼內(nèi)、外動(dòng)力條件強(qiáng)烈交織和轉(zhuǎn)化，使大型滑坡災(zāi)害發(fā)育（黃潤(rùn)秋，2007），具有災(zāi)害類型多、影響范圍廣、破壞程度大等特點(diǎn)，其中滑坡、泥石流、崩塌產(chǎn)生的直接危害最嚴(yán)重。青藏高原地區(qū)反傾層狀巖質(zhì)邊坡深層傾倒現(xiàn)象較常見(jiàn)，特別是金沙江、岷江等深切峽谷地區(qū)，傾倒深度可達(dá)200～300 m，對(duì)川藏鐵路、水資源開(kāi)發(fā)等重大工程建設(shè)產(chǎn)生影響。反傾層狀巖質(zhì)邊坡傾倒破壞嚴(yán)重影響山區(qū)經(jīng)濟(jì)社會(huì)發(fā)展，危及群眾生命財(cái)產(chǎn)安全，因此，需研究反傾層狀巖質(zhì)邊坡傾倒破壞特征。

目前，國(guó)內(nèi)外學(xué)者對(duì)反傾層狀巖質(zhì)邊坡的數(shù)值模擬分析可分為以下兩類：一類是將邊坡看作連續(xù)體進(jìn)行數(shù)值計(jì)算（如極限平衡法、有限元法等），是對(duì)真實(shí)情況的數(shù)值近似表征，分析結(jié)構(gòu)應(yīng)力與變形后，再進(jìn)行邊坡穩(wěn)定性判斷，但對(duì)裂縫等不連續(xù)問(wèn)題和結(jié)構(gòu)破壞后運(yùn)動(dòng)過(guò)程的處理較困難。韓貝傳等（1999）、王宇等（2013）、代仲海等（2018）、王霄等（2018）構(gòu)建了反傾層狀巖質(zhì)邊坡有限元模型，探討了坡體變形破壞機(jī)制及穩(wěn)定性影響因素。另一類是將巖石看作離散體進(jìn)行數(shù)值計(jì)算（如離散元法、非連續(xù)變形分析方法），以不連續(xù)力學(xué)為基礎(chǔ)，將巖體結(jié)構(gòu)看作離散的塊體，塊體之間采用接觸連接，不僅可分析塊體運(yùn)動(dòng)、變形、破壞等，還可對(duì)塊體間的相互作用進(jìn)行分析。然而，該方法采用全一階位移近似和小角度假定，且在塊體接觸、參數(shù)取值、邊界處理等方面存在不足，影響含裂縫的不連續(xù)結(jié)構(gòu)計(jì)算結(jié)果準(zhǔn)確性。劉云鵬等（2012）、王章瓊等（2014）、陶志剛等（2017）等使用UDEC 分析了層狀巖質(zhì)邊坡破壞模式。

1992 年，石根華首次提出數(shù)值流形法（NMM），充分發(fā)揮以上兩類方法的優(yōu)勢(shì)，解決了有限元、DDA和解析法存在的計(jì)算問(wèn)題，可同時(shí)處理連續(xù)和非連續(xù)問(wèn)題，具有廣泛應(yīng)用前景。NMM 以數(shù)值流形為核心，在非連續(xù)變形方法的塊體系統(tǒng)理論基礎(chǔ)上，融入有限元法和解析法，是連續(xù)分析方法。因此，該方法既可計(jì)算塊體內(nèi)部變形，也可模擬塊體系統(tǒng)力學(xué)行為（裴覺(jué)民，1997）。

針對(duì)反傾層狀巖質(zhì)邊坡的研究多采用有限元、離散元或非連續(xù)變形分析方法等,采用NMM 進(jìn)行地震力作用下的邊坡破壞研究較少。本文采用NMM 建立反傾邊坡模型，對(duì)地震作用下不同內(nèi)摩擦角的反傾層狀巖質(zhì)邊坡進(jìn)行傾倒破壞模擬，并分析傾倒破壞過(guò)程。

1 邊坡概況

滑坡區(qū)主要出露石英片巖、砂巖等硬巖及少量片麻巖和石灰?guī)r等軟巖，巖體結(jié)構(gòu)面發(fā)育。石英片巖是坡內(nèi)出露的主要巖性，成層性較好，表面風(fēng)化強(qiáng)烈。堆積體主要成分為石英片巖夾少量片麻巖，區(qū)域巖性變化較小。滑坡堆積體方量4 500 萬(wàn)m3，水平向長(zhǎng)度1 400 m，高度方向落差1 000 m，滑坡堆積體距坡趾最遠(yuǎn)達(dá)1 022 m?；聟^(qū)主要受3 組節(jié)理控制，傾向分別為375°∠40°、78°∠82°和185°∠78°，巖層表現(xiàn)為反傾向。通過(guò)調(diào)查分析，繪制貢扎滑坡體平面圖和主剖面圖，分別如圖1、2 所示。

圖1 貢扎滑坡平面圖Fig. 1 Plan of gong-za landslide

圖2 貢扎滑坡主剖面圖Fig. 2 Main section of gong-za landslide

2 室內(nèi)力學(xué)試驗(yàn)

點(diǎn)荷載試驗(yàn)是獲取巖石力學(xué)參數(shù)的主要方法之一，在進(jìn)行室內(nèi)點(diǎn)荷載試驗(yàn)的基礎(chǔ)上，通過(guò)儀器油表獲得巖石點(diǎn)荷載強(qiáng)度。為確定數(shù)值模擬主要物理力學(xué)參數(shù)，進(jìn)行室內(nèi)點(diǎn)荷載試驗(yàn)。由于對(duì)貢扎滑坡進(jìn)行了現(xiàn)場(chǎng)實(shí)地勘察，并從現(xiàn)場(chǎng)帶回了多塊堆積體巖石（圖3），因此開(kāi)展室內(nèi)點(diǎn)荷載試驗(yàn)較方便，試驗(yàn)數(shù)據(jù)具有可信度。

圖3 試驗(yàn)所用硬巖和軟巖Fig. 3 Hard and soft rocks used in the test

根據(jù)點(diǎn)荷載計(jì)算原理，進(jìn)行徑向測(cè)試時(shí)，等效巖芯直徑De即為巖芯直徑D。進(jìn)行不規(guī)則巖塊軸向測(cè)試時(shí)，等效巖芯直徑De計(jì)算如式（1）所示。當(dāng)De取50 mm 時(shí)，可換算出直徑為50 mm 標(biāo)準(zhǔn)試件點(diǎn)荷載強(qiáng)度，由于原型巖石形狀不規(guī)則，故采用式（2）、（3）進(jìn)行計(jì)算。通過(guò)點(diǎn)荷載強(qiáng)度Is（50）及式（4）可推導(dǎo)出巖石單軸抗壓強(qiáng)度值σc及單軸抗拉強(qiáng)度στ，進(jìn)而推導(dǎo)出黏聚力c及內(nèi)摩擦角φ。

式中：W為兩加荷點(diǎn)最小截面寬度（mm）；F為修正系數(shù)；Is為未修正點(diǎn)荷載強(qiáng)度（MPa）；m為修正指標(biāo)，一般取值為0.45；k為修正系數(shù)，當(dāng)巖石為硬巖時(shí)k=0.8～0.9，當(dāng)巖石為軟巖時(shí)k=0.6～0.7。

計(jì)算得到貢扎滑坡地區(qū)硬巖黏聚力c=17.23 MPa，內(nèi)摩擦角φ=47.39°；軟巖點(diǎn)荷載強(qiáng)度Is（50）=0.74 MPa，單軸抗壓強(qiáng)度σc=14.8 MPa，單軸抗拉強(qiáng)度στ=1.11 MPa，黏聚力c=2.02 MPa，內(nèi)摩擦角φ=41.57°。

3 算例分析

3.1 計(jì)算模型

為簡(jiǎn)化數(shù)值模擬計(jì)算過(guò)程，采用圖4 所示簡(jiǎn)化邊坡模型，其中坡角為30°、巖層傾角為65°、邊坡長(zhǎng)1 400 m、高1 200 m。模型邊界固定采用位移邊界固定法，對(duì)左右邊界及底部進(jìn)行約束。各巖層均采用莫爾-庫(kù)侖準(zhǔn)則巖石本構(gòu)模型，同時(shí)，在邊坡上布置3 組編號(hào)為1～6 的監(jiān)測(cè)點(diǎn)，以便分析內(nèi)摩擦角對(duì)監(jiān)測(cè)點(diǎn)位移的影響。

圖4 簡(jiǎn)化邊坡模型Fig. 4 simplified model of slope

3.2 模型參數(shù)選取

坡體主要包括石英片巖、砂巖及少量片麻巖和石灰?guī)r等，石英片巖是坡內(nèi)出露的主要巖體，成層性較好，表面風(fēng)化強(qiáng)烈。堆積體主要成分為石英片巖夾少量片麻巖，通過(guò)室內(nèi)試驗(yàn)得到相關(guān)力學(xué)參數(shù)如表1 所示。動(dòng)力系數(shù)為1，單步允許最大位移率為0.01，總時(shí)間步選為10 000 步，泊松比為0.25，超松弛系數(shù)為1.35。

表1 巖石力學(xué)參數(shù)取值Table 1 Rock mechanical parameters

3.3 數(shù)值計(jì)算過(guò)程

主要研究?jī)?nèi)摩擦角φ對(duì)巖質(zhì)邊坡傾倒破壞的影響，取計(jì)算總時(shí)間步數(shù)為10 000 步，時(shí)間步長(zhǎng)為0.008 s，最大位移比為0.001。計(jì)算模型如圖5 所示，在邊坡底部輸入日本巖手地震記錄到的地震波（圖6）。

圖5 邊坡網(wǎng)格劃分模型圖Fig. 5 Grid division model of slope

圖6 加速度時(shí)程曲線Fig. 6 Acceleration time history curves

NMM 計(jì)算步驟如下：（1）對(duì)含節(jié)理、裂縫的巖體進(jìn)行流形法網(wǎng)格剖分，定義物理和數(shù)學(xué)網(wǎng)格；（2）進(jìn)行第1 步加載，計(jì)算地震作用下的應(yīng)力分布及變形；（3）計(jì)算各裂縫尖端應(yīng)力強(qiáng)度因子；（4）讀取每個(gè)流形單元數(shù)據(jù)并計(jì)算其接觸方式；（5）判斷接觸方式并計(jì)算變形；（6）形成法向彈簧和剪切彈簧的加法矩陣，進(jìn)行迭代求解；（7）計(jì)算變形模式。

模擬結(jié)果如圖7 所示，由圖7 可知，隨著內(nèi)摩擦角的增大，塊體之間的摩擦系數(shù)增加，塊體抗滑力增加。因此坡體從開(kāi)始破壞到新的平衡狀態(tài)和達(dá)到最大位移所需的時(shí)間越短。隨著內(nèi)摩擦角的逐漸增大，摩擦系數(shù)增大，當(dāng)增大至一定程度時(shí)，阻止塊體滑動(dòng)的力將大于塊體滑動(dòng)力，這時(shí)斜坡上的巖石塊體處于穩(wěn)定狀態(tài)，如圖7(e)所示。

圖7 內(nèi)摩擦角不同時(shí)邊坡破壞特征Fig. 7 Failure characteristics of slope with different friction angles

3.4 模擬結(jié)果分析

通過(guò)觀察數(shù)值模擬結(jié)果可知，不同內(nèi)摩擦角下邊坡最終破壞狀態(tài)不同，但邊坡上部均開(kāi)始出現(xiàn)裂縫，然后上部破壞，繼而推動(dòng)邊坡下部巖體滑動(dòng)，導(dǎo)致邊坡傾倒破壞。使用Origin 軟件進(jìn)行數(shù)據(jù)處理，得到監(jiān)測(cè)點(diǎn)位移如圖8 所示。

圖8 監(jiān)測(cè)點(diǎn)位移圖Fig. 8 Displacement diagram of monitoring points

內(nèi)摩擦角<40°時(shí)，大部分巖塊產(chǎn)生整體滑移，靠近坡頂處的巖塊發(fā)生輕微轉(zhuǎn)動(dòng)，推動(dòng)前面的巖塊加速滑動(dòng)，呈傾倒-滑移模式。內(nèi)摩擦角＞40°時(shí)，坡頂巖塊首先產(chǎn)生滑動(dòng)并轉(zhuǎn)動(dòng)驅(qū)使前面的巖塊滑動(dòng)，推動(dòng)坡腳處巖塊產(chǎn)生滑動(dòng)，呈漸進(jìn)式傾倒破壞，產(chǎn)生整體性破壞的可能性較小。

由圖8 可知，前15 s 滑坡體整體移動(dòng)，說(shuō)明在水平地震作用下整體下滑，位移變化差異??；然后下部坡體在水平力作用下拋出，中間的傾倒區(qū)域塊體出現(xiàn)傾覆現(xiàn)象，最終坡體傾倒破壞。坡腳位置的監(jiān)測(cè)點(diǎn)5、6水平位移遠(yuǎn)大于其他位置，說(shuō)明在地震來(lái)臨時(shí)，邊坡下部的巖體對(duì)地震作用更敏感，移動(dòng)距離較大。

4 結(jié)語(yǔ)

采用二維數(shù)值流形法分析了在其他巖體物理力學(xué)參數(shù)保持常量時(shí)，僅改變內(nèi)摩擦角，地震作用下反傾層狀巖質(zhì)邊坡傾倒破壞特性。結(jié)果表明，內(nèi)摩擦角<40°時(shí)，坡體首先在水平地震作用下整體下滑，然后下部坡體在水平力作用下拋出，中間的傾倒區(qū)域塊體出現(xiàn)傾覆現(xiàn)象，呈傾倒-滑移模式；內(nèi)摩擦角>40°時(shí)，少量巖塊發(fā)生傾倒破壞，邊坡較穩(wěn)定。