地震作用下黃土滑坡離散元數值模擬研究

關曉迪

(西安理工大學土木建筑工程學院，陜西 西安 710048)

黃土因其形成環境的特殊性，表現出大孔隙、水敏性、弱膠結、力學強度低、濕陷性強等特點，使得黃土地區滑坡、崩塌等地質災害頻發，形成機理復雜，嚴重制約了黃土地區的社會經濟發展。據調查，自20世紀60年代以來，在陜西省范圍內發生的滑坡就多達1 131處[1]，甘肅東部新老滑坡有4萬余處[2]。作為黃土滑坡代表，甘肅黑方臺滑坡一直備受關注，對于黑方臺黃土滑坡的形成機理、滑坡類型、分布規律、監測預警及時間預測預報等方面，學者們做了大量的研究[3-6]，為該滑坡周圍的人們生產生活及區域規劃等提供理論依據和指導。

數值技術克服了傳統方法不能對幾何形態復雜的邊坡進行計算，不能體現材料的各向異性等缺陷，可更加全面地考慮影響滑坡運動破壞的各種因素。近些年來由于計算機技術的成熟及軟件的不斷完善，數值模擬方法在滑坡領域得到了廣泛應用，模擬效果逐步趨于實際，已成為滑坡運動過程分析中的重要手段。胡煒等[7]以黑方臺焦家崖頭黃土滑坡為研究對象，采用有限差分數值模型模擬了滑坡的運動過程，證明了此類灌溉滲透型黃土滑坡具有高速遠程滑動的特征。賈俊等[8]采用離散元數值方法，對黑方臺地區飽水黃土軟弱層誘發滑坡的形成過程和滑體的運動學特征進行了研究，分析了滑體高速運動及穩定堆積過程中的相關特征。王念秦等[9]采用PFC2D軟件模擬了涇陽南塬滑坡滑體的沖擊過程，分析了滑體的刮鏟堆積及其速度場的變化情況。通過掃描巖土體顆粒輪廓進行建模以及FLAC3D耦合PFC3D的真實三維地形建模已經成為了多數學者常用的研究方法[10-12]，滑坡模擬的研究方向也擴展到了滑體能量分析、滑帶摩擦熱研究等更多領域[13-16]。對于邊坡地震動力分析的離散元模擬工作，石崇等[17-18]、李祥龍等[19-20]、胡訓健等[21]、卞康等[22]學者做了較多研究。施鳳根[23]運用PFC3D對文家溝滑坡的高速遠程運動特征進行了研究，通過數值分析與現場調查結果進行對比，揭示了文家溝滑坡的運動機理。

當下，盡管國內外專家學者對黃土滑坡的破壞過程、運動機理等方面做了大量的工作，取得了豐富的研究成果，但仍有一些問題需要進一步研究，如地震荷載作用下黃土滑坡的運動滑移特征等。本文主要以甘肅黑方臺滑坡為例，對不同地震荷載下以及滑源區不同位置顆粒的運動特征進行分析，找出地震荷載作用下滑源區顆粒的運動規律，為黃土地區地震滑坡的防治和規劃提供一定的參考。

1 模型建立

黑方臺位于中國甘肅省中部的臨夏永靖縣，是祁呂賀山字型構造與隴西旋卷構造的復合部位，區內新構造運動強烈，以差異性上升為主要特征，由于受黃河侵蝕，形成了多級階地地貌。黑方臺屬黃河Ⅳ級基座階地。黑方臺發育的滑坡多達77處，結合黑方臺滑坡所處的地質環境條件，該地區的滑坡大致分為6個區段，分別為新源端、黨川段、黃茨段、焦家崖段、焦家段和陳家段。本文以2019年2月發生的甘肅黑方臺焦家6號滑坡為例展開研究，圖1為焦家6號滑坡的無人機影像，根據AA′剖面做出滑坡的縱剖面見圖2。然后利用PFC2D程序，建立滑坡的二維模型(圖3)，模型底部寬540 m，高140 m，模型中共包含顆粒8 005個，分為滑源區顆粒(黃色)和非滑源區顆粒(藍色)。由于滑坡處于黃土地區，黃土顆粒之間的接觸與PFC中自帶的接觸黏結接觸模型較為接近，故而選取接觸黏結模型作為本次模擬的接觸模型。通過試錯法以及大量的數值試驗來對細觀參數進行修改，最終得到模型的模擬參數，見表1。

圖2 焦家6號滑坡縱剖面(m)

圖3 焦家6號滑坡PFC2D模型

表1 模擬參數

2 模型邊界以及動荷載設置

