基于分形理論的非飽和土多相流相對滲透系數(shù)模型

曹 帥，寧金成，楊慶國

(1.河南交通職業(yè)技術(shù)學(xué)院 公路學(xué)院，鄭州 450015； 2.重慶交通大學(xué) 土木建筑學(xué)院，重慶 400074)

1 研究背景

非飽和土的滲透破壞與暴雨誘發(fā)滑坡、高庫大壩水力破壞、洪澇管涌破壞等工程問題[1-2]密切相關(guān)，而相對滲透系數(shù)則是表征非飽和土滲流規(guī)律最為關(guān)鍵的參數(shù)之一。因此，合理確定非飽和土的相對滲透系數(shù)對于準確把握非飽和土的滲流力學(xué)行為，進而為工程設(shè)計及災(zāi)害防治提供依據(jù)，具有重要的理論價值與工程意義。

眾所周知，非飽和土的相對滲透系數(shù)受土體細觀孔隙分布特征(如孔隙的孔徑大小、不同孔徑孔隙的體積分布、孔隙形狀等)及飽和度影響極大[3-5]。然而，由于非飽和土在較高飽和度時，相對滲透系數(shù)對飽和度的變化非常敏感，飽和度的微小變化會使相對滲透系數(shù)發(fā)生相當大的改變;而在較低飽和度時，非飽和土的相對滲透系數(shù)很小，造成該系數(shù)測量難度大、費時長，使得直接通過物理測定獲得該參數(shù)的準確值變得極為困難。因此，提出各種預(yù)測模型來預(yù)測非飽和土的相對滲透系數(shù)成為較為實際的方法。不少學(xué)者對該問題進行了深入研究[6-15]，如Burdine[6]和Mualem[7]通過提出理想的土壤孔隙模型建立了非飽和土液相滲流模型，該模型依賴于土水特征曲線試驗數(shù)據(jù)來計算液相的相對滲透系數(shù)。Van Genuchten[8]通過提出經(jīng)驗的土水特征曲線模型,并分別與Mualem[6]和Burdine[7]的非飽和多孔介質(zhì)廣義模型相結(jié)合，得到了三參數(shù)非飽和土的相對滲透系數(shù)模型，該模型一定程度上反映了土壤持水與孔隙分布規(guī)律的關(guān)系，但模型參數(shù)為經(jīng)驗參數(shù)，無具體物理含義。Yang和Mohanty[9]結(jié)合Kosugi[10]的土水特征曲線函數(shù)，推導(dǎo)了多孔介質(zhì)非濕潤相相對滲透率的一般表達式，并使用該廣義公式研究了Burdine[6]、Mualem[7]、Alexander等[11]、Skaggs[12]模型對非濕潤相滲透率的預(yù)測精度。

土壤孔隙形狀的不規(guī)則性和連通性使得準確估算非飽和土滲透系數(shù)仍是一個具有挑戰(zhàn)性的課題，尤其在土壤多相流運移模型方面。為研究土壤多相流運移模型，本文將土壤孔隙簡化為一簇服從分形定律的迂曲毛細管，通過建立液相與氣相在多孔介質(zhì)中的輸運特征模型，結(jié)合Young-Laplace方程與Buckingham-Darcy定律得到了非飽和土的液相與氣相的相對滲透系數(shù)模型，該模型僅包含2個參數(shù)，模型參數(shù)物理意義明確且均可通過試驗測得。最后，通過8組試驗數(shù)據(jù)驗證了本文模型的合理性與適用性。

2 Young-Laplace方程

毛細管模型為理解空氣-水-土交界面上的液體壓力與蒸汽壓之間的關(guān)系提供了物理解釋。一般情況下，可用圓柱面、橢球面和球面等簡單的幾何形狀來近似表示空氣-水-土交界面。Young-Laplace方程[16]用雙曲率模型(橢球面如圖1所示)建立了基質(zhì)吸力(ua-uw)與交界面幾何形狀關(guān)系的公式，其表達式為

圖1 橢球體狀氣-液-固交界面Fig.1 Ellipsoid gas-liquid-solid interface

(1)

式中：ua與uw分別為氣相和液相的壓力；Ts為水表面張力；θ為接觸角；r1、r2分別為橢圓的半長軸和半短軸。

為方便計算且不失一般性，本文對毛細管模型做出如下假設(shè)：

(1)土體孔隙可簡化為橫截面為圓形的迂曲毛細管。

(2)對于不同孔徑的毛細管，接觸角θ均相等且為常數(shù)。

(3)不考慮溫度變化對水表面張力的影響，即水表面張力Ts為常數(shù)。

由r1=r2=r，ua-uw=rwgh，代入式(1)可得

(2)

式中：h為基質(zhì)水頭；ρw為液相密度；g為重力加速度。

3 分形尺度分析及模型建立

3.1 分形尺度分析

