不同粗糙度裂隙滲流特性數(shù)值模擬研究

李治豪，陳世江

(內(nèi)蒙古科技大學(xué) 礦業(yè)研究院，內(nèi)蒙古 包頭 014000)

巖體中的裂隙錯(cuò)綜復(fù)雜，形態(tài)各異，流體在裂隙中的運(yùn)移滲流受多重因素的影響，而裂隙本身表面粗糙度就是其中一個(gè)十分重要的因素。張戈[1]基于Boltzmann方法模擬不同JRC(節(jié)理粗糙度系數(shù))粗糙單裂隙內(nèi)的流體滲流，通過(guò)考慮粗糙特性和迂曲度的雙重影響對(duì)粗糙裂隙模型進(jìn)行了計(jì)算，并對(duì)立方定律進(jìn)行了一定程度的修正；王鵬飛[2]利用3D打印技術(shù)制作裂隙然后澆筑出JRC=1～20的裂隙試件，對(duì)不同粗糙度、不同裂隙開(kāi)度的貫通裂隙的滲流特性進(jìn)行了研究，得出了在圍壓作用下粗糙度、裂隙開(kāi)度與裂隙滲透率之間的關(guān)系；潘汝江[3]采用數(shù)值模擬的方法研究了曲折度對(duì)單裂隙在不同粗糙度下滲流機(jī)理的影響，模擬結(jié)果顯示立方定律并不適用于壓差較大的情況，并結(jié)合曲折度對(duì)立方定律進(jìn)行了修訂；曲冠政[4]以Lattice Boltzmann方法研究了剪切滑移作用對(duì)裂縫中流體滲流的影響，研究了滲流能力和滲透率受裂縫面滑移作用影響的變化情況；王學(xué)敏[5]以慈林山礦為研究背景，采用COMSOL數(shù)值模擬軟件分析了水力壓裂過(guò)程中不同裂縫角和不同應(yīng)力差對(duì)巖石破裂的影響，研究結(jié)果表明應(yīng)力差和裂縫角的變化會(huì)改變巖石的起裂壓力；覃源[6]為研究節(jié)理粗糙度對(duì)水流在裂隙中滲流的影響，采用了分形維數(shù)的方法，通過(guò)研究發(fā)現(xiàn)，節(jié)理表面起伏角、伸長(zhǎng)率和分維數(shù)是影響水流流速的敏感因素;王珂[7]針對(duì)在應(yīng)力和滲流侵蝕耦合作用條件下巖石裂隙表面粗糙度的變化情況，探究了巖石裂隙的滲流機(jī)理演變規(guī)律，得到了應(yīng)力作用下裂隙表面粗糙度不同時(shí)其滲流量、裂隙寬度和滲透率的變化規(guī)律；馬德宜[8]對(duì)砂巖劈裂裂隙面JRC進(jìn)行了計(jì)算，并且以巖石滲透儀為平臺(tái)進(jìn)行了試驗(yàn)，得到了滲流量與圍壓、軸壓及JRC值之間的關(guān)系，并且以此建立了滲流量預(yù)測(cè)模型；劉樂(lè)[9]采用三軸水力壓裂系統(tǒng)對(duì)巖石試件進(jìn)行了水力壓裂試驗(yàn)，結(jié)果表明裂隙拓展面積會(huì)隨壓力增大而增大，且當(dāng)流量超過(guò)一定值時(shí)水力壓裂效果會(huì)明顯降低；劉基[10]對(duì)煤層頂板砂巖進(jìn)行了井下放水試驗(yàn)，以此來(lái)獲取滲透系數(shù)，并通過(guò)對(duì)比分析為礦井涌水量計(jì)算提供依據(jù)；蒙學(xué)禮[11]對(duì)裂隙巖體在不同圍壓下的滲流進(jìn)行了試驗(yàn)研究，結(jié)果表明裂隙的粗糙度越小，裂隙內(nèi)部的非線性滲流特征越明顯；高超[12]針對(duì)單裂隙巖體的滲流中裂隙尺寸效應(yīng)進(jìn)行了數(shù)值模擬，結(jié)果表明裂隙開(kāi)度對(duì)巖體裂隙滲透性的影響比粗糙度和水力梯度的影響更顯著；孫可明[13]將物理試驗(yàn)與數(shù)值模擬相結(jié)合，探究了應(yīng)力加載下不同粗糙度裂隙巖心的滲流特性，結(jié)果表明裂隙滲透率與粗糙度呈現(xiàn)正相關(guān)關(guān)系；張燁[14]通過(guò)人工制作粗糙裂隙巖體試件，改變裂隙開(kāi)度和水力梯度進(jìn)行了滲流試驗(yàn)，結(jié)果表明滲流速度較高時(shí)，水力梯度與流速呈現(xiàn)出線性關(guān)系。

