基于有限体积法的采空区CO浓度场模型构建

王月红，蒋冀萍，高萌，周宁

(华北理工大学 矿业工程学院，河北 唐山 063210；2.河北省矿山开发与安全技术重点实验室，河北 唐山 063210)

引言

采空区煤自燃发火是煤矿重大灾害之一，煤自燃往往会产生有害气体(CO)，CO是一种特殊的气体产物，伴随在整个煤氧复合反应过程中，被公认为是非常有效的预测预报自燃火灾的指标气体[1]。为了更好地了解煤自燃的发展过程，有效遏制火灾事故的发生，采空区自然发火数值模拟研究受到更多重视。20世纪80年代以来，国内外学者针对采空区或地面煤堆的自燃条件建立了多种煤自然发火数学模型，用来预测采空区或煤堆的自然发火危险性[2-8]，邵磊[9]通过建立采空区流场物理数学模型对采空区的氧浓度进行了数值模拟，Schmal D[10]给出了储煤煤堆的数学模型。该研究应用优化后的有限体积算法来构建采空区CO浓度场模型，能够在物理意义上表示控制体的通量平衡，为预测矿井火灾问题提供了有效的新思路。

1有限体积算法

基于物理量守恒原则的有限体积算法在流体工程领域中得到了广泛的应用，从整体意义上讲虽然选取不同的控制体是等价的，但在每个节点上的条件不是“全等”的，最终形成不同的控制方程，所以控制体的选取直接关系到有限体积法的计算精度，如图1所示为不同控制体的确定方法。

图1 控制体示意图

如图1所示，图1(a)是选取以P0为顶点的三角单元，六边形M1M2M3M4M5M6为控制体的确定方法；图1(b)是选取以P0为顶点的三角单元，六边形P01P02P03P04P05P06为控制体的确定方法；图1(c)是选取以zi为顶点的三角单元，多边形M1Q1M2Q2M3Q3M4Q4M5Q5M6Q6为控制体的确定方法。

本文根据图1(d)图的方法选择有限体积法的控制体，即该区域由平行于副线的折线圈出，每一个折线通过三角形元素的重心。采取这种新的控制体方法，不仅赋予了每部分明确的物理意义，同时反映了其科学性和合理性。因此基于有限体积法建立数学模型，既避免了研究泛函的问题、简化了建模过程，又能保证计算精度。

2采空区CO浓度场模型

2.1 CO浓度场方程

运用有限体积法建立采空区氧浓度场数值模型，对采空区CO模型进行研究。根据质量守恒定律：在单位时间内，漏风风流通过任意封闭区域，CO质量的变化量(WC)主要由以下几个方面引起：(1)由于渗透风流中的CO流出与流入任意封闭曲面之差(W1)；(2)CO气体生成量(W2)；(3)浓度差异引起的弥散进出任意封闭曲面的CO之差(W3)，即：

WC=W1+W2+W3

(1)

其中：

(2)

(3)

(4)

(5)

将式(2)、(3)、(4)以及(5)代入(6)得：

(6)

为了便于计算，将式(6)中质量浓度改为摩尔浓度，即为式(7)。

(7)

式中：Cco为氧气浓度，mol/m3；Kco为氧气扩散系数常数；v为通过采空区内某点气体的渗流速度，m/s；u(t)为松散煤体内的氧气消耗速度，m/s。

2.2 数值模型离散方程

移动坐标下采空区的CO浓度方程为：

(8)

对图1(d)中的内部节点进行质量守恒分析，多边形ABCDEF控制体为分析对象。此时，式(11)中的边界Г为图1(d)中的边界A-B-C-D-E-F-A，即6个小三角形拼接的多边形ABCDEF。于是式(9)可表示为：

(9)

式中：Cco为节点5、7的第k个小控制体的平均浓度，mol/m3；Δyk，Δxk为节点5、7的第k个小控制体的对边在x，y轴摄影长度，m；ΔSk为节点5、7的第k个小控制体的面积，m2。

