四川省漢源縣中海村滑坡動力學特征數(shù)值分析

徐奕梓， 樊曉一,2， 張友誼， 田述軍， 溫 翔， 劉浩南， 鄭榆楓， 廖洪陽

(1.西南石油大學土木工程與測繪學院，四川 成都 610500； 2.西南石油大學工程安全評估與防護研究院，四川 成都 610500； 3.西南科技大學土木工程與建筑學院，四川 綿陽 621010)

0 引言

四川省漢源縣地處橫斷山脈北段東緣，為川西高原與四川盆地之間的過渡地帶，區(qū)內(nèi)地形起伏大，地質(zhì)構(gòu)造復雜且地震活動強烈，是滑坡災害的頻發(fā)區(qū)[1-2]。近年來，發(fā)生了多起地質(zhì)災害，導致嚴重的人員傷亡和財產(chǎn)損失[3-5]。

災難性滑坡常具有隱蔽性、突發(fā)性特征，現(xiàn)場監(jiān)測和實時獲取災害參數(shù)存在較大的難度，模型試驗和數(shù)值模擬成為了理解和分析這類災害的有效方法。滑坡運動距離和堆積特征是滑坡防災減災和工程防護的主要技術(shù)參數(shù)，前人利用室內(nèi)實驗模型探究了滑坡運動距離與滑坡體積之間的關(guān)系以及滑坡運動距離的影響因素[6-7]，總結(jié)了體積、垂直高度、基底摩擦角、地形條件和不同偏轉(zhuǎn)角度等因素對堆積特性的影響[8-9]。計算機的發(fā)展為滑坡運動過程的模擬帶來更多的可能性，DAN-W模型、離散元模型、DAN3D模型和Massflow模型是4種基于不同理論方法和本構(gòu)模型建立起來的模型，均能對滑坡失穩(wěn)破壞和運動堆積過程進行演化分析。國內(nèi)學者分別采用這4種模型對舟曲泄流滑坡、三溪村滑坡、水城滑坡、甘肅黑方臺滑坡的運動特征、堆積特征、能量轉(zhuǎn)化特征和沖擊鏟刮效應展開數(shù)值模擬，取得了較好的研究成果[10-13]。

降雨、地震、水庫蓄水、人類工程活動、斜坡結(jié)構(gòu)、特殊地層結(jié)構(gòu)、地質(zhì)構(gòu)造等是滑坡發(fā)生的主要因素[3,13-15]，地形條件、巖土體組成則影響了滑坡的運動過程和致災區(qū)域。2020年8月21日發(fā)生的中海村滑坡在地形條件、巖土體組成上不同于該區(qū)域的其他滑坡，本文通過現(xiàn)場調(diào)查并結(jié)合滑坡區(qū)域地形地貌、地層巖性和氣候水文條件，分析滑坡災后特征，在此基礎(chǔ)上建立Massflow數(shù)值模型，反演整個滑坡動態(tài)過程并對滑坡動力學參數(shù)和運動堆積形態(tài)進行深入探討，旨在揭示影響滑坡啟動的關(guān)鍵因素和動力學特征，為該類滑坡運動機理的研究提供參考和借鑒。

1 滑坡區(qū)地質(zhì)環(huán)境條件

1.1 地理位置及地形地貌

滑坡位于距漢源縣城約5 km的富泉鎮(zhèn)中海村，地處大渡河左岸，地貌類型屬于侵蝕構(gòu)造中低山地貌，海拔高程820～1 060 m，斜坡坡向N15°E，滑坡區(qū)兩側(cè)各發(fā)育1條沖溝，沖溝在坡體下部省道S306處匯聚。斜坡平均坡度17°，滑坡區(qū)分別在海拔高程940 m、1 000 m處分布2級緩坡臺地(圖1)。

圖1 研究區(qū)地理位置及地形地貌

1.2 地層巖性

滑坡區(qū)地層巖性單一，主要以第四系昔格達組半成巖為主(圖2)。半成巖是介于軟巖和松散土體之間的特殊巖體，含泥量較高，具有明顯的水平層理構(gòu)造。根據(jù)現(xiàn)場鉆孔結(jié)果顯示，中海村滑坡山體巖性在垂直方向自上而下依次為強風化層—泥巖半成巖—泥巖(圖3)，強風化層和泥巖半成巖主要呈黃褐色，泥巖呈灰白色。出露在滑坡區(qū)域周圍的地層由老至新依次有寒武系川盆西南筇竹寺組—婁山關(guān)組、中上二疊統(tǒng)喀大崩組+奔子欄組、上二疊統(tǒng)揚子西南緣峨眉山玄武巖組+宣威組/黑泥哨組。

