大面積采空區(qū)煤自燃環(huán)境下氣體運(yùn)移規(guī)律研究

張鵬宇，陳曉坤，趙 亮，魏高明，王 奇

（1.西安科技大學(xué) 安全科學(xué)與工程學(xué)院，陜西 西安 710054；2.西安科技大學(xué) 煤火災(zāi)害防治陜西省重點(diǎn)實(shí)驗(yàn)室，陜西 西安 710054；3.兗礦能源集團(tuán)公司 鮑店煤礦，山東 鄒城 273500）

隨著現(xiàn)代化礦井開采強(qiáng)度與機(jī)械化程度不斷提升，礦井采空區(qū)面積和范圍隨之增加，形成大面積采空區(qū)[1]。在工作面正常推采期間，采掘活動(dòng)以及頂板周期來(lái)壓使得沿空側(cè)預(yù)留煤柱壓酥破碎，形成大量裂縫，導(dǎo)致采空區(qū)之間相互貫通，為煤的氧化升溫提供漏風(fēng)通道[2]。隨著大面積采空區(qū)持續(xù)供氧，采空區(qū)遺煤長(zhǎng)時(shí)間被氧化，增大了采空區(qū)煤自燃危險(xiǎn)性。采空區(qū)煤自燃產(chǎn)生如CH4、CO2、CO、C2H6、C2H4、C3H8、H2、SO2、H2S 等有害氣體，會(huì)在漏風(fēng)通道內(nèi)運(yùn)移并涌出[3-5]。因此，采空區(qū)煤自燃環(huán)境作用下氣體非控運(yùn)移有可能引發(fā)瓦斯、煤塵爆炸和礦井通風(fēng)系統(tǒng)紊亂等一系列嚴(yán)重的次生災(zāi)害[6-7]。采空區(qū)煤-氣-熱共生問(wèn)題日益嚴(yán)峻，大面積采空區(qū)流場(chǎng)特性和有害氣體運(yùn)移機(jī)制亟待進(jìn)一步深入研究。

國(guó)內(nèi)外學(xué)者對(duì)采空區(qū)流場(chǎng)特征進(jìn)行了大量研究，主要采用現(xiàn)場(chǎng)觀測(cè)、相似模擬實(shí)驗(yàn)和數(shù)值模擬方法[8-10]。由于采空區(qū)條件復(fù)雜，人員無(wú)法進(jìn)入其內(nèi)部，采空區(qū)隱蔽災(zāi)害信息判斷和氣體流場(chǎng)探測(cè)極為困難[11]。物理相似模擬實(shí)驗(yàn)雖然在一定程度上解決了現(xiàn)場(chǎng)對(duì)于采空區(qū)垮落環(huán)境下束管監(jiān)測(cè)系統(tǒng)監(jiān)測(cè)點(diǎn)位破壞和數(shù)據(jù)局限性等問(wèn)題，但在實(shí)驗(yàn)中不可避免產(chǎn)生的誤差，無(wú)法在全域條件下研究采空區(qū)時(shí)空演化規(guī)律。而數(shù)值模擬可以彌補(bǔ)這一缺陷，并進(jìn)一步為物理相似模擬實(shí)驗(yàn)的研究提供驗(yàn)證依據(jù)，為現(xiàn)場(chǎng)工程應(yīng)用的可行性提供指導(dǎo)[12-13]。許多研究人員已經(jīng)建立了多物理場(chǎng)耦合模型，旨在分析采空區(qū)中的氣體滲流[14-15]、氣體濃度場(chǎng)[16-17]、以及溫度場(chǎng)[18-19]變化規(guī)律，這些研究為采空區(qū)模擬復(fù)雜流場(chǎng)和氣體運(yùn)移提供了必要的理論基礎(chǔ)。由于已建立數(shù)值模型通常基于單一工作面采空區(qū)的漏風(fēng)規(guī)律進(jìn)行研究，缺乏考慮采空區(qū)連通后氣體流場(chǎng)的運(yùn)移規(guī)律；此外，忽略了由于采空區(qū)內(nèi)部氣體密度不同，煤自燃溫度導(dǎo)致采空區(qū)壓力以及氣體滲流和擴(kuò)散能力的改變。為此，針對(duì)大面積采空區(qū)建立了一個(gè)涉及多孔介質(zhì)的氣體運(yùn)輸、擴(kuò)散以及能量傳遞之間復(fù)雜相互作用的多物理場(chǎng)耦合模型。考慮在重力影響下，煤自燃溫度會(huì)使大面積采空區(qū)局部氣體運(yùn)動(dòng)變得更加復(fù)雜，根據(jù)大面積采空區(qū)壓力、溫度和氣體濃度分布規(guī)律，對(duì)采空區(qū)氣體流場(chǎng)形成過(guò)程進(jìn)行分析，為大面積采空區(qū)的漏風(fēng)控制和煤自燃防治提供科學(xué)指導(dǎo)。

