李 鵬，戴廣龍，葉慶樹，唐明云，聶士斌，宋小林

（1.神東煤炭集團(tuán)有限責(zé)任公司，內(nèi)蒙古 鄂爾多斯017000；2.安徽理工大學(xué) 能源與安全學(xué)院，安徽 淮南232001；3.中國(guó)礦業(yè)大學(xué) 安全工程學(xué)院，江蘇 徐州221116）

我國(guó)西北煤炭資源具有埋藏淺、基巖薄、煤層群間距近的特點(diǎn)，煤層開采時(shí)，極易形成貫通工作面、采空區(qū)和地表之間的漏風(fēng)通道，從而容易造成采空區(qū)內(nèi)遺煤自燃。因此，研究淺埋煤層采空區(qū)的漏風(fēng)規(guī)律具有重要意義。鄔劍明[1]等利用SF6瞬時(shí)釋放技術(shù)測(cè)定了大柳塔礦活雞兔井21304 工作面井上下漏風(fēng)規(guī)律；王建文[2]等通過示蹤氣體SF6測(cè)定了檸條塔煤礦N1201 工作面地表漏風(fēng)通道，得出采礦塌陷形成裂隙帶是漏風(fēng)的主要通道；徐會(huì)軍[3]等采用數(shù)值模擬方法模擬了封堵上、下隅角和注氮對(duì)試驗(yàn)工作面采空區(qū)漏風(fēng)和氧濃度的影響；趙啟峰[4]等通過采用“地表氣壓與井下采空區(qū)氣壓差、工作面進(jìn)風(fēng)巷與回風(fēng)巷風(fēng)量差、示蹤氣體地面瞬時(shí)釋放法”確定了平朔礦區(qū)淺埋綜放面地表漏風(fēng)的各參數(shù)；賀飛[5]等借助Fluent 數(shù)值模擬軟件對(duì)酸刺溝煤礦6下109工作面采空區(qū)流場(chǎng)進(jìn)行了模擬，得出地表漏風(fēng)逐漸向下部工作面回風(fēng)側(cè)采空區(qū)偏移。從以上文獻(xiàn)可知，眾多學(xué)者在采用數(shù)值模擬手段分析淺埋煤層工作面采空區(qū)漏風(fēng)規(guī)律時(shí)，大多是研究抽出式通風(fēng)方式下的漏風(fēng)規(guī)律，但對(duì)壓入式通風(fēng)方式下的漏風(fēng)規(guī)律研究較少，而且對(duì)于從空氣能位角度分析工作面與地表的空氣能位差對(duì)采空區(qū)漏風(fēng)的影響研究也較少。為此以采用壓入式通風(fēng)方式的活雞兔井12下206 工作面為研究對(duì)象，基于采空區(qū)“O”型垮落壓實(shí)[6]和Bachmat 提出的非線性滲流理論[7]，并結(jié)合采空區(qū)“豎三帶”計(jì)算及現(xiàn)場(chǎng)地面裂隙觀測(cè)結(jié)果，建立了三維滲透率分布模型，利用專業(yè)的流體力學(xué)軟件Fluent 對(duì)12下206 工作面采空區(qū)流場(chǎng)進(jìn)行模擬分析，并在此基礎(chǔ)上研究了工作面與地表之間能位差對(duì)采空區(qū)漏風(fēng)的影響。

1 工作面概況

12下206 工作面位于神東大柳塔煤礦活雞兔井的12 煤復(fù)合區(qū)下分層，12 煤采用分層開采。東側(cè)為12 煤復(fù)合區(qū)的火燒邊界，西北始于12 煤南側(cè)輔運(yùn)巷，南側(cè)為12下208 綜放工作面，北側(cè)為12下204采空區(qū)，正上方為12上206 采空區(qū)。該工作面采用放頂煤工藝開采，煤層層間距為直接頂，距離上覆采空區(qū)層間距為0.8～3.8 m，平均層間距2.5 m。正常回采期間工作面煤機(jī)采動(dòng)采高設(shè)計(jì)3.7 m（±0.2 m），放煤高度2.3 m，整體工作面采煤高度6 m。工作面采用“U”型壓入式通風(fēng)方式，供風(fēng)量為1 100 m3/min。

2 漏風(fēng)數(shù)值模擬

2.1 數(shù)學(xué)模型

2.1.1 采場(chǎng)流體守恒方程

根據(jù)流體的質(zhì)量和動(dòng)量守恒，建立的采場(chǎng)流體數(shù)學(xué)模型如式（1）[7-8]。

