趙紀(jì)飛，黃嘉悅，侯曉坤，李同錄，藺曉燕

(1.長(zhǎng)安大學(xué) 地質(zhì)工程與測(cè)繪學(xué)院，陜西 西安 710054； 2. 中國(guó)電建集團(tuán)西北勘測(cè)設(shè)計(jì)研究院有限公司 陜西 西安 710054； 3. 中國(guó)有色金屬工業(yè)西安勘察設(shè)計(jì)研究院，陜西 西安 710054； 4. 西安石油大學(xué) 地球科學(xué)與工程學(xué)院，陜西 西安 710065)

灌溉誘發(fā)的黑方臺(tái)黃土滑坡泥流機(jī)理分析*

趙紀(jì)飛1,2，黃嘉悅3，侯曉坤1，李同錄1，藺曉燕4

(1.長(zhǎng)安大學(xué) 地質(zhì)工程與測(cè)繪學(xué)院，陜西 西安 710054； 2. 中國(guó)電建集團(tuán)西北勘測(cè)設(shè)計(jì)研究院有限公司 陜西 西安 710054； 3. 中國(guó)有色金屬工業(yè)西安勘察設(shè)計(jì)研究院，陜西 西安 710054； 4. 西安石油大學(xué) 地球科學(xué)與工程學(xué)院，陜西 西安 710065)

針對(duì)黑方臺(tái)灌溉觸發(fā)的黃土滑坡泥流，通過(guò)野外測(cè)量查明坡體地質(zhì)結(jié)構(gòu)，室內(nèi)基本物理力學(xué)性質(zhì)測(cè)測(cè)試、三軸試驗(yàn)和土水特征曲線的測(cè)試等，建立二維有限元滲流模型和彈塑性模型，分析灌溉條件下地下水位的變化特征和坡體的穩(wěn)定性；并提取滑面上的應(yīng)力狀態(tài)分析其剪應(yīng)力和抗剪強(qiáng)度的變化特點(diǎn)和屈服特征，結(jié)果表明該滑坡為純黃土滑坡，黃土的飽和含水率高于液限，擾動(dòng)易發(fā)生液化；黃土下層的粉質(zhì)黏土層飽和滲透系數(shù)小，儲(chǔ)水能力強(qiáng)，易積水形成潛水位。灌溉條件下，地下水位逐漸升高，邊坡穩(wěn)定性開(kāi)始基本保持不變或緩慢降低，其快速降低始于水位開(kāi)始接觸滑面；灌溉前滑面上的抗剪強(qiáng)度均大于剪應(yīng)力，隨著灌溉進(jìn)行，坡腳先發(fā)生屈服，且范圍不斷擴(kuò)大，由于應(yīng)力重分布作用，滑坡后緣也逐漸發(fā)生屈服，最終沿整個(gè)滑面上土體的剪應(yīng)力均達(dá)到其抗剪強(qiáng)度。該滑坡為牽引式滑坡。

黃土滑坡；泥流；野外測(cè)量；室內(nèi)試驗(yàn)；有限元模型；牽引式滑坡；甘肅黑方臺(tái)

黃土地區(qū)為干旱半干旱氣候，為了農(nóng)業(yè)發(fā)展及人們生活需要，自1970年代以來(lái)，修建了大量的引水灌溉工程。這些工程對(duì)黃土地區(qū)的農(nóng)業(yè)生產(chǎn)、經(jīng)濟(jì)發(fā)展發(fā)揮了重要的作用。但大量的引水灌溉不同程度的改變了灌區(qū)的自然地質(zhì)條件，造成了灌區(qū)的生態(tài)破壞，從而引發(fā)了一系列滑坡災(zāi)害。甘肅黑方臺(tái)是灌溉誘發(fā)的黃土滑坡、滑坡泥流最典型最集中的區(qū)域。該地區(qū)自1968年提水灌溉到1980年代初期方臺(tái)村塬邊開(kāi)始出現(xiàn)滑坡以來(lái)，發(fā)生了許多大規(guī)模的滑坡和滑坡泥流，累積造成37人死亡，100多人受傷，并迫使4所小學(xué)搬遷等。對(duì)典型案例進(jìn)行專門研究，查明其形成機(jī)理,對(duì)災(zāi)害預(yù)測(cè)、減少人民生命和財(cái)產(chǎn)損失，具有重要的科學(xué)意義和社會(huì)意義。

