土體非飽和滲透特性對降雨入滲的影響

陳旺旺，劉 剛，劉 暢，鐘佳思，童富果

(三峽大學(xué) 水利與環(huán)境學(xué)院， 湖北 宜昌 443002)

連續(xù)性降雨、強降雨等是觸發(fā)泥石流、滑坡等地質(zhì)災(zāi)害的主要誘因[1-4]。非飽和土中雨水入滲，產(chǎn)生水氣二相流的運動，繼而改變邊坡土壤的非飽和特性，引發(fā)邊坡失穩(wěn)[5-6]。土體的非飽和滲流特性包含了土水特征關(guān)系、水-氣滲透關(guān)系，是土體微觀結(jié)構(gòu)在宏觀上的綜合反映。因此，進行不同非飽和滲流特性條件下的入滲過程進行數(shù)值模擬，對于揭示多個物理特性變量下土體的入滲規(guī)律具有重要意義。

非飽和滲流及其對邊坡穩(wěn)定性的影響有不少學(xué)者做了相關(guān)研究。朱偉等[7-8]通過設(shè)計一維降雨入滲試驗，研究了非飽和滲流特性對降雨入滲水量的影響。豐光亮等[9]在室內(nèi)模擬了人工降雨入滲非飽和土柱試驗，提出鄂西恩施地區(qū)入滲影響區(qū)入滲前鋒運移規(guī)律和降雨強度、降雨歷時對入滲率的影響規(guī)律。Neuman[10]將有限元方法應(yīng)用到求解飽和-非飽和滲流問題。陳善雄等[11]和陳守義[12]用積分有限差分方法模擬了降雨條件下土體中水分的運動情況，并對降雨條件下非飽和土坡等的穩(wěn)定性的分析方法做了研究。上述的試驗方法和數(shù)值模擬多建立在液相單相流理論上[13]。Lam等[14]假定土壤包含飽和和非飽和區(qū)域，在不飽和區(qū)域氣體為連續(xù)體且與大氣相通，利用二維有限元模型計算非飽和土體的應(yīng)力狀態(tài)。孫冬梅等[15]運用積分有限元差分方法模擬了土質(zhì)邊坡的降雨入滲過程，定量分析了空氣阻力對水流入滲的影響。降雨入滲過程是一個涉及水氣二相流耦合的非飽和滲流過程，以往涉及土體非飽和滲透特性的研究多集中于外界因素對降雨入滲規(guī)律的探討，而與非飽和滲透特性直接相關(guān)的土水特征關(guān)系、水-氣滲透關(guān)系對降雨入滲的影響也非常值得研究。

因此，本文基于非飽和土-水-氣二相流理論，采用有限單元法對不同非飽和滲透特性條件下的入滲過程進行計算，探究了土水特征關(guān)系和水相對滲透曲線對穩(wěn)定入滲強度的影響。對了解降雨誘發(fā)邊坡失穩(wěn)的機制，提高滑坡的預(yù)報和防治具有重要意義，也對城市防洪、農(nóng)業(yè)節(jié)水灌溉、污染物質(zhì)傳輸?shù)戎T多工程問題提供參考。

1 計算方程與數(shù)值模型

1.1 水氣兩相流控制方程

土體視為液、固、氣三相混合的連續(xù)介質(zhì)，各相物質(zhì)在土中的運動可描述為關(guān)于時間和空間的偏微分方程組。土壤內(nèi)部的液相流動主要受水壓力梯度的驅(qū)動，水分通過土壤中的水流通道入滲，土體的斷面面積大于實際水流的過水面積，因此Darcy定律計算的水流速度小于土中水流的真實流速。根據(jù)液相質(zhì)量守恒方程得到孔隙水的非飽和流動方程[16]：

(1)

式中：ρw為水的密度，kg/m3；μw為液體黏滯系數(shù)，Pa·s；Krw為水相對滲透系數(shù)；K為本征滲透系數(shù)，m2；Sr為水飽和度；φ為土壤孔隙率；Pw為孔隙水壓力；Qw為內(nèi)源項。

氣相在土壤中的流動受到重力和氣體壓力梯度的驅(qū)動，同時土中氣體僅通過氣流通道，土體的斷面面積大于實際氣流通過的面積，故由Darcy定律得來的氣流速度小于土中氣流的真實流速。根據(jù)氣相質(zhì)量守恒方程得到孔隙氣的非飽和流動方程[16]：

(2)

式中：ρg為氣的密度，kg/m3；μg為液體黏滯系數(shù)，Pa·s；Krg為氣相對滲透系數(shù)；Pg為孔隙水壓力。