圖2 昔格達組分布示意圖[16]

圖3 滑坡現(xiàn)場巖石特征

1.3 降雨特征

漢源縣屬于亞熱帶季風濕潤氣候，冬暖夏涼，四季分明。2010—2020年期間，研究區(qū)年平均氣溫16 ℃，年降雨量600 mm，其中2020年降雨量為858.9 mm (圖4a)。2020年全年降雨量不均勻分布，7—8月為雨季，降雨量分別為167.5 mm和351.5 mm，約占全年降雨量的60%(圖4b)。8月強降雨主要集中在滑坡發(fā)生的前7天，16日出現(xiàn)最大降雨，日降雨量達56.5 mm (圖4c)。經(jīng)實地調(diào)查，滑坡災后現(xiàn)場表面出現(xiàn)了流水特征，認為滑坡發(fā)生前受連續(xù)降雨的影響，雨水沿裂縫結(jié)構(gòu)面下滲使得動靜水壓力增加，飽和巖土體強度降低并引起應力的變化，導致天然巖體結(jié)構(gòu)遭到破壞產(chǎn)生變形，而半成巖因成巖度不高且具遇水易軟化膨脹的特性，對于強風化層來說滲透性較大，因此，強降雨和半成巖特性是導致滑坡災害發(fā)生的關(guān)鍵因素。

(a) 年平均氣溫及降雨量 (b) 月降雨量 (c) 日降雨量

2 滑坡基本特征

2.1 滑坡規(guī)模及運動演變特征

滑坡形態(tài)整體呈長條形狀，最大水平運動距離約600 m，最大寬度約為260 m，前后高差約180 m(后緣高程1 035 m，前緣高程855 m)，滑坡運動等效摩擦系數(shù)為0.3?；聠訒r，約15×104m3的山體在前期連續(xù)強降雨的影響下失穩(wěn)并快速滑出，隨即破碎解體，在上部緩坡平臺的下部，滑坡體向下加速運動并鏟刮坡體上表面的松散土體，使體積逐漸增加至20×104m3； 運動過程中受地形影響，運動方向多次發(fā)生改變?；麦w沖擊掩埋了山坡平臺上和山腳的民房，部分滑坡邊緣區(qū)房屋嚴重受損，最終導致S360省道道路中斷，造成7人死亡2人失聯(lián)，8棟民房被毀(圖5)。

(a) 滑坡前 (b) 滑坡后

2.2 滑坡分區(qū)及運動堆積特征

該滑坡發(fā)生在古滑體的殘坡積層[13]。對災后現(xiàn)場的勘察顯示，滑坡區(qū)域呈階梯狀分布，通過當?shù)鼐用瘾@悉，1989年滑坡源區(qū)后緣曾發(fā)生滑坡。在滑坡現(xiàn)場，利用無人機遙感技術(shù)生成滑后DEM數(shù)據(jù)，對滑坡前后高程變化進行了差分計算(圖6)。

圖6 滑坡前后的高程變化

結(jié)合實地調(diào)查結(jié)果和高程變化特征，將致災范圍沿滑坡運動方向分為滑源區(qū)、堆積區(qū)Ⅰ、鏟刮區(qū)和堆積區(qū)Ⅱ共4個部分。

2.2.1 滑源區(qū)及堆積區(qū)Ⅰ

滑源區(qū)位于S360省道北側(cè)，主滑方向為N8°E，據(jù)現(xiàn)場測量，滑源區(qū)寬度約為165 m，滑坡后壁的最大下錯距離為14 m，且壁面完整平滑，形成了典型的“圈椅狀”滑坡壁，具有明顯的剪切摩擦痕跡(圖7)。野外滑坡痕跡及特點顯示，滑源區(qū)土體失穩(wěn)向下滑動，滑源區(qū)后部瞬間缺少支撐，受到巖土體向下運動時的張拉和剪切作用，使得滑源區(qū)后緣部出現(xiàn)2條明顯的裂縫。第一條為直線形裂縫，裂縫高程為1 051 m，長度為25 m，裂縫發(fā)展方向為E40°S和W40°N，屬于張拉裂縫。第二條裂縫為弧形裂縫，裂縫高程為1 043 m，長度為11 m，裂縫發(fā)展方向為東南E30°S和W40°N。由于第二條裂縫的位置臨近滑源區(qū)，受到張拉力的同時也受到了滑源區(qū)坡體向下運動時的剪切力，使得第二條裂縫出現(xiàn)下錯，最大下錯距離約為1 m，屬于剪切裂縫。

