高瓦斯煤層采空區(qū)流場滲流模型研究

中國40%的厚煤層自然起火，其中大部分是空區(qū)煤自殺火災(zāi)，占煤自發(fā)火災(zāi)的60%。通過綜合開采和上煤技術(shù)，利用厚煤層，由于回收率相對(duì)較低，采空區(qū)頂板上的磚塊數(shù)量大，且上海角頂板的磚塊數(shù)量頻繁超標(biāo)。所以對(duì)于自燃厚煤層來說,有必要研究采空區(qū)的滲流特性,分析采空區(qū)的流場分布,為防治遺煤自燃和瓦斯超限制定相應(yīng)措施。筆者基于國內(nèi)外研究現(xiàn)狀討論了采空區(qū)多孔介質(zhì)的特性,結(jié)合O形圈理論將采空區(qū)劃分為不同區(qū)域,并研究了各個(gè)區(qū)域的孔隙率;利用連續(xù)性方程、動(dòng)量方程、瓦斯動(dòng)力彌散方程和滲流半經(jīng)驗(yàn)公式建立風(fēng)流在采空區(qū)流動(dòng)的基本方程組,為采空區(qū)風(fēng)流參數(shù)求解提供理論基礎(chǔ);結(jié)合多孔介質(zhì)以及風(fēng)流流態(tài)分析了數(shù)值模擬參數(shù)的設(shè)定;最后通過實(shí)例分析了采空區(qū)的風(fēng)流流場分布規(guī)律,為預(yù)防采空區(qū)遺煤自燃和治理上隅角瓦斯超限提供相應(yīng)措施。1流場模擬的發(fā)展一般將采空區(qū)視為是由多孔介質(zhì)充填的立體空間。采空區(qū)的滲流特性包括多孔介質(zhì)自身的特性(如孔隙率、滲透率)和風(fēng)流在多孔介質(zhì)中滲流特性(如層流、過渡流和紊流等)。多孔介質(zhì)的孔隙率和滲透率是影響整個(gè)采空區(qū)漏風(fēng)流場的重要參數(shù)?；诓蓜?dòng)垮落分布特征,孔隙率和滲透率的眾多理論模型和經(jīng)驗(yàn)公式相繼提出并被應(yīng)用到流場模擬中:錢鳴高等指出在采空區(qū)走向上,頂板由初始垮落到壓實(shí)近似服從負(fù)指數(shù)變化關(guān)系,文獻(xiàn)根據(jù)“砌體梁”結(jié)構(gòu)推導(dǎo)出了滲透率隨采空區(qū)深度變化的數(shù)學(xué)公式;宋振騏指出傾向上基本頂斷裂的發(fā)展進(jìn)程和出現(xiàn)時(shí)間與走向上是同步的,在采空區(qū)傾向上和走向上的壓力顯現(xiàn)規(guī)律基本相同;李樹剛認(rèn)為采空區(qū)可以細(xì)化為自然堆積區(qū)、承壓碎脹區(qū)和壓實(shí)區(qū),各個(gè)區(qū)的碎脹系數(shù)根據(jù)承載壓力的不同而不同(圖1)。1856年Darcy提出了單相流滲流線性公式,第一次描述了流速與壓降的關(guān)系;Carman對(duì)Darcy公式進(jìn)行了深入討論,在Kozeny關(guān)系式的基礎(chǔ)上運(yùn)用毛細(xì)管模型推導(dǎo)出了多孔介質(zhì)單相流在低雷諾數(shù)的壓降公式;Reynolds指出高流速下壓降與流速的關(guān)系偏離Darcy公式;Forchheimer指出多孔介質(zhì)中單相流的壓降由黏性力和慣性力2部分組成;Ergun在前人的基礎(chǔ)上,認(rèn)為壓降為流速一次方和二次方之和。2多孔介質(zhì)在采空區(qū)的特性2.1滲透率與孔隙率的關(guān)系描述多孔介質(zhì)的主要參數(shù)有孔隙率n、碎脹系數(shù)Kp、滲透率e。孔隙率與碎脹系數(shù)的關(guān)系為n=1?1/Kp(1)n=1-1/Κp(1)由Blake-Kozeny公式,滲透率與孔隙率的關(guān)系式為e=D2m150n3(1?n)2(2)e=Dm2150n3(1-n)2(2)其中,Dm為多孔介質(zhì)平均粒子直徑,m。2.2工作面平面范圍特征據(jù)巖層移動(dòng)規(guī)律和O形圈理論,采空區(qū)中央呈O形一樣的壓實(shí)區(qū),碎脹系數(shù)較小,越靠近采空區(qū)邊界越大;采空區(qū)走向和傾向的頂板垮落相似。不考慮煤層傾角影響,采空區(qū)冒落情況沿工作面中部對(duì)稱,采空區(qū)3個(gè)區(qū)的平面范圍如圖2所示。相似試驗(yàn)表明A1+B1≈(1～2)L1(L1為周期跨距),A2+B2≈0.8L2(L2為初次跨距),A3+B3≈L2。但由于靠近工作面處有支架支撐,靠近工作面的采空區(qū)煤巖所受壓力應(yīng)比靠近開切眼處的壓力要小,因此A3+B3值會(huì)比初始跨距小,而且隨著工作面的推移,壓實(shí)區(qū)會(huì)緩慢向走向兩側(cè)發(fā)展,A3+B3的值會(huì)越來越小。1管鄰近工作面區(qū)域碎脹系數(shù)A區(qū)巖層冒落頂成拱形,垮落煤巖塊未承受壓力,為自然膨脹區(qū),碎脹系數(shù)可取為定值。盡管靠近工作面或開切眼的區(qū)域和兩側(cè)區(qū)域垮落形式相似,但由于受力后的垮落程度不同,A1、A2、A3值不同,碎脹系數(shù)亦不同。在A區(qū)中,支架后方以及開切眼附近的傾向區(qū)域的平均碎脹系數(shù)設(shè)為Kpa1,采空區(qū)進(jìn)風(fēng)巷和回風(fēng)巷附近區(qū)域平均碎脹系數(shù)設(shè)為Kpa2。2“回”形采空區(qū)由圖2可知,B區(qū)呈“回”字形。該區(qū)域的煤巖平均碎脹系數(shù)為Kpb。很多文獻(xiàn)粗略地將采空區(qū)劃分為近采空區(qū)和遠(yuǎn)采空區(qū),而沒有考慮周邊冒落情況。由于采空區(qū)周邊是主要漏風(fēng)通道,結(jié)合自然堆積區(qū)和壓實(shí)區(qū)劃分成“回”形是合理的。該區(qū)域的煤巖碎脹系數(shù)受頂板壓力的影響變化,距離壓實(shí)區(qū)越近,碎脹系數(shù)越小。3不考慮采空區(qū)高度的二維平面根據(jù)O形圈理論,壓實(shí)區(qū)近似為一橢拋體,由于頂板周期冒落,壓實(shí)區(qū)包絡(luò)線(不考慮采空區(qū)高度的二維平面)近似為一矩形和2個(gè)橢圓形半圓合成的邊界。在開切眼和工作面支架附近為半橢圓形,在兩巷為一直線。一般采空區(qū)深部對(duì)所研究的問題影響不大,所以取值固定,該區(qū)域煤巖平均碎脹系數(shù)為Kpc。3在廣闊的空間里，流的流動(dòng)和數(shù)學(xué)模型3.1雷諾數(shù)計(jì)算re采空區(qū)風(fēng)流流動(dòng)屬于滲流流動(dòng),由于采空區(qū)各處煤與矸石壓實(shí)程度差異很大,各處風(fēng)壓及距進(jìn)回風(fēng)邊界距離很大,因此采空區(qū)各點(diǎn)的流速相差很大,各點(diǎn)風(fēng)流的雷諾數(shù)亦相差很大。采空區(qū)的多孔介質(zhì)內(nèi)風(fēng)流流態(tài)用滲流雷諾數(shù)Re表示:Re=ρve1/2/μ(3)Re=ρve1/2/μ(3)式中:v為當(dāng)?shù)貪B流風(fēng)速,m/s;μ為混合氣體絕對(duì)黏度,一般井下風(fēng)流取1.8×10-5Pa·s;ρ為混合氣體密度,井下一般取1.2kg/m3。當(dāng)Re≤10為層流狀態(tài);10<Re≤100為過度流態(tài);Re>100為紊流狀態(tài)。若取v=2m/s,n=0.2(自然堆積區(qū)孔隙率),μ=1.8×10-5Pa·s,ρ=1.2kg/m3,Dm=0.2m,計(jì)算得Re=240,為紊流狀態(tài)。隨風(fēng)流流向采空區(qū)深部,采空區(qū)逐漸壓實(shí),風(fēng)速逐漸減小,最終轉(zhuǎn)為層流直至風(fēng)流靜止。因此采空區(qū)風(fēng)流是從紊流過度到層流的非線性滲流。3.2流暢流的基本公式1采空區(qū)瓦斯源實(shí)際狀況一般認(rèn)為在亞音速下的風(fēng)流為不可壓氣體,那么井下采空區(qū)風(fēng)流穩(wěn)定的連續(xù)性方程為divV?=0(4)divV→=0(4)其中,V?V→為滲流速度矢量。井下風(fēng)流一般含有瓦斯,采空區(qū)、煤壁也有瓦斯涌出源,所以如果要分析采空區(qū)瓦斯流動(dòng)的特點(diǎn),要加上瓦斯動(dòng)力彌散方程。?c?t=??xi(Dij?c?xj)???xi(cui)+ICH4(5)?c?t=??xi(Dij?c?xj)-??xi(cui)+ΙCΗ4(5)式中:c為采空區(qū)某點(diǎn)的瓦斯質(zhì)量濃度,kg/m3;t為瓦斯源涌出時(shí)間,s;Dij為瓦斯動(dòng)力彌散系數(shù)的9個(gè)分量(i,j=1,2,3),m2/s;ui為i方向上采空區(qū)某點(diǎn)風(fēng)流平均流速,m/s;xi、xj分別為i、j方向上采空區(qū)某點(diǎn)距源相距離,m;ICH4為瓦斯源項(xiàng),kg/(m3·s)。2采空區(qū)內(nèi)的動(dòng)力壓降ρ??xj(uiuj)=??p′?xi+?τij?xj+ρgi+Si(6)ρ??xj(uiuj)=-?p′?xi+?τij?xj+ρgi+Si(6)式中:uj為j方向上采空區(qū)某點(diǎn)風(fēng)流平均流速,m/s;p′為混合氣體的靜壓;τij為切應(yīng)力張量;gi為i方向的重力體積力(可以忽略);Si為i方向的外部體積力源項(xiàng)。采空區(qū)由于流體黏性以及顆粒存在,會(huì)造成流體動(dòng)量的損失。眾所周知,層流狀態(tài)(即小雷諾數(shù))下,多孔介質(zhì)內(nèi)單位長度壓降遵從達(dá)西滲流定律:V?=?Eμ?p(7)V→=-Eμ?p(7)式中:E為滲透率,二階張量,m2,對(duì)于各向同性介質(zhì),E蛻化為標(biāo)量e;?p為單位長度流體總壓力損失,Pa/m。Ergun總結(jié)了大雷諾數(shù)下,多孔介質(zhì)中單相流半經(jīng)驗(yàn)非線性滲流公式,對(duì)該方程略作改動(dòng),擴(kuò)展為三維空間,對(duì)于各向同性的多孔介質(zhì)單相流壓降方程可表達(dá)為??p=[150μD2m(1?n)2n3+1.75ρDm(1?n)n3V]V?(8)-?p=[150μDm2(1-n)2n3+1.75ρDm(1-n)n3V]V→(8)其中,V為的模量,m/s,當(dāng)V趨近于0時(shí),也趨近于0,此時(shí)式(8)即轉(zhuǎn)變?yōu)槭?7)。由于采空區(qū)同時(shí)存在紊流、過渡流及層流,所以在采空區(qū)運(yùn)用式(8)概括是合適的。實(shí)際上,式(8)即為風(fēng)流在多孔介質(zhì)中的動(dòng)量守恒方程,表示風(fēng)流受多孔介質(zhì)的影響而造成的動(dòng)量損失與風(fēng)速的關(guān)系。結(jié)合式(4)—式(8)可以求得風(fēng)流的運(yùn)動(dòng)狀態(tài)。4多孔介質(zhì)模型采用Fluent軟件對(duì)采空區(qū)風(fēng)流流場進(jìn)行數(shù)值模擬時(shí),由于采空區(qū)存在多孔介質(zhì),所以多孔介質(zhì)的參數(shù)設(shè)定是數(shù)值模擬的重點(diǎn)。多孔介質(zhì)需輸入的主要參數(shù)有孔隙率、黏性阻力系數(shù)或慣性阻力系數(shù)。1區(qū)域合作擬合基于上述論述,孔隙率可以用碎脹系數(shù)表示,而碎脹系數(shù)可以利用試驗(yàn)擬合。所設(shè)定的值因區(qū)域而異。由于取值不固定,需要利用udf(usersdefinedfunction,用戶自定義函數(shù))導(dǎo)入Fluent。2多孔介質(zhì)模型在利用Fluent進(jìn)行數(shù)值模擬時(shí),需要設(shè)置多孔介質(zhì)的參數(shù)。根據(jù)Fluent的規(guī)定,多孔介質(zhì)相當(dāng)于一個(gè)動(dòng)量源,對(duì)于簡單的多孔介質(zhì)有:S?=?p=?(C1iμ+C2iρV/2)V?(9)S→=?p=-(C1iμ+C2iρV/2)V→(9)其中:S?S→為以Si為標(biāo)量的多孔介質(zhì)動(dòng)量源矢量,N/m3;C1i為黏性阻力系數(shù);C2i為慣性阻力系數(shù)。對(duì)于各向同性的多孔介質(zhì),各個(gè)方向阻力系數(shù)相同,可將其下標(biāo)去掉。對(duì)比式(8)、式(9)可知,C1=150(1-n)2/(D2mn3)=1/e,C2=3.5(1-n)/(Dmn3)。阻力系數(shù)同樣需要利用udf函數(shù)輸入。5使用實(shí)例5.1工作面動(dòng)煤巖碎脹系數(shù)實(shí)測分析陽泉三礦主采15號(hào)煤為高瓦斯易自燃煤層。該煤層總厚7.08m,煤凈厚6.71m,煤層傾角1～7°,平均5°,頂板為砂質(zhì)頁巖。其中K8206工作面傾向長252m,走向長1579m;采用綜采放頂煤的采煤方法,采高2.8m,工作面平均推進(jìn)速度為2～4m/d;采用U型通風(fēng)方式,工作面配風(fēng)量1650m3/min。由于未考察該工作面的采空區(qū)滲流特性,筆者借鑒該礦80606工作面采空區(qū)現(xiàn)場實(shí)測及相似試驗(yàn)結(jié)果。采空區(qū)自然堆積區(qū)和O形壓實(shí)區(qū)碎脹系數(shù)Kpa1=1.58,Kpa2=1.25,Kpc=1.11。承壓碎脹區(qū)的碎脹系數(shù)介于自然堆積區(qū)和壓實(shí)區(qū)的碎脹系數(shù)之間,煤巖碎脹系數(shù)見表1。據(jù)考察,碎脹系數(shù)在由自然堆積區(qū)到壓實(shí)區(qū)是逐漸變小的,因此可將承壓碎脹區(qū)劃分為多個(gè)環(huán)形區(qū)域,每個(gè)區(qū)域的碎脹系數(shù)近似認(rèn)為是相同的,由采空區(qū)邊至采空區(qū)中部遞減。區(qū)域劃分及參數(shù)設(shè)置利用udf編寫。5.2采空區(qū)注意事項(xiàng)如圖3所示,模擬回采距開切眼150m的情況,坐標(biāo)原點(diǎn)為開切眼與回風(fēng)巷底板交點(diǎn)處。模型沿采空區(qū)中部對(duì)稱。進(jìn)口風(fēng)速3.05m/s,氣體成分不含甲烷。幾何邊界參數(shù)根據(jù)相似試驗(yàn)?zāi)M的距離確定,見表2。其他參數(shù)見表3。根據(jù)采空區(qū)防滅火“三帶”的標(biāo)準(zhǔn)之一,利用風(fēng)速不同,可將采空區(qū)劃分為冷卻帶、氧化帶、窒息帶,對(duì)應(yīng)的“三帶”的風(fēng)速范圍分別為>0.24、0.10～0.24、<0.10m/min。只有在氧化帶內(nèi)的煤體可能發(fā)生氧化自燃,所以確定“三帶”范圍是非常必要的。由于K8206遺煤厚約0.5m,所以截取高度0.5m的風(fēng)速面,查看三帶分布是合理。根據(jù)圖4可知,采空區(qū)氧化帶分布幾乎對(duì)稱,進(jìn)、回風(fēng)巷氧化帶范圍距支架后方距離分別為30～90、20～88m,采空區(qū)中部的氧化帶距支架8～40m。由于O形圈的存在,采空區(qū)中部煤巖較壓實(shí),孔隙率較小,阻力較大,風(fēng)速梯度大,而且該位置的工作面和采空區(qū)壓差很小,造成距支架很近的地方風(fēng)速急劇減小至0.24m/min。因此形成了采空區(qū)中部氣流較緩慢,寬度較小,而兩端較大的趨勢。綜合起來看采空區(qū)氧化帶為8～90m。距底板高1m時(shí)瓦斯體積分?jǐn)?shù)分布如圖5所示。由圖5可知,采空區(qū)回風(fēng)巷瓦斯?jié)舛茸罡?尤其是末端,達(dá)到9.5%以上,從采空區(qū)回風(fēng)巷末端至上隅角瓦斯體積分?jǐn)?shù)降至5%,瓦斯?jié)舛热匀缓艽?由于通風(fēng)已經(jīng)很大(規(guī)定小于4m/s)仍滿足不了治理瓦斯要求,因此必須采取相應(yīng)措施。5.3帶覆蓋區(qū)域防滅火措施根據(jù)以上分析,針對(duì)自然發(fā)火,在初采期就必須建立預(yù)測預(yù)報(bào)系統(tǒng),采取均壓防滅措施或設(shè)置堵漏風(fēng)裝置;加快工作面的推進(jìn);由于三帶范圍為8～90m,在距支架后方8m處就要設(shè)置防滅火管口,例如注氮防火或灌漿滅火。對(duì)于瓦斯治理,由于風(fēng)排瓦斯量局限,應(yīng)向采空區(qū)及工作打鉆孔進(jìn)行采空區(qū)抽采及預(yù)抽,或者掘巷進(jìn)行立體抽采。采取