張玉柱

（中煤科工集團(tuán)重慶研究院有限公司，重慶 400037）

煤作為一種非均質(zhì)多孔介質(zhì)材料，其內(nèi)部包含大量的孔隙、裂隙結(jié)構(gòu)，這些結(jié)構(gòu)對(duì)煤層氣的存儲(chǔ)以及運(yùn)移有著重要影響[1]。一般認(rèn)為，煤基質(zhì)內(nèi)的孔隙主要作為存儲(chǔ)介質(zhì)，在外部壓力作用下，孔隙內(nèi)自由氣體會(huì)流向裂隙介質(zhì)。裂隙網(wǎng)絡(luò)為氣體的主要滲流通道，氣體的流動(dòng)主要發(fā)生在裂隙介質(zhì)中，煤層中的自然發(fā)育的裂隙控制著煤儲(chǔ)層的滲透性[2-3]。然而煤中裂隙的分布往往是無(wú)序且不規(guī)則的，其特征難以詳細(xì)描述，但其空間分布規(guī)律可能遵循一些宏觀統(tǒng)計(jì)規(guī)律。此外，由于煤的非均質(zhì)性較強(qiáng)，孔隙結(jié)構(gòu)和裂隙結(jié)構(gòu)對(duì)煤層氣的貢獻(xiàn)也不相同，如何定義孔隙、裂隙的物理屬性對(duì)研究煤層氣流動(dòng)特性至關(guān)重要。近年來(lái)，隨著圖像處理技術(shù)的發(fā)展，采用數(shù)字圖像研究煤巖體材料的非均質(zhì)特性成為一種簡(jiǎn)單有效的方法，可以通過(guò)不同的像素值來(lái)反映不同結(jié)構(gòu)的特征[4-6]。因此，基于Monte Carlo方法生成隨機(jī)裂隙網(wǎng)絡(luò)圖像，采用圖像處理軟件將其轉(zhuǎn)化為二值化圖像，再導(dǎo)入COMSOL中，通過(guò)圖像函數(shù)功能再將其轉(zhuǎn)化為計(jì)算域，并定義了孔、裂隙結(jié)構(gòu)的孔隙度和滲透率函數(shù)，建立了煤層瓦斯流動(dòng)的數(shù)值計(jì)算模型，進(jìn)一步分析了不同裂隙結(jié)構(gòu)特征對(duì)煤層滲透率的影響。

1 建立隨機(jī)裂隙網(wǎng)絡(luò)模型

實(shí)際煤層中含有大量幾何結(jié)構(gòu)復(fù)雜的裂隙網(wǎng)絡(luò)，詳細(xì)描述實(shí)際煤中裂隙的幾何特征是極為困難的，但它們的空間分布形態(tài)往往具有宏觀統(tǒng)計(jì)特征。因此，可以采用某種概率統(tǒng)計(jì)密度函數(shù)來(lái)描述其幾何特征參數(shù)（裂隙產(chǎn)狀、密度、裂隙跡長(zhǎng)和開(kāi)度）。目前，采用Monte Carlo法建立離散裂隙網(wǎng)絡(luò)模型是應(yīng)用最為廣泛、也是公認(rèn)最為準(zhǔn)確的一種隨機(jī)建模方法[7]，它能夠較好的模擬煤層的非均質(zhì)特性。因此，采用該方法建立煤層離散裂隙網(wǎng)絡(luò)模型。裂隙幾何參數(shù)見(jiàn)表1。隨機(jī)裂隙網(wǎng)絡(luò)模型如圖1，模型尺寸為30 m×30 m。

表1 裂隙幾何參數(shù)Table 1 Fracture geometric parameters

圖1 隨機(jī)裂隙網(wǎng)絡(luò)模型Fig.1 Stochastic fracture network model

2 煤層瓦斯流動(dòng)的數(shù)值模擬

煤體為雙重孔隙介質(zhì)，主要由孔隙和裂隙構(gòu)成，氣體在煤層中的流動(dòng)是1個(gè)極其復(fù)雜的過(guò)程，研究假設(shè)氣體在煤層中的流動(dòng)為單向流動(dòng)，符合達(dá)西定律[8-12]；不考慮氣體的吸附和解吸特性，僅研究氣體在孔裂隙中的流動(dòng)特性。

