沿空留巷采空区煤自燃堵漏控风机制数值模拟研究

司俊鸿，王乙桥，程根银，李 林，高世杰

(1.华北科技学院 应急技术与管理学院，北京 101601；2.内蒙古伊泰煤炭股份有限公司宏景塔一矿，内蒙古 鄂尔多斯 017000)

沿空留巷是2个相邻工作面共享1条回采巷道、采空区不留设保护煤柱的煤炭地下无煤柱开采技术之一。受采煤工艺、机械化程度、支护材料的限制，初期的沿空留巷技术仅适用于薄煤层。沿空留巷的巷旁支护方式有木垛、矸石袋、单体支架、钢筋混凝土墙等，切顶卸压自动成巷无煤柱开采技术基于“长壁开采切顶短臂梁”理论形成，通过在采空区侧定向切顶卸压，可消除煤体上方集中应力[1-2]，具有合理开发煤炭资源、提高煤炭回采率、降低巷道掘进率、改善矿井采掘接替紧张的特点，能够提高煤矿的经济效益[3]。

沿空留巷工作面通常采用“两进风、尾巷回风”的Y型通风方式，该方式下采空区密封性较差、漏风情况严重[4-6]，尾巷侧采空区的开放程度会影响采空区的风流运移规律[7-8]。由于采空区位置特殊、人员无法进入、监测设备布置困难，故对采空区风流场的研究主要采用现场实测、实验室相似模拟与数值模拟有限元分析相结合的方法[9-12]，分析不同工作面通风参数、不同采空区漏风强度条件下采空区内的气体运移规律[13-14]。

根据煤自燃氧化理论，可知漏风向采空区供氧是遗煤自燃的主要原因之一。研究采空区漏风场、氧气浓度场、温度场等多场耦合作用下的采空区气体流动及传热传质规律[15-16]，探究工作面回采速度、漏风速率、遗煤放热量与温度场、氧气浓度场、速度场等参数之间的关系[17-18]，是研究采空区煤自燃问题的关键。沿空留巷Y型通风改变了采区开拓布置，解决了上隅角瓦斯浓度超限问题[19]、缩短了通风路线、降低了矿井通风阻力，且由于材料巷和运输巷均在进风流中，改善了煤矿井下作业环境。

随着煤炭开采工艺的革新，准确掌握沿空留巷采空区风流运移规律、确定合理的尾巷侧采空区开放程度、分析采空区的风流控制机制，是研究开放式采空区煤自燃氧化生热、聚热及风流散热的平衡关系中亟待解决的主要问题，对预防和治理采空区内遗煤自燃具有重要意义。

1 矿井概况

图1 316综采工作面通风系统示意图

2 沿空留巷采空区控风机制

2.1 风流运移影响因素

沿空留巷采空区风流运移，属于多孔介质固—气耦合问题，其影响因素主要包括采空区结构因素、流体因素、物理化学因素，以及采矿方法和支护方式等[20]。开放式沿空留巷采空区风流运移的主要影响因素指标如图2所示。

图2 风流运移影响因素指标框图

在开采工艺、支护手段及通风方式确定的情况下，主要分析多孔介质特征(渗透率)、气体成分、压力损失对采空区风流运移的影响规律。

2.2 采空区控风机制

根据多孔介质固—气耦合气体流动理论，沿空留巷采空区漏风源为工作面进风和相邻采空区漏风，漏风汇合为尾巷侧回风，采空区的漏风风量与漏风风阻、漏风源汇间的阻力损失满足阻力定律，漏风源汇间的阻力损失越大，采空区的漏风风量越大；采空区漏风风阻与多孔介质漏风通道的导通孔隙尺寸大小呈负相关关系，导通孔隙越大，漏风风阻越小，漏风风量越大。因此，增大漏风风阻，或降低漏风源汇间的阻力损失，可降低采空区的漏风风量。

3 沿空留巷采空区渗透率分布

渗透率分布是影响采空区风流场模拟结果的重要指标。由于沿空留巷采空区渗透率分布具有不对称的特点，以下隅角(主要进风巷和316综采工作面交汇处)作为物理模型的坐标原点，定义x轴方向为采空区深度方向，y轴方向为工作面倾斜方向，z轴方向为高度方向(见图1)。

