基于降雨入滲的Newmark模型改進(jìn)及地震滑坡危險(xiǎn)性預(yù)測(cè)研究

雷 真,李林銳,隆交鳳,陳敬男,楊 洋

(云南大學(xué) 建筑與規(guī)劃學(xué)院,云南 昆明 650500)

0 引言

滑坡是地震與降雨誘發(fā)的一種主要次生災(zāi)害,造成的經(jīng)濟(jì)損失往往比地震或強(qiáng)降雨造成的直接損失大得多。導(dǎo)致滑坡發(fā)生的因素有多種,其中降雨與地震兩種因素尤為重要。降雨滲透將導(dǎo)致斜坡土體含水量與重度增加,下滑力增大,抗剪強(qiáng)度降低。尤其在土體與巖體的交界面上,由于土、巖交界地層形成的年代普遍較晚。上部覆蓋的土層孔隙較大,也較松散,下部巖層又表現(xiàn)出較強(qiáng)的力學(xué)性質(zhì),在降水滲透的情況下,更容易導(dǎo)致滑坡失穩(wěn)。Green-Ampt(GA)降雨入滲模型是目前降雨滑坡分析中應(yīng)用最為廣泛的模型之一[1]。傳統(tǒng)模型計(jì)算簡(jiǎn)單、易于求解,但是該模型求解的是水平面,而未考慮斜面。目前基于該模型,大量學(xué)者對(duì)其進(jìn)行改進(jìn)修正。針對(duì)初始含水率因素,王千等[2]提出了初始含水率為線性分布時(shí)的GA入滲模型,Chen等[3]提出了初始含水率均勻的斜面降雨入滲GA模型。彭振陽等[4]基于分層假設(shè)的GA模型,引入濕潤(rùn)層等效導(dǎo)水率改進(jìn)了傳統(tǒng)GA入滲模型。劉衛(wèi)濤等[5]考慮土體非飽和特性,對(duì)傳統(tǒng)GA模型進(jìn)行修正,并應(yīng)用到斜坡穩(wěn)定性研究中。李秀珍等[6]結(jié)合Mein-Larson入滲模型進(jìn)行了飽和土、非飽和土與近似飽和土三種情況下的斜坡穩(wěn)定性分析,對(duì)Green-Ampt入滲模型應(yīng)用與飽和土與近似非飽和土斜坡穩(wěn)定性模型進(jìn)行了改進(jìn)。

Newmark累積位移模型[7]是Newmark于1965年首次提出且隨著GIS技術(shù)的不斷創(chuàng)新和發(fā)展,該方法被國(guó)內(nèi)外學(xué)者大量的應(yīng)用于地震滑坡危險(xiǎn)性評(píng)估研究中[8-12]。Newmark累積位移回歸方程可基于強(qiáng)震記錄數(shù)據(jù)進(jìn)行經(jīng)驗(yàn)擬合,Jibson[8]根據(jù)世界范圍內(nèi)30次地震的2 270條強(qiáng)震記錄得到的數(shù)據(jù)建立了多種位移回歸方程;徐光興等[9]采用Newmark方法基于強(qiáng)震動(dòng)記錄提出斜坡永久位移預(yù)測(cè)模型,驗(yàn)證了采取CAR臨界加速度比為參數(shù)的模型更加適用。許多學(xué)者采用Newmark方法對(duì)“5·12”大地震受災(zāi)區(qū)域進(jìn)行了地震滑坡危險(xiǎn)性評(píng)估分析與區(qū)劃圖的編制[10-11]。劉甲美等[12]基于Newmark模型開展了九寨溝地震滑坡快速評(píng)估,并得出該次地震滑坡主要以淺表型碎屑流和小規(guī)模崩塌居多的結(jié)論。金凱平等[13]提出對(duì)Newmark累積位移的修正,并將該模型應(yīng)用于蘆山地震進(jìn)行地震滑坡危險(xiǎn)性區(qū)劃。