以上研究成果都在一定程度上探究了粗糙單裂隙的滲流特性，但對(duì)于不同JRC值裂隙的滲流特性研究還不夠全面。因此，開(kāi)展不同JRC值下粗糙裂隙的滲流研究仍具有十分重要的意義。筆者在前人的研究基礎(chǔ)上，建立起不同JRC值的裂隙模型，并且用COMSOL軟件進(jìn)行數(shù)值模擬，探究不同JRC值下裂隙的滲流速度、滲透率和應(yīng)力變化規(guī)律，并針對(duì)不同裂隙開(kāi)度模型進(jìn)行比較，以確保模擬結(jié)果的準(zhǔn)確性。

1 數(shù)值模型的建立

不同粗糙度的單裂隙微觀滲流模型的構(gòu)建是基于確定的JRC值的Barton十條曲線進(jìn)行的，如圖1所示。

圖1 具有確定JRC值的Barton標(biāo)準(zhǔn)曲線

JRC為裂隙結(jié)構(gòu)面的粗糙度系數(shù)，通過(guò)對(duì)標(biāo)準(zhǔn)Barton曲線在高精度分辨率下的提取后導(dǎo)入CAD中進(jìn)行格式轉(zhuǎn)換，然后利用COMSOL Multiphysics建立數(shù)值模型。裂隙的開(kāi)度是Barton曲線在豎直方向平移的距離，由于自然界巖體中的裂隙開(kāi)度相對(duì)較小，為使模擬與實(shí)際裂隙開(kāi)度一致，故以微米為模擬裂隙開(kāi)度的單位建立微觀裂隙模型。結(jié)合實(shí)際裂隙情況并保證模擬計(jì)算過(guò)程收斂，構(gòu)建的初始裂隙通道的開(kāi)度為10 μm，裂隙長(zhǎng)度為200 μm，長(zhǎng)度與開(kāi)度的比值為20，具有一定的代表性。裂隙通道入口邊界處流體流入的滲透壓為100 Pa。模擬網(wǎng)格為物理場(chǎng)控制網(wǎng)格，裂隙通道內(nèi)部為自由三角形網(wǎng)格，為確保模擬的精度，將邊界處設(shè)為矩形四邊形網(wǎng)格。

2 數(shù)值模擬方法及條件

基于不可壓縮的Navier-Stokes方程，考慮重力加速度的影響，在遵循質(zhì)量守恒、動(dòng)量守恒和能量守恒的基礎(chǔ)上，對(duì)流體場(chǎng)中的各個(gè)宏觀量進(jìn)行模擬計(jì)算，得到滲流場(chǎng)中的速度、壓力場(chǎng)中的壓力及渦流場(chǎng)中的渦流形式。其方程形式與模型信息如下[9]:

(1)

(2)

式中:ρ為流體的密度，kg/m3；t為時(shí)間步長(zhǎng)，s;u為速度矢量，u=(u,v)，m/s；p為標(biāo)準(zhǔn)壓力，Pa；ν為動(dòng)力黏度，Pa·s。

