破碎花崗巖非達西滲流的試驗研究

賀香蘭，周佳慶，魏 凱，王 敏

(1.武漢大學水資源與水電工程科學國家重點實驗室，武漢 430072；2. 武漢大學水工巖石力學教育部重點實驗室，武漢 430072)

0 引 言

裂隙是巖體滲流的主要通道，控制著巖體的滲流特性，隨著高壩樞紐工程和抽水蓄能電站的興建，巖體可能承受高滲壓、高水力梯度作用，此時節(jié)理裂隙等結構面易張開擴展，甚至發(fā)生水力劈裂，巖體滲流往往表現(xiàn)出非達西特征，特別是對于裂隙充分發(fā)育的破碎巖體，滲流特性更為復雜。黃先伍等[1]、李順才等[2]在破碎巖石構成的多孔介質(zhì)滲流試驗的基礎上，研究了非達西滲流特性與孔隙率的關系；Thauvin和Mohanty[3]利用數(shù)值模擬的方法從細觀結構的角度研究了多孔介質(zhì)非達西滲流的成因；胡大偉等[4]研究了單軸壓縮破壞的大理巖峰后非線性滲流特征演化規(guī)律。

對于非達西流，傳統(tǒng)的達西滲流理論已不能準確表征巖體的滲透特性，因此滲流流態(tài)的判斷和分析理論的選擇是研究介質(zhì)滲流的關鍵。為此，國內(nèi)外學者對達西流與非達西流的臨界條件開展了大量研究，主要采用雷諾數(shù)Re和福希海默數(shù)F0作為臨界參數(shù)[5]。孟如真等[6]、李健等[7]分別通過現(xiàn)場壓水試驗和室內(nèi)滲流試驗，得到了量級為10～100的臨界雷諾數(shù)；Zimmerman等[8]通過對單裂隙砂巖進行滲流試驗和數(shù)值模擬，發(fā)現(xiàn)Re大于20時，滲流滿足Forchheimer非達西方程；Javadi等[9]開展了粗糙裂隙面滲流試驗，得到的臨界雷諾數(shù)介于0.001～25；同時還有文獻提出了高達2 300的臨界雷諾數(shù)[10]。另一方面，Zeng和Grigg等[11]根據(jù)福希海默數(shù)F0的定義及其物理意義，建議取F0=0.11作為達西流和非達西流的臨界值；隨后，Macini等[12]通過研究天然砂粒的滲流試驗，得到值為0.4左右的臨界福希海默數(shù)，該值對應28%的非達西效應；文獻[13]在進行數(shù)值模擬時采用了F0=0.01作為臨界值。通過眾多研究可以發(fā)現(xiàn)，雷諾數(shù)Re和福希海默數(shù)F0臨界值的取值范圍變化較大，且目前尚未得到臨界值的變化規(guī)律，因此采用參數(shù)Re和F0判斷滲流流態(tài)存在一定的主觀性和不確定性。

綜上所述，盡管滲流的非達西性得到了大量的研究，但涉及破碎巖體的非達西滲流的研究成果相對較少，并且達西流和非達西流的臨界問題仍未得到有效地解決?；诖耍疚耐ㄟ^開展破碎花崗巖在不同圍壓條件下的高水力梯度非達西滲流試驗，研究破碎巖體復雜裂隙介質(zhì)的滲流特性，并在試驗結果基礎上，提出一種簡單直觀方便的非達西流判定方法。該方法可將臨界值的不確定性問題轉化為實際工程目標需求的容許誤差，具有較好的實用價值。

1 試驗系統(tǒng)與步驟

1.1 試驗系統(tǒng)介紹

破碎巖體非達西滲流試驗在武漢大學與法國里爾科技大學(USTL)聯(lián)合研發(fā)的Triaxial Cell三軸試驗系統(tǒng)上進行，見圖1。該系統(tǒng)由圍壓、軸壓和孔隙水壓3套獨立加載部分構成，具有應力、應變、位移和流量4種控制加載方式，可施加最大圍壓60 MPa，最大軸壓375 MPa，最大滲透水壓力50 MPa。其中滲透水壓力由升級后的高精度液壓泵PHMP 50-1000型施加，該流量泵可容納100 mL體積水，流量范圍為0.000 41～66.8 mL/min。系統(tǒng)出口端接入測量精度為0.01 g的SA2202S-CW 型高精度自動電子天平，以測量出口水流流量。

