基于Scoops 3D模型的庫(kù)岸邊坡穩(wěn)定性研究

王志浩，黨進(jìn)謙，郭 釗

（西北農(nóng)林科技大學(xué)水利與建筑工程學(xué)院，陜西楊凌 712100）

水庫(kù)在運(yùn)行期間，為了滿足防洪、灌溉、發(fā)電等需要，庫(kù)水位需要進(jìn)行周期性的漲落，在這過(guò)程中岸坡的穩(wěn)定性會(huì)出現(xiàn)不同程度的下降，進(jìn)而發(fā)生庫(kù)岸滑坡災(zāi)害。國(guó)內(nèi)外庫(kù)岸滑坡災(zāi)害屢見(jiàn)不鮮，如湖南省資水河柘溪水庫(kù)的塘巖光滑坡［1］、三峽庫(kù)區(qū)的千江坪滑坡［2］、意大利Vajont 水庫(kù)的左岸滑坡［3］等均對(duì)人民生命財(cái)產(chǎn)安全造成重大損失。

目前，國(guó)內(nèi)外針對(duì)庫(kù)岸邊坡有著大量的研究成果。鄭穎人［4］、吳瓊［5］通過(guò)建立岸坡飽和與非飽和非穩(wěn)定滲流模型，探討了庫(kù)水位波動(dòng)情況下岸坡內(nèi)地下水滲流場(chǎng)的變化情況；江強(qiáng)強(qiáng)［6］、肖詩(shī)榮［7］采用物理模型試驗(yàn)研究岸坡在不同水位波動(dòng)、降雨作用下的穩(wěn)定性演化規(guī)律，揭示滑坡的變形過(guò)程；仉文崗［8］等使用FLAC 進(jìn)行數(shù)值模擬，表明庫(kù)水位下降對(duì)岸坡穩(wěn)定性影響較大，且水位下降過(guò)程存在最危險(xiǎn)水位；楊敏［9］等采用離散元軟件FINAL 探索了錨索加固后岸坡穩(wěn)定性分析方法；Xia［10］等通過(guò)對(duì)三峽區(qū)域石榴樹(shù)包滑坡多年的監(jiān)測(cè)研究，表明岸坡淺層滑坡主要受降雨影響，而庫(kù)水變動(dòng)是造成岸坡深層滑坡的主要原因，D R Montgomery 和W E Dietrich［11］建立了可考慮降雨的SHALSTAB 穩(wěn)定滲流模型；J W Godt 等［12］建立了TRIGRS 模型。史笑凡［13］等提出了不確定性邊坡穩(wěn)定性Bootstrap-BCa估計(jì)方法；王禹璽［14］等基于Copula 理論和Bootstrap 方法提出了邊坡系統(tǒng)可靠度區(qū)間的高效估計(jì)方法；但是，對(duì)三維真實(shí)岸坡整體區(qū)域穩(wěn)定性進(jìn)行的動(dòng)態(tài)研究較少。

所以，本文以山美水庫(kù)為研究對(duì)象，首先采用解析法即包辛涅斯克潛水非穩(wěn)定性方程計(jì)算庫(kù)水位升降過(guò)程中岸坡地下水位，然后采用三維斜坡穩(wěn)定性計(jì)算模型Scoops 3D 進(jìn)行庫(kù)岸邊坡穩(wěn)定性的計(jì)算，并就巖土力力學(xué)參數(shù)對(duì)水庫(kù)岸坡的穩(wěn)定性影響程度進(jìn)行分析。此論文主要目的是研究庫(kù)水位的升降速率和巖土體力學(xué)參數(shù)對(duì)岸坡穩(wěn)定性的影響，為水庫(kù)運(yùn)行和庫(kù)岸滑坡評(píng)價(jià)提供參考。

1 山美水庫(kù)區(qū)域概況

1.1 地質(zhì)、水文概況

山美水庫(kù)位于福建省泉州市西北部，地處晉江上游的東溪，是晉江流域唯一一座綜合利用的大型水利樞紐工程，水庫(kù)集雨面積1 023 km2，總庫(kù)容6.55 億m3。研究區(qū)雨量年內(nèi)分布極不均勻，多集中在4-9月，占全年雨量的70%～80%，而10-12月僅占6%～10%，且夏季多暴雨、大雨。研究區(qū)域地層發(fā)育不完全，出露的地層主要為上三疊統(tǒng)-侏羅系、第四紀(jì)。其流域地勢(shì)由西北向東南顯著傾斜，西部山嶺縱橫，地形多以低山與丘陵為主，東部呈階梯狀，散布著山間小盆地。