式中：ρ 為流體的密度，kg/m3；τ 為時(shí)間，s； xi為i 方向上坐標(biāo)變量，m；ui、uj為i、j 方向上流體的速度，m/s；qm為流體質(zhì)量源項(xiàng)，kg/（m3·s）；p 為流體的壓力，Pa；τij為流體的剪應(yīng)力，Pa；ρgi為流體的重力，N/m3；Fi為流體的動(dòng)量源項(xiàng)，N/m3；在采空區(qū)內(nèi)，流體屬于層流或者過渡流，F(xiàn)i可用式（2）進(jìn)行求解。

式中：μ 為流體的動(dòng)力黏度，Pa·s；k 為采空區(qū)內(nèi)滲透率，m2；C2為多孔介質(zhì)內(nèi)部阻力系數(shù)，1/m；｜u｜為i 方向速度u 的標(biāo)量值，m/s。

2.1.2 采空區(qū)漏風(fēng)阻力系數(shù)模型

多孔介質(zhì)非線性滲流運(yùn)動(dòng)方程[9]：

式中：▽h 為壓力梯度，無(wú)因次；a 為黏性力系數(shù)，s/m；b 為慣性力系數(shù)，s2/m2；v 為運(yùn)動(dòng)黏度，m2/s；g為重力加速度9.81 m/s2；n 為孔隙率，無(wú)因次；β 為介質(zhì)顆粒形狀系數(shù)，無(wú)因次，取1.5；Dm為平均調(diào)和粒徑，m，需要根據(jù)工作面實(shí)際開采條件進(jìn)行反求，在本研究中通過模擬反求[7]。

滲透率k 是孔隙率n 的函數(shù)，可用Blake-Kozeny 公式計(jì)算：

孔隙率的大小與采空區(qū)垮落巖石的碎脹系數(shù)有關(guān)，因此可通過采空區(qū)垮落巖石的碎脹系數(shù)Kp得出：

根據(jù)一般的礦壓觀測(cè)規(guī)律，在采空區(qū)內(nèi)垮落和壓實(shí)之間的碎脹系數(shù)Kp近似呈負(fù)指數(shù)衰減規(guī)律變化。為了使得采空區(qū)滲透率模型更能反映實(shí)際，碎脹系數(shù)分布模型除了表征碎脹系數(shù)在x，y 方向上的不同，還應(yīng)體現(xiàn)其在垂直方向上的不同，本次建立的碎脹系數(shù)模型如式（6）[10]。

式中：Kp0為初始垮落碎脹系數(shù)，無(wú)因次；Kp1為壓實(shí)后的碎脹系數(shù)，無(wú)因次；Kp，i（x，y）為采空區(qū)平面上的碎脹系數(shù)分布，無(wú)因次，i 取0 和1；a1、a2、a3分別為x、y、z 方向上的衰減率，m－1；d1、d2、d3分別為采空區(qū)坐標(biāo)（x、y、z）處距x、y、z 邊界的距離，m。

為了較為準(zhǔn)確地描述采空區(qū)在垂直方向上的碎脹系數(shù)分布規(guī)律，主要以采空區(qū)“豎三帶”為依據(jù)來建立12下206 工作面采空區(qū)碎脹系數(shù)分布模型。根據(jù)現(xiàn)場(chǎng)資料，12 煤頂板屬于中硬巖，結(jié)合相關(guān)中硬巖下垮落帶[11]與導(dǎo)水?dāng)嗔褞в?jì)算公式，可以計(jì)算得出12下206 工作面后方采空區(qū)垮落帶和導(dǎo)水?dāng)嗔褞У母叨?，具體的“豎三帶”[12]分布示意如圖1。

圖1 12 下206 采空區(qū)“豎三帶”分布圖

按照綜放工作面采空區(qū)礦壓顯現(xiàn)的一般規(guī)律，并結(jié)合現(xiàn)場(chǎng)對(duì)地面裂隙觀測(cè)[13]情況分析，對(duì)采空區(qū)內(nèi)不同區(qū)域的初始碎脹系數(shù)、壓實(shí)后的碎脹系數(shù)及衰減率取值（表1）。將表1 取值代入式（6）中，即可獲得采空區(qū)不同位置的碎脹系數(shù)分布模型[7]。

表1 采空區(qū)內(nèi)不同區(qū)域的碎脹系數(shù)模型參數(shù)取值