根据式(9)，任取一内部节点l，设其与n个单元相关联，则可建立如下质量方程：

(10)

设第k个三角形单元的3个节点分别为i、j、m，下面分析第k个三角形单元对m点质量方程的贡献Wmk：

(11)

同理可得出第k个三角形单元对i、j点能量方程的贡献Wik、Wjk：单元k对其3个节点质量守恒方程贡献，写成矩阵表达式：

(12)

式中：

(l,n=i,j,m)

(13)

式中单元的孔隙率(n)可根据单元的3个节点孔隙率确定大小，系数Kco根据实际情况而定，同样办法对单元速度(vk)进行处理。

3 CO浓度场数学模型中参数的确定及模拟

3.1 实验装置及过程

该实验采用煤缓慢升温氧化实验，实验主要仪器由煤升温氧化炉、实验温度控制仪、气相色谱仪、计算机采集信息系统、数据处理系统等部分组成，其余实验装置有锤子、筛子和天平。

实验选取某煤矿4个地点的煤样，将选取的大块煤样捣碎筛分出相同粒径的实验煤样，用天平称取130 g后放入煤样罐，用于实验。煤样准备好后，将煤样罐装入煤样氧化炉，实验过程中通入的气体为空气。通过温度控制器预先设定好恒温箱内温度和采样温度，当恒温箱内温度达到预定温度时，测量煤样温度，同时利用气相色谱仪取其气样进行气体分析。

3.2 基础数据分析

在煤升温缓慢氧化实验中，主要测得的是氧化炉出口处的气体成分及浓度百分比。图2所示为通过origin软件将所得到的实验数据绘制的温度与CO、CO2和O2消耗量的关系曲线。

图2 温度与气体消耗量变化关系

由图2(b)知，当煤样温度较低时，没有CO气体产生，CO气体的百分比为0；煤样温度大约在90 ℃时，CO气体开始产生；温度为160 ℃时，氧气的消耗量突然增大，同时CO和CO2气体产生量突然增加，且由图中数据得出，CO2气体产生量要大于CO气体产生量。

3.3 放热强度数据分析

煤体在氧化和自燃过程中持续放热，采空区温度上升所需能量主要是来自于这一过程。煤体的放热强度可以表征其自燃剧烈程度。根据化学动力学和化学平衡理论，由实验得出所需数据，计算出煤体的放热强度q(t)，计算公式如式(14)。

(14)

将实验所测得煤样的放热强度转换到氧气浓度为21%新鲜风流中的放热强度，计算公式见式(15)所示。

q0(t)=q(t)×(Cin/Cave)

(15)

式中：Cave为进出口的氧气平均浓度，%。Cave=(Cin+Cout)/2

计算结果如图3所示。

图3 放热强度q0(t)与温度t(℃)之间的关系图

由图3分析发现，放热强度q0(t)与温度t(℃)两者之间的关系近似指数关系，但经计算得，两者的相关系数较小，小于0.9。若将其设为多项式关系时，虽然相关系数比较大，但考虑到计算简便问题，还是将其看做指数关系。故设其方程式如式(16)所示。

q0(t)=m×en×t

(16)

其回归方程相关系数如表1所示。

表1 放热强度与温度之间的相关系数

由图3和表1可知，在环境温度较低时，煤样的放热强度低；当温度上升到一定数值时，放热强度增大，且增大的速率也越来越大。这正符合放热强度与温度之间的指数关系。

建立采空区自然发火模型所确定的参数，如表2所示。

表2 基本数据

3.4 采空区CO浓度场数值模拟

程序设计界面基于VB编程软件，能够简单有效地将模型方程解算出结果。根据得出的CO浓度场数学模型，利用VB软件编写程序解算方程。设计输入参数页面如图4所示。

图4 程序设计界面