堆積區(qū)Ⅰ位于斜坡3級平臺上，堆積方向保持不變，區(qū)域內(nèi)最大堆積厚度約8 m。在滑坡發(fā)生后平臺上的一棟居民房屋完全被掩埋，導致2位老人失聯(lián)。

圖7 滑坡裂縫

2.2.2 鏟刮區(qū)

鏟刮區(qū)位于2級平臺上。部分滑坡體在堆積區(qū)Ⅰ開始堆積后，另一部分滑坡體沿著下墊面向下運動，運動方向發(fā)生改變，主滑方向為N24°E。該段滑坡下墊面場地為梯田，表層土體為強風化層，由于常年受到人類活動的干預，土體較為松軟，且殘積坡積物較厚，滑坡體向下運動時不斷鏟刮表面的松散土體使滑坡體體積增加，鏟刮表層土體的同時攜卷滑坡區(qū)兩側(cè)土體?；麦w向下鏟刮一段距離后由于地形的阻擋作用運動方向再次偏轉(zhuǎn)，發(fā)生第二次堆積。

2.2.3 堆積區(qū)Ⅱ

第二次運動方向偏轉(zhuǎn)后，即滑坡體進入最終堆積階段(1級平臺)。主堆積方向為N22°W，沿堆積方向的最大水平距離約為120 m，最大堆積厚度為5 m。在堆積區(qū)Ⅱ存在坡度較大的陡坎(圖8)，滑坡體拋射加速后沖毀并掩埋了公路內(nèi)外側(cè)的房屋。根據(jù)滑坡的原始地形和堆積特征，中海村滑坡運動地形為階梯型地形。階梯型滑坡相對于凹面型滑坡和坡腳型滑坡而言，最大水平運動距離最小，但在坡腳位置具有最大沖擊力[17]。因此，雖然滑坡在堆積區(qū)Ⅱ運動距離短，但導致了嚴重的房屋損壞和掩埋。

圖8 滑坡現(xiàn)場

3 滑坡數(shù)值模擬

3.1 Massflow介紹

動態(tài)反演分析是在獲取滑坡基本參數(shù)和地形數(shù)據(jù)條件下，借助數(shù)值分析工具進行的動態(tài)特征分析。Massflow從流體力學的N-S方程出發(fā)，基于深度積分的連續(xù)介質(zhì)力學理論，提出具有二階精度的高效能數(shù)值仿真模型，其高效且擴展自如的特點能讓軟件使用者根據(jù)自身需求使用Fortran語言對源代碼進行添加和修改[18]。該軟件可用來模擬如山體滑坡、碎屑流、泥石流、堰塞湖、潰壩、雪崩等多種山區(qū)地質(zhì)災害[19]，不僅可有效反演災害動態(tài)過程，模擬結(jié)果還能為工程實踐提供參考。近年來，該方法在地質(zhì)災害研究中取得了豐碩的成果，已應用在如香港滑坡[20]、北京房山滑坡[21]、國外的諾拉泥石流[20]、國內(nèi)的紅椿溝泥石流[22]以及冷漬溝泥石流[23]的工程治理效果評價中。

3.2 模型參數(shù)選取

在Massflow已構(gòu)建的Coulomb-viscous模型、Voellmy模型、Manning模型這3種物理摩擦模型中，選取適用于滑坡、碎屑流和巖崩等地質(zhì)災害的Coulomb-viscous模型，該摩擦模型需考慮的物理力學參數(shù)有： 密度ρ、黏聚力c、基底摩擦系數(shù)μ、孔隙水壓力系數(shù)ru、基底摩擦角δ和內(nèi)摩擦角φ。以上模擬參數(shù)主要根據(jù)文獻調(diào)研[24-25]和模擬校準確定： 泥巖半成巖的密度范圍為1.8～2.3 g/cm3，取2.05 g/cm3，黏聚力c為28 kPa，內(nèi)摩擦角φ為25°，基底摩擦系數(shù)采用等效摩擦系數(shù)(f=0.3)，基底摩擦角則利用等效摩擦系數(shù)的反正切函數(shù)計算(δ約為17°)，孔隙水壓力系數(shù)控制在0.6～0.9范圍?？紤]到體積膨脹，將初始體積乘以系數(shù)1.2。