眾所周知，非飽和土土體中分布著大量復(fù)雜的孔隙，其分布及自身尺寸則直接影響了土體的水力特性。為描述該特性，本文采用一簇具有不同半徑的毛細管來代替這些孔隙，如圖2所示。

圖2 分形多孔介質(zhì)孔隙簡化模型Fig.2 Simplified pore modelof fractal porous media

在土中取一個具有代表性的圓柱體單元(半徑為R0，長度為L)，假設(shè)土體孔隙可用一簇具有不同半徑的毛細管代替，其半徑r在一定范圍內(nèi)，即rminr的毛細管的累積數(shù)量近似服從以下分形定律[17-18]，即

(3)

式中：ε為分形標度，此處指土壤孔隙半徑；Df為分形維數(shù)。

Yu[18]提出分形維數(shù)Df可由式(4)計算。

(4)

式中：dE為歐式空間維度數(shù)，在二維歐式空間里dE=2，在三維歐式空間里dE=3；φ為孔隙率。在二維歐式空間里，多孔介質(zhì)的孔隙分形維數(shù)一般有1

將式(4)對r微分，可得到毛細管半徑從r到r+dr的毛細管累積面積，即

(5)

由式(5)可知代表毛細管的累積數(shù)量隨著孔徑的增大而減少。

假設(shè)每根毛細管中的水都處于可流動狀態(tài)，根據(jù)熱力學(xué)局部平衡理論，在具有一定含水量的土體中，當土中水運動時(一般分為干燥過程和濕潤過程)，存在一個臨界半徑R(R為毛細管恰好飽和狀態(tài)時的臨界孔隙半徑)，使得R為控制著毛細管飽和狀態(tài)和非飽和狀態(tài)的臨界半徑，即毛細管半徑r=R時恰好為飽和狀態(tài)，其水頭高度恰好為L，則當毛細管半徑為rmin≤r≤R時，毛細管為飽和狀態(tài)；當R

3.2 非飽和土多相流相對滲透系數(shù)模型建立

非飽和多孔介質(zhì)中的多相流問題可以分為多個單相流問題。在具有代表性的圓柱體單元中，某一特定含水量下的非飽和多孔介質(zhì)的孔隙氣壓力ua和孔隙水壓力uw為常數(shù)。因此，根據(jù)Bear[19]提出的Navier-Stokes方程特解，單根孔隙的氣相和液相體積流率可分別表示為：

(6)

(7)

式中：μa和μw分別為氣相和液相的流體動力黏度；ua和uw分別為非飽和土多孔介質(zhì)的孔隙氣壓力和孔隙水壓力；ra和rw分別為液相和氣相密度；z為垂直坐標；g為重力加速度。

結(jié)合式(5)和式(7)可得具有代表性的圓柱體單元橫截面的液相總體積流率Qw，即

(8)

結(jié)合式(2)和式(8)，液相的總體積流率可用基質(zhì)水頭表示為

(9)

由Buckingham-Darcy定律[20]可知，非飽和多孔介質(zhì)的液相總流率還可表示為

(10)

式中kw為非飽和多孔介質(zhì)液相的滲透率。結(jié)合式(9)與式(10)可得

(11)

同理可得氣相(R～rmax)的滲透率ka可表示為

(12)

由式(11)和式(12)可知，當h=hmax時，土壤為干燥的，此時土壤孔隙被氣相填滿；當h=hmin時，土壤孔隙被水充滿，土壤是飽和的。因此，飽和土壤的滲透率ksat可表示為

(13)

結(jié)合式(11)—式(13)可得非飽和土的液相和氣相的相對滲透系數(shù)，即：

(14)

(15)

式中：kr,w和kr,a分別為非飽和土液相和氣相的相對滲透系數(shù)；hmin和hmax分別為最小和最大基質(zhì)水頭。土壤的飽和度Se可由飽和土定義求得，即

式中:Vw為土中孔隙水的體積；Vv為孔隙總體積。

因此，聯(lián)立式(2)與式(16)可得飽和度與基質(zhì)水頭關(guān)系的表達式，即

結(jié)合式(14)、式(15)和式(17)可得非飽和多孔介質(zhì)液相與氣相的相對滲透系數(shù)與飽和度關(guān)系的表達式，即：

式中ξ=hmax/hmin。

式(18)和式(19)即為本文提出的非飽和土多相流(氣相與液相)相對滲透系數(shù)模型，該模型僅包含Df和ξ2個參數(shù)，且Df和ξ可分別通過計盒法和軸平移技術(shù)測得。

4 模型驗證與討論

本部分采用不同土類的相對滲透系數(shù)試驗數(shù)據(jù)驗證所提出模型的可行性與合理性，這些土包括Sable de Riviere土、Gllat 砂質(zhì)壤土和Rubicon 砂質(zhì)壤土[21]；Grimsel Test Site土[22]；Mixed 砂和Oakley 砂[23]。采用MatLab全局搜索功能擬合試驗數(shù)據(jù)并將得到的模型參數(shù)列于表1。

