部分充填裂隙網(wǎng)絡(luò)巖體的滲透特性分析

鄭 雙，張 斌，符文熹

(1. 四川大學(xué)水利水電學(xué)院，成都 610065；2. 四川大學(xué)水力學(xué)與山區(qū)河流保護(hù)國(guó)家重點(diǎn)實(shí)驗(yàn)室，成都 610065; 3. 四川省公路規(guī)劃勘察設(shè)計(jì)研究院有限公司，成都 610041)

0 引 言

西南地區(qū)是我國(guó)水利能源開(kāi)發(fā)的主要地區(qū)，目前已修建和在建大量的水利工程，大壩特別是高壩多位于基巖上，基巖的滲透性往往決定了整個(gè)壩基的滲漏[1]。然而，西南地區(qū)地質(zhì)條件十分復(fù)雜，巖體在長(zhǎng)期內(nèi)外動(dòng)力地質(zhì)作用下，巖體內(nèi)部常常發(fā)育大量的裂隙，形成裂隙網(wǎng)絡(luò)巖體[2]。因此，研究裂隙巖體的滲透特性具有重大工程意義。

對(duì)于裂隙巖體，可采用裂隙網(wǎng)絡(luò)方法進(jìn)行滲流分析。Wittke[3]不考慮巖塊的滲流，以真實(shí)的裂隙網(wǎng)絡(luò)為基礎(chǔ)，利用線單元法，創(chuàng)建了線素法，最早提出了巖體裂隙網(wǎng)絡(luò)滲流模型。Wilson和Witherspoon[4]將有限單元法引入到線素法，計(jì)算了壩下理想裂隙網(wǎng)絡(luò)巖體的滲流特性。王恩志[5]利用圖論的思想對(duì)巖體裂隙網(wǎng)絡(luò)進(jìn)行了數(shù)學(xué)描述，給出了描述裂隙網(wǎng)絡(luò)連通關(guān)系的數(shù)學(xué)方法，提出了一個(gè)新的地下水滲流數(shù)值模型。莫海鴻[6]基于拓?fù)鋵W(xué)的理論，提出了裂隙網(wǎng)絡(luò)的離散模型。趙紅亮[7]通過(guò)計(jì)算機(jī)模擬構(gòu)筑了三維裂隙網(wǎng)絡(luò)系統(tǒng)。經(jīng)過(guò)不斷的發(fā)展，巖體裂隙網(wǎng)絡(luò)滲流模型逐步走向成熟，得到了廣泛的應(yīng)用。

對(duì)于裂隙網(wǎng)絡(luò)巖體的滲透特性，通常基于立方定理計(jì)算單裂隙的滲流量，通過(guò)裂隙組合進(jìn)行分析[8]。Zimmerman[9]通過(guò)研究單相流在粗糙巖石裂隙中的運(yùn)動(dòng)，基于立方定理給出了Navier-Stokes方程的解，從理論上說(shuō)明了立方定理的適用條件(雷諾數(shù)Re<<1)。柴軍瑞[10]研究了裂隙水流為紊流且不滿足立方定理?xiàng)l件下的裂隙滲透特性，給出了巖體裂隙網(wǎng)絡(luò)非穩(wěn)定滲流的控制方程。朱紅光[11]在立方定理的基礎(chǔ)上，推導(dǎo)了交叉裂隙的滲流計(jì)算公式，并將其應(yīng)用到巖體裂隙網(wǎng)絡(luò)分析中。劉日成[12]為了研究裂隙網(wǎng)絡(luò)的非線性滲流，建立了兩種離散裂隙網(wǎng)絡(luò)模型，根據(jù)相應(yīng)邊界條件，對(duì)Navier-Stokes方程進(jìn)行求解，并且通過(guò)室內(nèi)透水試驗(yàn)進(jìn)行了驗(yàn)證，提出粗糙的裂隙表面會(huì)導(dǎo)致滲流量降低，但是誤差可以忽略不計(jì)。

綜上所述，巖體裂隙的粗糙度、開(kāi)度、交叉情況以及裂隙流的穩(wěn)態(tài)和非穩(wěn)態(tài)問(wèn)題都得到了廣泛細(xì)致的研究，然而，由于風(fēng)化和侵蝕作用，天然巖體裂隙中常常存在一定厚度的充填物，此時(shí)直接采用立方定理進(jìn)行計(jì)算，將產(chǎn)生較大誤差。為了研究含充填介質(zhì)裂隙網(wǎng)絡(luò)巖體的滲透特性，本文建立了部分充填條件下的單裂隙滲透模型，結(jié)合Navier-Stokes方程和Brinkman-extended Darcy方程，得到了單裂隙流場(chǎng)分布特征，推求出單裂隙滲流量的解析解；在無(wú)充填的情況下，該滲流量表達(dá)式可以簡(jiǎn)化為經(jīng)典的立方定理。后基于裂隙網(wǎng)絡(luò)滲流原理，構(gòu)建了含充填介質(zhì)裂隙巖體滲流網(wǎng)絡(luò)的數(shù)學(xué)模型，可以得到裂隙網(wǎng)絡(luò)節(jié)點(diǎn)水頭。結(jié)合案例分析，探討了充填介質(zhì)對(duì)裂隙網(wǎng)絡(luò)巖體滲透特性的影響。