然而目前采用Newmark模型進(jìn)行地震滑坡危險(xiǎn)預(yù)測(cè)中考慮降雨因素的研究較少,且降雨和地震耦合作用必將提高滑坡發(fā)生的概率,因此分析兩因素同時(shí)誘導(dǎo)滑坡發(fā)生的情況是非常有必要的。本文基于GA降雨入滲模型提出了修正的Newmark累積位移模型,對(duì)斜坡安全系數(shù)Fs進(jìn)行推導(dǎo)。以云南省魯?shù)榭h某區(qū)域?yàn)槔_展降雨和地震耦合作用下的滑坡危險(xiǎn)性預(yù)測(cè),并對(duì)預(yù)測(cè)結(jié)果進(jìn)行分析,為今后的滑坡危險(xiǎn)性預(yù)測(cè)提供一定的理論依據(jù)。

1 GA模型修正及地下水位確定

1.1 GA模型修正

降雨入滲后,地下水位可以通過一個(gè)水文模型來描述,如圖1所示。假設(shè)該水文模型的初始含水率均勻分布,圖中α為斜坡坡腳,q為降雨強(qiáng)度,則可通過Green-Ampt降雨入滲模型并簡(jiǎn)化為一維入滲問題推導(dǎo)濕潤(rùn)鋒深度h。

降雨累積入滲量與濕潤(rùn)鋒深度可通過質(zhì)量守恒理論表示為:

(1)

式中:F為降雨入滲總量;θs為土體飽和含水率;θi為土體初始含水率;h為濕潤(rùn)鋒深度。

在降雨強(qiáng)度q小于土體飽和導(dǎo)水系數(shù)KS時(shí),可認(rèn)為降雨全部被土體吸收,濕潤(rùn)鋒深度h可表示為:

(2)

如果降雨強(qiáng)度q大于土體飽和滲透系數(shù)KS,隨著降雨時(shí)間的增長(zhǎng),地面將在某個(gè)時(shí)刻開始積水,假設(shè)該臨界時(shí)刻為降雨歷時(shí)t=tp且濕潤(rùn)鋒深度以上土體均為飽和土體,則在積水前入滲深度可由式 (2)計(jì)算;令Sf為濕潤(rùn)峰深度h處的基質(zhì)吸力,則積水后垂直于斜面的降雨入滲深度可表示為:

(3)

Sf可由王紅閃等[14]提出的方法計(jì)算,即:

(4)

式中:S為土壤吸濕率。因此降雨入滲深度h隨時(shí)間t增加可表示為:

(5)

根據(jù)式(5)可以得到在臨界時(shí)刻tp時(shí),由于降雨入滲深度為關(guān)于時(shí)間t的連續(xù)型函數(shù),將tp帶入式(5)得到。

(6)

將式(3)代入整理可以得到臨界時(shí)刻tp的關(guān)系式。

(7)

1.2 地下水位的確定

隨著降雨入滲過程的發(fā)展,飽和土體內(nèi)地下水將沿著滲流管網(wǎng)平行于坡面進(jìn)行流動(dòng),因此濕潤(rùn)鋒深度也隨之發(fā)生變化,此過程可假設(shè)平行坡面的水流由一個(gè)體積為V的容器控制,則對(duì)于該體積,假定土體性質(zhì)是恒定的,則在時(shí)間t流出體積V的水流量可以用高度變化Δh(圖1)來表示:

(8)

式中:x為元素單位大小;s是元素單位間的距離。

根據(jù)Beven提出的飽和帶內(nèi)任何點(diǎn)的水力梯度等于模型坡度的假設(shè),則V可表示為:

V=KShstsinα

(9)

則修正后的入滲深度:

Hm=h±Δh

(10)