表1 不同土類的相對滲透系數(shù)模型參數(shù)Table 1 Parameters of relative permeability coefficientmodel for different soil types

圖3為Sable de Riviere土、Gllat 砂質(zhì)壤土和Rubicon 砂質(zhì)壤土液相的相對滲透系數(shù)試驗數(shù)據(jù)與本文模型模擬曲線的對比。由圖3可知，本文模型模擬曲線與試驗數(shù)據(jù)吻合良好。其中，對于Sable de Riviere土的液相相對滲透系數(shù)試驗數(shù)據(jù)，當土體飽和時土壤的相對滲透系數(shù)為1，即土壤滲透率為飽和滲透率。當土壤的飽和度<1時，土壤的相對滲透率也隨之不斷減小，直至飽和度為0時，相對滲透系數(shù)也變?yōu)?，此時土壤氣相的滲透率達到最大值。Gllat 砂質(zhì)壤土和Rubicon 砂質(zhì)壤土的相對滲透系數(shù)試驗數(shù)據(jù)也表現(xiàn)出了類似的特征。

圖3 不同非飽和土液相的相對滲透率試驗數(shù)據(jù)與本文模型模擬曲線對比Fig.3 Comparison between test data of liquid-phaserelative permeability coefficient of different unsaturatedsoils and the simulated curve by the present model

圖4為Grimsel Test Site土、Mixed砂和Oakley砂的氣相相對滲透系數(shù)與本文模型模擬曲線對比。

圖4 不同非飽和土氣相的相對滲透率試驗數(shù)據(jù)與本文模型模擬曲線對比Fig.4 Comparison between test data of gas-phaserelative permeability of different unsaturated soils andthe simulated curve by the present model

由圖4可知本文模型可較好地反映非飽和多孔介質(zhì)的氣相滲流規(guī)律。對于Grimsel Test Site土氣相的相對試驗數(shù)據(jù)，當土壤完全干燥時，土壤的氣相滲透率最大且其值隨著飽和度的增大而不斷減小。值得注意的是，當土壤的氣相傳導(dǎo)能力接近0時，土壤的飽和度并未到達1，這主要是因為當土壤含水量達到一定值以后會將土壤孔隙堵住,從而限制了氣相的輸運。Mixed砂和Oakley砂的氣相相對滲透系數(shù)試驗數(shù)據(jù)也表現(xiàn)出了與Grimsel Test Site土類似的滲流特征。

孔隙分布規(guī)律對非飽和土的滲流規(guī)律有著不可忽視的影響。圖5表明了分形維數(shù)Df對非飽和土多相滲流的影響。由圖5可知，在同一飽和度下分形維數(shù)越大，土壤液相的相對滲透系數(shù)越小，氣相的相對滲透系數(shù)越大。這是因為分形維數(shù)越大表明土壤的較小孔隙占比越大，由于孔隙水表面張力的存在，土壤的導(dǎo)水能力越差，而氣相的輸運能力相對就越強。

圖5 孔隙尺寸分布對非飽和土相對滲透系數(shù)的影響Fig.5 Influence of pore size distribution on relativepermeability coefficient of unsaturated soil

由上述分析可知，本文提出的多相相對滲透系數(shù)模型可較好地模擬非飽和土的多相滲流規(guī)律。由于本文在重構(gòu)土壤孔隙系統(tǒng)時是將土壤孔隙簡化為一簇光滑的毛細管，因此忽略了孔隙的粗糙度、連通性及不規(guī)則性對土壤滲流規(guī)律的影響，造成本文模型在預(yù)測土壤的相對滲透系數(shù)時仍存在一些誤差，這些因素對土壤的相對滲透系數(shù)的影響也是筆者進一步研究的目標。

5 結(jié) 論

本文在可用毛細管系統(tǒng)近似代替土孔隙的假設(shè)上，基于Young-Laplace方程和分形理論，探討了非飽和土的多相滲流特征，并得到了以下結(jié)論：

(1)將土壤孔隙簡化為等截面迂曲毛細管，假設(shè)土壤孔隙尺寸分布服從分形定律并重構(gòu)了土壤孔隙系統(tǒng)模型。

(2)結(jié)合Young-Laplace方程與Buckingham- Darcy定律得到了非飽和土的液相與氣相的相對滲透系數(shù)模型，該模型僅包含2個參數(shù)且模型參數(shù)物理意義明確。采用6組不同土類的試驗數(shù)據(jù)對本文提出的模型進行了驗證，結(jié)果顯示本文提出的多相相對滲透系數(shù)模型與試驗數(shù)據(jù)吻合良好，表明了本文模型可較好地模擬非飽和土的多相相對流滲透系數(shù)。

(3)討論了分形維數(shù)Df對非飽和土液相與氣相的滲流特性的影響。由分析可知同一飽和度下，分形維數(shù)Df越大，土壤液相的相對滲透系數(shù)越小，氣相的相對滲透系數(shù)越大。