1.2 試驗步驟

試驗巖樣取自廣東省陽江抽水蓄能電站高壓隧洞區(qū)，為中?；◢弾r，經(jīng)取芯機、剛性切割機和打磨機加工后制備得到Φ50 mm×100 mm的圓柱狀標準巖石試樣，見圖2。破碎巖體滲流試驗過程如下。

(1)將飽和的標準巖石試樣在預定圍壓條件下施加軸向壓力，直至試樣破壞，得到破碎巖樣。

(2)以位移加載模式控制巖樣保持在破碎初始狀態(tài)，并以0.1 MPa的滲透水壓持續(xù)進行12 h的過流，保證試樣充分飽和。

(3)在設定的圍壓條件下，逐級施加滲透水壓(見圖3)，測量并記錄每級水壓穩(wěn)定后的出口水流流量，完成該圍壓下的滲流試驗后，改變圍壓，重復上述過程，直至圍壓達到30 MPa。

圖1 非達西滲流試驗裝置Fig.1 Sketch of non-Darcy flow test equipment

圖2 試驗前后花崗巖試樣Fig.2 Granite samples before and after test

圖3 滲透水壓加載過程曲線Fig.3 Curves of the applied water pressure

2 試驗結果與分析

試驗采用陽江破碎花崗巖標準巖樣，每個試樣開展了不少于10級圍壓下的不同水力梯度的滲流試驗。試驗后的3個破碎巖樣見圖2(b)。由圖2(b)可見，試樣破碎程度較高，縱向裂隙發(fā)育，部分貫通試樣，同時存在少量橫向裂紋。本文給出了試樣1～3巖石試樣在不同圍壓條件下兩端水力梯度▽P與滲流流量Q之間的試驗曲線，見圖4。由圖4可見，▽P～Q關系曲線皆呈下凸狀，表明隨著水力梯度的增大，滲流流量的增長幅度逐漸減小。 換言之，隨著流量的增大，增加同等的流量所需消耗的水頭也越來越大。另外，對于同一試樣，圍壓增大可能引起巖樣壓縮，裂隙壓密，過水通道減小，導水性降低，因此對于相同的流量，圍壓越大，所對應的水力梯度也越大，故在▽P～Q關系曲線圖中表現(xiàn)為圍壓越大，曲線越偏離橫軸，偏離程度與圍壓大小的變化基本保持一致。但也存在例外，例如試樣3，在圍壓σ3達到25.0 MPa后，繼續(xù)增大σ3對▽P～Q曲線的偏離作用略顯微弱，因為此時裂隙的壓密程度接近最大，巖樣透水性基本不再變化。

圖4 試樣在不同圍壓下的▽P～Q關系曲線Fig.4 Variation of pressure gradient ▽P with flow rate Q under various confining stresses

圖4所示的▽P～Q非線性關系主要由水流慣性阻力作用的逐漸增強產(chǎn)生[14]，可采用Forchheimer(1901)提出的二項式方程來描述其滲流規(guī)律：

-▽P=AQ+BQ2

(2)

式中：▽P為壓力梯度，Pa/m；Q為滲流流量，m3/s；A為線性項系數(shù)，(Pa·s)/m4；B為非線性項系數(shù)，(Pa·s2)/m7；Ah為過流面積，m2；μ為流體的動力黏度，Pa·s；ρ為密度，kg/m3；k為過流介質(zhì)的固有滲透率，m2；β為非達西系數(shù)，m-1，又被稱為紊流因子、非達西滲流慣性系數(shù)，是研究滲流非線性的重要參數(shù)[15-17]。

由式(1)可知，F(xiàn)orchheimer方程。在形式上僅僅是對達西方程增加了一個二次項(慣性項)，意義上則表示非達西滲流的能量損失不再僅由黏滯阻力控制，而是由黏滯阻力和慣性阻力2者共同決定[18]，換言之，滲流非達西效應取決于慣性阻力作用的相對大小。