對(duì)黃土邊坡破壞機(jī)理的研究有定性法和定量法。定性法是根據(jù)邊坡破壞時(shí)的工程地質(zhì)條件推斷邊坡破壞形式，如地下水位、水源、邊坡幾何形狀等[1-2]；或根據(jù)室內(nèi)三軸或直剪試驗(yàn)得到的土樣的應(yīng)力應(yīng)變特征進(jìn)而推測(cè)邊坡的破壞模式[3-4]。但這只是對(duì)邊坡破壞機(jī)理的推測(cè)。定量法是采用極限平衡法、有限差分法等計(jì)算邊坡穩(wěn)定性。傳統(tǒng)計(jì)算方法是建立在飽和土力學(xué)基礎(chǔ)上的，采用飽和土的抗剪強(qiáng)度理論和強(qiáng)度參數(shù)。隨著非飽和土力學(xué)的發(fā)展，非飽和土抗剪強(qiáng)度的建立，人們逐漸意識(shí)到水分入滲會(huì)引起基質(zhì)吸力降低，使得土體的抗剪強(qiáng)度減小，邊坡穩(wěn)定性降低[5-8]。近年來(lái)，不斷有學(xué)者采用考慮飽和非飽和滲透過(guò)程的數(shù)值分析方法研究入滲條件下邊坡穩(wěn)定性[9-12]，該法考慮了水分入滲過(guò)程中坡體內(nèi)基質(zhì)吸力的降低及其對(duì)邊坡穩(wěn)定性的影響，一定程度上反映了邊坡的破壞過(guò)程。但對(duì)邊坡穩(wěn)定性分析仍是采用的極限平衡法，只能得到邊坡整體的穩(wěn)定性，并未能反映滑面上抗剪強(qiáng)度和剪應(yīng)力的變化情況。

綜上所述，本文以甘肅黑方臺(tái)一典型滑坡泥流為例，在查明坡體地質(zhì)結(jié)構(gòu)的基礎(chǔ)上，取原狀樣進(jìn)行基本物理力學(xué)指標(biāo)的測(cè)試、三軸試驗(yàn)和土水特征曲線測(cè)試，并建立數(shù)值模型，進(jìn)行飽和非飽和滲流分析和彈塑性有限元分析得到灌溉條件下邊坡的水分場(chǎng)和應(yīng)力場(chǎng)，在此基礎(chǔ)上結(jié)合非飽和土的抗剪強(qiáng)度理論計(jì)算邊坡的穩(wěn)定系數(shù)，通過(guò)分析滑面上抗剪強(qiáng)度和剪應(yīng)力的變化情況揭示邊坡的變形破壞機(jī)理. 所用軟件為Geo-studio。

1 理論介紹

非飽和滲流包括穩(wěn)態(tài)滲流和瞬態(tài)滲流。利用質(zhì)量守恒定律和達(dá)西定律，可得到滲流的控制方程(二維)：

(1)

式中：kx(hm)、kz(hm)分別為x方向和豎直方向滲透系數(shù)；hm為基質(zhì)吸力水頭,其大小為ψ/g,ψ為基質(zhì)吸力，g為重力加速度，ρ為水密度，θ為體積含水率，t為時(shí)間。

若式(1)右側(cè)為0，為穩(wěn)態(tài)滲流；若等式右側(cè)不為0，則為瞬態(tài)滲流。非飽和滲流的控制方程中并未考慮土骨架變形，故該模型為非水力耦合模型。公式的求解需土水特征曲線、非飽和滲透系數(shù)、邊界條件以及初始含水狀態(tài)等。

非飽和滲流的不同時(shí)間步，通過(guò)彈塑性有限元分析可計(jì)算得到不同狀態(tài)下坡體的應(yīng)力場(chǎng)和應(yīng)變場(chǎng)，進(jìn)而分析其破壞機(jī)理。土體的本構(gòu)模型可采用理想彈塑性模型，即破壞前土體為線彈性體，當(dāng)應(yīng)力達(dá)到其強(qiáng)度時(shí)，土體破壞。其中飽和區(qū)土體的強(qiáng)度由常用的莫爾庫(kù)強(qiáng)度公式確定，見(jiàn)式(2)；而非飽和土的抗剪強(qiáng)度公式目前主要有單變量強(qiáng)度理論(Bishop)和雙變量強(qiáng)度理論，考慮到抗剪強(qiáng)度隨飽和度的變化，Vanapalli 等[13], Oberg等[14], Khalili等[15]等提出利用土水特征曲線預(yù)測(cè)非飽和土的抗剪強(qiáng)度。由于本文的有限元計(jì)算是建立在非飽和滲流分析的基礎(chǔ)上，開(kāi)始計(jì)算時(shí)邊坡的水分場(chǎng)已知，利用土水特征曲線預(yù)測(cè)非飽和土的抗剪強(qiáng)度簡(jiǎn)單易行，因此本文采用Vanapalli[13]提出的非飽和土強(qiáng)度公式進(jìn)行計(jì)算，見(jiàn)式(3)。值得注意的是由于軟件并沒(méi)有專門計(jì)算非飽和區(qū)土體變形的本構(gòu)模型，非飽和區(qū)土體也采用彈性理想塑形本構(gòu)模型，但其彈性模量取值較飽和區(qū)彈性模量大。盡管這種計(jì)算不如采用非飽和本構(gòu)模型計(jì)算結(jié)果準(zhǔn)確，但其仍能反映出邊坡應(yīng)力場(chǎng)變化的一些重要特點(diǎn)(后面將會(huì)描述)。

τf=c’+σ’tanφ’；

(2)

(3)