1.2 庫(kù)水位變動(dòng)特征

近年來(lái)，山美水庫(kù)庫(kù)水位在70～95 m 之間呈現(xiàn)規(guī)律性地波動(dòng)，可以將庫(kù)水位的波動(dòng)過(guò)程分為4 個(gè)階段（見(jiàn)圖1），具體劃分如表1所示。本論文主要針對(duì)庫(kù)水位蓄水、下降階段進(jìn)行庫(kù)岸邊坡的穩(wěn)定性研究。

表1 山美水庫(kù)運(yùn)行階段劃分表Tab.1 Shanmei Reservoir Operation Stage Division Table

2 庫(kù)岸邊坡穩(wěn)定性研究方法

2.1 庫(kù)岸邊坡地下水浸潤(rùn)線計(jì)算

地下水是影響庫(kù)岸邊坡穩(wěn)定性的重要因素，所以，在庫(kù)岸邊坡穩(wěn)定性研究中準(zhǔn)確、合理地確定庫(kù)水位升降作用下的浸潤(rùn)線分布具有重要的意義。在進(jìn)行庫(kù)水位變動(dòng)引起的岸坡地下水位變化計(jì)算時(shí)，解析法具有方便快捷的優(yōu)勢(shì)。目前常用的解析解法有線性解、近似解、簡(jiǎn)化解，此論文計(jì)算庫(kù)岸邊坡浸潤(rùn)線變化時(shí)使用Boussinesq 一維非穩(wěn)定滲流方程。計(jì)算中的基本假設(shè)為：采用裘布依假定、潛水層為均質(zhì)含水層。由這些假設(shè)條件下的潛水非穩(wěn)定性運(yùn)動(dòng)微分方程可由包辛涅斯克（Boussinesq）方程得到，即：

式中：h是指岸坡內(nèi)部浸潤(rùn)線高程，m；K表示巖土體滲透系數(shù)，m/d；μ為給水度；t表示水位變動(dòng)時(shí)間，d；H表示巖土體含水層厚度。

計(jì)算簡(jiǎn)圖如圖2所示。依據(jù)基本假定可以得到水位升降的半無(wú)限含水層中岸坡地下水非穩(wěn)定滲流計(jì)算滲流模型：

式中，u表示土體中水位變化高度；用水位升降始末高程的平均值hm代替巖土體含水層厚度，簡(jiǎn)化方法將其看為常量，即取hm等于H。

將上述的滲流計(jì)算模型利用拉普拉斯積分變換和逆變換即可得到微分方程如下：

為方便計(jì)算對(duì)η-M(η)函數(shù)曲線進(jìn)行擬合，得到的多項(xiàng)式擬合公式為：

于是就可以得到庫(kù)水位升降時(shí)浸潤(rùn)線的簡(jiǎn)化計(jì)算公式，其表達(dá)式為：

式中：h0，0為庫(kù)水升降前的水位，d；v0為庫(kù)水位升降速率，m/d，水位下降時(shí)取負(fù)值，水位上升時(shí)取為正值。

2.2 Scoops 3D模型

Scoops 3D 是在2015 由美國(guó)地質(zhì)調(diào)查局的M E Reid 等［15］開(kāi)發(fā)的一個(gè)三維斜坡穩(wěn)定性計(jì)算模型，此模型以數(shù)字高程模型為基礎(chǔ)，采用球形滑面，可系統(tǒng)全面的搜索大量的潛在滑坡并計(jì)算其穩(wěn)定性。在進(jìn)行穩(wěn)定性計(jì)算時(shí)，此模型可以耦合多種不同環(huán)境、力學(xué)等條件參數(shù)進(jìn)行分析，常見(jiàn)的有：斜坡體的孔隙水壓力、坡體的地下水位高度、巖土體力學(xué)參數(shù)、地震荷載、真實(shí)復(fù)雜地形等。

Scoops 3D 模型基于一個(gè)橢球形的潛在滑面利用三維柱體極限平衡分析法計(jì)算潛在滑坡的穩(wěn)定系數(shù)［16，17］，進(jìn)行穩(wěn)定性計(jì)算的基本單元為三維柱體，主要過(guò)程為，首先在設(shè)定的數(shù)字高程模型范圍上方形成大量的球心，然后生成若干的搜索球體和潛在滑坡三維柱體。球面和數(shù)字高程模型相交形成的不規(guī)則坡體內(nèi)所有三維柱體構(gòu)成潛在的滑動(dòng)體，模型中潛在三維滑坡體如圖3所示。