根据Ergun方程，采空区渗透率与孔隙率满足：

(1)

式中：k为渗透率，m2；Dm为采空区破碎煤岩平均粒径，m；ε(x,y,z)为采空区(x,y,z)坐标点的孔隙率。

根据砌体梁及关键层理论，建立工作面采空区三维孔隙率分布数学模型[21]：

(2)

式中：m为煤层厚度，m；W为上覆岩层下沉量，m；H为采空区岩块垮落后的高度，m；a、b、c和d为待定系数。

经测定，工作面岩块碎胀系数Kp=1.29，关键层的厚度∑h=6.28 m，关键层的拉伸强度αT=1.3 MPa，岩块自身质量及其上部岩层质量q=2.36×109kg，沿空留巷侧顶板垮落后尾巷附近煤岩体孔隙率小于留煤柱采空区。

计算得到关键层最大下沉量Wmax=0.97 m，Dm=0.4 m，H=10.59 m，y方向上与坐标原点的距离Ly=246 m。根据沿空留巷采空区的边界条件：

(3)

由式(3)求解得a=ln[m/(m-W(x))]/Ly，b=0，c=[1-ε(x,y)]/[Hε(x,y)]，d=ε(x,y)。

推导得到316综采工作面底板处(z=0 m)的采空区渗透率分布函数公式：

(4)

采空区渗透率分布云图如图3所示。

图3 采空区渗透率分布云图(z=0 m)

4 沿空留巷采空区物理模型的建立

4.1 物理模型和网格划分

采空区物理模型可分为综采工作面、主要进风巷、辅助进风巷、尾巷和采空区，沿空留巷采空区物理模型网格划分如图4所示。

图4 沿空留巷采空区物理模型网格划分

采用Tet/Hybrid网格划分方法，在综采工作面、采空区靠近进风巷、尾巷区域采用了网格局部加密方法，具体尺寸及网格划分情况如表1所示。

表1 沿空留巷采空区物理模型参数

4.2 边界条件设置

各区域均为流体类型(Fluid Type)，将采空区设置为多孔介质，渗透率表达式见式(4)。采用自定义函数(UDF)编写代码，加载到采空区多孔介质特征中。模型边界设置如表2所示。

表2 模型边界条件

5 沿空留巷采空区风流运移规律

5.1 采空区风流场分布

采用Fluent数值模拟软件，计算宏景塔一矿316综采工作面的风流分布，绘制出距底板1 m水平位置(z=1 m)采空区风速及O2体积分数云图，如图5 所示。

(a)风速

由图5(a)可以看出，沿空留巷侧的采空区漏风风速大于非留巷侧。受风流的惯性影响，在下隅角处向采空区深部方向的漏风风速较大，形成了风流冲击区；从工作面中部111 m开始，风速等值线间距变宽，此区域的风速逐渐增大，大于3.57×10-5m/s。

开放式沿空留巷尾巷毗邻采空区，属于采场低压能区，根据采空区“O”形圈理论，采空区漏风由工作面进入采空区后，受中间压实区的影响，气体向尾巷侧运移，氧化带呈现L型分布，如图5(b)所示。

5.2 沿空留巷采空区“三带”划分

为了清楚划分Y型通风方式下采空区“三带”分布，将沿空留巷采空区的氧化带按工作面倾斜方向和尾巷平行方向分为横向氧化带和竖向氧化带。横向氧化带位于采空区高压能区域，漏风强度较大，氧化带较宽；竖向氧化带位于采空区低压能区域，氧化带宽度较小。

现实中，判定采空区O2体积分数分布比较困难。基于数值模拟结果，提出采空区横向氧化带平行四边形面积判定法，采用横向氧化带与工作面倾斜方向的夹角α和非沿空留巷侧采空区氧化带宽度Hp构成的平行四边形面积作为横向氧化带范围划分依据，为现场采空区“三带”划分提供参考。采空区氧化带定位如图6 所示。

图6 采空区氧化带定位图

根据图5(b)中O2体积分数划分采空区“三带”，316综采工作面的横向散热带在靠近工作面 21 m 内，纵向散热带在靠近尾巷8 m内；横向氧化带平均宽度为100 m，纵向氧化带平均宽度为4 m。根据平行四边行面积判定法可知α=6°、Hp=100 m、氧化带区域面积为2.34×104m2。