1 大面積采空區(qū)多物理場(chǎng)耦合模型

1.1 基本守恒方程

流體運(yùn)動(dòng)必須遵循3 個(gè)規(guī)則，即質(zhì)量守恒、動(dòng)量守恒和能量守恒[20]。

質(zhì)量守恒方程又稱為連續(xù)性方程，是指單位時(shí)間內(nèi)微元體中流體質(zhì)量的增加等于同一時(shí)間間隔內(nèi)流入該微元體的凈質(zhì)量，如式（1）：

式中：ρ 為混合氣體的密度，kg/m3；ε 為采空區(qū)孔隙度；u 為氣體滲流速度，m/s；Sm為采空區(qū)遺煤耗氧源項(xiàng)，kg/（m3·s）；t 為時(shí)間，s，。

動(dòng)量守恒方程描述了流體系統(tǒng)所受外力與流體系統(tǒng)動(dòng)量變化之間的關(guān)系，即微元體中流體動(dòng)量的增加率等于作用在微元體上各種力之和，如式（2）[21]：

式中：p 為氣體靜壓，Pa；τ 為黏性應(yīng)力張量；ρgn為重力體積力，kN/m3；α 為多孔介質(zhì)滲透率，m2；μ 為氣體動(dòng)力黏度，Pa·s。

能量守恒方程建立在物質(zhì)運(yùn)動(dòng)變化過(guò)程中控制體內(nèi)能量轉(zhuǎn)換的等量關(guān)系，決定了采空區(qū)高溫點(diǎn)的溫度分布和變化[22]。采空區(qū)氣體在運(yùn)移過(guò)程中，溫度的變化同時(shí)也會(huì)引起采空區(qū)氣體動(dòng)能和內(nèi)能改變，如式（3）：

式中：ρf、ρs為氣體和固體介質(zhì)的密度，kg/m3；Cpf、Cps為氣體和固體介質(zhì)的定壓比熱容，J/（kg·K）；kf、ks為氣體和固體介質(zhì)的導(dǎo)熱系數(shù)，W/（m·K）；T 為溫度，K；hi為組分i 的焓，kJ/mol；Ji為組分i 的擴(kuò)散通量，kg/（m2·s）；Sf為遺煤的耗氧放熱源項(xiàng)，W/m3。

煤氧化反應(yīng)速率r 與氧氣濃度和環(huán)境溫度有關(guān)，可根據(jù)Arrhenius 方程定義為式（4）：

式中：A 為指前因子，s-1；E 反應(yīng)活化能，kJ/mol；R 為普適氣體常數(shù)，kJ/（mol·K）；CO2為氧氣濃度，kmol/m3；n 為反應(yīng)級(jí)數(shù)。

指前因子A 和反應(yīng)活化能E 可根據(jù)煤低溫氧化實(shí)驗(yàn)獲得，利用Fluent 二次開發(fā)工具UDF 編寫遺煤耗氧及放熱源項(xiàng)。

1.2 粒子輸送方程

粒子輸送方程是指任一瞬時(shí)系統(tǒng)內(nèi)物理量（質(zhì)量、動(dòng)量和能量）隨時(shí)間的變化率等于該時(shí)間控制體內(nèi)物理量的變化率與通過(guò)控制體表面凈通量之和，是流體質(zhì)點(diǎn)物理量在流場(chǎng)中在時(shí)間空間上輸運(yùn)傳送變化過(guò)程[23]。粒子輸送方程如式（5）：

式中：ci為氣體組分i 的質(zhì)量分?jǐn)?shù)；Qi為組分i的增減源項(xiàng)，kg/（m3·s）。

1.3 大面積采空區(qū)孔隙率和滲透率方程

在采空區(qū)形成過(guò)程中，采空區(qū)上覆巖層垮落破碎，垮落巖石填充采空區(qū)形成多孔介質(zhì)空間。在采空區(qū)邊界煤巖層的支撐作用和地層的垮落壓實(shí)作用下，上覆巖層的變形呈現(xiàn)出一定垮落規(guī)律，用碎脹系數(shù)Kp表征煤巖塊堆積與壓實(shí)狀況，采空區(qū)走向上碎脹系數(shù)Kp（x）如式（6）：