模型使用力矩平衡法即簡(jiǎn)化的Bishop 法進(jìn)行三維斜坡穩(wěn)定性系數(shù)計(jì)算，Scoops 3D中基本計(jì)算公式為［18］：

式中：Am，n為岸坡潛在滑面的表面積；αm，n為岸坡潛在滑面的視傾角；βm，n為岸坡潛在滑面的傾角；Wm，n為潛在滑動(dòng)三維柱體的重量；um，n為孔隙水壓力，作用在三維柱體潛在滑面上；km，n為三維柱體中心受到的水平的振動(dòng)荷載；cm，n和φm，n滑面處的巖土體抗剪強(qiáng)度參數(shù)；Rm，n為轉(zhuǎn)動(dòng)軸到潛在三維柱體滑面的距離，在三維中隨著柱體的位置等改變Rm，n也不斷發(fā)生變化。

3 研究區(qū)數(shù)據(jù)參數(shù)

3.1 DEM數(shù)據(jù)

此次研究中，山美水庫(kù)的表面地形DEM數(shù)據(jù)由地理信息數(shù)據(jù)庫(kù)得到，使用ArcGIS 軟件進(jìn)行處理后得到分辨率為5 m 的高精度DEM，將其作為模型的基礎(chǔ)地形高程數(shù)據(jù)。庫(kù)岸邊坡高程以及坡度數(shù)據(jù)分別如圖4、5所示。研究區(qū)庫(kù)岸邊坡的面積約為16.1 km2。

3.2 巖土體參數(shù)

模型計(jì)算所需的主要物理力學(xué)參數(shù)為：庫(kù)岸邊坡巖土體的土體重度、飽和抗剪強(qiáng)度以及土體的滲透系數(shù)。在Scoops 3D模型計(jì)算時(shí)，根據(jù)山美水庫(kù)地質(zhì)勘察資料選取巖土體物理力學(xué)參數(shù)，其中黏聚力為20 kPa，內(nèi)摩擦角為16°，滲透系數(shù)為0.05 m/d，土體部分飽和重度為18 kN/m3，飽和重度為20 kN/m3。

3.3 計(jì)算工況的確定

3.3.1 不同庫(kù)水位升降速率下庫(kù)岸邊破穩(wěn)定性

近年來(lái)，山美水庫(kù)水位在90～75 m之間漲落，水位下降期的平均水位降速為0.8 m/d。為了研究庫(kù)水位升降速率對(duì)庫(kù)岸邊坡穩(wěn)定性的影響，設(shè)置水位變動(dòng)范圍為90 m 至75 m，選取庫(kù)水位升降速率分別為0.5、0.8、2 和3 m/d，計(jì)算不同水位升降速率下水位變動(dòng)最后時(shí)刻瞬時(shí)狀態(tài)下山美水庫(kù)庫(kù)岸邊坡穩(wěn)定性的變化情況，據(jù)此評(píng)價(jià)水位升降速率對(duì)庫(kù)岸邊坡穩(wěn)定性的影響。

3.3.2 庫(kù)岸邊坡穩(wěn)定性對(duì)巖土參數(shù)的敏感性

模擬計(jì)算過(guò)程中使用畢肖普強(qiáng)度理論進(jìn)行岸坡穩(wěn)定性計(jì)算，朱銳［19］等通過(guò)室內(nèi)直剪試驗(yàn)，揭示了干濕循環(huán)作用下膨脹土邊坡土體黏聚力出現(xiàn)明顯降低、內(nèi)摩擦角呈小幅度減小的衰減規(guī)律，據(jù)此確定黏聚力c、內(nèi)摩擦角φ、土體重度γ為岸坡體穩(wěn)定性主要影響因素。對(duì)庫(kù)岸邊坡穩(wěn)定性進(jìn)行大量的數(shù)值模擬，來(lái)獲得庫(kù)岸邊坡在主要影響因素變化時(shí)穩(wěn)定性的變化規(guī)律，進(jìn)而進(jìn)行敏感性分析。模擬中水位波動(dòng)條件設(shè)定為由90 m 以0.8 m/d 的速度降至75 m。通過(guò)控制其中兩個(gè)影響因素，另外一個(gè)因素作為變量進(jìn)行計(jì)算，具體模擬方案見(jiàn)表2。

表2 模型巖土力學(xué)參數(shù)變化表Tab.2 Model rock and soil mechanics parameter change table

3.4 模型參數(shù)