式中：τf,為抗剪強(qiáng)度；c’,φ’為有效強(qiáng)度參數(shù)；σ’為有效應(yīng)力；σ為總應(yīng)力；ua為孔隙氣壓；uw為孔隙水壓；(ua-uw)為基質(zhì)吸力；θ為體積含水率；θs為飽和體積含水率；θr為殘余體積含水率。

若破壞面已知，則根據(jù)滑面上的主應(yīng)力和剪應(yīng)力通過(guò)式(4)和式(5)可得到下滑力，跟據(jù)式(2)和式(3)強(qiáng)度理論可計(jì)算得到其抗剪強(qiáng)度[16]。以破壞面長(zhǎng)度L為x軸，沿滑動(dòng)面上的下滑力、抗剪強(qiáng)度為Y軸，繪制τf,τa曲線，通過(guò)分析曲線的變化特征，可闡明邊坡的破壞過(guò)程；此外滑面上剪應(yīng)力和抗剪強(qiáng)度積分的比值即為邊坡的穩(wěn)定系數(shù)[17-18]，見(jiàn)式(6)。該計(jì)算方法無(wú)假設(shè)條件，其計(jì)算結(jié)果取決于有限元計(jì)算精度。將此計(jì)算結(jié)果同常用的極限平衡法的計(jì)算結(jié)果進(jìn)行對(duì)比分析，以檢驗(yàn)計(jì)算結(jié)果的準(zhǔn)確性。

(4)

(5)

(6)

式中：σα，τα分別為滑動(dòng)面上的正應(yīng)力和剪應(yīng)力；a為滑面傾角；xA和xB分別為剪入口和剪出口的橫坐標(biāo)；Fs為穩(wěn)定系數(shù)。

2 室內(nèi)試驗(yàn)及結(jié)果

本文所選案例為黑方臺(tái)焦家一典型滑坡泥流，實(shí)測(cè)地質(zhì)剖面如圖1所示。該剖面土層厚度共52 m，其中0～24 m為Q3黃土，呈淡黃色，結(jié)構(gòu)疏松多孔，具有垂直節(jié)理，即馬蘭黃土；24～44 m為Q2黃土，呈淡黃色，較致密，具有垂直節(jié)理；44～52 m為沖積粉質(zhì)粘土，具有水平層理，呈棕紅色，粘粒含量高，致密堅(jiān)硬；該層之下為厚5 m的卵石層，磨圓度好，為黃河四級(jí)階地堆積；卵石層之下為緩傾的白堊系紫紅色砂泥巖，基巖產(chǎn)狀為172°∠27°。滑坡發(fā)生時(shí)的地下水位埋深約34 m，位于Q2黃土層。而滑坡剪出口位于Q2黃土層中，距坡頂約45 m，為純黃土滑坡?；掳l(fā)生后，隨即轉(zhuǎn)化為泥流，順開(kāi)闊的坡面流下，形成坡面泥流，至坡腳處才停留堆積下來(lái)，水平滑距約280 m。分別在Q3黃土、Q2黃土和粉質(zhì)黏土層中取土進(jìn)行基本物理參數(shù)和土水特征曲線的測(cè)試。由于邊坡破壞發(fā)生在黃土層中，故只對(duì)Q3和Q2黃土進(jìn)行飽和三軸試驗(yàn)測(cè)得其強(qiáng)度參數(shù)，每層土的取樣位置見(jiàn)圖1。基本物理參數(shù)見(jiàn)表1，其中顆粒分析采用 Battersize 激光粒度儀測(cè)定，飽和滲透系數(shù)采用常規(guī)變水頭滲透試驗(yàn)測(cè)定。由表1可知，Q3和Q2黃土的飽和含水率均大于液限，說(shuō)明當(dāng)土樣受擾動(dòng)時(shí)，亦發(fā)生液化。黃土的飽和滲透系數(shù)遠(yuǎn)遠(yuǎn)大于粉質(zhì)粘土。

表1 Q3黃土、Q2黃土和粉質(zhì)黏土層的基本物理性質(zhì)參數(shù)

圖1 黑方臺(tái)焦家滑坡地質(zhì)剖面圖

2.1 三軸試驗(yàn)

本次試驗(yàn)所用試驗(yàn)儀器為南京土壤儀器廠生產(chǎn)的SLB-1型應(yīng)力應(yīng)變控制式三軸儀。試樣直徑6.18 cm， 高12.5 cm。按《土工試驗(yàn)規(guī)范》操作步驟將其試樣安置在三軸儀的壓力室后，先采用水頭飽和法將土體內(nèi)的大部分氣體排出，此過(guò)程需要1～2 d，然后采用反壓飽和的方法將土體飽和，反壓為200 kPa，待孔壓系數(shù)B>95%時(shí)飽和結(jié)束。之后打開(kāi)排水閥門，分別將試樣在100 kPa、200 kPa、300 kPa、400 kPa和500 kPa的圍壓下固結(jié)，待超孔隙水壓力完全消散或消散95%時(shí)，固結(jié)結(jié)束。關(guān)閉排水閥門，開(kāi)始剪切，剪切速率按規(guī)范要求設(shè)置為0.1mm/min。共10組試驗(yàn)。