采空區(qū)內(nèi)部的孔隙度分布與巖石的破碎程度相關(guān)，當(dāng)巖石破碎程度較小時(shí)，巖塊尺寸較大，孔隙度也變大；在采空區(qū)中心區(qū)域，上覆巖層在垮塌過(guò)程中破裂較為劇烈，在上覆巖層的自重壓實(shí)作用下，巖塊尺寸變小，孔隙度也越小。因此，有研究表明采空區(qū)孔隙度分布與邊界距離呈現(xiàn)“O”型分布形狀[24]。

滲透率是在壓差作用下，介質(zhì)中允許流體通過(guò)的能力，是用于描述流體在介質(zhì)內(nèi)傳導(dǎo)液體能力的參數(shù)。其中，滲透率越大，流體通過(guò)多孔介質(zhì)的速度越快。采空區(qū)孔隙率ε 與碎脹系數(shù)Kp如式（7），采空區(qū)滲透率α 根據(jù)Blake-Kozeny 公式計(jì)算得到[25]：

式中：Dp為垮落帶煤巖塊平均直徑，m。

2 大面積采空區(qū)數(shù)值模型

2.1 工程實(shí)例

鮑店煤礦位于山東省濟(jì)寧市，隸屬兗州煤業(yè)股份有限公司，鮑店煤礦煤樣的基礎(chǔ)參數(shù)為：水分1.910%，灰分12.790%，揮發(fā)分13.040%，固定碳52.760%，碳元素72.660%，氫元素4.274%，氮元素1.276%，臨界溫度70 °C。該地區(qū)煤樣為氣煤，易發(fā)生自燃。5316 工作面鄰近多個(gè)采空區(qū)，煤層埋藏深、通風(fēng)系統(tǒng)和地質(zhì)條件復(fù)雜、采空區(qū)范圍大、內(nèi)部遺煤多、漏風(fēng)通道多。5316 工作面推采過(guò)程中鄰近采空區(qū)漏風(fēng)流場(chǎng)極為復(fù)雜，存在煤自燃隱患。以U 型通風(fēng)5316 工作面采空區(qū)和鄰近五采區(qū)的5309、5310、5311 和5312 串聯(lián)采空區(qū)為物理模型進(jìn)行3D 數(shù)值模擬研究。大面積采空區(qū)模型如圖1。

圖1 大面積采空區(qū)模型Fig. 1 Large area goaf model

圖1（a）是工作面的二維示意圖，工作面回采率約90%，兩側(cè)留有5 m 煤柱。5309、5310、5311 和5312 采空區(qū)已封閉，5316 為生產(chǎn)面。CFD 模型如圖1（b），大面積采空區(qū)采用六面體網(wǎng)格，為提高計(jì)算精度，對(duì)邊界層網(wǎng)格進(jìn)行加密處理。

2.2 模型驗(yàn)證

5309 和5310 工作面采空區(qū)預(yù)留鉆孔，通過(guò)氣相色譜監(jiān)測(cè)采空區(qū)內(nèi)部氣體體積分?jǐn)?shù)，反映采空區(qū)氣體滲流情況，大面積采空區(qū)氧氣體積分?jǐn)?shù)的實(shí)測(cè)值與模擬值的對(duì)比驗(yàn)證如圖2。

圖2 大面積采空區(qū)氧氣體積分?jǐn)?shù)的實(shí)測(cè)值與模擬值的對(duì)比驗(yàn)證（風(fēng)速v=1.5 m/s，高度切片Z=1 m）Fig. 2 Verification of the measured and simulated values of oxygen concentration in large area goaf（ventilation flux v=1.5 m/s, extracted slices Z=1 m）

數(shù)值模擬結(jié)果與現(xiàn)場(chǎng)實(shí)測(cè)O2體積分?jǐn)?shù)變化規(guī)律基本一致，由于煤柱漏風(fēng)，隨著埋深的增加，采空區(qū)內(nèi)O2體積分?jǐn)?shù)迅速下降。風(fēng)流從進(jìn)風(fēng)口到回風(fēng)口的過(guò)程，5309 和5310 工作面O2體積分?jǐn)?shù)也有所降低。對(duì)于5316 工作面采空區(qū)，當(dāng)進(jìn)風(fēng)側(cè)采空區(qū)距工作面72 m 時(shí)，從散熱區(qū)到氧化升溫區(qū)，O2體積分?jǐn)?shù)降低到18%，距工作面150 m 時(shí)，從氧化區(qū)到窒息區(qū)，O2體積分?jǐn)?shù)降低至8%，氧化帶寬為78 m；回風(fēng)側(cè)距工作面82 m 時(shí)，O2體積分?jǐn)?shù)降低至18%，從散熱區(qū)到氧化升溫區(qū)，距工作面109 m 時(shí)，O2體積分?jǐn)?shù)降低8%，從氧化區(qū)到窒息區(qū)，氧化帶寬度為27 m。