1 含充填介質(zhì)裂隙滲透模型

1.1 理論模型

本文為了研究充填物對(duì)裂隙巖體滲透特性的影響，建立如圖1所示的概化模型，裂隙寬度為e，填充物厚度為b1，空裂隙的厚度為b2，模型長(zhǎng)度為L(zhǎng)。假設(shè)裂隙中滲流流速為u，充填介質(zhì)中滲流流速為w。模型存在下述的假設(shè)，①流體不可壓縮，為層流；②滲流只沿著x軸的流動(dòng)；③裂隙中滲流u符合Navier-Stokes方程[13]，充填物中的滲流w符合Brinkman-extended Darcy方程[14]；④充填物與裂隙滲流的交界面符合流速相等和剪應(yīng)力連續(xù)的邊界條件[15]。

圖1 巖體裂隙的概化模型Fig.1 Generalized model of rock mass fracture

根據(jù)假設(shè)③，裂隙中流速u(mài)和充填介質(zhì)中滲透流速w分別滿足Navier-Stokes方程和Brinkman-extended Darcy方程，同時(shí)它們還滿足連續(xù)性方程。

(1)

(2)

式中：p為流體壓強(qiáng)，M/(L·T2)；ρ為水流密度，M/L3；n為孔隙率；η為水的動(dòng)力黏度，M/(L·T)；K為充填介質(zhì)的滲透率，L2；▽為Hamilton算子。

又由模型假設(shè)b，可得：

(3)

分別代入式(1)、(2)，可得

(4)

求解式(4)中的微分方程，可得：

(5)

(6)

式中：ΔP分別為裂隙兩端的壓強(qiáng)差。

根據(jù)假設(shè)，式(5)和(6)滿足下述邊界條件：

(2)裂隙滲流與巖石交界面處流速為0，即y=b1+b2時(shí)，ux= 0。

(3)充填物中滲流與巖石交界面處流速為0，即y= 0時(shí)，wx= 0。

聯(lián)立上述邊界條件即可求出參數(shù)B1、B2、C1、C2。

(7)

根據(jù)流速分布，可以求得整個(gè)裂隙斷面的流量為：

(8)

可以發(fā)現(xiàn)當(dāng)填充比(β=b1/b2)為0，即b1= 0時(shí)，式(8)可以化簡(jiǎn)得到：

(9)

式(9)即為經(jīng)典的立方定理計(jì)算的裂隙流量。式(8)更具有普遍性，立方定理是它的一種特殊情況。

1.2 模擬試驗(yàn)

本文為了驗(yàn)證理論模型的正確性，開(kāi)展物理模型試驗(yàn)進(jìn)行驗(yàn)證。試驗(yàn)裝置設(shè)計(jì)圖如圖2所示，整個(gè)模型長(zhǎng)為30 cm，寬為20 cm，高為20 cm。其中的不透水巖塊利用亞克力版制作的長(zhǎng)方體模擬，長(zhǎng)為20 cm，寬為20 cm，高為9.5 cm，可以達(dá)到不透水的效果，充填物采用河砂模擬，經(jīng)過(guò)試驗(yàn)測(cè)定，滲透率K為3.54×10-8cm2，孔隙率n為0.39。為了防止填充物發(fā)生侵蝕破壞，采用鐵絲網(wǎng)進(jìn)行固定。根據(jù)設(shè)計(jì)進(jìn)行組裝，裝填完的試驗(yàn)裝置如圖3所示。

圖2 試驗(yàn)裝置設(shè)計(jì)圖Fig.2 Design drawing of test device

圖3 試驗(yàn)裝置Fig.3 Test device

參照土工試驗(yàn)方法標(biāo)準(zhǔn)進(jìn)行滲透試驗(yàn)，分別在填充比為25%、50%、75%，在不同水頭差ΔH條件下進(jìn)行試驗(yàn)。將試驗(yàn)結(jié)果記錄在表1中。表1中試驗(yàn)值均為多次試驗(yàn)的平均值。

將試驗(yàn)值與理論值進(jìn)行對(duì)比，如圖4所示。從圖4中可以發(fā)現(xiàn)，理論值與試驗(yàn)值較接近，誤差相對(duì)較小(小于10%)。在一定程度上可以驗(yàn)證本文提出含充填介質(zhì)的巖體裂隙的滲流模型。