圖2和圖3分別為Q3黃土和Q2黃土的三軸試驗(yàn)結(jié)果，從中可知兩種黃土的應(yīng)力應(yīng)變關(guān)系曲線和總(有效)應(yīng)力路徑曲線等基本相似。從圖2a和圖3a中可知，土體的應(yīng)力隨應(yīng)變的增大逐漸增大，基本表現(xiàn)出應(yīng)變硬化或理想塑形的特征，僅圖2a中400 kPa圍壓下的應(yīng)力應(yīng)變曲線有稍許軟化特征。圖2b和圖3b分別為Q3和Q2黃土的(有效)應(yīng)力路徑曲線，取應(yīng)變?yōu)?0%的狀態(tài)點(diǎn)為破壞點(diǎn)，繪制其(有效)破壞主應(yīng)力跡線(Kf′)Kf線, 根據(jù)其傾角a和截距b，利用公式(7)，可求得土體的有效抗剪強(qiáng)度參數(shù)和總強(qiáng)度參數(shù)。其中Q3黃土有效強(qiáng)度參數(shù)為c′=0.12 kPa,φ′=32.0°；總抗剪強(qiáng)度參數(shù)為c=0.15 kPa，φ=19.3°；Q2黃土的有效抗剪強(qiáng)度參數(shù)為c′=1.9 kPa,φ′=32.7°；總抗剪強(qiáng)度參數(shù)為c=0.7 kPa，φ=21.1°。對(duì)比可知，Q2和Q3黃土的抗剪強(qiáng)度參數(shù)大小接近。

(7)

2.2 土水特征曲線及非飽和滲透系數(shù)

采用張力計(jì)法測(cè)定Q2、Q3原狀黃土和粉質(zhì)粘土的土水特征曲線。所用張力計(jì)為浙江托普儀器公司生產(chǎn)的TEN型張力計(jì)，可測(cè)得的基質(zhì)吸力的范圍為0～100 kPa， 精度為2 kPa。 張力計(jì)使用前，采用真空飽和法將陶瓷頭飽和。試驗(yàn)所用土樣為大塊正方體型原狀樣，邊長(zhǎng)約為30 cm，周圍用塑料膜包裹嚴(yán)實(shí)防治水分散失。試驗(yàn)時(shí)將土樣用同張力計(jì)直徑相近的鉆頭打孔，孔底取樣測(cè)定其質(zhì)量含水率，清理孔隙中的散土后將充滿水的張力計(jì)插入孔中，并用細(xì)沙填充孔隙后，將土樣表層和張力計(jì)的接觸部位密封，記錄張力計(jì)讀數(shù)直到讀數(shù)基本穩(wěn)定。此穩(wěn)定值即為該土樣在此含水率下的基質(zhì)吸力。此時(shí)將張力計(jì)取出，土樣灑水，密封至少10 d，之后重復(fù)上述操作步驟進(jìn)行下一級(jí)基質(zhì)吸力的測(cè)量。

上述測(cè)試方法只能得到一些離散的數(shù)據(jù)點(diǎn)，國(guó)內(nèi)外學(xué)者提出多種數(shù)學(xué)模型反映非飽和土體基質(zhì)吸力和體積含水率的關(guān)系曲線，常見(jiàn)的模型有Gardenr模型、Van Genuchten模型、Fredlund and Xing模型。

圖2 Q3黃土三軸試驗(yàn)結(jié)果

圖3 Q2黃土三軸試驗(yàn)結(jié)果

本文采用VG 模型[19](式8)對(duì)試驗(yàn)所測(cè)數(shù)據(jù)進(jìn)行擬合，結(jié)果見(jiàn)圖4，從中可知，相關(guān)系數(shù)均大于0.99,擬合效果較好。

(8)

式中：θw為體積含水率；θr為殘余體積含水率；θs為飽和體積含水率；a、m、n為擬合參數(shù)，其中m=1-1/n。

圖4 Q3黃土、Q2黃土和粉質(zhì)黏土層的土水特征曲線

土體非飽和滲透系數(shù)的獲得方法有直接法和間接法[20]。直接法是采用瞬態(tài)剖面法、穩(wěn)態(tài)法等測(cè)試手段直接測(cè)量非飽和滲透系數(shù)。該法操作復(fù)雜，耗時(shí)長(zhǎng)。間接法是土體的土水特征曲線法和相應(yīng)的經(jīng)驗(yàn)公式求取非飽和滲透系數(shù)。盡管預(yù)測(cè)結(jié)果和實(shí)測(cè)值稍有差異，但數(shù)量級(jí)大體相同，且可基本反映出非飽和滲透系數(shù)隨飽和度的變化情況。考慮到試驗(yàn)成本問(wèn)題，本文采用后者得到各層土體的非飽和滲透系數(shù)，采用VG模型[19](式9)計(jì)算得到的Q2、Q3黃土和粉質(zhì)粘土的非飽和滲透系數(shù)見(jiàn)圖5，從中可知Q3黃土的非飽和滲透系數(shù)最大，且黃土的非飽和滲透系數(shù)遠(yuǎn)遠(yuǎn)大于粉質(zhì)粘土，這有利于水分在黃土中積聚。

