基于MATLAB的氧浓度对采空区“三带”划分的数值模拟

于 敏，栗继祖，康立勋，冯国瑞，张晓亮

(1.太原理工大学矿业工程学院，山西太原030024;2.太原理工大学科技处，山西太原030024; 3.国电永寿煤业有限责任公司，陕西咸阳713401)

由于现代煤矿有逐步走向大型化的趋势，使得工作面的数量和搬家的次数逐渐减少，巷道的布置变得更加简化，从而减少了煤柱的预留，但也使矿压突显的情况变得非常严重，从而造成采空区漏风严重。当采空区的遗煤达到自燃所需的3个条件时(煤具有自燃倾向性，且以破碎形式存在;有持续供氧的条件;有积聚氧化热的环境，且上述3个条件需持续时间不得小于煤炭的自然发火期)［1］，就容易发生采空区遗煤的自燃。因此对采空区自燃“三带”，即散热带、氧化自燃带和窒熄带［2］的正确划分对预防采空区遗煤自燃是非常必要的。

由于采空区里面的实际情况复杂，采空区“三带”划分至今没有一个统一的标准，且现有的划分方法又具有一定的局限性，但大多数专家认可的划分方法有以下3种:

(1)采空区漏风风速V(散热带:V＞0.02m/ s;氧化带:0.02m/s≥V≥0.001m/s;窒熄带:V＜0.001m/s)。

(2)采空区氧气浓度C(窒熄带:C≤7%;氧化带:18%＞C＞7%;散热带:C≥18%)。

(3)采空区温度升高的速率 (目前对此划分指标的研究尚不成熟，很少作为“三带”划分的指标)［3－5］。

在工程实践中应用最为广泛的是采空区氧气浓度的划分标准，这也是目前最有效的划分方法。本文以采空区氧气浓度划分指标为依据，在数据分析部分，利用现场实测的数据，运用MATLAB强大的三维数据处理能力和ORIGIN的数值拟合，对采空区自燃“三带”进行划分。再利用GAMBIT和FLUENT来进行建模和模拟采空区的内部环境，所得到的采空区氧气浓度的模拟结果与采用MATLAB和ORIGIN数值拟合的采空区“三带”的划分结果进行比较，验证MATLAB和ORIGIN的数据分析方法对“三带”划分的准确性。

1 采空区氧气分布规律MATLAB数值模拟

山西鼎盛煤业位于太原西山煤田复式向斜的东翼，煤层基本上走向为北东－南西，向南东倾伏的单斜构造，地层倾角3°左右。1202工作面开采煤层为二叠系山西组2号煤层，煤厚3.0m左右，厚度稳定，全区可采结构简单，普含一层夹石，夹石厚度0.2～0.5m左右，煤层倾角平均5°左右。工作面走向长度850m，工作面采长175m。2号煤层顶板为泥岩或砂质泥岩，厚7m左右，煤层底板为细砂岩，厚度3m左右。2号煤层为自然发火煤层，发火期为50～109d;煤尘具有强爆炸性，其爆炸指数为40.32%。

1.1 采空区束管布置及数据的收集

沿1202工作面倾向布置5个测点，1号测点距工作面距离为20m，测点间距为30m，铺设单芯束管，束管用直径133mm(4寸)钢管保护，防止压坏束管，采样器探头部位用坑木加以保护。埋入采空区的采样器探头处与套管密封连接，密封段采用聚氨酯材料密封，附近采用大块矸石防护，以防止浮煤、水、泥浆堵管和抽取到套管内气体［6］。可实现同时对5个监测地点取气进行分析，及时发现发火预兆，束管布置见图1。1202工作面采空区氧气浓度的现场观测数据见表1。

图1 束管布置

1.2 对收集气体数据的处理和分析

在MATLAB中绘制立体浓度分布图时，工具库中的meshgrid函数用于生成平面网格点，利用mesh函数绘制由线条框构成的网格曲面，由于绘制的图形网格线条框分明，不能很好地体现采空区内气体浓度的渐变过程，因此用surf和shading interp函数来对生成的曲面进行网格内填充颜色和颜色插值处理，使得表面显得光滑。在绘制等值曲线图时，使用griddata函数对数据进行平面差分，利用contour函数进行等值线图绘制［7－9］。通过对实测数据的处理，得到了采空区氧气浓度立体分布图 (图2)和采空区氧气浓度等值曲线图 (图3)。

表1 实测氧气浓度值 %

图2 采空区氧气浓度立体分布

由图2和图3可以得出以下结论:

(1)在采空区后端的等高线比前端的等高线密集，这是由于采空区末端岩石处于重新压实区，空隙率低，供氧条件变差，致使氧气浓度迅速降低。

图3 采空区氧气浓度等值曲线