圖4 試驗(yàn)值與理論值對(duì)比Fig.4 comparison between experimental and theoretical values

2 巖體裂隙網(wǎng)絡(luò)模型

裂隙網(wǎng)絡(luò)滲流模型首先由Wittke提出，后經(jīng)過(guò)Louis、Wilson和王恩志等人的發(fā)展，現(xiàn)可以用于二維穩(wěn)定和非穩(wěn)定裂隙流的計(jì)算[16]。裂隙網(wǎng)絡(luò)模型如圖5所示，以巖體裂隙的交點(diǎn)為節(jié)點(diǎn)，兩個(gè)節(jié)點(diǎn)之間的裂隙為一個(gè)線元。假設(shè)滲流域內(nèi)有n個(gè)節(jié)點(diǎn)，與第i個(gè)節(jié)點(diǎn)相連的線元有m個(gè)，在節(jié)點(diǎn)i處，根據(jù)水流均衡原理[17]，可得節(jié)點(diǎn)i處的水流控制方程為：

表1 試驗(yàn)記錄表Tab.1 Test record table

圖5 裂隙網(wǎng)絡(luò)模型Fig.5 Fracture network model

(10)

式中：Qi為節(jié)點(diǎn)i的源項(xiàng)，i取值為：1，2，3,…,n。

若裂隙中存在充填物，則裂隙中的滲流量可按式(8)進(jìn)行計(jì)算。將式(8)帶入式(10)可得：

(11)

式中：bj-1為線元j中填充物的厚度；bj-2為線元j中空裂隙的厚度；wj-x為填充物中滲流的流速；uj-x為空裂隙中的滲流流速。

式(11)為裂隙中存在充填物條件下，節(jié)點(diǎn)i處的流量控制方程。由于整個(gè)滲流域內(nèi)總共有n個(gè)節(jié)點(diǎn)，則整體的滲流控制方程為：

(12)

式中：A為裂隙網(wǎng)絡(luò)的銜接矩陣，其為n×m階矩陣。

它是用來(lái)描述滲流域中線元與節(jié)點(diǎn)的鏈接關(guān)系的矩陣，根據(jù)文獻(xiàn)[5]可以求得。Q= {Q1,Q2, …，Qn}T。結(jié)合模型的邊界條件，可以得到裂隙中存在填充物條件下，巖體裂隙網(wǎng)絡(luò)滲流的數(shù)學(xué)模型為：

(13)

式中：Γ1和Γ2分別為已知水頭邊界；h1和h2分別為相應(yīng)邊界上的水頭。

通過(guò)Matlab編程，求解式(13)，可以得到各個(gè)節(jié)點(diǎn)的水頭。

3 算 例

如圖6所示，大壩下的巖體中存在兩組相互垂直的理想裂隙，巖塊不透水，邊界AC、BD、CD、EF分別為不透水邊界。滲流沿著裂隙網(wǎng)絡(luò)流動(dòng)。邊界AE的水頭為100 m，邊界BF的水頭為20 m，AB長(zhǎng)度為100 m，AC長(zhǎng)度為40 m。裂隙的寬度為0.01 m，深度為0.1 m。水平向的裂隙含有厚度為0.005 m的充填介質(zhì)，充填比β=50%，滲透率K=3.5×10-8m2。根據(jù)裂隙網(wǎng)絡(luò)滲流原理，對(duì)節(jié)點(diǎn)進(jìn)行編號(hào)，整個(gè)滲流域內(nèi)共有51個(gè)節(jié)點(diǎn)，66個(gè)線元。

圖6 模型與節(jié)點(diǎn)編號(hào)和節(jié)點(diǎn)詳圖Fig.6 Model node numbered and node detail

在裂隙分別為充填和不充填的條件下，根據(jù)本文的數(shù)學(xué)模型(含充填介質(zhì))進(jìn)行計(jì)算，得到各節(jié)點(diǎn)的水頭，如圖7所示。根據(jù)節(jié)點(diǎn)編號(hào)，可以分為5個(gè)區(qū)間(1～9、10～20、21～31、32～42、43～51)，分別表示2條水平邊界(AB、CD)和3條水平向裂隙上的節(jié)點(diǎn)水頭。可以發(fā)現(xiàn)，在每個(gè)區(qū)間內(nèi)，兩種方式計(jì)算的水頭變化趨勢(shì)保持一致，都是逐漸減小，但是具有一定的差距，明顯在裂隙含充填介質(zhì)條件下，計(jì)算得到的水頭變化幅度更大，具有更大的最大值和更小的最小值。

圖7 各節(jié)點(diǎn)水頭Fig.7 Water head of each node