在孔隙水壓力系數(shù)分別為0.6、0.7、0.8和0.9情況下進行試錯試驗，通過反向分析得到在孔隙水壓力系數(shù)為0.8時，模擬結(jié)果與實際情況更加吻合。為了量化分析數(shù)值模擬堆積厚度分布，將模擬結(jié)果在滑坡主滑方向、滑源區(qū)、堆積區(qū)Ⅰ和堆積區(qū)Ⅱ進行分析(圖6)，使其與航空圖像和野外測量數(shù)據(jù)進行對比驗證，對比結(jié)果如圖9所示。由于模型中沒有考慮表面微觀地貌的變化和房屋對滑坡體運動的影響，剖面3-3′和4-4′的研究結(jié)果存在一定差異，但出現(xiàn)了較為一致的變化特色，因此，數(shù)值模擬結(jié)果具有可行性。

(a) 剖面1-1′ (b) 剖面2-2′

(c) 剖面3-3′ (d) 剖面4-4′

3.3 滑坡全運動過程模擬

3.3.1 滑坡體運動形態(tài)變化特征

采用Massflow模型的仿真模擬結(jié)果如圖10所示，以10 s間隔繪制不同時間步長的滑坡運動形態(tài)，滑坡運動過程共持續(xù)約70 s。為了更直觀地觀察不同時間厚度的變化特征，在主滑方向進行剖切，從整個堆積過程可見，滑坡最大運動距離約600 m，滑坡體主要堆積在3個平臺上，最大堆積厚度為8 m(圖11)。

圖10 中海村滑坡不同時刻運動厚度分布

圖11 中海村滑坡不同時刻厚度分布

3.3.2 運動速度特征

滑坡運動速度的大小對整個滑坡致災強度有著至關(guān)重要的影響。圖12中顯示了不同時間段滑坡體的速度分布特征。10 s時最大運動速度位于滑坡體前部，20 s時最大運動速度位于滑坡體中部，30 s時最大滑坡運動速度位于滑坡體中后部； 在此期間最大運動速度為16.2 m/s，均出現(xiàn)在3級平臺上。40～60 s階段，當滑坡體運動至1級平臺時滑坡體前緣速度開始減小，由于慣性作用，后緣顆粒的能量傳遞到前緣顆粒，使得前緣顆粒繼續(xù)向前運動，掩埋了對側(cè)房屋，而此時滑源區(qū)物質(zhì)向3級平臺上繼續(xù)堆積，直到70 s時完全停止運動。

圖12 中海村滑坡不同時刻運動速度分布

在目前已有的眾多滑坡各區(qū)段運動速度估算方法中，最為常見且有效的方法是謝德格爾(A. E. Scheidegger)法[26]，其計算公式如下：

。

(1)

式中：V為運動速度，m/s；g為重力加速度,通常取9.8 m/s2；H為研究點距滑坡最高點的垂直距離,m；L為研究點距滑坡后緣的水平距離,m；f為等效摩擦系數(shù)。

經(jīng)計算，整個滑坡運動速度曲線如圖13所示。對比2種方法分析的滑坡運動速度，在整個運動速度的演化過程中，滑源區(qū)以外的滑坡體運動速度規(guī)律基本一致且出現(xiàn)3個峰值，峰值速度均出現(xiàn)在3個緩坡平臺的后部，其原因在于坡度由小到大的顯著變化致使運動速度增加，從大到小的變化因滑體撞擊地面而速度減小。在滑坡整體的速度演化分布上，由于沒有考慮地形因素對速度的影響，導致計算值大于數(shù)值模擬值。

圖13 運動速度對比

4 結(jié)論

(1)野外調(diào)查和地形高程變化特征顯示，中海村滑坡前后緣高差180 m，最大水平運動距離600 m，總堆積體積為20×104m3?；聟^(qū)域總體可分為滑源區(qū)、堆積區(qū)Ⅰ、鏟刮區(qū)和堆積區(qū)Ⅱ這4個區(qū)域。

(2)降雨和半成巖特性是導致滑坡災害發(fā)生的關(guān)鍵因素。滑坡所在區(qū)域地層巖性為第四系昔格達組泥巖半成巖，該巖層具有明顯的水平層理構(gòu)造，透水性強，遇水易軟化膨脹； 滑坡發(fā)生前的連續(xù)強降雨增加了地下水的入滲，誘發(fā)了滑坡。

(3)數(shù)值模擬結(jié)果表明，孔隙水壓力系數(shù)為0.8時模擬結(jié)果與實際情況吻合?；逻\動持續(xù)時間約為70 s，最大堆積厚度為8 m，最大運動速度為16.2 m/s。

(4)中海村滑坡速度在運程上的分布演化存在3個峰值，峰值速度和主堆積區(qū)均分布在3個緩坡平臺上。