圖5 Forchheimer方程系數(shù)A和B與圍壓的關系曲線Fig.5 Variation of coefficients A and B with confining stress

采用Forchheimer非達西滲流方程對▽P～Q數(shù)據(jù)進行擬合，得到各曲線的擬合系數(shù)見表1。由表1可知，F(xiàn)orchheimer方程對非線性滲流的擬合效果很好，相關系數(shù)R2值均介于0.983 6和0.999 9之間，接近于1。為了更為直觀地分析，圖5給出了擬合系數(shù)A、B與圍壓的關系曲線。由圖5可知，線性項系數(shù)A隨圍壓嚴格單調(diào)遞增，而非線性項系數(shù)B則隨著圍壓的增大呈總體上升趨勢。由式(2)可知，系數(shù)A反映的是破碎巖體的過流能力，圍壓增大，巖樣內(nèi)部的裂隙被壓密，過流能力減弱，導致系數(shù)A增大；系數(shù)B表征流動的慣性效應，破碎巖樣的裂隙分布復雜，形態(tài)各異，在圍壓增大的過程中，縱向裂隙閉合，橫向裂隙擴展，2者對滲流的非線性具有相反的作用機制，且作用強度具有波動性，但總體上縱向裂隙作用強于橫向裂隙作用，因而系數(shù)B在局部存在起伏波動，但總體上增大。

由Forchheimer方程的式(1)和式(2)可知，當忽略滲流流體的物理特性變化時，線性項系數(shù)A取決于巖石固有滲透率k，且與k成反比，非線性項系數(shù)B則由非達西系數(shù)β確定，并與之成正比。因此圖5所示的系數(shù)A、B的變化關系在一定程度上反映了非達西系數(shù)β與固有滲透率k之間存在的相關關系；另一方面，在裂隙和多孔介質(zhì)的滲流研究中，非達西系數(shù)β與固有滲透率k之間廣泛采用冪函數(shù)關系描述[16,17,19,20]?；诒疚牡臐B流試驗結果，計算出3個破碎花崗巖巖樣在各圍壓下的固有滲透率k與非達西系數(shù)β，并將2者關系曲線繪制于圖6中，采用冪函數(shù)關系式進行擬合，相關系數(shù)R2=0.916 7，擬合程度較高，表達式如下：

β=6.962×10-13k-1.709

(3)

3 非達西滲流評估

如前所述，雷諾數(shù)Re和福希海默數(shù)F0是劃分達西流和非達西流常用的指標，但是確定臨界值卻十分困難，因此本文在試驗結果分析得到的非達西系數(shù)與固有滲透率關系式的基礎上，建立另一類簡單直觀的非達西滲流評估方法。

表1 Forchheimer方程擬合系數(shù)Tab.1 The best-fitted coefficients of Forchheimer Equation

圖6 非達西系數(shù)β與固有滲透率k關系曲線Fig.6 Variation of non-Darcy coefficient β with intrinsic permeability k

3.1 水力梯度比ΦP

假設破碎巖體的滲流流量Q已知，則由式(1)計算得到的非達西水力梯度|-▽P|F與由達西定律計算所得的達西水力梯度|-▽P|D2者的比值ΦP滿足[14]：

(4)

上式表明非達西與達西水力梯度的比值ΦP由福希海默數(shù)F0決定，因此參數(shù)ΦP也具有表征滲流非達西程度的意義。將式(2)代入式(4)，并考慮式(3)所得的β～k冪函數(shù)關系，則有：

(5)

可知，ΦP是固有滲透率k和滲流流速v的函數(shù)。

根據(jù)式(4)和式(5)，圖7(a)和圖7(b)分別做出了ΦP～F0關系曲線和不同ΦP值下的v～k關系曲線。由圖7可知，ΦP大于1，表明滲流慣性效應的存在將引起額外的壓力損失，因此產(chǎn)生與達西流相等的流量需要更大的水力梯度。在評估滲流的非線性時，只要確定了水力梯度比ΦP的值，就可從圖7(b)找到對應的臨界v～k曲線，該曲線上部區(qū)域的滲流為非達西流，下部區(qū)域滲流的非線性則忽略不計，可采用達西定律進行滲流分析。雖然忽略非達西效應的滲流計算結果必然存在偏差，但只要在容許范圍內(nèi)則影響不大。假設某工程要求水力梯度的容許相對誤差為10%，即對應臨界水力梯度比ΦP=1.11(該值對應Zeng等[11]提出的達西流與非達西流的臨界值F0=0.11)，圖7(b)中ΦP=1.11對應v～k曲線即是相應的臨界曲線。