圖8為滲流域內(nèi)的水頭等值線圖，可以發(fā)現(xiàn)，兩種方式計(jì)算的結(jié)果具有明顯的差別。在靠近滲流入口和溢出面的位置，裂隙含充填介質(zhì)條件下，計(jì)算的水頭等值線更加稀疏，水力比降較小，說(shuō)明考慮充填介質(zhì)會(huì)導(dǎo)致裂隙巖體的滲流量減小。另一方面，在滲流域的中部，含有充填介質(zhì)裂隙巖體的等水頭線較密，水力比降較大，滲流的流速增大，可能會(huì)導(dǎo)致局部的侵蝕破壞，影響裂隙巖體的強(qiáng)度。

在裂隙無(wú)充填的條件下，基于立方定理計(jì)算的滲流量為2.25 m3/s，在裂隙存在充填物的條件下，根據(jù)本文數(shù)學(xué)模型計(jì)算的滲流量為1.43 m3/s，相比下降了36.44%。可見(jiàn)，巖體裂隙填充物對(duì)裂隙巖體滲流特性的影響較大，不能忽略不計(jì)。

圖8 滲流域內(nèi)水頭等值線圖Fig.8 Isometric map of water head in seepage basin

為了研究充填介質(zhì)的滲透率K對(duì)滲流量Q的影響，在不同填充比β的工況下，計(jì)算得到滲透率K與流域滲流量Q的關(guān)系，從圖9中可以看出，隨著滲透率下降，滲流量逐漸減小，當(dāng)滲透率較大時(shí)，滲透率變化對(duì)流量的影響較大。當(dāng)滲透率低于1×10-7cm2，滲流量隨著滲透率的變化率將趨于0。填充度(β)越高，滲流量越小，當(dāng)填充度為0時(shí)，計(jì)算的滲流量和文獻(xiàn)基于立方定理計(jì)算的滲流量相等，這個(gè)也印證了式(9)的推導(dǎo)。

圖9 不同填充度下滲透率與流量的關(guān)系Fig.9 Relationship between permeability and seepage volume with different filling degree

含填充介質(zhì)裂隙巖體的滲透特性不僅與填充介質(zhì)的滲透系數(shù)有關(guān)，還和裂隙的填充比β有關(guān)，在實(shí)際的工程中，它可以通過(guò)野外調(diào)查獲得。為了研究填充比對(duì)含充填介質(zhì)裂隙巖體滲透特性的影響，在充填介質(zhì)具有不同滲透率的情況下，計(jì)算填充比與滲流量的關(guān)系，獲得圖10中的曲線，發(fā)現(xiàn)滲流量與填充比具有非線性的關(guān)系，曲線具有兩個(gè)拐點(diǎn)。當(dāng)填充比較小，低于30%時(shí)，曲線斜率較小，滲流量隨填充比的變化率較小。當(dāng)填充比較高，位于30%到70%之間時(shí)，曲線斜率較陡，滲透量隨填充比的變化率較大。當(dāng)填充比高于70%，曲線再次變緩，滲透量隨填充比的變化率再次減小。

從不同滲透系數(shù)的關(guān)系曲線上看，填充比越高，充填介質(zhì)滲透系數(shù)對(duì)滲流量的影響也越大，這是由于當(dāng)填充比較低時(shí)(β<50%)，巖體裂隙的滲透量主要由空裂隙中的滲流量控制，滲流量受填充介質(zhì)的影響較小；當(dāng)填充比較高時(shí)(β>50%)，巖體裂隙的滲透量主要由填充介質(zhì)中的滲流量控制，滲流量受填充介質(zhì)的影響較大。

圖10 不同滲透率下填充比與流量的關(guān)系Fig.10 relationship between filling ratio and seepage volume under different permeability

4 結(jié) 論

通過(guò)對(duì)部分充填條件下的單裂隙滲透模型的理論分析，推導(dǎo)出了含充填介質(zhì)的巖體裂隙的流速分布，獲得了單裂隙的流量公式。當(dāng)填充度為0時(shí)，該滲流量公式可以化簡(jiǎn)到立方定理，室內(nèi)滲流試驗(yàn)的結(jié)果驗(yàn)證了理論分析的正確性。后結(jié)合裂隙網(wǎng)絡(luò)的滲流原理，建立了含充填介質(zhì)的裂隙網(wǎng)絡(luò)的數(shù)學(xué)模型，可以求解裂隙巖體的滲流場(chǎng)。

在算例中，以填充比和填充介質(zhì)的滲透 系數(shù)為指標(biāo)，討論了填充介質(zhì)性質(zhì)對(duì)裂隙巖體滲透性質(zhì)的影響。充填介質(zhì)滲透率越小、填充比越高，裂隙巖體的滲流量越低。充填介質(zhì)的存在會(huì)還導(dǎo)致裂隙兩端的水力比降增加，提高裂隙的滲流流速，增加局部發(fā)生沖刷侵蝕破壞的可能性，進(jìn)而影響裂隙巖體的整體強(qiáng)度。