自然界中的邊坡均為半無限介質，但受限于計算負荷，數值模型的計算區域只能是有限的，因此對滑坡模型的底部、左側和右側的邊界顆粒進行設定，對邊界顆粒固定位移，圖4所示，在滑源區對12個點的顆粒進行監測，從左向右、從上向下依次編號為1—12號，監測顆粒的位移、速度以及運行軌跡，其中1、2、5、8、11、12號點為近坡面點，其中1號點位于滑坡后緣，11、12號點位于滑坡前緣，12號點相比于11號點位置更靠前，滑源區中部有三列點，每一列又有3個監測點，每一列點的x坐標是一致的，主要目的是監測滑源區顆粒隨深度增加，對地震動荷載的響應規律。綠色顆粒為邊界顆粒，可傳遞地震荷載，地震荷載為汶川地震時的加速度時程曲線，峰值加速度1.0g，振動持時50 s，見圖5。

圖4 模型邊界以及監測顆粒設置

圖5 地震波加速度時程曲線

3 初始加速度下模擬結果

3.1 滑坡運動過程力鏈分布

滑坡滑動過程主要是滑源區顆粒的運動，也體現在滑源區顆粒間接觸的破壞。圖6所示，不同時刻滑坡的運動狀態，0 s時，滑坡處于穩定狀態；30 s時，滑坡已經開始滑動，在地震波的影響下，滑源區顆粒經過流通區，向下滑移并開始堆積，可以看到，在滑坡滑動過程中，滑源區、流通區和堆積區都有顆粒，60 s時，滑坡運動已經結束，可以看到，此時流通區基本上沒有顆粒，顆粒主要在堆積區，滑源區也有一些殘積土，但是殘留在滑源區的土體從力鏈結構分布上看，較為松散，堆積區的土體也重新形成新的結構連接，但是力鏈結構也較為松散。

a)0 s時滑坡處于穩定狀態

3.2 滑源區顆粒位移分析

通過對滑源區不同位置顆粒位移的分析可以獲得整體滑源區顆粒的位移情況，見圖7。由圖7a可以看出，滑源區最外層顆粒的位移也是有很大差異性的，12號點顆粒，即滑坡體最前緣顆粒的位移量是最小的，位移量基本為0；而1號點顆粒作為滑坡體后緣顆粒，其位移量也是非常小的。可以看出，在滑坡滑動過程中并非是整個滑坡體整體下滑，滑坡體前緣和后緣都會有部分土體存在殘留的情況。而11號點顆粒，作為滑坡體前緣顆粒，在整個滑坡體的推動下，位移量最大為392 m，其次是位于滑坡體中部位置的5號點顆粒位移量為345 m，兩者相差13.6%，而且在最開始的滑坡過程中，5號點顆粒的位移量要高于11號點顆粒的位移量，說明在地震荷載加載過程中，滑坡體中部位置的顆粒最先滑動，隨即帶動滑坡前緣的顆粒進行滑動，進而整個滑坡體開始滑動。由圖7b可以看出，位于滑坡體中部位置的5、6、7號點的初始位移是不一樣的，說明滑源區不同深度的顆粒對地震荷載的響應程度也不一樣，近滑坡面的顆粒對于地震荷載的響應程度更大，而且隨著地震持時的增加，7號點的位移量基本不變，6號點顆粒的位移量是緩慢增加，然后逐漸平穩，而5號點顆粒的位移量是急劇增加，然后在30 s后開始逐漸趨于平穩。

a)滑源區近坡面顆粒

3.3 滑源區顆粒速度分析

對滑源區不同位置顆粒速度進行分析可以看到整個滑坡體在滑動過程中各個位置顆粒的運動狀態，見圖8。由圖8a可以看出，12號點顆粒(滑坡體最前緣的顆粒)以及1號點顆粒(滑坡體最后緣的顆粒)的速度都是較小的，在滑坡運動過程中，速度也變化不大，而其他近坡面位置顆粒的速度較大，而且速度變化也比較大，和位移變化規律類似，11號點顆粒的速度變化幅度是最大的，速度峰值達到28 m/s，其次是5、8號點。由圖8b可以看出，滑源區不同深度的顆粒速度時程曲線差異性較大，7號點的速度時程曲線在0附近波動，而6號點的速度變化不大，其速度峰值為6.5 m/s，但是要大于7號點顆粒的速度，而位于近坡面位置的5號點顆粒速度最大為25 m/s，其速度峰值是6號點速度峰值的3.8倍，表明隨著深度的增加，顆粒所受到的限制和約束越多，進而導致其速度變化越小，而近坡面的顆粒由于所受的限制和約束較少，速度波動較大。

a)滑源區近坡面顆粒

4 結論

采用離散單元法對地震動荷載下黃土滑坡進行數值模擬，分析地震荷載作用下滑源區黃土運動特征，主要結論如下：①地震荷載作用下的黃土滑坡存在滑源區、流通區以及堆積區3個部分，滑源區仍然會有一定黃土滑坡的殘積土存在，但是結構相較于未滑動前更加松散；②滑源區顆粒的位移差異性較大，處于滑源區前端的顆粒的位移量最大，隨著深度的增加，顆粒的位移量逐漸減小，且差異性較大。