基于分形理論的非達西滲流參數研究

馮雪冬，張毅，由晗楊，李文英，馬艷飛

(1. 山東理工大學 資源與環境工程學院, 山東 淄博 255049；2.山東理工大學 農業工程與食品科學學院, 山東 淄博 255049)

石油在開采、煉制、貯運及使用過程中往往伴隨著嚴重的環境問題，其中地下水污染是環境問題之一[1-2]。研究地下水污染，需要從石油污染物滲流方程的角度出發，達西定律作為多孔介質中飽和滲流的基本方程，描述了流體滲流速度與壓力梯度之間的線性關系，目前被廣泛應用于地下水、石油工程等領域[3]。然而，理論研究與室內試驗發現，當流體黏度或滲流速度增大時，滲流速度與壓力梯度之間呈非線性關系，達西定律不再成立。因此，對非達西滲流的研究成為當前滲流領域的研究熱點。

從現有研究成果來看，描述非達西滲流的方程有很多，Forchheimer方程因其具有明確的物理意義，被廣泛應用于非達西滲流的計算中[4]。但目前國內對Forchheimer方程中滲流參數的研究主要集中在數值模擬方面[5-6]，多孔介質的粒徑構成及孔隙結構等特征與非達西滲流參數的關系研究卻相對較少。同時，由于多孔介質的顆粒和孔隙具有不規則性、自相似性和非線性等特點[7]，使得多孔介質孔隙結構與非達西滲流參數之間的關系很難用傳統線性分析方法進行定量化描述。

分形理論作為定量描述多孔介質的新方法，其分形維數值既可以準確反映多孔介質的粒徑構成，又能真實描述多孔介質的孔隙結構[8]。因此，本文利用分形理論，探討了多孔介質分形特征與Forchheimer型非達西滲流參數之間的關系，以期為石油工程領域的非達西滲流研究提供一定的理論依據。

1 材料與方法

1.1 試驗材料

試驗所用多孔介質來源于山東省東營市廣饒縣和淄博市臨淄區的農田，按0～20 cm，20～40 cm和40～60 cm深度由上向下分層采集各土層的土壤樣品，每層重復取樣3次，其中廣饒縣0～20 cm，20～40 cm和40～60 cm的土壤編號依次為GR1、GR2和GR3，臨淄區0～20 cm，20～40 cm和40～60 cm的土壤編號依次為LZ1、LZ2和LZ3；所采多孔介質樣品帶回實驗室自然風干、除雜，過20目篩后，儲存在密封容器中；采用Mastersizer 2000激光粒度分析儀測定多孔介質樣品的顆粒組成，并按照國際制分級標準來劃分土壤質地類型。

試驗所用流體為0#柴油，購自中國石化齊魯石化股份有限公司，密度為0.82 g·cm-3，動力黏度為3.11 mPa·s，油水界面張力為37.30 mN·m-1，表面張力為23.90 mN·m-1。

1.2 滲流試驗

試驗采用自制的一維土柱滲流試驗裝置，如圖1所示。為保證試樣的均勻性，試驗時按標準化分層方法進行裝樣，每次2 cm厚，裝樣高度為20 cm；進油端采用蠕動泵進樣，并且在供油管上設有分流球閥，以控制進油端油壓；土柱出口端連接恒定液面高度的油箱，以保持土柱出流端液面高度恒定；土柱進出油端分別裝有測壓管，用于測定進出油端的油壓；柴油的滲流速度由量筒和秒表進行測定。

圖1 試驗裝置示意圖Fig.1 The schematic diagram of the experiment equipment

1.3 滲流理論方程

1856年法國工程師達西通過大量試驗研究，發現了滲流速度與壓力梯度之間的線性關系[9]，達西滲流方程可表示為：

(1)

式中：J為壓力梯度，MPa·m-1；ρ為流體密度，kg·L-1；g為重力加速度，m·s-2；K為滲透系數，m·s-1；v為滲流速度，m·s-1。

Forchheimer型非達西滲流方程與達西滲流方程相比，其主要區別在于增加了一個二次項，意義上則表示壓力梯度的損失不再完全由黏滯力所決定，而是由黏滯力和慣性力共同決定[6]，其具體表達形式為：

(2)

式中：β為非達西滲流因子，m-1；Kf為非達西滲透系數，m·s-1。

1.4 分形維數計算方法

本試驗用粒徑的重量分布對多孔介質分形維數進行計算[10-11]：

(3)

