采空區(qū)自然發(fā)火多場(chǎng)耦合數(shù)值模擬研究

黎 力，李治剛，奚 弦，賀 偉

(1.太原理工大學(xué)礦業(yè)工程學(xué)院，山西 太原 030024； 2.太原理工大學(xué) 原位改性采礦教育部重點(diǎn)實(shí)驗(yàn)室，山西 太原 030024；3.中國(guó)礦業(yè)大學(xué)安全工程學(xué)院,江蘇 徐州 221008；4.中國(guó)礦業(yè)大學(xué)煤炭資源與安全開(kāi)采國(guó)家重點(diǎn)實(shí)驗(yàn)室，江蘇 徐州 221008)

采空區(qū)自然發(fā)火多場(chǎng)耦合數(shù)值模擬研究

黎 力1,2，李治剛1，奚 弦3,4，賀 偉1,2

(1.太原理工大學(xué)礦業(yè)工程學(xué)院，山西 太原 030024； 2.太原理工大學(xué) 原位改性采礦教育部重點(diǎn)實(shí)驗(yàn)室，山西 太原 030024；3.中國(guó)礦業(yè)大學(xué)安全工程學(xué)院,江蘇 徐州 221008；4.中國(guó)礦業(yè)大學(xué)煤炭資源與安全開(kāi)采國(guó)家重點(diǎn)實(shí)驗(yàn)室，江蘇 徐州 221008)

采空區(qū)遺煤自然發(fā)火嚴(yán)重制約著煤礦的安全生產(chǎn)，其受到滲流場(chǎng)、擴(kuò)散場(chǎng)、化學(xué)反應(yīng)場(chǎng)以及多孔介質(zhì)傳熱場(chǎng)共同作用影響，基于這一特點(diǎn)，本文采用有限元數(shù)值模擬手段對(duì)采空區(qū)內(nèi)部自然發(fā)火多場(chǎng)耦合作用進(jìn)行細(xì)致研究，并通過(guò)模擬結(jié)果與現(xiàn)場(chǎng)實(shí)測(cè)數(shù)據(jù)比對(duì)驗(yàn)證了模擬的可行性以及模擬結(jié)果的可靠性。研究結(jié)果表明：依靠傳統(tǒng)的三帶劃分方式并不能很好的確定采空區(qū)可能自然發(fā)火位置，此外，工作面漏風(fēng)入口后部煤體溫度變化明顯，在實(shí)際工作中該區(qū)域范圍應(yīng)作為重點(diǎn)監(jiān)控對(duì)象。

采空區(qū)；自然發(fā)火；數(shù)值模擬；多場(chǎng)耦合

煤炭資源是我國(guó)經(jīng)濟(jì)發(fā)展的主體能源之一，然而煤炭自然發(fā)火嚴(yán)重制約著我國(guó)煤礦的安全高效生產(chǎn)。據(jù)統(tǒng)計(jì)，我國(guó)自然發(fā)火礦井占全國(guó)重點(diǎn)煤礦的50%以上，由于煤炭自然引起的礦井火災(zāi)占總數(shù)的90%以上[1]，由此可見(jiàn)，煤炭自然發(fā)火的防治是煤礦安全生產(chǎn)的重要任務(wù),其中采空區(qū)自然發(fā)火防治為重中之重[2]。近年來(lái)，使用數(shù)值模擬手段研究采空區(qū)內(nèi)部煤炭自然發(fā)火受到越來(lái)越多的研究學(xué)者關(guān)注[3-5]。針對(duì)采空區(qū)自然發(fā)火這一包含滲流場(chǎng)、擴(kuò)散場(chǎng)、化學(xué)反應(yīng)場(chǎng)以及多孔介質(zhì)熱傳遞場(chǎng)這一多場(chǎng)耦合特性，本文使用COMSOL Multiphysics軟件對(duì)采空區(qū)自然發(fā)火進(jìn)行了瞬態(tài)模型模擬求解。并將實(shí)測(cè)數(shù)據(jù)與模擬數(shù)據(jù)進(jìn)行比對(duì)驗(yàn)證，為研究采空區(qū)自然發(fā)火提供可靠的數(shù)據(jù)支持。

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

1.1 采空區(qū)流體運(yùn)移控制方程