在求解水、氣二相流控制方程時，孔隙水壓力和孔隙氣壓力之間恒有Pc≡Pg-Pw，式中Pc為基質(zhì)吸力。方程(1)、方程(2)構(gòu)成時間和空間的非線性偏微分方程，求解時對方程組的空間離散化采用Garlerkin有限單元法，時間域離散采用一維差分方法，即假設(shè)在時間域的時間段內(nèi)變量隨時間的線性變化。

1.2 本構(gòu)模型關(guān)系與參數(shù)

方程(1)和(2)有Pw，Pg，Sr，Krw，Krg5個未知量，求解時還需要引入3個本構(gòu)模型關(guān)系，即土水特征關(guān)系、水相對滲透關(guān)系、氣相對滲透關(guān)系，反映了土體的含水率與基質(zhì)吸力、滲水能力、滲氣能力的關(guān)系。

土水特征關(guān)系選取適用于大多數(shù)土體的Van Genuchten模型[17]：

(3)

式中：Sr為水飽和度；Srw為殘余水飽和度；P0為進氣值；m為與土體孔徑分布有關(guān)的參數(shù)。

水相對滲透關(guān)系選取為Van Genten-Mualem模型[18]：

(4)

式中：k為與材料特性有關(guān)的參數(shù)。

氣相對滲透關(guān)系選取Brooks & Corey模型[19]：

(5)

式中：Srg為殘余氣飽和度。

方程涉及的其他計算參數(shù)為：ρw=1 t/m3，ρg=1.29 kg/m3，g=9.8 N/kg，φ=0.4，μw=1.0×10-3N·s/m2，μg=1.84×10-5N·s/m2，k=2.88×10-13m2。

1.3 幾何模型

對于土質(zhì)邊坡而言，自然降雨的入滲方向主要為垂向，水分入滲也主要集中在坡體淺層區(qū)域，故本文在數(shù)值模擬時用到的幾何模型為淺層均質(zhì)土體。幾何模型的大小由降雨入滲過程中水分入滲深度決定，而水分的入滲深度與土體特性、邊界條件等因素相關(guān)的同時，還與降雨入滲時間相關(guān)?，F(xiàn)選定幾何模型高3 m，頂部單元厚度為0.1 m，底部單元厚度為0.2 m，共計42個節(jié)點20個單元。為研究土體非飽和滲透特性對入滲強度的影響，假定入滲強度遠小于降雨強度，地表徑流對入滲的影響可忽略不計，頂部及底部均為透水透氣邊界，且其邊界氣壓力為大氣壓，側(cè)邊界為不透水不透氣邊界。

2 土體非飽和滲流特性與入滲強度的關(guān)系

2.1 土水特征關(guān)系對入滲強度的影響

土水特征關(guān)系表述了土壤基質(zhì)吸力和土壤飽和度的關(guān)系，本文采用的Van Genuchten模型中，參數(shù)m和參數(shù)P0值的變化會直接影響關(guān)系曲線的結(jié)果形式，本節(jié)將對這兩個參數(shù)分別研究。為便于分析比較，采用相對入滲強度之比來衡量滲流特性變化對入滲強度的影響強弱。入滲強度之比是一個無量綱變量，定義為不同滲透特性下的相對入滲強度與某一滲透特性下的相對入滲強度的比值。

(1) 參數(shù)m變化對入滲強度的影響。參數(shù)m是反映土體孔隙結(jié)構(gòu)特征的重要參數(shù)之一，與孔徑分布指數(shù)相關(guān)。通過模擬不同參數(shù)m條件下的入滲過程，探究孔徑分布對土壤入滲強度的影響。根據(jù)相關(guān)資料可知參數(shù)m的常見取值范圍[20]為0.3～0.5，計算過程中，參數(shù)m取0.30、0.35、0.40、0.45、0.50一共5組(見圖1)，飽和度取值從0.15～1.00共18組，參數(shù)P0=0.278，k=0.9，相互組合共計模擬90組不同入滲過程?；|(zhì)吸力隨飽和度的變化率因參數(shù)m的改變而改變，參數(shù)m越小，基質(zhì)吸力在飽和度越低時的影響越大，而當飽和度較高時，參數(shù)m對基質(zhì)吸力的影響較小。

圖1不同參數(shù)m條件下土水特征曲線

如圖2所示，參數(shù)m的改變只影響飽和度較低時的基質(zhì)吸力，隨著飽和度的增大，參數(shù)m的改變對基質(zhì)吸力影響不再明顯，故參數(shù)m主要影響飽和度較低時的土體穩(wěn)定入滲強度，當飽和度較高趨近于1時，土壤內(nèi)部的入滲通道是水分入滲的主要通道，幾乎不受基質(zhì)吸力的影響，故土體穩(wěn)定入滲強度幾乎不受參數(shù)m影響，且其極小值的大小和位置不隨參數(shù)m的改變而改變。