為提高模型計(jì)算效率和精度，需預(yù)先設(shè)定搜索范圍即三維搜索矩陣，在此范圍內(nèi)Scoops 3D 模型將生成大量的球形滑面球心。研究的區(qū)域分辨率為5 m，岸坡最高點(diǎn)高程為261 m，共有1 464×1 956 個(gè)柵格。據(jù)此定義三維搜索矩陣的起點(diǎn)為1×1號(hào)柵格，搜索矩陣的最低高程設(shè)定為300 m，最高高程為500 m，搜索半徑增量為0.5 m。庫(kù)岸區(qū)域模型參數(shù)見(jiàn)表3。

表3 山美水庫(kù)庫(kù)岸區(qū)域模型參數(shù)Tab.3 Model parameters of Shanmei reservoir bank area

4 結(jié)果分析

根據(jù)節(jié)2 的岸坡穩(wěn)定性研究方法和節(jié)3 的參數(shù)，并按照3.3節(jié)中設(shè)定工況進(jìn)行山美水庫(kù)庫(kù)岸邊坡穩(wěn)定性分析。結(jié)果分析中參考相關(guān)標(biāo)準(zhǔn)［20］依據(jù)穩(wěn)定性系數(shù)對(duì)岸破穩(wěn)定性進(jìn)行分級(jí)。分級(jí)標(biāo)準(zhǔn)見(jiàn)表4。

表4 岸坡穩(wěn)定性劃分標(biāo)準(zhǔn)Tab.4 Bank slope stability classification standard

4.1 不同水位升降速率下岸坡穩(wěn)定性計(jì)算

4.1.1 不同水位降速下岸坡穩(wěn)定性計(jì)算結(jié)果

將不同水位降速條件下，庫(kù)水位降至75 m時(shí)岸坡體瞬時(shí)浸潤(rùn)線分布狀態(tài)作為地下水位參數(shù)，并將巖土體力學(xué)等參數(shù)一并輸入到Scoops 3D 模型中進(jìn)行計(jì)算，得到山美水庫(kù)庫(kù)岸邊坡在90 m 常水位和不同水位下降速率下庫(kù)水位為75 m 時(shí)的岸坡安全系數(shù)，并進(jìn)行穩(wěn)定性分級(jí)，結(jié)果如圖6所示。其中水位下降速率為0.5 和0.8 m/d 時(shí)岸坡穩(wěn)定性變化極小，穩(wěn)定性分級(jí)圖基本一致，故在圖6中只繪制了v=0.8 m/d 降速下的穩(wěn)定性分級(jí)圖。在各工況間的趨勢(shì)對(duì)比如圖7所示。

由圖6、7 可以得到：①庫(kù)水位在90 m 高水位運(yùn)行期間，庫(kù)岸邊坡大部分具有較高的穩(wěn)定性，穩(wěn)定性等級(jí)在較高以上的岸坡占全庫(kù)岸段的68.77%，低穩(wěn)定性等級(jí)的岸坡占全庫(kù)岸段的14.37%。②庫(kù)水位從90 m 下降至75 m 過(guò)程中，當(dāng)下降速度為0.8 m/d 時(shí)，低穩(wěn)定性岸坡占全庫(kù)岸段的19.44%，較90 m 常水位時(shí)增加了5.07%。岸坡出現(xiàn)了新的零星分布的低穩(wěn)定性區(qū)域，右側(cè)庫(kù)岸段低穩(wěn)定和較低穩(wěn)定性區(qū)域較之前明顯擴(kuò)大，部分低穩(wěn)定性區(qū)域有點(diǎn)狀分布連通成片。③庫(kù)水位下降速度由0.5 m/d 增加到0.8 m/d 時(shí)，庫(kù)岸邊坡穩(wěn)定性基本不變。④庫(kù)岸邊坡的穩(wěn)定性隨著水位下降速度的增大逐漸較小，下降速度增加到2、3 m/d 時(shí)，低穩(wěn)定性岸坡占全庫(kù)岸段的百分比分別為21.82%、24.84%，相對(duì)0.8 m/d面積分別增加了0.38、0.87 km2。因此山美水庫(kù)在泄洪等排水過(guò)程中，庫(kù)水位下降速率應(yīng)控制在0.8 m/d之內(nèi)，且在此范圍內(nèi)根據(jù)需要可適當(dāng)增加減速，但應(yīng)減少2 m/d、3 m/d高降速的情況。

4.1.2 不同水位升速下岸坡穩(wěn)定性計(jì)算結(jié)果