一般認(rèn)為采空區(qū)內(nèi)部密布著塊徑不同的破碎巖體，同時(shí)破碎巖體之間所形成的裂隙結(jié)構(gòu)在整個(gè)采空區(qū)空間均存在[4]，因此多孔介質(zhì)這一材料界定可適用采空區(qū)巖體[6]。采空區(qū)滲流場(chǎng)可用達(dá)西定律表示，見(jiàn)式(1)。

(1)

采空區(qū)內(nèi)部孔隙率與滲透率之間相互關(guān)聯(lián)，采空區(qū)遺煤滲透率可由多孔介質(zhì)的Kozeny Carman公式得出其經(jīng)驗(yàn)計(jì)算公式，見(jiàn)式(2)。

(2)

采空區(qū)氣體濃度的變化規(guī)律符合菲克定律，所以，氧氣氣體濃度擴(kuò)散規(guī)律符合對(duì)流擴(kuò)散方程[7]，見(jiàn)式(3)。

▽(-D▽C)+U▽C=R′

(3)

式中：D為擴(kuò)散系數(shù)，m2/s；C為氣體濃度，mol/ m3；U為漏風(fēng)風(fēng)流速度，m/s；R′為反應(yīng)率mol/(m3·s)。

1.2 采空區(qū)遺煤氧化反應(yīng)化學(xué)反應(yīng)動(dòng)力原理

采空區(qū)具有自燃傾向性的煤層想要發(fā)生自燃，必須有一個(gè)溫度升高的過(guò)程，而溫度的升高主要和采空區(qū)煤與氧氣之間的物理化學(xué)復(fù)合作用以及散熱條件有關(guān)。采空區(qū)遺留的松散煤體與氧氣發(fā)生氧化反應(yīng)時(shí)，它的主要?dú)怏w為CO和CO2等物質(zhì)[8]，根據(jù)化學(xué)反應(yīng)動(dòng)力學(xué)，其反應(yīng)過(guò)程參數(shù)關(guān)系可表示為式(4)。其化學(xué)表達(dá)式可以表示為式(5)[9]。

f(CO2,CCO,CCO2;t)=0

(4)

煤+O2→產(chǎn)物+Q

(5)

式中，Q為單位時(shí)間化學(xué)反應(yīng)產(chǎn)熱，J/mol。

式(5)化學(xué)反應(yīng)反應(yīng)速率見(jiàn)式(6)[10]。

R′=C(O2)×r

(6)

式中：R′化學(xué)反應(yīng)速率， mol/(m3·s)；C(O2)氧氣濃度，mol/m3；r化學(xué)反應(yīng)速率常數(shù)，1/s。1889年Arrhenius提出了化學(xué)反應(yīng)速率的求解公式，并提出反應(yīng)速率為溫度的函數(shù)，見(jiàn)式(7)。

(7)

式中：r為化學(xué)反應(yīng)速率常數(shù)，該參數(shù)在不同的環(huán)境之下具有不同的量綱，本文為1/s；r0為頻率因子，量綱與化學(xué)反應(yīng)速率常數(shù)相同，1/s；E為活化能，J/mol；R為氣體常數(shù)，8 135 J/(mol·K)；T為溫度，K。

化學(xué)反應(yīng)產(chǎn)熱計(jì)算見(jiàn)式(8)。

(8)

式中：Qr為化學(xué)反應(yīng)產(chǎn)熱，W/m3；QP為生成物質(zhì)數(shù)量：Qre為反應(yīng)物質(zhì)數(shù)量；νij為化學(xué)計(jì)量數(shù)；hi為標(biāo)準(zhǔn)摩爾生成焓，J/mol。

1.3 采空區(qū)多孔介質(zhì)傳熱

采空區(qū)內(nèi)部熱量傳遞極為復(fù)雜，根據(jù)前述采空區(qū)巖體為多孔介質(zhì)這一特點(diǎn)，結(jié)合傳熱理論，認(rèn)為熱量在采空區(qū)內(nèi)部傳遞方式包括固體熱傳導(dǎo)、流體與固體、流體與流體對(duì)流傳熱，其傳熱控制方程見(jiàn)式(9)。

(9)

式中：ρp為采空區(qū)巖體密度，m3/kg；Cp,p為巖體比熱容，J/(kg·K)；Cp為氣體常壓比熱容，J/(kg·K)；βp為巖體導(dǎo)熱系數(shù)，W/(m·K)；β為氣體導(dǎo)熱率，W/(m·K)。

2 模型描述

2.1 模型描述與參數(shù)選取