本文基于ArcGIS平臺(tái)采用D8單流向算法按上述式子計(jì)算Hm。D8流向算法即8個(gè)方向建模,是ArcGIS平臺(tái)水文分析重要的工具,其工作原理如圖2所示。本文假設(shè)雨水入滲至土體與巖石交界面處不再下滲,最后將修正后的入滲深度視為地下水位深度Hw,即Hw=Hm。

圖2 D8算法工作原理Fig.2 Working principle of D8 algorithm

2 Newmark模型的修正

Newmark累積邊坡位移模型是Newmark于1965年提出的一種極限平衡理論模型,雖然該模型未考慮豎向地震作用,但仍在地震滑坡危險(xiǎn)性預(yù)測(cè)分析中被廣大學(xué)者使用。其原理是將邊坡視為沿斜面滑動(dòng)的剛塊,當(dāng)外力達(dá)到使邊坡由穩(wěn)定到失穩(wěn)的最小力時(shí),邊坡開始滑動(dòng),在地震發(fā)生時(shí)可將地震引起的外力用mac表示且方向平行于坡面。在降雨期間,隨著雨水入滲,土體含水量逐漸增加,斜坡抗剪強(qiáng)度降低,孔隙水壓力不斷增大,致使滑坡發(fā)生的概率增加。本文假設(shè)在地震發(fā)生時(shí)孔隙水壓力不隨著地震作用發(fā)生變化,不考慮地震所產(chǎn)生的超孔隙水壓力,并將Green-Ampt模型與Newmark模型相結(jié)合,對(duì)Newmark模型進(jìn)行改進(jìn)(圖3)。

圖3 改進(jìn)的Newmark模型Fig.3 Improved Newmark model

在無降雨期間,坡體可視為干燥條件不考慮孔隙水壓力的作用,則靜態(tài)安全系數(shù)FS可表示為:

(11)

(12)

在降雨期間,隨著降雨時(shí)間的增加,孔隙水壓力不斷增大,此時(shí)安全系數(shù)FS采用如下式計(jì)算:

FZ=csatb+(P-U)tanφsatcosα

(13)

FX=Psinα

(14)

P=b[γ(H-Hw)+Hwγsat]

(15)

U=γwbHw

(16)

(17)

地震作用下,當(dāng)坡體達(dá)到極限平衡時(shí),即FS=1,則臨界加速度ac可表示為:

ac=(FS-1)gsinα

(18)

式中:P為坡體重度(kN);H為坡體厚度(m);U為孔隙水壓力(kN);c和csat分別為土體的天然和飽和黏聚力(MPa);φ和φsat分別為土體的天然和飽和內(nèi)摩擦角(°);m為坡體質(zhì)量;g為重力加速度;ac為致使坡體下滑失穩(wěn)的臨界加速度;γ、γsat和γw分別表示土體天然重度、土體飽和重度與水的重度(kN/m3);b為坡體長(zhǎng)度(m)。

目前計(jì)算滑坡累積位移DN的模型有很多,代表性的模型主要分為三種,第一種是積分模型:

DN=?[a(t)-ac]dtdt

(19)

式中:a(t)為加速度時(shí)程,由于中強(qiáng)震動(dòng)記錄比較缺乏且大量精確積分難度較高,因此該模型應(yīng)用較為不便;第二種模型則是以臨界加速度比(CAR)為參數(shù)計(jì)算的模型;第三種則是以Arias強(qiáng)度與CAR為參數(shù)的計(jì)算模型。通過事實(shí)分析,Arias強(qiáng)度雖能較好地反映出地震動(dòng)特性,但其獲得的過程復(fù)雜且精度未見提高,因此不建議使用[9]。

lg(DN)=0.546MW-0.425CAR-1.086±0.263

(20)

3 模型應(yīng)用

3.1 試驗(yàn)區(qū)基本地質(zhì)特征及計(jì)算參數(shù)