3 結(jié)果與討論

3.1 煤自燃環(huán)境下大面積采空區(qū)溫度場(chǎng)演化規(guī)律

大面積采空區(qū)溫度隨時(shí)間的演化云圖如圖3，大面積采空區(qū)溫度隨高度變化云圖如圖4，大面積采空區(qū)溫度分布曲線如圖5。

圖3 大面積采空區(qū)溫度隨時(shí)間的演化云圖（風(fēng)速v=1.5 m/s，高度切片z=1 m）Fig.3 Cloud images of temperature evolution with time in large area goaf（ventilation flux v=1.5 m/s,extracted slices z=1 m）

圖4 大面積采空區(qū)溫度隨高度變化云圖（風(fēng)速v=1.5 m/s，第30 d）Fig.4 Cloud images of temperature variation with height in large area goaf（ventilation flux v=1.5 m/s,30th day）

圖5 大面積采空區(qū)溫度分布曲線（風(fēng)速v=1.5 m/s，第30 d，切片y=320 m）Fig.5 Temperature distribution curves of large area goaf（ventilation flux v=1.5 m/s, 30th day,extracted slices y=320 m）

由圖3 可知：隨著時(shí)間的推移，采空區(qū)高溫區(qū)域向四周擴(kuò)散，并隨風(fēng)流運(yùn)動(dòng)方向遷移，通過(guò)溫度的等值線可以看出溫度向外擴(kuò)散和蔓延的幅度較小，與環(huán)境溫度（300 K）形成較大的溫度梯度；由于風(fēng)流流經(jīng)隅角，在隅角處積聚，風(fēng)流對(duì)于熱量耗散作用較小，高溫點(diǎn)逐漸趨向隅角處，形成高溫點(diǎn)向隅角偏移現(xiàn)象。

由圖4 和圖5 可知：受傳導(dǎo)和對(duì)流作用影響，溫度在采空區(qū)豎直方向上進(jìn)行傳遞，高溫區(qū)域面積沒有明顯改變；在切片位置為z=1、5、10 m 處的最高溫度分別為346、339、336 K，溫度上升梯度分別為1.75 K/m 和0.6 K/m；在離開遺煤區(qū)域后，距火源距離越遠(yuǎn)，溫度上升梯度越小，這是由于采空區(qū)屬于多孔介質(zhì)空間，內(nèi)部的巖石和空氣的導(dǎo)熱能力較差，導(dǎo)致熱量不斷積聚，從而形成高溫點(diǎn)并能夠持續(xù)的維持其影響范圍。

3.2 煤自燃環(huán)境下大面積采空區(qū)氣體運(yùn)移規(guī)律

壓力是促使大面積采空區(qū)內(nèi)部氣體運(yùn)移的直接動(dòng)力，掌握采空區(qū)氧氣體積分?jǐn)?shù)分布和流動(dòng)規(guī)律是防治自然發(fā)火的關(guān)鍵技術(shù)基礎(chǔ)。大面積采空區(qū)不同高度下靜壓分布云圖如圖6，大面積采空區(qū)靜壓分布曲線如圖7，大面積采空區(qū)氧氣濃度分布云圖如圖8，大面積采空區(qū)氧氣濃度分布曲線如圖9，不同風(fēng)速條件下大面積采空區(qū)氧氣云圖如圖10。

圖6 大面積采空區(qū)不同高度下靜壓分布云圖（風(fēng)速v=1.5 m/s，第30 d）Fig.6 Cloud images of static pressure distribution at different heights in large area goaf（ventilation flux v=1.5 m/s, 30th day）

圖8 大面積采空區(qū)氧氣體積分?jǐn)?shù)分布云圖（風(fēng)速v=1.5 m/s，第30 d）Fig.8 Cloud images of oxygen concentration distribution in large area goaf（ventilation flux v=1.5 m/s, 30th day）

圖9 大面積采空區(qū)氧氣體積分?jǐn)?shù)分布曲線（風(fēng)速v=1.5 m/s，第30 d）Fig.9 Oxygen concentration distribution curves of large area goaf（ventilation flux v=1.5 m/s, 30th day）

圖10 不同風(fēng)速條件下大面積采空區(qū)氧氣云圖（切片高度z=10 m）Fig.10 Oxygen cloud images of large area goaf under different wind speeds（extracted slices z=10 m）