結(jié)合式（5）和式（4），可得出采空區(qū)不同位置的滲透率分布模型，綜合式（3），推導(dǎo)得到Fluent 的黏性阻力系數(shù)C1和慣性阻力系數(shù)C2：

最后，利用自定義函數(shù)UDF 進(jìn)行編譯實(shí)現(xiàn)[14-15]。

2.2 物理模型和基本假設(shè)

1）物理模型。12下206 工作面距離地表為90 m，工作面進(jìn)、回巷道寬×高為5.4 m×3.5 m。結(jié)合12下206 工作面的實(shí)際情況，對(duì)物理模型進(jìn)行簡(jiǎn)化，建立了三維采場(chǎng)物理模型，模型參數(shù)見表2。

2）基本假設(shè)。①采場(chǎng)流體為不可壓流體；②采空區(qū)內(nèi)介質(zhì)具有各向同性。

表2 12 下206 采場(chǎng)物理模型參數(shù)

2.3 漏風(fēng)模擬結(jié)果及分析

2.3.1 模擬參數(shù)設(shè)置

由于地表裂隙多，分布規(guī)律不明顯，且漏風(fēng)速率無(wú)法確定，因此將地表設(shè)置為壓力出口，工作面進(jìn)風(fēng)口和回風(fēng)口分別設(shè)置為壓力入口和壓力出口。采場(chǎng)網(wǎng)格模型及邊界條件如圖2。

圖2 12 下206 采場(chǎng)網(wǎng)格模型及邊界條件

根據(jù)邊界條件設(shè)置，需要現(xiàn)場(chǎng)測(cè)試工作面與地表之間的能位值。由于工作面進(jìn)風(fēng)隅角和回風(fēng)隅角的風(fēng)流速度為0.9 m/s，則計(jì)算得到的空氣動(dòng)能小于1 Pa，因此在計(jì)算工作面進(jìn)、回風(fēng)隅角的能位時(shí)，忽略空氣動(dòng)能的影響。根據(jù)井下對(duì)工作面進(jìn)、回風(fēng)巷和采空區(qū)對(duì)應(yīng)地表進(jìn)行的能位測(cè)試，12下206 工作面能位測(cè)試見表3。

表3 12 下206 工作面能位測(cè)試

2.3.2 模擬結(jié)果及分析

按照上述模擬參數(shù)，對(duì)采空區(qū)進(jìn)行模擬。進(jìn)、回風(fēng)巷與地表的漏風(fēng)量大小可通過改變式（4）中平均粒徑Dm的大小來調(diào)整，以此模擬得到與現(xiàn)場(chǎng)實(shí)際基本吻合的漏風(fēng)量。通過不斷調(diào)整，當(dāng)Dm=0.03 m 時(shí)，模擬得到工作面的進(jìn)風(fēng)為1 092 m3/min，回風(fēng)為966m3/min，從工作面漏至地表的風(fēng)量為126 m3/min，與現(xiàn)場(chǎng)實(shí)測(cè)的采空區(qū)漏風(fēng)量119 m3/min 基本吻合，模擬得到的采場(chǎng)壓力分布如圖3。由圖3 可得，采空區(qū)垂直方向上壓力呈層狀分布，最底部壓力最大，隨著高度增加，壓力逐漸減少，對(duì)應(yīng)地面氣壓為最小，這與地表大氣壓分布趨勢(shì)一致。

圖3 12 下206 采場(chǎng)壓力分布圖

為直觀分析12下206 工作面采空區(qū)平面上的壓力分布，繪制了z=2 m 剖面的壓力分布圖（圖4）。從圖4 可以看出，由于工作面通風(fēng)阻力的影響，使得在工作面進(jìn)風(fēng)側(cè)往回風(fēng)側(cè)方向上的壓力逐漸減小。由于該工作面屬于正壓通風(fēng)，工作面能位高于采空區(qū)，使得在采空區(qū)內(nèi)x 方向上，隨著x 增加，壓力逐漸減小，因此，在氣體壓力的作用下，工作面部分風(fēng)量逐漸流入采空區(qū)深部，從圖中的流線圖可以明顯的看出采空區(qū)內(nèi)的氣體流動(dòng)方向。

圖4 12 下206 采場(chǎng)z=2 m 剖面壓力分布圖