圖1 采空區(qū)模擬幾何模型

圖1中模型被用作研究采空區(qū)內(nèi)部的氣體運(yùn)移以及三帶劃分，采空區(qū)區(qū)域?yàn)?20 m×300 m，G1為液壓支架，左側(cè)上下S1與S2為采面漏風(fēng)進(jìn)口與出口，寬度均為5 m，求解域單元數(shù)為3 924，求解自由度為15 518。模型中各計(jì)算參數(shù)設(shè)置為：空氣摩爾質(zhì)量M29 g/mol；采空區(qū)孔隙率φ0.29；煤體在自然堆積狀態(tài)下的滲透率5.38×10-6m2，在受周?chē)鷰r體載荷影響的情況下滲透率2.6×10-6m2，在圍巖壓實(shí)穩(wěn)定地區(qū)煤體的滲透率情況下滲透率1.3×10-6m2；動(dòng)力黏度μ2.01×10-5Pa·s；氣體擴(kuò)散系數(shù)D2.88×10-2m2/s，假設(shè)采空區(qū)煤體氧化反應(yīng)主要產(chǎn)物為CO2與CO，其反應(yīng)速率分別為1.026×10-3mol/(m3·s)和3.4×10-4mol/(m3·s)；氣體常壓比熱容Cp29.099 J/(kg·K)；巖體比熱容Cp,p1.26kJ/(kg·K)。

適當(dāng)?shù)倪吔缗c初始條件是確保模擬準(zhǔn)確性的關(guān)鍵。該模擬模型中，滲流場(chǎng)采空區(qū)內(nèi)部初始?jí)毫0為0.1 MPa，漏風(fēng)入口端氣體壓力Pi為0.5 MPa，漏風(fēng)風(fēng)速為0.55 m/s，出口端壓力為0.1 MPa，其余邊界為零流量邊界；擴(kuò)散場(chǎng)采空區(qū)內(nèi)氧氣初始濃度為0 mol/m3，入口邊界條件為9.286 mol/m3；溫度場(chǎng)初始條件為298.15 K，除工作面與采空區(qū)頂端邊界條件為定溫邊界條件298.15 K，其余邊界條件為絕熱邊界條件。

2.2 模型數(shù)據(jù)與實(shí)際數(shù)據(jù)對(duì)比

為了驗(yàn)證該數(shù)據(jù)擬合模型的可靠性，將模擬得到的采空區(qū)溫度數(shù)據(jù)與實(shí)際監(jiān)測(cè)數(shù)據(jù)進(jìn)行對(duì)比驗(yàn)證，如圖2所示。圖2中模擬數(shù)據(jù)與實(shí)測(cè)數(shù)據(jù)趨勢(shì)一致，數(shù)據(jù)基本吻合，模擬模型對(duì)采空區(qū)內(nèi)部溫度分布復(fù)現(xiàn)性較好，模擬結(jié)果能夠有效的反應(yīng)采空區(qū)真實(shí)狀況。

3 模擬結(jié)果分析與討論

3.1 采空區(qū)內(nèi)部氣壓與氣體流速分布

圖2 采空區(qū)溫度模擬數(shù)據(jù)與實(shí)測(cè)數(shù)據(jù)對(duì)比

由圖3可見(jiàn)，采空區(qū)內(nèi)部的氣體壓力分布呈“對(duì)稱(chēng)”分布，液壓支架后部自漏風(fēng)進(jìn)口至漏風(fēng)出口壓力降低明顯，隨著距工作面距離的變大，壓力降低幅度減小，至100 m左右以外壓力趨于一致。此外，壓力變化較大區(qū)域分布在工作面漏風(fēng)進(jìn)出口位置，其余部分變化很小。

圖3 采空區(qū)內(nèi)氣體壓力分布

圖4 采空區(qū)內(nèi)氣流速度分布

圖4為采空區(qū)內(nèi)部氣體流速與流向情況分布，其中箭頭方向?yàn)闅怏w流動(dòng)方向。采空區(qū)內(nèi)部氣流速度在進(jìn)出口處變化梯度較大，氣流速度最大點(diǎn)處位于液壓支架邊緣地區(qū)。氣體流向呈環(huán)形。氣體流速較高區(qū)域主要集中在采空區(qū)的煤體自然堆積區(qū)域，其他區(qū)域幾乎沒(méi)有可監(jiān)測(cè)的氣體流動(dòng)。由此可見(jiàn)，如用漏風(fēng)風(fēng)速劃定采空區(qū)三帶并不合理。