圖5 Q3 、Q2黃土和粉質(zhì)粘土滲透系數(shù)函數(shù)曲線

(9)

式中：kw為非飽和土的滲透系數(shù)；ks為飽和滲透系數(shù)；a、m、n為式9擬合參數(shù)。

3 滑坡破壞機(jī)理分析

圖6 數(shù)值模型圖

根據(jù)實(shí)測(cè)剖面建立計(jì)算模型(圖6)，土層自下而上分別為Q3黃土、Q2黃土，粉質(zhì)粘土和第三系泥巖，其土層厚度依次為24 m、20 m和8 m。值得注意的是在邊坡地質(zhì)剖面(圖1)中，粉質(zhì)粘土下層為砂卵石層，考慮到其對(duì)上部邊坡穩(wěn)定性影響小，且厚度較薄，故將其省去。邊坡的破壞面根據(jù)實(shí)測(cè)值確定(圖1)，剪出口設(shè)定在Q2黃土層底部，剪入口在坡頂。有限元網(wǎng)格尺寸設(shè)置為1 m，程序根據(jù)需要自動(dòng)劃分為三角形網(wǎng)格和四邊形網(wǎng)格，共10 993個(gè)單元。由于Geo-studio軟件的有限元求解采用增量法，故對(duì)滲流分析和彈塑性有限元分析均需設(shè)定其初始條件，即初始水分場(chǎng)和應(yīng)力場(chǎng)。初始水分場(chǎng)通過(guò)年均降雨條件下的穩(wěn)態(tài)滲流分析得到；初始應(yīng)力場(chǎng)可由軟件可自行計(jì)算獲得。此外有限元的求解須知其邊界條件，滲流分析的邊界包括水頭邊界和流量邊界，具體為地下水位的位置和地表灌溉的流量，此外還包括坡面為潛在滲流面(圖6)。黑方臺(tái)地區(qū)灌溉的地下水位位于粉質(zhì)粘土層內(nèi)，灌溉流量根據(jù)當(dāng)?shù)氐膶?shí)際灌溉量確定，據(jù)調(diào)查該區(qū)的灌溉主要分為春灌和冬灌，分別集中在5月份和11月份，灌溉流量為設(shè)置為16 mm/d，其他月份為停灌階段,流量邊界設(shè)置為0 mm/d。該滲流過(guò)程模擬了黑方臺(tái)地區(qū)自1968年灌溉以來(lái)到1980年該區(qū)開(kāi)始頻繁出現(xiàn)滑坡期間的滲流過(guò)程，共12年。

彈塑性有限元計(jì)算所需的邊界條件為位移邊界條件，即邊坡底部為固定端，各網(wǎng)格單元節(jié)點(diǎn)處沒(méi)有橫向位移和豎向位移，邊坡左端各網(wǎng)格節(jié)點(diǎn)控制水平位移[16]。上述實(shí)驗(yàn)已獲得非飽和滲流計(jì)算和彈塑性有限元計(jì)算的所需的土層參數(shù)，其中Q2，Q3黃土的彈性模量根據(jù)三軸試驗(yàn)的應(yīng)力應(yīng)變關(guān)系曲線，采用鄧肯-張模型計(jì)算所得。此外由于Q3黃土區(qū)多為非飽和區(qū)，其彈性模量取值較飽和彈性模量大，根據(jù)滲流模擬結(jié)果計(jì)算得到該區(qū)的平均基質(zhì)吸力，采用Oh & Vanapalli提出的非飽和區(qū)土彈性模量經(jīng)驗(yàn)公式計(jì)算得到[21]，見(jiàn)式(10)。由于粉質(zhì)粘土和泥巖的性質(zhì)參數(shù)對(duì)邊坡應(yīng)力場(chǎng)的計(jì)算結(jié)果影響不大，故其彈性模量和泊松比取經(jīng)驗(yàn)值。具體參數(shù)見(jiàn)表2。

(10)

式中：E為彈性模量；Es為飽和彈性模量；(ua-uw)為基質(zhì)吸力；sr為飽和度；a和k為擬合參數(shù)；對(duì)于IP為9.1的粉土其取值可分別為0.1，2。