將不同水位上升速率條件下庫(kù)水位升至90 m 瞬時(shí)地下水位狀態(tài)等條件輸入到Scoops 3D 模型中進(jìn)行庫(kù)岸穩(wěn)定性分析，計(jì)算得到庫(kù)水位上升速率分別為0.8、2、3 m/d 時(shí)庫(kù)庫(kù)邊坡穩(wěn)定性系數(shù)并進(jìn)行分級(jí)，各工況間的變化趨勢(shì)對(duì)比如圖8所示。

由統(tǒng)計(jì)可知：①庫(kù)水位上升工況下，庫(kù)岸邊坡整體穩(wěn)定性系數(shù)較庫(kù)水位下降工況下大，處于高穩(wěn)定區(qū)域占全庫(kù)岸段的百分比增加了5%～10%，低穩(wěn)定性區(qū)域也較水位下降時(shí)有所降低，降低約3%～6%。②庫(kù)水位上升過(guò)程中，隨著庫(kù)水位上升速度的增加岸坡穩(wěn)定性逐漸增大，上升速度為2、3 m/d 相對(duì)0.8 m/d 低穩(wěn)定性區(qū)域分別降低了1.17%、1.57%。所以在水庫(kù)運(yùn)行期間，水位上升速度在3 m/d范圍內(nèi)可適當(dāng)增大。

據(jù)此可得到，研究區(qū)庫(kù)岸滑坡穩(wěn)定性對(duì)庫(kù)水位升降速很敏感。庫(kù)水位下降高度相同條件下，隨著庫(kù)水位下降速率的增大庫(kù)岸穩(wěn)定性顯著降低；而庫(kù)水位上升高度保持不變時(shí)，庫(kù)水位上升速度的增大對(duì)庫(kù)岸邊坡的穩(wěn)定性有利。

4.2 巖土體參數(shù)敏感性分析

采用單因素敏感性分析方法按照3.3.2 節(jié)設(shè)計(jì)工況對(duì)庫(kù)岸邊坡滑坡穩(wěn)定性進(jìn)行研究，庫(kù)岸滑坡穩(wěn)定性分級(jí)計(jì)算結(jié)果如圖9所示。

由圖9可知，巖土體力學(xué)參數(shù)的變化對(duì)庫(kù)岸邊坡穩(wěn)定性影響很顯著，將這些指標(biāo)按照標(biāo)準(zhǔn)歸一化處理進(jìn)而進(jìn)行分析。采用敏感度計(jì)算公式［21，22］：

式中：L(xk)為巖土體參數(shù)xk的敏感度；F*為xk取力學(xué)參數(shù)標(biāo)準(zhǔn)值xk*時(shí)對(duì)應(yīng)的穩(wěn)定性區(qū)域占比；Δxk、ΔF為巖土體力學(xué)參數(shù)x*k和穩(wěn)定性占比F的變動(dòng)值。

計(jì)算過(guò)程中，對(duì)每個(gè)力學(xué)參數(shù)，選取庫(kù)岸邊坡低穩(wěn)定性區(qū)域占岸坡總庫(kù)岸段的百分比進(jìn)行敏感度計(jì)算。計(jì)算得到，改變黏聚力c、內(nèi)摩擦角φ、土體重度γ，敏感度分別為1.14、3.77、1.26。由此可知，各個(gè)巖土體力學(xué)參數(shù)對(duì)岸坡穩(wěn)定性的影響程度大小為：內(nèi)摩擦角φ＞土體重度γ＞黏聚力c。

5 結(jié) 論

（1）運(yùn)用Scoops 3D 模型計(jì)算不同工況下岸坡穩(wěn)定性的結(jié)果表明，庫(kù)水位下降速率為0.5、0.8、2、3 m/d 時(shí)，低穩(wěn)定性庫(kù)岸區(qū)域占去全庫(kù)岸段的百分比為：18.95%、19.44%、21.82%、24.84%；庫(kù)水位上升速率為0.8、2、3 m/d 時(shí)，低穩(wěn)定性庫(kù)岸區(qū)域占去全庫(kù)岸段的百分比為：18.59%、17.42%、17.02%。即庫(kù)岸邊坡穩(wěn)定性隨庫(kù)水位下降速率的增大而降低，而較大的水位上升速率對(duì)岸坡穩(wěn)定性有利。

（2）采用單因素敏感性分析方法研究了庫(kù)岸邊坡穩(wěn)定性對(duì)巖土體力學(xué)參數(shù)的敏感性，結(jié)果表明，敏感度由大到小依次是內(nèi)摩擦角φ、土體重度γ、黏聚力c?！?/p>