根據(jù)模擬結(jié)果，得到12下206 采場(chǎng)不同剖面的流場(chǎng)分布圖，如圖5 和圖6。根據(jù)圖5 可知，在12下206 工作面上的風(fēng)速最大，最大風(fēng)速達(dá)到1.0 m/s，采空區(qū)內(nèi)漏風(fēng)流場(chǎng)近似呈“O”型分布，越靠近采空區(qū)邊界，流速越大，越靠近采空區(qū)中部，流速越小?？拷ぷ髅?zhèn)鹊牟煽諈^(qū)風(fēng)流速度大的區(qū)域明顯大于采空區(qū)上下“兩道”，它也是由于工作面?zhèn)鹊乃槊浵禂?shù)的衰減率更小所致。從圖6 可知，垂直方向上不同高度的采空區(qū)流場(chǎng)分布圖也是基本一致。工作面后方采空區(qū)，漏風(fēng)流速大，且區(qū)域?qū)?，說明靠近工作面附近的采空區(qū)是漏風(fēng)的主要區(qū)域。由x=70 m 和y=70 m 剖面流速等值線看出，隨著垂直高度的增加，采空區(qū)漏風(fēng)的主要區(qū)域逐漸變窄，而且均出現(xiàn)了突變現(xiàn)象，在垂直方向上，整體漏風(fēng)流場(chǎng)近似呈現(xiàn)“浴盆”狀，這主要是因?yàn)?2下206 采空區(qū)的垮落帶的滲透率比裂隙帶更大所致。從圖中x=70 m 和y=70 m 剖面的流線方向可知，采空區(qū)內(nèi)的風(fēng)流方向是由采空區(qū)底部流向地表的，這與地表漏風(fēng)測(cè)試的結(jié)果也是基本吻合。

圖5 采場(chǎng)z=2 m 剖面流場(chǎng)分布云圖

圖6 采場(chǎng)x=70 m、y=70 m 和z=2 m 剖面流場(chǎng)分布圖

3 工作面與地表能位差對(duì)采空區(qū)漏風(fēng)影響

為了研究空氣不同能位差對(duì)工作面采空區(qū)漏風(fēng)的影響，本次模擬假定采空區(qū)內(nèi)空氣位能不變，即為0.917 kPa，通過改變工作面與地表氣壓差來獲得工作面與地表之間的能位差。擬定了4 種不同的工作面與地表之間的能位差，模擬結(jié)果見表4。根據(jù)表4得到12下206 工作面與地表之間不同能位差下的漏風(fēng)量變化規(guī)律（圖7）。

同樣根據(jù)圖7，擬合得到能位差與漏風(fēng)量之間的關(guān)系函數(shù)，如式（8）：

式中：QL為12下206 工作面采空區(qū)漏風(fēng)量，m3/min；△pn為12下206 工作面與地表之間的能位差，kPa。

表4 不同能位差下的地表漏風(fēng)量分布

圖7 能位差與采空區(qū)漏風(fēng)量變化規(guī)律

從圖7 可知，隨著12下206 工作面與地表的能位差增大，工作面流至采空區(qū)的漏風(fēng)量也增大。當(dāng)工作面與地表壓差為0.937 kPa 時(shí)，該壓差與采空區(qū)內(nèi)的空氣柱位能近似相等，此時(shí)工作面與地表之間的能位差為0.02 2 kPa，該條件下的漏風(fēng)量為44 m3/min，漏風(fēng)量為最小。根據(jù)圖7，可以擬合得到工作面采空區(qū)漏風(fēng)量與能位差之間呈指數(shù)函數(shù)關(guān)系。

從以上分析可知，除了及時(shí)對(duì)地表裂隙進(jìn)行回填以增大采空區(qū)漏風(fēng)風(fēng)阻外，還可以采用均壓措施實(shí)時(shí)調(diào)整井下能位，盡量使得它與地表能位相等，以減少采空區(qū)漏風(fēng)。

4 結(jié) 論

1）壓入式通風(fēng)下淺埋煤層工作面采空區(qū)內(nèi)漏風(fēng)流場(chǎng)近似呈“O”型分布，靠近工作面附近的采空區(qū)是漏風(fēng)的主要區(qū)域。

2）采空區(qū)主要漏風(fēng)區(qū)域隨垂直高度增加而逐漸變窄，整體漏風(fēng)流場(chǎng)近似呈“浴盆”狀。

3）壓入式通風(fēng)方式下，工作面與地表之間的能位差決定了工作面采空區(qū)與地表之間的漏風(fēng)量大小，能位差越大，漏風(fēng)量也越大；工作面與地表之間能位差與采空區(qū)漏風(fēng)量呈指數(shù)函數(shù)關(guān)系。