圖7為非飽和滲流分析計(jì)算結(jié)果，從中可知隨著灌溉的進(jìn)行，地下水位不斷上升初始水位線位于粉質(zhì)粘土層，經(jīng)過(guò)2年的灌溉和停灌，水位線由初始的 26 m 上升到42 m，第6年水位達(dá)到Q2黃土上部，約50 m的位置，第12年的水位上升到52 m，達(dá)到了Q3黃土，穿越了滑帶，在坡腳以泉水的形式溢出。圖中 AB 線為從坡頂?shù)椒圪|(zhì)粘土底層的直線，自上而下每隔2 m設(shè)定一個(gè)點(diǎn)，共25個(gè)，讀取每個(gè)點(diǎn)上的孔隙水壓力，繪制成曲線如圖8所示，從圖8可知，隨著水位線的上升，坡體飽和區(qū)土體的孔隙水壓力壓不斷增大，而非飽和土土體的孔隙水壓力在停灌階段基本保持不變；在灌溉階段，非飽和區(qū)上部負(fù)孔隙水壓力會(huì)減小?？傮w上，隨著灌溉過(guò)程的進(jìn)行，飽和區(qū)不斷增大，非飽和區(qū)出頂部5 m范圍內(nèi)的負(fù)孔隙水壓不斷增大縮小的波動(dòng)外，其他部位的孔隙水壓力基本維持恒定。根據(jù)式(3)非飽和土強(qiáng)度理論可知，土體的抗剪強(qiáng)度隨基質(zhì)吸力的減小而減小，當(dāng)土體飽和時(shí)，若應(yīng)力狀態(tài)不變其抗剪強(qiáng)度最小。根據(jù)滲流分析結(jié)果和非飽和土抗剪強(qiáng)度理論可初步得到邊坡的失穩(wěn)主要同地下水位的上升有關(guān)，而這一現(xiàn)象將在后面穩(wěn)定系數(shù)分析的結(jié)果中具體體現(xiàn)出來(lái)。

彈塑性有限元分析可得到邊坡的應(yīng)力場(chǎng)，通過(guò)提取不同灌溉年份時(shí)滑面上不同節(jié)點(diǎn)的應(yīng)力狀態(tài)(σx,σy,τxy)，通過(guò)式(4)、(5)和式(2)或式(3)計(jì)算得到其剪應(yīng)力和抗剪強(qiáng)度，以潛在滑面為x軸，以剪應(yīng)力和抗剪強(qiáng)度為y軸，分析不同階段剪應(yīng)力和抗剪強(qiáng)度的變化特征。其中非飽和區(qū)的抗剪強(qiáng)度還需提取滑帶上土體的基質(zhì)吸力(ua-uw)。圖9為不同階段滑帶上各節(jié)點(diǎn)的應(yīng)力分布曲線，分別為灌溉前(a)、第3年10月停灌(b)、第8年5月灌溉(c)和第11年5月灌溉后剪應(yīng)力和抗剪強(qiáng)度沿滑面上的分布。由圖9a可以看出，初始狀態(tài)下整個(gè)滑面上土體的抗剪強(qiáng)度均大于剪應(yīng)力，邊坡處于穩(wěn)定狀態(tài)。隨著灌溉時(shí)間的增長(zhǎng)以及土體中孔隙水壓力的增大，土體的抗剪強(qiáng)度和剪應(yīng)力都在不斷地發(fā)生變化，抗剪強(qiáng)度不斷地降低，由圖9b可知灌溉第3年時(shí)，小范圍的坡腳部位剪應(yīng)力達(dá)到其抗剪強(qiáng)度，土體屈服.由圖9c可以看出，當(dāng)持續(xù)灌溉到第8年時(shí)，盡管邊坡處于穩(wěn)定狀態(tài)，其安全系數(shù)為1.11(見(jiàn)圖10)，但坡腳處的土體仍然發(fā)生屈服，且屈服范圍較第三年灌溉時(shí)大；此外滑坡后壁處也有部分土體屈服，這同應(yīng)力重分布相關(guān)；當(dāng)灌溉到第12年時(shí)，屈服范圍又?jǐn)U大，幾乎整個(gè)滑面上土體的剪切應(yīng)力到達(dá)了其抗剪強(qiáng)度，坡體失穩(wěn)。從上述分析中不難看出，坡腳處的土體先屈服，坡腳失穩(wěn)，該滑坡為牽引式滑坡。

表2 模型中所需各土層參數(shù)

圖7 地下水位隨灌溉的變化圖

圖8 參照線剖面上的孔隙水壓力輪廓線

圖9 不同階段滑面節(jié)點(diǎn)抗剪強(qiáng)度和剪應(yīng)力曲線

圖10 穩(wěn)定系數(shù)隨灌溉時(shí)間變化曲線

將有限元滲流計(jì)算的滲流場(chǎng)耦合用于邊坡穩(wěn)定性計(jì)算，利用式(6)計(jì)算穩(wěn)定系數(shù)，最終得出穩(wěn)定系數(shù)隨灌溉時(shí)間的關(guān)系曲線如圖10所示。計(jì)算結(jié)果表明，隨著地下水位的上升，滑坡的穩(wěn)定系數(shù)降低，1970年11月灌溉后穩(wěn)定系數(shù)突降，由初始1.34降到1.22，分析其原因，此時(shí)的地下水位由初始水位上升到了與粉質(zhì)粘土層接觸的水平滑帶位置，粉質(zhì)粘土層為相對(duì)隔水層，滑帶土孔隙水壓力開(kāi)始上升，有效應(yīng)力逐漸減小，土的抗剪強(qiáng)度降低。自1971年到1977年，穩(wěn)定系數(shù)逐年降低，直到1980年5月灌溉之后，穩(wěn)定系數(shù)達(dá)到1.01，該滑坡處于臨界狀態(tài)。計(jì)算結(jié)果同黑方臺(tái)地區(qū)1980年開(kāi)始出現(xiàn)大規(guī)模滑坡的實(shí)際狀況相一致。