本文選取魯?shù)榭h某個(gè)區(qū)域作為試驗(yàn)區(qū)進(jìn)行該方法的滑坡預(yù)測(cè)應(yīng)用。魯?shù)榭h地處云貴高原西北部,滇東北高原南部(103°09′～103°40′E,26°59′～27°32′N)[15],為昭通市一縣城。南北向構(gòu)造與東北向構(gòu)造共同組成了魯?shù)榭h的主干構(gòu)造。其中,南北向構(gòu)造的主要斷層與褶皺有野牛塘斷層、五里牌斷層、韋家渡背斜、新街子向斜等,東北向構(gòu)造主要為騾馬口構(gòu)造帶[16]。魯?shù)榭h區(qū)巖性主要為灰?guī)r、白云巖、砂巖、泥巖、頁巖和玄武巖,其中玄武巖大片分布。泥巖與頁巖巖性較為軟弱,加之玄武巖風(fēng)化嚴(yán)重,導(dǎo)致研究區(qū)地質(zhì)環(huán)境整體脆弱,為滑坡發(fā)育提供基礎(chǔ)[17]。根據(jù)當(dāng)?shù)氐刭|(zhì)特征及降雨情況資料,該試驗(yàn)區(qū)相關(guān)計(jì)算參數(shù)選取列于表1。

通過該研究區(qū)Dem(30 m)數(shù)據(jù)生成坡度圖,圖4所示。根據(jù)魯?shù)榭h以往的地震記錄,本次模擬取MW=6.1。根據(jù)魯?shù)榈卣饡r(shí),美國(guó)地質(zhì)勘探局(USGS)所公布的PGA(PGA，地震動(dòng)峰值加速度)數(shù)據(jù)如圖5所示。由于地震誘發(fā)的滑坡大多數(shù)為淺層滑坡,根據(jù)前人經(jīng)驗(yàn)取滑坡厚度H=3 m[10-12]。

圖4 坡度圖[單位:(°)]Fig.4 Slope map [Unit:(°)]

圖5 PGA圖(單位:g)Fig.5 PGA map (Unit:g)

3.2 滑坡位移計(jì)算及危險(xiǎn)區(qū)劃

通過假設(shè)該試驗(yàn)區(qū)在均衡降雨強(qiáng)度q作用下,由表1參數(shù)可知假設(shè)的降雨強(qiáng)度q大于土體飽和滲透系數(shù)KS,則存在臨界降雨時(shí)刻tp使得地面開始積水,通過1.1節(jié)公式可推導(dǎo)tp:

表1 研究區(qū)參數(shù)表Table 1 Parameters of study area

(21)

則通過ArcGIS平臺(tái)計(jì)算得出研究區(qū)開始積水的臨界時(shí)刻tp(圖6),研究區(qū)中白色區(qū)域?yàn)槠露刃∮?0°的區(qū)域,依據(jù)經(jīng)驗(yàn)[13],該部分區(qū)域不發(fā)生滑坡,因此不參與計(jì)算。

圖6 臨界時(shí)刻tp圖(單位:h)Fig.6 Critical time tp map (Unit:h)

為方便計(jì)算,本文取該研究區(qū)tp數(shù)據(jù)的中值tp=5 h,分別計(jì)算降雨歷時(shí)5 h(此時(shí)地面還未積水)與降雨歷時(shí)10 h地下水位深度Hw(圖7)。

圖7 地下水位深度圖(單位:m)Fig.7 Depth of groundwater level (Unit:m)

通過計(jì)算分別得到無降雨時(shí)期(t=0 h)、降雨歷時(shí)t=5 h(地面無積水)以及降雨歷時(shí)t=10 h(地面積水)的地震滑坡累計(jì)位移Dn。Dn值的大小與滑坡發(fā)生的危險(xiǎn)性呈正相關(guān)的趨勢(shì),參考前人經(jīng)驗(yàn)并根據(jù)以往實(shí)際滑坡的觀測(cè)數(shù)據(jù)可以發(fā)現(xiàn),在安全區(qū)域或少量發(fā)生滑坡的區(qū)域,其累積位移計(jì)算值往往小于0.5 cm,危險(xiǎn)程度較低;而發(fā)生大規(guī)模地震滑坡災(zāi)害的區(qū)域,其累積位移值均處于2 cm以上[12,18],所以將累積位移值分為三類“DN<0.5 cm(危險(xiǎn)程度低)、0.5 cm2 cm(危險(xiǎn)程度高)”。最終結(jié)果如圖8所示。