由圖6 和圖7 可知：①無(wú)煤自燃的情況下，通風(fēng)作用會(huì)在進(jìn)風(fēng)側(cè)形成較大壓力，在回風(fēng)側(cè)形成負(fù)壓。靠近進(jìn)風(fēng)側(cè)采空區(qū)深部壓力隨著高度增加而增大，而靠近回風(fēng)側(cè)采空區(qū)壓力隨高度增大而急劇減小，這是由于進(jìn)風(fēng)側(cè)風(fēng)流隨漏風(fēng)通道進(jìn)入采空區(qū)，向采空區(qū)上部運(yùn)移，導(dǎo)致采空區(qū)上部壓力增大，最終風(fēng)流匯聚在回風(fēng)側(cè)流出，致使壓力驟降；②煤自燃產(chǎn)生高溫會(huì)使氣體膨脹，密度減小，從而在采空區(qū)底部形成負(fù)壓區(qū)域，相應(yīng)的受浮力作用影響，氣體受熱膨脹向上運(yùn)動(dòng)，在采空區(qū)頂部積聚形成高壓區(qū)域。

由圖8 和圖9 可知：①氧氣進(jìn)入采空區(qū)，受孔隙率分布影響，在向采空區(qū)深部擴(kuò)散過(guò)程中呈現(xiàn)出不同的濃度分布，且隨著采空區(qū)豎直高度越高氧氣體積分?jǐn)?shù)越低，這是由于氧氣密度大于空氣，受重力影響，導(dǎo)致氧氣在采空區(qū)底部積聚；②在煤自燃作用下，可以看到溫度產(chǎn)生的火風(fēng)壓作用使局部風(fēng)流方向發(fā)生改變，使得采空區(qū)煤自燃區(qū)域不斷卷吸周邊氧氣，產(chǎn)生強(qiáng)對(duì)流作用涌入采空區(qū)上部區(qū)域并在隅角處積聚。

加強(qiáng)通風(fēng)會(huì)導(dǎo)致采空區(qū)漏風(fēng)強(qiáng)度增大，進(jìn)而影響采空區(qū)氣體的運(yùn)移和積聚。

由圖10 可知：隨著進(jìn)風(fēng)風(fēng)速的增加，進(jìn)風(fēng)側(cè)的壓力明顯增大，采空區(qū)漏風(fēng)更加嚴(yán)重；因此，進(jìn)風(fēng)風(fēng)速越大，氧化升溫帶越寬，并且明顯從工作面向采空區(qū)深部移動(dòng)；隨著進(jìn)風(fēng)風(fēng)速?gòu)?.5 m/s 增加到2.0 m/s，5309 采空區(qū)進(jìn)風(fēng)側(cè)氧化升溫帶從距工作面55 m移動(dòng)到65 m；由于氧化升溫帶向深部轉(zhuǎn)移，采空區(qū)煤自燃的風(fēng)險(xiǎn)會(huì)更加嚴(yán)重。

4 結(jié) 語(yǔ)

建立了壓力場(chǎng)、溫度場(chǎng)、氣體濃度場(chǎng)的多場(chǎng)耦合模型，并將重力作用納入模擬煤自燃環(huán)境下大面積采空區(qū)模型，研究煤自燃環(huán)境下大面積采空區(qū)氣體運(yùn)移規(guī)律。同時(shí)，討論了不同風(fēng)速條件下O2濃度場(chǎng)的變化規(guī)律。

1）大面積采空區(qū)的漏風(fēng)和煤自燃高溫引起的熱浮力效應(yīng)是影響采空區(qū)氣體運(yùn)移主要因素。采空區(qū)發(fā)生煤自燃時(shí)，由于內(nèi)部巖石和空氣導(dǎo)熱系數(shù)較低，高度方向上溫度梯度從1.75 K/m 下降至0.6 K/m，熱量在遺煤區(qū)域積聚，引起周圍氣體的漏風(fēng)和運(yùn)移。

2）煤自燃高溫產(chǎn)生的火風(fēng)壓作用使局部風(fēng)流方向發(fā)生改變，使得采空區(qū)煤自燃區(qū)域不斷卷吸周邊氧氣，產(chǎn)生強(qiáng)對(duì)流作用涌入采空區(qū)上部區(qū)域并在隅角處積聚。

3）隨著漏風(fēng)源風(fēng)速的增大，大面積采空區(qū)進(jìn)風(fēng)側(cè)氧化升溫帶從距工作面55 m 移動(dòng)到65 m，增大了大面積采空區(qū)煤自燃的風(fēng)險(xiǎn)，同時(shí)加快了采空區(qū)深部有害氣體向工作面的運(yùn)移。