將此計(jì)算結(jié)果同常用的極限平衡法計(jì)算結(jié)果進(jìn)行對(duì)比，見(jiàn)圖11。其中極限平衡法所用軟件為Geo-studio中的Slope模塊，計(jì)算方法為M-P法，計(jì)算不需邊界條件，所需材料參數(shù)為SWCC、容重、有效強(qiáng)度參數(shù)c’和φ’，具體見(jiàn)表2。通過(guò)對(duì)比分析可知，極限平衡法計(jì)算結(jié)果和本文有限元法的計(jì)算結(jié)果不盡相同，但差異較小，其相關(guān)系數(shù)為0.99。

圖11 穩(wěn)定系數(shù)計(jì)算結(jié)果對(duì)比圖

4 結(jié)論

通過(guò)野外剖面測(cè)量，室內(nèi)試驗(yàn)，和數(shù)值模擬等手段，結(jié)合非飽和土力學(xué)理論研究了甘肅黑方臺(tái)一滑坡泥流的失穩(wěn)機(jī)制，并得到以下結(jié)論：

(1) 該滑坡為純黃土滑坡，坡體前緣有廣闊的臨空面。黃土的飽和含水率大于液限，擾動(dòng)時(shí)易發(fā)生液化；黃土層下的粉質(zhì)粘土層飽和滲透系數(shù)小，儲(chǔ)水能力強(qiáng)，有利于水分積聚，形成潛水位。

(2) 該地區(qū)自1967年灌溉以來(lái)，地下水位不斷上升，非飽和區(qū)的基質(zhì)吸力變化不大。結(jié)合非飽和土強(qiáng)度理論，邊坡穩(wěn)定性的降低主要同地下水位上升相關(guān)，其穩(wěn)定系數(shù)的快速下降始于地下水位開(kāi)始接觸潛在滑面。

(3) 剪應(yīng)力和抗剪強(qiáng)度的對(duì)比分析結(jié)果表明，滑面上初始邊坡的抗剪強(qiáng)度均高于剪應(yīng)力，隨著地下水位的上升，坡腳土體先屈服，且屈服范圍不斷擴(kuò)大；由于應(yīng)力重分布，滑坡后緣的土體也逐漸屈服。該滑坡為牽引式滑坡。

(4) 野外調(diào)查中若在黃土層中發(fā)現(xiàn)泉流等地下水排泄口，表明地下水位已上升到黃土層內(nèi)部，邊坡處于欠穩(wěn)定狀態(tài)，應(yīng)及時(shí)做好防范措施，如距離坡面一定范圍內(nèi)設(shè)置為危險(xiǎn)區(qū)、及時(shí)遷移相關(guān)農(nóng)戶等。

[1] 雷祥義. 陜西涇陽(yáng)南緣黃土滑坡災(zāi)害與引水灌溉的關(guān)系[J]. 工程地質(zhì)學(xué)報(bào), 1994, 3(1): 56- 64.

[2] 王家鼎. 高速黃土滑坡一種機(jī)理—飽和黃土蠕動(dòng)液化[J].地質(zhì)評(píng)論, 1992, 38(6): 532-539.

[3] 戴福初，陳守義，李焯芬. 從土的應(yīng)力應(yīng)變特性探討滑坡發(fā)生機(jī)理[J]. 巖土工程學(xué)報(bào), 2000, 22(1):127-130.

[4] 金艷麗, 戴福初. 灌溉誘發(fā)黃土滑坡機(jī)理研究[J]. 巖土工程學(xué)報(bào)， 2007，29(10): 1493-1499.

[5] Fredlund D G and Rahardjo H. Soil Mechanics for Unsaturated Soils[M].Wiley, New York, NY, USA,1993.

[6] Brand E W. Some thoughts on rainfall-induced slope failures[C]// Proc of 10th Int Conf on Soil Mech and Found Eng, Brookfield, 1981:373-376.

[7] Lee L M, Gofar N and Rahardjo H. A simple model for preliminary evaluation of rainfall-induced slope instability[J]. Engineering Geology ,2009,108(3/4): 272-285.

[8] 簡(jiǎn)文星，許強(qiáng)，童龍?jiān)?，?三峽庫(kù)區(qū)黃土坡滑坡降雨入滲模型研究[J]. 巖土力學(xué), 2013, 34(1): 3527-3535.

[9] 金艷麗, 戴福初. 地下水位上升下黃土斜坡穩(wěn)定性分析[J]. 工程地質(zhì)學(xué)報(bào), 2007,15(5): 599-606.

[10] 朱立鋒, 胡煒, 賈俊，等.甘肅永靖黑方臺(tái)地區(qū)灌溉誘發(fā)型滑坡發(fā)育特征及力學(xué)機(jī)制[J].地質(zhì)通報(bào), 2013, 32(6): 840-846.

[11] 谷天峰, 朱立峰, 胡煒，等.灌溉引起地下水位上升對(duì)斜坡穩(wěn)定性的影響[J]. 現(xiàn)代地質(zhì), 2015, 29(2): 408-415.