圖7 水力梯度比值ΦP評估圖Fig.7 Evaluation of pressure gradient ratio ΦP

3.2 流量比ΦQ

假設破碎巖體的水力梯度|-▽P|已知，則由式(2)解得的非達西流量QF為[14]：

(6)

定義流量比ΦQ為非達西流量QF與由達西定律計算得到的達西流量QD的比值。為求得ΦQ，聯(lián)立式(2)和達西方程，則：

AQD=AQF+BQ2F

(7)

故：

(8)

聯(lián)立式(2)、式(3)和式(8)，得：

(9)

顯然，ΦQ是固有滲透率k和水力梯度|-▽P|的函數(shù)。

同樣地，根據(jù)式(8)和式(9)分別得出了ΦQ～F0關系曲線和不同ΦQ值下的|-▽P|～k關系曲線，見圖8。由圖8(a)可知，流量比ΦQ的值為0～1，因為在相同的水力梯度作用下，非達西流動的慣性項不可忽略，致使非達西流量小于達西流量。在評估滲流的非線性時，只要流量比ΦQ的值被確定后，就可以從圖8(b)找到對應的臨界|-▽P|～k曲線，該曲線下部區(qū)域的滲流可忽略非線性，采用達西定律進行滲流分析，其他區(qū)域的滲流作為非達西流處理。若流量的容許相對誤差定為20%，則臨界流量比ΦQ=0.83，圖8(b)中ΦQ=0.83的|-▽P|～k曲線即是相應的臨界曲線。

圖8 流量比值ΦQ評估圖Fig.8 Evaluation of flow rate ratio ΦQ

與雷諾數(shù)Re和福希海默數(shù)F0相比，水力梯度比ΦP和流量比ΦQ2個評估參數(shù)的物理意義更為清晰直觀，通過限制非達西水力梯度或者非達西流量相對于達西流情況的偏離程度來判定滲流的非線性，對工程應用具有重要的意義。例如油氣開發(fā)更關注開采量，因此可選擇表示流量偏離的參數(shù)ΦQ來考慮滲流非線性，其臨界值大小由工程實踐中可容許的流量偏差決定；當分析與巖體滲透壓力有關的穩(wěn)定性問題時，選擇水力梯度比ΦP來判定滲流非線性更為合適。但值得指出的是，式(5)和式(9)是在式(3)的基礎上推導得到的，冪函數(shù)關系式(3)的擬合系數(shù)值與巖體介質(zhì)的滲流特性有關，因此對于不同的透水介質(zhì)，式(5)和式(9)中各系數(shù)值也會有所差異。

4 結 語

本文對3個陽江花崗巖破碎后，開展了不同圍壓條件下的高水力梯度滲流試驗，并對之進行了分析討論，在試驗結果的基礎上提出了一種簡單直觀的非達西滲流評估方法。具體得到如下結論。

(1)試驗結果表明，破碎巖體中的滲流在高水力梯度下表現(xiàn)出顯著的非線性，采用Forchheimer方程擬合▽P～Q關系的相關系數(shù)為0.983 6～0.999 9。

(2)隨著圍壓的上升，與巖樣過水能力有關的線性項系數(shù)A嚴格單調(diào)遞增，反映滲流慣性阻力作用的非線性系數(shù)B呈總體增大趨勢；非達西系數(shù)β與固有滲透率k具有相關性較好的冪函數(shù)關系。

(3)評估滲流的非線性時，采用水力梯度比ΦP和流量比ΦQ更加簡單直觀，且物理意義清晰，可根據(jù)工程具體情況，通過限制非達西流與達西流情況下2者的水力梯度比或者流量比來確定滲流非線性的臨界條件，決定是否采用達西定律進行滲流分析。