模擬中流體為不可壓縮介質(zhì)流動(dòng)，流體介質(zhì)的動(dòng)力黏度為1×10-3Pa·s，密度為1×103kg/m3，重力加速度為9.8 m/s2。

以光滑單裂隙通道滲流過(guò)程的廣義立方定律為基礎(chǔ)，對(duì)粗糙裂隙進(jìn)行立方定律修正后，滲透率k的計(jì)算表達(dá)式為[16]:

(3)

式中:n為裂隙通道的數(shù)量；b為裂隙的開(kāi)度，μm；l為裂隙通道的長(zhǎng)度，μm；A為過(guò)水通道橫截面面積，μm2；ξ為相對(duì)粗糙度系數(shù)。

3 數(shù)值模擬結(jié)果及分析

3.1 不同JRC值時(shí)裂隙通道滲流速度場(chǎng)變化情況

根據(jù)每條不同JRC值的裂隙通道滲流模擬結(jié)果的速度場(chǎng)可知，裂隙通道表面的粗糙度對(duì)滲流速度的影響是比較明顯的，如圖2所示?？梢钥闯?，從裂隙通道的邊界處至裂隙通道的中心處，滲流速度由外向內(nèi)是在逐漸增大的，在裂隙通道的中心處速度達(dá)到最大。保持裂隙入口處的滲透壓不變，增大裂隙通道表面的粗糙度系數(shù)JRC值，結(jié)果表明，隨著JRC值的逐漸增大，裂隙通道內(nèi)部滲流過(guò)程的速度場(chǎng)會(huì)發(fā)生明顯的減小，且通過(guò)對(duì)比10條裂隙的速度場(chǎng)發(fā)現(xiàn)，裂隙表面粗糙度越大，通道內(nèi)部高流速速度場(chǎng)的寬度越窄。

圖2 不同JRC值下數(shù)值模擬宏觀速度場(chǎng)

3.1.1 裂隙通道粗糙度與滲流過(guò)程最大流速和平均流速之間的關(guān)系

JRC值對(duì)裂隙通道內(nèi)最大流速和平均流速的影響情況如圖3所示。

圖3 不同JRC值時(shí)裂隙滲流過(guò)程最大流速和

由圖3可見(jiàn)，在滲透壓為100 Pa時(shí)，隨著裂隙通道粗糙度的逐級(jí)增大，在裂隙內(nèi)部流體的最大流速、裂隙出口邊界處的最大流速和裂隙通道流出的平均流速基本都呈現(xiàn)逐漸減小的趨勢(shì)。流速關(guān)系為:裂隙內(nèi)部最大流速>裂隙出口最大流速>裂隙出口處的平均流速。其中，裂隙內(nèi)部最大流速與裂隙出口處的最大流速較為接近，但是內(nèi)部的最大流速要大于出口邊界處的最大流速。在前期JRC值較小時(shí)，裂隙內(nèi)部與出口邊界處的流速變化較為平緩，且兩者之間差距較小，在JRC值達(dá)到10～12之后，內(nèi)部的最大流速變化有一定程度的起伏，但下降程度較小，而出口邊界的最大流速則開(kāi)始隨著JRC值的增大下降幅度變得較大。裂隙出口處的平均流速與出口邊界處的最大流速變化規(guī)律相似，在JRC值較小時(shí)速度變化較為平緩，在JRC值達(dá)到 10～12之后開(kāi)始下降，與出口最大流速變化情況不同的是，出口平均流速呈階段性的減小趨勢(shì)。

3.1.2 裂隙通道粗糙度與滲流過(guò)程中渦流場(chǎng)之間的關(guān)系

不同JRC值時(shí)裂隙通道內(nèi)速度場(chǎng)等值線分布圖如圖4所示。

圖4 不同JRC值時(shí)裂隙通道內(nèi)速度場(chǎng)等值線分布圖