2.1 COMSOL中的圖像函數(shù)功能

在COMSOL中可以使用圖像數(shù)據(jù)來(lái)表示二維材質(zhì)分布情況，通過(guò)顏色或灰度來(lái)標(biāo)識(shí)具有不同材質(zhì)的區(qū)域。以這種方式使用的圖像可以有許多來(lái)源，例如掃描電子顯微鏡（SEM）、計(jì)算機(jī)斷層掃描（CT）或磁共振成像等，導(dǎo)入的圖像作為通用的COMSOL插值函數(shù)，可用于任何建模目的。

研究利用COMSOL中的圖像函數(shù)功能，將處理過(guò)后的二值化裂隙網(wǎng)絡(luò)圖像（圖1）導(dǎo)入COMSOL中進(jìn)行數(shù)值建模。將導(dǎo)入的圖像函數(shù)定義為im（x，y），像素值范圍從0到1。通過(guò)COMSOL中圖像函數(shù)功能導(dǎo)入的裂隙網(wǎng)絡(luò)圖像如圖2，其中紅色部分代表煤基質(zhì)，藍(lán)色部分為裂隙。

圖2 通過(guò)COMSOL中圖像函數(shù)功能導(dǎo)入的裂隙網(wǎng)絡(luò)圖像Fig.2 Fracture network image imported by image function in COMSOL

由于煤的非均質(zhì)特性，其煤基質(zhì)與裂隙物理特性差異較大。因此，為了從圖像顏色代碼定義煤體物理特性，對(duì)孔隙度φ和滲透率K實(shí)現(xiàn)以下關(guān)系：

式中：λm為煤基質(zhì)所對(duì)應(yīng)的像素值；λf為裂隙所對(duì)應(yīng)的像素值；φm為煤基質(zhì)的孔隙度，φf(shuō)為裂隙的孔隙度；Km為煤基質(zhì)的孔隙滲透率；Kf為裂隙滲透率。

2.2 煤層氣流動(dòng)控制方程

假定氣體在煤層中的流動(dòng)是定場(chǎng)和等溫的，忽略氣體的壓縮特性，則氣體在孔、裂隙介質(zhì)中的流動(dòng)應(yīng)滿足Navier-Stokes方程[7]，即：

式中：μ為流體動(dòng)力黏度，瓦斯動(dòng)力黏度取11.067μPa·s；犖P為氣體壓力梯度；V為流動(dòng)速度。

煤體為雙重孔隙介質(zhì)模型，氣體主要在煤體中孔隙以及裂隙中流動(dòng)，因此煤層滲透率K可描述為：

流體在裂隙中的流動(dòng)方程常采用立方定量描述：

式中：b為裂隙開(kāi)度。

2.3 模型邊界條件

研究設(shè)定的模型邊界條件如圖3。

圖3 模型邊界條件Fig.3 Model boundary conditions

模型尺寸為30 mm×30 mm，模型左側(cè)為進(jìn)氣口壓力p1=0.1 MPa；右側(cè)出口處壓力為p0=0 MPa，上下邊界為不滲透邊界。同時(shí)，為了保證計(jì)算精度和節(jié)約計(jì)算時(shí)間，采用自適應(yīng)三角形網(wǎng)格對(duì)裂隙區(qū)域進(jìn)行細(xì)化處理。

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

結(jié)合邊界條件以及相關(guān)數(shù)值模擬參數(shù)，通過(guò)COMSOL中的達(dá)西穩(wěn)態(tài)流求解方程（式（3）），得到了裂隙網(wǎng)絡(luò)中氣體壓力和速度的分布特征，數(shù)值模擬結(jié)果如圖4。

圖4 數(shù)值模擬結(jié)果Fig.4 Numerical simulation results

從圖4可明顯看出，氣體壓力從入口到出口逐漸降低，氣體在裂隙中的流動(dòng)速度要遠(yuǎn)大于基質(zhì)滲流速度。此結(jié)果也表明煤層滲透率主要由裂隙網(wǎng)絡(luò)所控制，裂隙對(duì)煤層滲透率的貢獻(xiàn)要遠(yuǎn)遠(yuǎn)大于基質(zhì)孔隙的滲透率，這也與實(shí)際情況一致。從模擬結(jié)果來(lái)看，基于裂隙網(wǎng)絡(luò)圖像來(lái)模擬煤層氣流動(dòng)是一種可行的方法。