6 沿空留巷采空区堵漏控风试验

控制漏风源是解决漏风问题的重要途径，提出封堵下隅角(方案1)、封堵上隅角(方案2)、同时封堵上下隅角(方案3)共3种漏风源堵漏方法，每个位置封堵幕墙的长度均为10 m，采用数值模拟方法研究沿空留巷采空区的风流变化规律。

绘制出3种控风条件下z=1 m水平的O2体积分数分布云图，如图7所示。

(a)封堵下隅角

漏风源封堵后的横向、纵向氧化带分布情况如表3所示。

表3 漏风源封堵后横向、纵向氧化带分布情况

由表3可知，采用漏风源堵漏方法后，3种封堵方案的横向氧化带平均宽度分别减小了29、28、26 m，氧化带平行四边形夹角分别增加了9°、6°、4°，纵向氧化带平均宽度均为4 m。因此，漏风源堵漏方法减小了横向氧化带的平均宽度，降低了煤氧化自燃的概率，但对纵向氧化带的影响较小。

由图7可以看出，采空区两侧O2体积分数的边界，能直观表现氧化带区域变化规律，在316综采工作面走向236、1 m处，分别沿采空区深部方向间隔10 m提取未堵漏控风(用0表示)及3种封堵方案的O2体积分数数据，绘制O2体积分数随采空区深度变化曲线，如图8所示。图中“1-0”代表在x=1 m处未采用堵漏控风措施时的O2体积分数，“1-1”代表在x=1 m处采取堵漏控风方案1时的O2体积分数，其他以此类推。

图8 O2体积分数随采空区深度变化曲线

由图8可知，采用漏风源堵漏控风方法后，采空区各采样点的O2体积分数均减小，3种封堵方案在236 m处O2体积分数的最大变化量分别位于采空区105.26、126.32、105.26 m处，其值分别为-0.06、-0.04、-0.05；在1 m处O2体积分数的最大变化量分别位于采空区157.90、168.42、157.90 m处，其值分别为-0.035 9、-0.033 6、-0.047 8。

沿采空区深度方向，分别在316综采工作面走向236、1 m处提取O2体积分数为18%和8%的点的y轴坐标，整理得到不同方案采空区“三带”的划分范围，结果如表4所示。

表4 不同方案采空区“三带”分布

由表4可以看出，位于316综采工作面1 m处，采取封堵下隅角(方案1)时，能使1 m处O2体积分数为18%的点靠近工作面；采取封堵上隅角(方案2)时，能使1 m处O2体积分数为18%的点远离工作面；采取3种方案封堵后O2体积分数为8%的点均靠近工作面。在工作面1 m处，方案2的横向氧化带宽度变化最大，减小到了60.30 m。

位于316综采工作面236 m处，采取3种方案后O2体积分数18%的点均远离工作面；采取3种封堵方案后O2体积分数为8%的点均靠近工作面。在工作面236 m处，方案1的横向氧化带平均宽度变化最大，减小到了80.40 m。

根据图5和图7计算氧化带区域面积，结果见表5。

表5 不同方案氧化带区域面积

由表5可以看出，氧化带区域面积相对变化量最大的是封堵下隅角(方案1)，减小了13.80%，故采用方案1的采空区氧化带区域面积变化最大。

7 结论

1)建立了开放式沿空留巷采空区风流运移影响因素指标，其影响因素主要包括采空区结构因素、流体因素、物理化学因素，以及采矿方法和支护方式等。

2)将沿空留巷采空区的氧化带按工作面倾斜方向和尾巷平行方向分为横向氧化带和竖向氧化带。提出了采空区横向氧化带平行四边形判定法，采用横向氧化带与工作面倾斜方向的夹角α和非沿空留巷侧采空区氧化带宽度Hp构成的平行四边形面积作为横向氧化带范围划分依据。

3)采用堵漏控风方法研究了沿空留巷采空区风流运移规律，结果表明堵漏措施减小了横向氧化带的宽度、增大了横向氧化带与工作面的夹角，不影响纵向氧化带的宽度。封堵下隅角时，采空区氧化带区域面积变化最大，减小了13.80%。