由圖4可以看出，在流體滲流的過(guò)程中，滲流通道內(nèi)有渦流場(chǎng)產(chǎn)生，而渦流場(chǎng)的大小及分布情況則取決于粗糙度系數(shù)JRC值的大小。在JRC值較小時(shí)，滲流通道中的渦流旋渦較為整齊，但隨著JRC值的增大，滲流通道中的渦流場(chǎng)明顯變得紊亂，且隨著粗糙度的增大，產(chǎn)生的渦流旋渦也在變大，渦流旋渦逐漸變得不規(guī)則，不再呈現(xiàn)出近似圓形的形狀。在JRC值增大的過(guò)程中，渦流場(chǎng)由通道中間逐漸向兩邊延伸，不僅僅在通道中間產(chǎn)生，在滲流通道的入口處及出口處也產(chǎn)生了渦流場(chǎng)。渦流場(chǎng)會(huì)改變滲流過(guò)程中流體運(yùn)動(dòng)方向，并且會(huì)使能量脈動(dòng)，產(chǎn)生壓力損失。模擬結(jié)果表明，流體在裂隙內(nèi)的滲流過(guò)程中，由于粗糙度的改變，在流動(dòng)過(guò)程中因流體的壓差改變而產(chǎn)生渦流，而渦流場(chǎng)的存在進(jìn)而影響了裂隙通道內(nèi)的流體速度場(chǎng)。

3.2 裂隙開(kāi)度對(duì)滲流通道各流速的影響

為了分析裂隙開(kāi)度對(duì)滲流通道各流速的影響規(guī)律，分別計(jì)算裂隙開(kāi)度在10、12、20 μm時(shí)各JRC值下的內(nèi)部最大流速、出口最大流速和出口平均流速，結(jié)果如圖5所示。

(a)裂隙開(kāi)度為10 μm

由圖5總體來(lái)看，裂隙通道內(nèi)部最大流速、出口最大流速和出口平均流速的曲線增長(zhǎng)幅度是基本接近的，在裂隙開(kāi)度增大幅度為20%時(shí)，最大流速和平均流速的增長(zhǎng)幅度接近60%；在裂隙開(kāi)度增大幅度為100%時(shí)，最大流速和平均流速的增長(zhǎng)幅度接近300%。由此可知，粗糙裂隙的流速變化程度和裂隙開(kāi)度的變化程度為3倍的正比例關(guān)系，但變化趨勢(shì)與裂隙開(kāi)度無(wú)關(guān)。

為了保證模擬試驗(yàn)的準(zhǔn)確性，以200 Pa的滲透壓用相同步驟進(jìn)行二次模擬，再將模擬得出的內(nèi)部最大流速、出口最大流速和出口平均流速與100 Pa滲透壓時(shí)的數(shù)據(jù)進(jìn)行對(duì)比，結(jié)果如圖6所示。

圖6 滲透壓分別為100 Pa和200 Pa時(shí)滲流速度對(duì)比曲線

由圖6可以看出，100 Pa和200 Pa的滲透壓下的滲透速度的變化趨勢(shì)接近，不過(guò)隨著滲透壓的增大，裂隙通道內(nèi)的流速也在增大，流速曲線所反映的趨勢(shì)拐點(diǎn)也更為明顯，但增長(zhǎng)趨勢(shì)基本不變。由此可見(jiàn)，滲透壓會(huì)改變滲流速度的大小，但不會(huì)改變滲流速度的變化趨勢(shì)。

3.3 粗糙裂隙在不同裂隙開(kāi)度和不同JRC條件下的滲透率變化情況

分別計(jì)算在裂隙開(kāi)度為10、12、20 μm時(shí)裂隙通道滲流過(guò)程中的滲透率，并且將計(jì)算結(jié)果繪制成曲線圖進(jìn)行分析，如圖7所示。

圖7 不同裂隙開(kāi)度下不同JRC值的滲透率