□

[1] 黃先伍, 唐 平, 繆協(xié)興, 等. 破碎砂巖滲透特性與孔隙率關系的試驗研究[J]. 巖土力學, 2005,26(9):1 385-1 388.

[2] 李順才, 繆協(xié)興, 陳占清, 等. 承壓破碎巖石非 Darcy 滲流的滲透特性試驗研究[J]. 工程力學, 2008,25(4):85-92.

[3] Thauvin F, Mohanty K K. Network modeling of non-Darcy flow through porous media[J]. Transport in Porous Media, 1998,31(1):19-37.

[4] 胡大偉, 周 輝, 謝守益, 等. 峰后大理巖非線性滲流特征及機制研究[J]. 巖石力學與工程學報, 2009,28(3):451-458.

[5] Zhou J Q, Hu S H, Fang S, et al. Nonlinear flow behavior at low Reynolds numbers through rough-walled fractures subjected to normal compressive loading[J]. International Journal of Rock Mechanics and Mining Sciences, 2015,80:202-218.

[6] 孟如真, 胡少華, 陳益峰, 等. 高滲壓條件下基于非達西流的裂隙巖體滲透特性研究[J]. 巖石力學與工程學報, 2014,33(9):1 756-1 764.

[7] 李 健, 黃冠華, 文 章, 等. 兩種不同粒徑石英砂中非達西流動的實驗研究[J]. 水利學報, 2008,39(6):726-732.

[8] Zimmerman R W, Al-yaarubi A, Pain C C, et al. Non-linear regimes of fluid flow in rock fractures[J]. International Journal of Rock Mechanics and Mining Sciences, 2004,41:163-169.

[9] Javadi M, Sharifzadeh M, Shahriar K, et al. Critical Reynolds number for nonlinear flow through rough-walled fractures: The role of shear processes[J]. Water Resources Research, 2014,50(2):1 789-1 804.

[10] 周創(chuàng)兵, 陳益峰, 姜清輝, 等. 復雜巖體多場廣義耦合分析導論[M]. 北京: 中國水利水電出版社, 2008.

[11] Zeng Z, Grigg R. A criterion for non-Darcy flow in porous media[J]. Transport in Porous Media, 2006,63(1):57-69.

[12] Macici P, Mesini E, Viola R. Laboratory measurements of non-Darcy flow coefficients in natural and artificial unconsolidated porous media[J]. Journal of Petroleum Science and Engineering, 2011,77(3):365-374.

[13] Liu R, Li B, Jiang Y. Critical hydraulic gradient for nonlinear flow through rock fracture networks: the roles of aperture, surface roughness, and number of intersections[J]. Advances in Water Resources, 2016, 88:53-65.

[14] Chen Y F, Zhou J Q, Hu S H, et al. Evaluation of Forchheimer equation coefficients for non-Darcy flow in deformable rough-walled fractures[J]. Journal of Hydrology, 2015,529:993-1 006.

[15] Firdaouss M, Guermond J L, Quéré L P. Nonlinear corrections to Darcy's law at low Reynolds numbers[J]. Journal of Fluid Mechanics, 1997,343:331-350.

[16] 許 凱, 雷學文, 孟慶山, 等. 非達西滲流慣性系數(shù)研究[J]. 巖石力學與工程學報, 2012,31(1):164-170.

[17] Geertsma J. Estimating the coefficient of inertial resistance in fluid flow through porous media[J]. Society of Petroleum Engineers Journal, 1974,14(5):445-450.

[18] Zhou J Q, Hu S H, Chen Y F, et al. The friction factor in the forchheimer equation for rock fractures[J]. Rock Mechanics and Rock Engineering, 2016,49(8):3 055-3 068.

[19] Li D, Engler T W. Literature review on correlations of the non-Darcy coefficient[C]∥SPE Permian Basin Oil and Gas Recovery Conference. Society of Petroleum Engineers, 2001.

[20] 王 媛, 秦 峰, 夏志皓, 等. 深埋隧洞涌水預測非達西流模型及數(shù)值模擬[J]. 巖石力學與工程學報, 2012,31(9):1 862-1 868.