式中：M(r

對該式兩邊取對數，以lgR，lg (M(r

D=3-b

(4)

2 結果與分析

2.1 壓力梯度與滲流速度的關系

調節分流球閥，測定不同壓力梯度下柴油的滲流速度，將試驗結果繪制成圖，如圖2所示。由圖2可知，多孔介質中柴油的滲流曲線均呈類拋物線，表明滲流速度與壓力梯度不再滿足線性關系，很可能發生了非達西滲流。為進一步證實上述分析結果的正確性，采用Forchheimer方程對試驗結果進行擬合，擬合結果見表1。由表1可知，Forchheimer方程的擬合效果較好，相關系數均大于0.96，說明柴油在滲流過程中確實發生了非達西滲流。

表1 滲流速度與壓力梯度的曲線擬合方程及相關系數Tab.1 The fitting equations of seepage velocity and pressure gradient

2.2 多孔介質的分形特征

圖3為lgR與lg (M(r

與廣饒地區相比，臨淄地區多孔介質的分形維數值相對較大，但兩地區分形維數值的變化趨勢也存在相似之處：隨著土壤深度增加，分形維數值也相應增大，在40～60 cm土層都達到最大值。從表2中還可以發現，多孔介質的分形維數值與粉粒、粘粒含量的變化趨勢基本一致，即多孔介質中粉粒、粘粒含量越高，土壤顆粒分形維數值就越大。

2.3 分形維數與粒徑含量的關系

為進一步探討分形維數值與粒徑含量的關系，對粒徑含量與分形維數值進行回歸分析，擬合結果如圖4所示。由回歸分析結果可以看出，分形維數值與砂粒含量之間的相關性達到極顯著水平(P<0.01)，相關關系為負相關(R=-0.985)；分形維數值與粉粒含量之間呈極顯著正相關關系(P<0.01，R=0.970)，與粘粒含量也呈極顯著正相關關系(P<0.01，R=0.995)。分形維數值與多孔介質粒徑含量的相關性較高，這表明分形維數值能夠在一定程度上反映多孔介質的粒徑組成，即分形維數值的變化與粉粒、粘粒含量密切相關，分形維數值越大，粉粒、粘粒含量越高，質地越細，微小孔隙越多，結構越復雜，這與張世熔等[14]的研究結果一致。

表2 不同多孔介質的粒徑組成及分形維數Tab.2 Particle size composition and fractal dimension of different porous media

(a)分形維數與砂粒含量的關系 (b)分形維數與粉粒含量的關系 (c)分形維數與黏粒含量的關系圖4 分形維數與粒徑含量的關系Fig.4 Relationship between fractal dimension and particle size content of soils

2.4 分形維數與非達西滲流參數的關系

從Forchheimer方程表達式可以看出，非達西滲透系數和非達西滲流因子分別代表著滲流過程中粘滯力對柴油的影響和慣性力對柴油的影響。當多孔介質孔隙較小、流速較低時，黏滯力占主導地位，慣性力可忽略，這時壓力梯度與滲流速度滿足線性關系，服從達西定律；當多孔介質孔隙增大、流速加快時，慣性力的影響逐漸增大，壓力梯度與滲流速度的線性關系發生偏移，即發生非達西滲流。因此，Forchheimer方程滲流參數的變化與多孔介質性質(尤其是孔隙結構)密切相關。

本研究中分形維數值能反映多孔介質的粒徑構成，而粒徑構成又影響著多孔介質的孔隙結構，因此通過建立分形維數值與非達西滲流參數之間的回歸方程，探討它們之間的相關關系，可以為進一步研究Forchheimer型非達西滲流參數特征提供一定的研究思路和方法。以非達西滲流參數，分形維數值為橫、縱坐標，進行曲線擬合，回歸方程及相關系數如圖5、圖6所示。

由圖5可以看出，分形維數值與非達西滲透系數之間呈線性負相關關系，且二者的相關程度較高。從微觀角度來看，分形維數值越大，多孔介質中粉粒、粘粒含量越高，多孔介質越緊實，微小孔隙越多，從而使得多孔介質的滲透性變差。與非達西滲透系數的變化趨勢相反，非達西滲流因子隨分形維數值的增大而增大(圖6)，反之，當分形維數值變小，多孔介質顆粒變粗，壓力梯度與滲流速度之間的非線性關系越不明顯。

綜上所述，多孔介質的粒徑大小、孔隙結構等在很大程度上決定著非達西滲流參數的大小。分形維數值較直觀地反映了多孔介質的顆粒組成，通過建立分形維數值與非達西滲流參數之間的相關關系，從而可以定量地描述Forchheimer型非達西滲流的參數特征。

3 結論

1) 柴油在多孔介質中的滲流為非達西滲流，Forchheimer方程能較為準確地描述壓力梯度與滲流速度之間的關系。

2) 多孔介質具有統計意義上的自相似性，即分形性，其中分形維數值與砂粒含量呈極顯著負相關，

圖5 分形維數與非達西滲透系數的關系Fig.5 Relationship between fractal dimension and non-Darcy permeability coefficient

圖6 分形維數與非達西滲流因子的關系Fig.6 Relationship between fractal dimension and turbulence factor

與粉粒、粘粒含量呈極顯著正相關。

3) 分形維數值與非達西滲透系數之間呈線性負相關關系，與非達西滲流因子之間呈線性正相關關系，因此，分形維數可以作為描述非達西滲流參數變化趨勢的重要指標。