3.1 裂隙長(zhǎng)度

為探究裂隙長(zhǎng)度對(duì)煤層滲透率的影響，通過(guò)Monte Carlo法建立不同裂隙長(zhǎng)度離散裂隙網(wǎng)絡(luò)模型，其它幾何特征參數(shù)不變，裂隙長(zhǎng)度變化范圍2～13 m。裂隙長(zhǎng)度與煤層滲透率的關(guān)系如圖5。

由圖5可知，隨著裂隙長(zhǎng)度的增加，滲透率也隨之增大，且裂隙長(zhǎng)度與滲透率之間呈現(xiàn)較好的冪律關(guān)系。這是因?yàn)榱严堕L(zhǎng)度越大，更容易連接其它裂隙，煤層內(nèi)裂隙網(wǎng)絡(luò)的連通性越好，煤層滲透性越好。

圖5 裂隙長(zhǎng)度與滲透率的關(guān)系Fig.5 Relationship between fracture length and permeability

3.2 裂隙密度

在二維空間中，裂隙密度為研究區(qū)域內(nèi)裂隙數(shù)量與面積的比值，它可以反映裂隙分布的復(fù)雜程度。為探究裂隙密度的改變對(duì)煤層滲透性的影響，固定其它幾何參數(shù)，設(shè)定不同的密度值，采用Monte Carlo方法采用生成裂隙網(wǎng)絡(luò)模型，密度變化范圍為0.4～2.0條/m2。裂隙密度與煤層滲透率的關(guān)系如圖6。

圖6 裂隙密度與滲透率的關(guān)系Fig.6 Relationship between fracture density and permeability

由圖6可知，隨著裂隙密度的增加，煤層滲透率也隨之增大，且裂隙密度與滲透率之間也呈現(xiàn)較好的冪律關(guān)系。其原因在于，裂隙密度的增加意味著煤層內(nèi)裂隙數(shù)量增加，越容易形成更多的滲流通道，煤層的滲透性顯著提高。

3.3 裂隙開(kāi)度

由式（6）可知，裂隙滲透率主要與裂隙開(kāi)度有關(guān)。因此，為分析裂隙開(kāi)度變化對(duì)煤層滲透率的影響，設(shè)定裂隙開(kāi)度的變化范圍為0.1～3 mm。裂隙開(kāi)度與滲透率的關(guān)系如圖7。

圖7 裂隙開(kāi)度與滲透率的關(guān)系Fig.7 Relationship between fracture opening and permeability

由圖7可知，煤層滲透率隨裂隙開(kāi)度的增大而增大，且裂隙開(kāi)度與滲透率之間同樣符合冪律增長(zhǎng)關(guān)系。對(duì)比裂隙長(zhǎng)度、密度與滲透率的關(guān)系（圖5、圖6），裂隙開(kāi)度的增加導(dǎo)致滲透率值迅速增大，其影響要遠(yuǎn)大于裂隙長(zhǎng)度以及密度對(duì)煤層滲透率的影響。

3.4 裂隙傾角

裂隙傾角為煤層氣流動(dòng)方向與裂隙之間的傾角。為分析裂隙傾角的變化對(duì)煤層滲透性的影響，研究設(shè)定傾角的變化范圍為10°～90°。裂隙傾角與煤層滲透率的關(guān)系如圖8。從圖8可以看出，當(dāng)裂隙傾角較小時(shí)，煤層滲透率值最大，隨著裂隙傾角的增大，滲透率值隨之降低，由圖中擬合曲線可知，裂隙傾角與滲透之間呈現(xiàn)較好冪律增長(zhǎng)關(guān)系。其原因在于，裂隙是氣體運(yùn)移的主要通道，裂隙傾角越大，氣體流動(dòng)阻力也越大，因此煤層滲透率值也越小。

圖8 裂隙傾角與滲透率的關(guān)系Fig.8 Relationship between fracture dip angle and permeability

4 煤層氣抽采模擬