圖2相對入滲強度與飽和度的關(guān)系曲線

(2) 參數(shù)P0變化對入滲強度的影響。參數(shù)P0為進氣值，通過模擬不同參數(shù)P0條件下的入滲過程，反映土體的最大孔隙尺寸對入滲強度的影響。參數(shù)P0的取值范圍為0.2～2.0，計算過程中，參數(shù)P0為0.2、0.6、1.0、1.4、1.8一共5組(見圖3)，飽和度取值從0.15～1.00共18組，參數(shù)m=0.359，k=0.9，相互組合共計模擬90組不同入滲過程?；|(zhì)吸力隨飽和度的變化率因參數(shù)P0的改變而改變，參數(shù)P0越大，基質(zhì)吸力在飽和度越低時的影響越大，隨著飽和度的增大，水相對滲透系數(shù)對入滲起主導(dǎo)作用，參數(shù)P0對基質(zhì)吸力的影響逐漸減小直至為0。

圖3不同參數(shù)P0條件下土水特征曲線

當飽和度增大時，入滲強度先減小后增大，相同飽和度的條件下，入滲強度隨參數(shù)P0的增大而增大(見圖4)。因為P0增大使得飽和度較高時的基質(zhì)吸力增大，入滲強度受基質(zhì)吸力的影響明顯，故入滲強度極小值點向飽和度較高處偏移。當飽和度趨近于1時，水分主要是通過土壤孔隙通道流動，水相對滲透系數(shù)是主要影響因素，基質(zhì)吸力此時的作用并不明顯。

圖4相對入滲強度與飽和度的關(guān)系曲線

2.2 水相對滲透曲線對入滲強度的影響

水相對滲透曲線反映了水相對滲透系數(shù)與飽和度的關(guān)系，是影響入滲強度的關(guān)鍵因素。本文采用Van Genten-Mualem模型，模型中參數(shù)k反映了水相對滲透系數(shù)隨飽和度的變化規(guī)律。本文參數(shù)k的取值范圍為0.5～0.9，計算過程中，參數(shù)k為0.5、0.6、0.7、0.8、0.9一共5組(見圖5)，飽和度取值從0.15～1.00共18組，參數(shù)P0=1.33，m=0.4，相互組合共計模擬90組不同入滲過程。水相對滲透系數(shù)隨著飽和度的增大而增大，土壤殘余含水和完全飽和時的水相對滲透系數(shù)分別為0和1，相同飽和度下水相對滲透系數(shù)隨著參數(shù)k的減小而減小。

圖5不同參數(shù)k條件下的水相對滲透曲線

由圖6可以看出穩(wěn)定入滲強度是取決于基質(zhì)吸力和水相對滲透系數(shù)兩者結(jié)果，水相對滲透系數(shù)減小時入滲強度減小，則參數(shù)k越小穩(wěn)定入滲強度也越小。隨著飽和度的增大，水相滲透系數(shù)的變化會影響孔隙通道的過水能力，所以飽和度較高時入滲強度受參數(shù)k的變化影響明顯。土體穩(wěn)定入滲強度的極小值位于飽和度較高時，參數(shù)k的減小會導(dǎo)致水相對滲透系數(shù)的降低進而導(dǎo)致入滲強度的極小值減小。

圖6相對入滲強度與飽和度的關(guān)系曲線

3 結(jié) 論

(1) 在非飽和入滲過程中，基質(zhì)吸力的變化直接影響水分入滲過程中的水壓力梯度。當飽和度較低、參數(shù)m較小時，基質(zhì)吸力對土體的入滲強度影響比較明顯；隨著飽和度的增大并趨近于1時，土體入滲強度幾乎不受參數(shù)m的影響。

(2) 當飽和度較低時，基質(zhì)吸力驅(qū)動水分入滲，入滲強度隨著參數(shù)P0的增大而變化明顯；當飽和度較大時，水相對滲透系數(shù)較大，孔隙通道作為水分主要下滲通道，參數(shù)P0增大對入滲強度的影響微弱。

(3) 土體穩(wěn)定入滲強度的極值點不僅受基質(zhì)吸力和相對滲透系數(shù)的影響，還與土體本征滲透性關(guān)系密切，土體本征滲透性決定了土體穩(wěn)定入滲強度極小值的位置。隨著本征滲透性的增大，土體滲水能力增強，極小值的位置向低飽和度運動，此時基質(zhì)吸力的影響可忽略不計；當本征滲透性足夠小時，入滲受基質(zhì)吸力的影響，此時穩(wěn)定入滲強度隨飽和度的增大而減小。