(2)采空区氧气浓度随工作面的距离增加而减小，进风处的氧气浓度降低速度相对于出风口的降低速度小，这是因为采空区进风处的顶板处于自然堆积区，垮落后具有相对较好的漏风通道，且在进风巷处具有很好的供氧条件。

(3)在采空区中部，相对于同一距离的氧气浓度中部小于两侧，这是由于在采空区中部处于水平方向的重新压实带，使得岩石被压实导致中部的空隙率和黏性阻力系数相对变大，氧气浓度降低。

(4)工作面推进到36～59m后，氧浓度降低到18%，从而进入到氧化带，在此之前是在散热带内，当工作面推进到163m以后，5号至1号的氧浓度逐渐下降至7%左右，随后进入到窒息带。

为了得到各个束管测点在某个氧气浓度时的确切的位置，本文利用ORIGIN软件将各个束管测点监测到的氧气浓度的实际数据进行了数值拟合，通过拟合曲线方程来确定各束管测点在采空区氧气浓度的详细分布情况，由于篇幅限制，本文只显示1号、2号测点的拟合曲线，见图4、图5。

图4 1202工作面1号测点数据拟合曲线

根据各测点氧气浓度拟合曲线得到采空区内氧浓度分布，见表2，表中的数据为距工作面的距离，单位为m。

图5 1202工作面2号测点数据拟合曲线

1.3 对采空区“三带”的划分

表2 采空区内各测点氧浓度分布

根据表2可以得出依据采空区氧气浓度所划分的采空区“三带”等值线图，如图6。由图6中采空区自燃“三带”的划分可以看出，采空区散热带宽度范围为36～59m，窒熄带是在距工作面200m之后的范围，氧化带介于散热带和窒息带之间，宽度约为130～150m。

图6 采空区“三带”分布

2 采空区氧气分布规律CFD模拟［11－12］

2.1 建立1202工作面模型

根据1202综放面的特点，利用GAMBIT建立工作面模型，模型参数为:采空区走向260m，倾向185m;巷道断面宽4m，高3m;工作面断面宽8m，高3m;模型高15m，U型负压通风。工作面的风量为 1440m3/min，采空区瓦斯涌出量为32m3/min，工作面煤壁瓦斯涌出量为9m3/min。各气体成分的体积分数为甲烷0.4%，氧气21%，其余为氮气。一般数据参数设置见表3。

2.2 采空区氧气浓度的数值模拟及分析

经过FLUENT模拟，得到采空区氧气浓度模拟分布图(图7)。根据采空区氧气浓度C(窒熄带:C≤7%;氧化带: 18%＞C＞7%;散热带:C≥18%)的划分标准，由图7可知“三带”的划分范围。散热带的范围＜57m;氧化带的范围是57～195m;窒息带的范围则是＞195m。与MATLAB进行的数值模拟划分“三带”结果相符合，由此可以证明MATLAB划分采空区“三带”的数值模拟是具有可行性的。

表3 采空区模型参数

图7 采空区氧气浓度模拟分布

3 结论

(1)利用MATLAB对数据的强大处理功能得到该工作面采空区氧浓度的立体分布和等值线图，利用ORIGIN得到束管测点氧浓度的拟合曲线方程，得到一种基于MATLAB的采空区“三带”划分的新方法。

(2)通过与CFD模拟结果相比较，二者以氧浓度为标准的三带划分结果大致接近，说明基于MATLAB的三带划分结果的准确性，且相对于CFD模拟此方法具有简洁、明了、易观察和不受地质条件影响的特点，显示出此方法的优越性。

(3)以氧浓度为采空区“三带”划分指标，利用MATLAB得到采空区“三带”和各测点“三带”的划分范围，判定结果可以为采空区防灭火提供直观、形象的指导。

［1］王省身，张国枢.矿井火灾防治［M］.徐州:中国矿业大学出版社，1990.

［2］王红刚.采空区漏风流场与瓦斯运移的叠加发放研究［D］.西安:西安科技大学，2009.

［3］杨胜强，张人伟，邸志前，等.综采面采空区自燃“三带”的分布规律［J］.中国矿业大学学报，2000，29(1).

［4］李宗翔，吴 强，王志清.自燃采空区耗氧－升温的区域分布特征［J］.煤炭学报，2009，34(5):667－672.

［5］徐精彩，余 锋，李树刚，等.石嘴山二矿2268综放采空区自燃危险区域划分研究［J］.煤炭学报，2003，28(3).

［6］张 勋.塔山煤矿综放面采空区自燃“三带”研究［D］.阜新:辽宁工程技术大学，2011.

［7］余明高，常绪华，贾海林，等.基于Matlab采空区自燃“三带”的分析［J］.煤炭学报，2010，35(4):600－604.

［8］常绪华.基于MATLAB综采面采空区自燃特性分析研究［D］.焦作:河南理工大学，2010.

［9］欧阳联华，王家林.一种在等值线图上任意截取剖面的MATLAB语言算法［J］.物探化探计算技术，2003，25(3).