圖8 滑坡預(yù)測(cè)區(qū)劃圖Fig.8 Landslide prediction zoning map

3.3 滑坡危險(xiǎn)區(qū)分布特征及結(jié)果分析

依據(jù)滑坡預(yù)測(cè)區(qū)劃圖可知,滑坡高危險(xiǎn)區(qū)域主要集中在PGA與坡度值高的地區(qū),表明地震滑坡危險(xiǎn)程度與坡度和PGA具有較高的相關(guān)性。隨著降雨時(shí)間的增加,危險(xiǎn)程度低的部分區(qū)域逐步轉(zhuǎn)變?yōu)槲ｋU(xiǎn)程度中的區(qū)域,部分中危險(xiǎn)程度區(qū)域轉(zhuǎn)變?yōu)楦呶ｋU(xiǎn)區(qū)域。在降雨期間,隨著降雨時(shí)間的增加(t=0 h→10 h),區(qū)域危險(xiǎn)性程度不斷上升,危險(xiǎn)面積占比同步提高。與未降雨情況相比,降雨歷時(shí)5 h、10 h情況下地震滑坡低危險(xiǎn)區(qū)域占比面積從51%分別下降至35%、33%,下降幅度分別為31%、35%。高危險(xiǎn)區(qū)域面積從1%分別提高至9%、12%,提高幅度分別為800%、1 100%。中危險(xiǎn)區(qū)域面積占比從48%分別增加至56%、55%,可以看出降雨達(dá)到5 h時(shí),低危險(xiǎn)區(qū)域面積占比下降幅度較大,高危險(xiǎn)區(qū)域面積占比顯著增加并隨著降雨時(shí)間增加而逐漸增大。中危險(xiǎn)區(qū)域面積占比出現(xiàn)了先增加后降低1%的情況,歸結(jié)原因是降雨5 h后到t=10 h時(shí)僅有極少部分的滑坡低危險(xiǎn)程度區(qū)域轉(zhuǎn)為了滑坡中危險(xiǎn)程度區(qū)域,而中危險(xiǎn)程度區(qū)域轉(zhuǎn)為高危險(xiǎn)程度區(qū)域的面積相對(duì)較多。

本文利用既定的研究區(qū)參數(shù)(黏聚力、內(nèi)摩擦角、土體重度),分析了坡度、地下水位深度(入滲深度)及PGA三種因素對(duì)滑坡位移值的影響(圖9)。

圖9(a)中為無降雨情況下,不同坡度(15°、20°、30°、35°)對(duì)應(yīng)計(jì)算出的累積位移曲線圖,在坡度值為15°、PGA為0.1g時(shí)所計(jì)算的位移值為負(fù)數(shù),表明坡體不會(huì)發(fā)生滑坡,本文中取位移值為-0.01,即坡體穩(wěn)定。由圖9(a)可知:(1)在PGA值為0.2g時(shí)坡度值15°對(duì)應(yīng)的位移值最小(0.44 cm),坡度值35°對(duì)應(yīng)的位移值最大(1.63 cm);(2)PGA為0.2g時(shí),不同坡度所對(duì)應(yīng)的累計(jì)位移值相差較大,而PGA為0.5g時(shí)不同坡度所對(duì)應(yīng)的累計(jì)位移值差距較小,主要原因是當(dāng)坡度值較小時(shí),PGA值越大對(duì)累計(jì)位移的影響度越大,坡度影響度隨之減小。