[12] 董英, 賈俊, 張茂省，等.甘肅永靖黑方臺(tái)地區(qū)灌溉誘發(fā)作用與黃土滑坡響應(yīng)[J]. 地質(zhì)通報(bào), 2013, 32(6): 893-898.

[13] Vanapalli S K,D G Fredlund,D E Pufahl,et al. Model for the prediction of shear strength with respect to soil suction[J] .Canadian Geotechnical Journal ,1996,33(3): 379-392.

[14] Khalili N,Khabbaz M H. A unique relationship for the determination of the shear strength of unsaturated soils[J]. Gotechnique,1998,48(5):681-687.

[15] Oberg A, Sallfors G. Determination of shear strength parameters of unsaturated silts and sands based on the water retention curve[J]. Geotech Test J GTJODJ ,1997,20(1):40-48.

[16] GeoSlope International Ltd. Sigma/W user’s guide for stress-deformation analysis[R]. Calgary, Alta: GEO-SLOPE International Ltd; 2007.

[17] 張常亮, 王阿丹, 邢鮮麗，等.侵蝕作用誘發(fā)黃土滑坡的機(jī)制研究[J].巖土力學(xué),2012, 33(5): 1585-1592.

[18] 侯曉坤,李同錄,李萍. 開(kāi)挖黃土高邊坡的應(yīng)力路徑及變形破壞機(jī)制分析[J]. 巖土力學(xué) ,2014，35(S2):558-556

[19] Van Genuchten,M.T. A closed form equation for predicting the hydraulic conductivity of unsaturated soils[J]. Soil Science Society of America Journal, 1980,44(5):892-898.

[20] 王紅. 非飽和黃土滲透性函數(shù)的試驗(yàn)研究[D]. 西安:長(zhǎng)安大學(xué)，2015.

[21] Vanapalli Sai and Won Oh.A model for predicting the modulus of elasticity of unsaturated soils using the soil-water characteristic curve[J]. International Journal of Geotechnical Engineering, 2010, 4(4): 425-433.

Abstract:For a loess flow slide in Heifangtai Gansu Province, field investigation for the geological structure, laboratory tests including basic physical parameters tests, traixial tests, SWCC tests, and finite seep and elastic-plastic numerical model are conducted to study the failure mechanism.The changing characters of groundwater table and the related slope stability are analyzed. The yield characteristics associated by comparing the shear stress and shear strength along the slide surface are studied by extracting the stress state along the failure surface. Results show that the flow-slide is a loess landslide, whose failure surface lies in the loess layer. The saturated water content of loess is higher than its liquid limit, meaning the loess is easy to liquefy when it is disturbed; the silt clay below the loess layer has a low saturated permeability and high water storage capacity, making it easy to form a water table above it. Under the irrigation conditions, the ground water table increases slowly, the factor of safety of the slope fluctuates at first and begin to decrease greatly when the groundwater table researches the slide surface. Before irrigation, the shear strength is higher than the shear strength along the slide surface. With the irrigation continuing, the slope toe yields first, and the yield area expands upwards. The following stress redistribution leads to the yield of soil at the crest of slope. At the end the whole soil along the slide surface yields with the stress researches the shear strength. The landslide is a retrogressive landslide.

Key words:loess flow-slide; field investigation; laboratory tests; finite numerical model; retrogressive landslide

Analysis of a Flow-slide in Heifangtai Induced by Irrigation

ZHAO Jifei1, 2, HUANG Jiayue3, HOU Xiaokun1, LI Tonglu1and LIN Xiaoyan4

(1.SchoolofGeologicalEngineeringandSurveying,Chang’anUniversity,Xi’an710054,China; 2.NorthwestEngineeringCorporationLimited,Xi’an710054,China; 3.Xi’anEngineeringandDesignResearchofChina,Xi’an710054,China;4.SchoolofEarthScienceandEngineering,Xi’anShiyouUniversity,Xi’an710065,China)

趙紀(jì)飛，黃嘉悅，侯曉坤，等. 灌溉誘發(fā)的黑方臺(tái)黃土滑坡泥流機(jī)理分析[J]. 災(zāi)害學(xué)，2017，32(4)：60-66. [ZHAO Jifei, HUANG Jiayue, HOU Xiaokun，et al. Analysis of a Flow-slide in Heifangtai Induced by Irrigation[J]. Journal of Catastrophology，2017，32(4)：60-66.

10.3969/j.issn.1000-811X.2017.04.010.]

X43；P642.131

A

1000-811X(2017)04-0060-07

2017-02-07

2017-04-17

國(guó)家自然科學(xué)基金項(xiàng)目(40772181, 41372329)；國(guó)家重點(diǎn)基礎(chǔ)研究發(fā)展計(jì)劃(973計(jì)劃)資助項(xiàng)目(2014CB744701)

趙紀(jì)飛(1987-)，男，山東日照人，博士研究生，工程師，主要從事邊坡工程方面的研究. E-mail：523327087@qq.com

10.3969/j.issn.1000-811X.2017.04.010