3.2 采空區(qū)內(nèi)部氧氣濃度分布

圖5 1 d時(shí)采空區(qū)內(nèi)部氧氣濃度分布

圖6 50 d時(shí)采空區(qū)內(nèi)部氧氣濃度分布

圖7 100 d時(shí)采空區(qū)內(nèi)部氧氣濃度分布

圖5～7為采空區(qū)內(nèi)部第1 d、50 d、100 d內(nèi)的氧氣濃度分布?？梢悦黠@看出，在第1 d時(shí)，采空區(qū)內(nèi)部氧氣主要分布于100 m以?xún)?nèi)，漏風(fēng)進(jìn)口處的氧氣濃度擴(kuò)散范圍要大于漏風(fēng)出口處，這是由于空氣的達(dá)西滲流在出口處較為劇烈，從而攜帶了更多的氧氣進(jìn)入采空區(qū)內(nèi)部。第1 d至第50 d之間，氧氣運(yùn)移效果十分明顯。至第50 d時(shí)，氧氣已經(jīng)擴(kuò)散到了采空區(qū)中部。而到第100 d時(shí)，氧氣擴(kuò)散范圍與第50 d時(shí)差異并不大。根據(jù)采空區(qū)三帶氧氣濃度劃分準(zhǔn)則，三帶劃分線應(yīng)為曲線形，大致范圍為0～100 m為散熱帶，100～180 m為氧化升溫帶，之后為窒息帶，這恰好大致對(duì)應(yīng)于采空區(qū)煤體自然堆積區(qū)、載荷影響區(qū)與壓實(shí)穩(wěn)定區(qū)。

3.2 采空區(qū)內(nèi)部溫度分布

圖8 1 d時(shí)采空區(qū)內(nèi)部溫度分布

圖9 50 d時(shí)采空區(qū)內(nèi)部溫度分布

圖10 100 d時(shí)采空區(qū)內(nèi)部溫度分布

圖8～10顯示了第1 d、50 d和100 d的采空區(qū)溫度變化。第1 d時(shí)，采空區(qū)溫度升高區(qū)域主要集中在液壓支架內(nèi)方100 m范圍內(nèi)，升溫幅度并不明顯，最高溫度僅為25.3 ℃。第50 d時(shí)，溫度升高區(qū)域移動(dòng)到100～200 m范圍，最高溫度上升至32.6℃，其前端由于達(dá)西滲流影響，大量熱量被帶至漏風(fēng)出口處，因此，在漏風(fēng)出口處有一個(gè)較緩的溫度變化帶，同時(shí)，漏風(fēng)出口處溫度較低的另一原因?yàn)樵撎幯鯕夂枯^低，煤與氧氣氧化產(chǎn)熱較少。第100 d時(shí)，采空區(qū)內(nèi)高溫區(qū)域較50 d時(shí)并沒(méi)有明顯的變化，但其最高溫度上升至37.3 ℃。值得注意的是，采空區(qū)高溫區(qū)域無(wú)論在哪個(gè)時(shí)刻都是大致處于漏風(fēng)入口后部，其可能導(dǎo)致原因?yàn)楦吆趿考暗蜌怏w流速造成的熱量的疊加累積。

6 結(jié)論

1)采空區(qū)自然發(fā)火是氣體滲流、擴(kuò)散，氧氣與采空區(qū)內(nèi)部煤體氧化反應(yīng)，以及采空區(qū)內(nèi)部巖體傳熱等效應(yīng)共同作用的結(jié)果，若只考慮其中一個(gè)或幾個(gè)因素并不能得到良好的研究效果。

2)基于COMSOL進(jìn)行計(jì)算的多場(chǎng)耦合采空區(qū)模擬模型能夠很好的匹配實(shí)際情況，為實(shí)際采空區(qū)自然發(fā)火防治提供指導(dǎo)。

3)采用工作面漏風(fēng)風(fēng)速劃定采空區(qū)三帶的指標(biāo)并不合理，采空區(qū)內(nèi)部自然三帶的劃分應(yīng)該綜合考慮氧氣濃度與溫度條件。

4)采空區(qū)內(nèi)部氧氣濃度與溫度影響區(qū)域在初期變化明顯，在后期則不然。同時(shí)，工作面漏風(fēng)入口處后部為煤體溫度升高最為顯著區(qū)域，即其自然發(fā)火可能性最大，在實(shí)際生產(chǎn)中需要重點(diǎn)監(jiān)控。