基于裂隙網(wǎng)絡(luò)圖像不僅能研究煤層氣流動(dòng)特性，同樣也可擴(kuò)展到工程應(yīng)用方面，例如圖像函數(shù)在瓦斯抽采方面的應(yīng)用?；贛onte Carlo方法建立了隨機(jī)裂隙網(wǎng)絡(luò)模型，模型尺寸為30 m×30 m，然后將建立的裂隙網(wǎng)絡(luò)圖像導(dǎo)入COMSOL中建立計(jì)算域，瓦斯抽采數(shù)值模型如圖9。

圖9 瓦斯抽采數(shù)值模型Fig.9 Numerical model of gas extraction

采用自由三角形網(wǎng)格劃分計(jì)算域，抽采孔直徑0.1 m，初始煤層瓦斯壓力為1 MPa，抽采孔壓力為標(biāo)準(zhǔn)大氣壓力0.1 MPa。

實(shí)際煤儲(chǔ)層中氣體質(zhì)量含量包括裂隙系統(tǒng)中的自由相氣體、吸收氣體以及基質(zhì)孔隙中的吸附氣體，氣體質(zhì)量含量mg可表示為[13-14]：

式中：Mg為甲烷氣體分子量，取16.04 g/mol；p為氣體壓力，MPa；T為煤層溫度，取310 K；R為氣體常數(shù)，取52.9 J/（mol·K）；ρgs為標(biāo)準(zhǔn)狀態(tài)下氣體密度，取0.716 kg/m3；ρs為煤體骨架密度，取1 350 kg/m3；VL朗繆爾體積常數(shù)，取0.036 m3/kg；pL是朗繆爾壓力常數(shù)，取3.034 MPa。

煤層瓦斯抽采過(guò)程中，氣體運(yùn)移擴(kuò)散是1個(gè)極為緩慢的過(guò)程。在氣體壓力影響下，基質(zhì)中的氣體逐漸向裂隙內(nèi)擴(kuò)散，裂隙網(wǎng)絡(luò)中各氣體組分的流動(dòng)應(yīng)滿足質(zhì)量守恒定律[13-14]：

式中：ρg為氣體密度，取0.648 kg/m3；t為時(shí)間。

不同時(shí)刻抽采孔周圍氣體壓力分布如圖10。從圖10可以看出，隨著抽采時(shí)間的推移，抽采孔周圍的氣體壓力逐漸降低。壓力在鉆孔附近迅速下降，距離抽采孔較遠(yuǎn)的地方，壓力下降速度較為緩慢。

圖10 不同時(shí)刻抽采孔周圍氣體壓力分布Fig.10 Gas pressure distribution around the extraction hole at different times

不同時(shí)間抽采孔周圍氣體速度分布如圖11。從圖11可以看出，裂隙中的氣體速度遠(yuǎn)大于基質(zhì)中的氣體速度，流速隨著氣體壓力的減小而降低。該模擬結(jié)果顯示了瓦斯抽采過(guò)程中煤層氣的流動(dòng)特征。

圖11 不同時(shí)間抽采孔周圍氣體速度分布Fig.11 Gas velocity distribution around the extraction hole at different times

5 結(jié)語(yǔ)

1）煤基質(zhì)滲透率遠(yuǎn)小于裂隙網(wǎng)絡(luò)滲透率，煤層滲透率主要為裂隙網(wǎng)絡(luò)所控制。

2）煤層滲透率隨著裂隙長(zhǎng)度、密度以及開(kāi)度的增加而增大，隨裂隙傾角的減小而增大；相比于裂隙長(zhǎng)度、密度，裂隙開(kāi)度對(duì)滲透率的影響最為顯著。

3）基于裂隙網(wǎng)絡(luò)圖像模擬瓦斯抽采是一種簡(jiǎn)單、可行的方法，可以反映整個(gè)抽采過(guò)程中各階段的煤層氣流動(dòng)特性，包括瓦斯流動(dòng)速度和壓力的時(shí)空演化、瓦斯流動(dòng)速度和壓力在裂隙與煤基質(zhì)之間的差異。

4）研究同樣存在一些局限，只是對(duì)非均質(zhì)煤體的滲流過(guò)程進(jìn)行示例研究，并未考慮應(yīng)力場(chǎng)、溫度場(chǎng)的對(duì)煤層瓦斯?jié)B流的影響，多場(chǎng)作用下煤層瓦斯的流動(dòng)特性的研究將是下一步工作的重點(diǎn)。