圖9 影響因子對(duì)滑坡位移影響圖Fig.9 Influence diagram of influencing factors on landslide displacement

圖9(b)為降雨情況下,降雨入滲深度Hw=0.6 m時(shí)不同坡度下計(jì)算出的累計(jì)位移,結(jié)果表明:(1) 降雨情況下,PGA=0.2g時(shí),坡度為15°時(shí)對(duì)應(yīng)的累計(jì)位移同樣為最小值(1.31 cm),坡度值35°對(duì)應(yīng)的位移值同樣為最大值(2.12 cm);(2)在坡度為15°、PGA分別為0.2g、0.5g時(shí)的累計(jì)位移相比于圖9(a)中無降雨情況下分別增加0.87 cm、0.84 cm,增加幅度分別為198%、30%,表明當(dāng)坡度與PGA等其他參數(shù)值相同時(shí),降雨入滲作用將增加滑坡位移值,當(dāng)PGA值與坡度值較小時(shí),增加幅度大;當(dāng)PGA值與坡度值較大時(shí),增加幅度小。

圖9(c)列出了當(dāng)坡度值30°,不同入滲深度(0.2 m、0.4 m、0.6 m、0.8 m)下的累計(jì)位移值曲線。結(jié)果表明:(1)PGA為0.1g時(shí),入滲深度為0.2 m所對(duì)應(yīng)的累計(jì)位移值最小(1.65 cm),入滲深度為0.8 m所對(duì)應(yīng)的位移值最大(1.82 cm),增加幅度為10.3%;(2)當(dāng)PGA值為0.5g時(shí),四種入滲深度所計(jì)算出的累積位移值幾乎相同。主要原因是坡度與PGA達(dá)到較高數(shù)值時(shí),降雨入滲因素對(duì)累積位移的影響程度較小。

4 結(jié)論

本文結(jié)合GA降雨模型與Newmark累積位移模型在降雨地震耦合作用下對(duì)安全系數(shù)Fs進(jìn)行了推導(dǎo),對(duì)魯?shù)榭h某區(qū)域開展無降雨情況、降雨無積水與降雨積水三種情況地震滑坡危險(xiǎn)性預(yù)測(cè)分析,并對(duì)滑坡位移值進(jìn)行坡度、降雨入滲及PGA的因素分析。

(1) 相對(duì)無降雨情況,加入降雨因素后,隨著降雨時(shí)間的增加,滑坡低危險(xiǎn)區(qū)域逐漸減少轉(zhuǎn)變?yōu)榛轮形ｋU(xiǎn)區(qū)域,地震滑坡高危險(xiǎn)程度區(qū)域增加,計(jì)算區(qū)域中滑坡高危險(xiǎn)程度區(qū)域占比上升幅度為1 100%,滑坡低危險(xiǎn)程度區(qū)域占比下降幅度為35%,說明降雨與地震耦合作用對(duì)誘發(fā)滑坡的影響不容忽視,降雨較大的增加了地震滑坡的可能性。

(2) 當(dāng)坡度較小時(shí),PGA值對(duì)滑坡累積位移影響相比坡度影響較大,當(dāng)坡度增加至一定數(shù)值時(shí),PGA對(duì)滑坡位移影響程度遠(yuǎn)小于坡度;與無降雨情況相比,降雨入滲將使滑坡位移值增大,入滲深度在PGA值較小時(shí)對(duì)滑坡位移值影響較大,而在PGA值較大時(shí)對(duì)滑坡位移值相對(duì)影響較小。

(3) 降雨與地震耦合作用下的計(jì)算模型與降雨強(qiáng)度q、降雨時(shí)間t、巖石物理力學(xué)參數(shù)以及Dem數(shù)據(jù)的精度有較大關(guān)系,在進(jìn)行耦合作用下滑坡預(yù)測(cè)分析時(shí),獲取精確的計(jì)算參數(shù)將極大提高該方法預(yù)測(cè)的精確度。