[1] 郭金偉. 采空區(qū)煤炭自燃發(fā)火防治技術(shù)研究[J]. 中國(guó)化工貿(mào)易, 2015,7(26):69.

[2] 嚴(yán)民杰. 納林廟煤礦采空區(qū)自然發(fā)火的綜合防治技術(shù)研究[J]. 中國(guó)礦業(yè), 2010,19(8):95-96.

[3] 李宗翔, 許端平, 劉立群. 采空區(qū)自然發(fā)火“三帶”劃分的數(shù)值模擬[J]. 遼寧工程技術(shù)大學(xué)學(xué)報(bào), 2002,21(5):545-548.

[4] 周西華, 郭梁輝, 孟樂(lè). 易自燃煤層綜放工作面采空區(qū)自然發(fā)火防治數(shù)值模擬[J]. 中國(guó)地質(zhì)災(zāi)害與防治學(xué)報(bào), 2012,23(1):83-87.

[5] 邵磊. 采空區(qū)自燃“三帶”的實(shí)測(cè)分析及數(shù)值模擬[D]. 淮南：安徽理工大學(xué), 2014.

[6] 孔祥言. 高等滲流力學(xué)[M]. 合肥: 中國(guó)科學(xué)技術(shù)大學(xué)出版社, 1999.

[7] 周春山. 煤自然實(shí)驗(yàn)過(guò)程中溫度場(chǎng)的數(shù)值模擬研究[D]. 太原：太原理工大學(xué), 2007.

[8] Deng J, Zhao J, Zhang Y, et al. Study on Coal Spontaneous Combustion Characteristic Temperature of Growth Rate Analysis[J]. Procedia Engineering, 2014,84:796-805.

[9] Qin B, Sun Q, Wang D, et al. Analysis and key control technologies to prevent spontaneous coal combustion occurring at a fully mechanized caving face with large obliquity in deep mines[J]. Mining Science and Technology, 2009,19(4):446-451.

[10] 許越. 化學(xué)反應(yīng)動(dòng)力學(xué)[M]. 北京: 化學(xué)工業(yè)出版社, 2005.

Multi-filed coupling numerical simulation of spontaneous combustion in goaf

LI Li1,2, LI Zhigang1, XI Xian3,4, HE Wei1,2

(1. College of Mining Engineering, Taiyuan University of Technology, Taiyuan 030024, China; 2. Key Laboratory of In-situ Property-improving Mining of Ministry of Education, Taiyuan University of Technology, Taiyuan 030024, China; 3.College of Safety Engineering, China University of Mining and Technology, Xuzhou 221000, China; 4.State Key Laboratory of Coal Resources and Safe Mining, China University of Mining and Technology, Xuzhou 221000, China)

Coal spontaneous combustion is one of the factors that endangers the produce in coal mine. It is affected by the combined action of seepage and diffusion of gas, chemical reaction of coal and oxygen as well as heat transfer in porous media. Based on this, a series of finite element numerical simulations were conducted to investigate the spontaneous combustion under the effect of different physical field, and the simulation results were compared with the data obtained by in-suit measurement in order to test the accuracy of the simulation data and the feasibility of the simulation model. It is found out that the numerical simulation model acts well in prediction of gas transport and temperature variation, besides it is also discovered that the traditional partition way of the three zone of spontaneous combustion in goaf can not precisely locate the position where spontaneous combustion might occur. In addition, the area behind the entrance of leaked air from working face to the goaf is the place where temperature changes significantly, and it is the most dangerous place that may take place coal fire, so during practical work of coal mining, more attentions should be attached

goaf; spontaneous combustion; numerical simulation; multi-field coupling

2016-06-05

太原理工大學(xué)2014年?；痦?xiàng)目資助(編號(hào):1205-04020102)；山西省研究生教育創(chuàng)新項(xiàng)目資助(編號(hào)：2016SY014)

黎力(1991-)，男，江蘇南京人，碩士研究生，主要從事煤礦數(shù)值模擬方面的研究工作,E-mail：lili0661@link.tyut.edu.cn。

李治剛(1978-)，男，遼寧大連人，博士，講師，主要從事煤礦安全方面研究工作,E-mail：zhiganglee2009@hotmail.com。

TD353

A

1004-4051(2017)03-0157-04