由圖7可以看出，當(dāng)裂隙開(kāi)度為10 μm時(shí)，在粗糙度逐漸增大的過(guò)程中，裂隙滲透率呈現(xiàn)出遞增的規(guī)律，且為標(biāo)準(zhǔn)的對(duì)數(shù)形式，擬合曲線的趨勢(shì)線方程為:y=9×10-19lnx+4×10-16，R2=1。改變裂隙的開(kāi)度至12 μm和20 μm，計(jì)算所得的滲透率變化與粗糙度的變化關(guān)系仍然保持對(duì)數(shù)函數(shù)關(guān)系，當(dāng)裂隙開(kāi)度為12 μm時(shí)，擬合曲線的趨勢(shì)線方程為:y=1×10-18lnx+6×10-16，R2=1；當(dāng)裂隙開(kāi)度為 20 μm 時(shí)，擬合曲線的趨勢(shì)線方程為:y=4×10-18lnx+2×10-15，R2=1。

3.4 裂隙通道不同JRC值及不同裂隙開(kāi)度下應(yīng)力分析

對(duì)于10級(jí)不同JRC值的粗糙度貫通裂隙，分別計(jì)算裂隙內(nèi)部的最大壓力和最小壓力，結(jié)果如圖8所示。

(a)裂隙通道內(nèi)最大壓力

由圖8可知，裂隙內(nèi)部的最大壓力隨著JRC值的增大波動(dòng)較大，并沒(méi)有呈現(xiàn)出較為明顯的遞增或遞減的規(guī)律性變化，反而裂隙內(nèi)部最小壓力隨著JRC值的增大而呈現(xiàn)出逐漸減小的趨勢(shì)，在JRC值為4～6時(shí)為最小壓力由正變負(fù)的轉(zhuǎn)折點(diǎn)，即在JRC小于4～6時(shí)，通道內(nèi)部的壓力最小值大于0 Pa，與滲流方向一致，而JRC值大于4～6時(shí)，通道內(nèi)部的壓力最小值小于0 Pa，與滲流方向相反。通過(guò)對(duì)比可以發(fā)現(xiàn)，最大壓力在JRC值為8～10時(shí)都達(dá)到了壓力值的最低點(diǎn)，最小壓力則隨著JRC值的增大持續(xù)減小。

考慮裂隙滲流過(guò)程中裂隙開(kāi)度對(duì)裂隙通道內(nèi)壓應(yīng)力的影響，從曲線圖總體趨勢(shì)來(lái)看，隨著裂隙開(kāi)度的增大，裂隙通道內(nèi)的最大壓力在增大，而裂隙通道內(nèi)的最小壓力卻是在減小。

4 結(jié)論

1)流體在粗糙裂隙通道內(nèi)的滲流過(guò)程中，裂隙內(nèi)的最大流速和出口處平均流速都隨著JRC值的增大而逐漸減小。最大流速受粗糙度的影響不明顯，而平均流速受粗糙度的影響較為明顯，平均流速下降較快且呈階段性減小趨勢(shì)。在裂隙開(kāi)度改變的同時(shí)，滲流速度的變化率是裂隙開(kāi)度的變化率的 3倍。

2)粗糙裂隙內(nèi)的滲流過(guò)程有渦流場(chǎng)的產(chǎn)生，JRC值越大，渦流圈的形狀越不規(guī)則，越來(lái)越紊亂，且渦流場(chǎng)的大小會(huì)影響滲流速度場(chǎng)的大小。

3)隨著JRC值的增大，裂隙內(nèi)的滲透率也在增大，呈現(xiàn)出典型的對(duì)數(shù)遞增形式，并且經(jīng)過(guò)對(duì)比，滲透率的增長(zhǎng)只與粗糙度有關(guān)，與裂隙開(kāi)度無(wú)關(guān)。

4)滲流過(guò)程中裂隙內(nèi)最大壓力受粗糙度的影響起伏波動(dòng)較大，但在JRC值為8～10時(shí)會(huì)達(dá)到最小值；裂隙內(nèi)的最小壓力會(huì)隨JRC值的增大而減小。裂隙開(kāi)度增大時(shí)，裂隙內(nèi)最大壓力會(huì)增大，而最小壓力卻會(huì)減小。