在实验的基础上，Tecplot9.0软件能够将解算的结果模拟出图，可以直观清晰地看出模拟结果，建立的CO浓度场模型如图5(c)所示，相同参数条件下的采空区氧浓度场和温度场如图5(a)、图5(b)所示。

图5 工作面长度为100m时氧浓度场、温度场、CO浓度场

4各因素对采空区CO浓度场的影响

在采空区内，由于影响CO浓度的因素较多，在这里选取工作面长度(m)和推进速度(m/d)2种因素进行对比，观察CO浓度场的变化情况。

4.1 工作面长度改变

在采空区其他参数不变的情况下，改变工作面长度，设置数值为50 m，工作面长度为50 m时CO浓度场、氧浓度场、温度场模拟结果如图6所示；设置数值为150 m时CO浓度场、氧浓度场、温度场模拟结果如图7所示进行方程解算。利用tecplot9.0模拟软件进行模拟。

图6 工作面长度为50 m时模拟结果

图7 工作面长度为150 m时模拟结果

4.2 推进速度改变

在采空区其他参数不变的情况下，改变推进速度，分别设置数值为2.0 m/d和2.8 m/d，进行方程解算。利用tecplot9.0模拟软件进行CO浓度场、氧浓度场和温度场模拟。模拟结果如图8、图9所示。

图8 推进速度为2.0 m/d时模拟结果

图9 推进速度为2.8 m/d时模拟结果

4.3 结果对比及分析

(1)根据图5、图6和图10，对比CO浓度场模型，当工作面长度为50 m时，CO浓度最高部位集中于采空区中后部，且CO浓度要高于工作面长度为100 m时浓度上隅角部位的最高点浓度，但上隅角部位的浓度不高。对比氧浓度场模型和温度场模型可知，工作面长度为50 m时，采空区温度降低，氧气浓度升高，氧气的消耗变慢，CO生成速率降低。分析出现该种情况的原因为：a. 当工作面长度较小时，生成的CO浓度增大；b. 新鲜风流从进风巷进入采空区后不利于向回风巷一侧流动，导致CO积聚在采空区图中部位。

图10 各因素变化时的CO最高浓度点值变化图

(2)根据图5、图7和图10，对比CO浓度场模型，当工作面长度为150 m时，采空区CO浓度最高部位同样出现在采空区上隅角且最高点的浓度变化不大，但采空区回风巷一侧的中部位置浓度也相对变高。对比氧浓度场模型和温度场模型，温度变化不大，但氧浓度相对较高。增大工作面长度，氧气的消耗速率降低导致CO的生成量降低，进而导致CO产生量降低，但由于通风原因，采空区回风巷一侧的浓度偏高。

(3)根据图6、图8，图9和图10，对比CO浓度场模型，当推进速度降低或增高时，CO的浓度变化不大，浓度最高点同样出现在采空区上隅角部位，最高点的CO浓度同样基本不变。不同的是当推进速度降低时，采空区中部位置的浓度会偏高一些，原因可能是当推进速度较慢时，采空区氧浓度相对较低，煤更易于发生不完全燃烧，此处产生相对较多的CO。

5结论

(1)基于有限体积法，得出采空区CO浓度场的数值模型。通过煤的缓慢氧化实验研究不同温度下煤的耗氧速度和放热强度，分析、整理实验数据，发现放热强度与温度之间的关系成指数关系。

(2)当工作面长度较小时，采空区的中后部易积聚CO且浓度最高；当工作面的长度增加时，CO的浓度逐渐积聚在采空区上隅角部位；采空区工作面长度继续增大，超过一定数值时，采空区的最高浓度点位置不变，最高浓度值不变，但回风巷一侧的浓度逐渐增高。

(3)改变推进速度时，CO浓度场分布变化不大，且浓度值变化不大，但推进速度较小时，也可能会造成CO在采空区回风巷一侧的积聚。