动态回采工作面煤壁瓦斯涌出数值计算及现场应用

刘 彦 青

(1.煤炭科学技术研究院有限公司 安全分院，北京 100013；2. 煤炭资源高效开采与洁净利用国家重点实验室,北京 100013)

0 引 言

回采工作面瓦斯涌出量预测是回采工作面配风量确定以及回采煤层瓦斯抽采设计最重要依据之一，回采工作面煤壁为工作面最主要瓦斯涌出源，工作面煤壁瓦斯涌出预测对于实现工作面瓦斯涌出准确预测至关重要，目前所使用的矿井瓦斯涌出量预测方法中开采层瓦斯涌出量计算公式为半经验公式，无法实现精准预测。

为适应现阶段煤炭行业精准化开采发展趋势，以实现工作面瓦斯涌出量的准确预测为目标，众多研究人员采用数学建模、现场实测、数值模拟等方法对回采工作面煤壁瓦斯涌出规律进行研究：文献[1-2]基于气固耦合理论对工作面前面煤层应力场和瓦斯流场进行了模拟分析；程国强等[3]采用区别于传统网格数值模拟方法的SPH无网格模拟技术对煤层瓦斯流场进行模拟分析；刘世奇等[4]构建了CO2注入煤储层结构演化-流体运移-储存能力THM全耦合数学模型，并采用COMSOL Multiphysics仿真软件进行了模拟分析；文献[5]基于实测的单位面积工作面煤壁瓦斯涌出强度函数推导建立了整个回采工作面煤壁瓦斯涌出量计算公式；文献[6]基于回采工作面瓦斯、风量实测参数构建工作面瓦斯分源计算数学模型，计算回采工作面煤壁瓦斯涌出量。

在前人研究成果的基础上，笔者重点研究工作面推进过程对工作面煤壁瓦斯涌出影响规律，工作面推进速度的数量级为10-3～10-4m/s，与靠近工作面煤壁的煤层区域瓦斯渗流速度较为接近，而且工作面煤壁涌出瓦斯来源于靠近工作面煤壁的煤层区域，工作面推进速度必然影响工作面煤壁瓦斯涌出，当工作面推进速度较大时，靠近工作面煤壁处的煤层瓦斯不会通过煤层渗流涌出到工作面，而随采煤机采落形成块状落煤一起被运出工作面，工作面推进速度较大导致煤壁暴露时间短，煤壁瓦斯涌出量会增大，前人的研究中未考虑过该问题。以往的煤层瓦斯渗流计算模型中将煤层瓦斯近似为单一的甲烷气体成分，而实际煤层瓦斯中往往还含有一定量的氮气与二氧化碳成分，基于上述工作面煤壁前方煤层瓦斯渗流计算模型存在的问题，笔者构建考虑工作面推进动态过程的回采工作面煤层多组分瓦斯渗流计算模型，对回采工作面煤壁瓦斯涌出量进行准确预测。

1 回采工作面概况

以郑州煤电股份有限公司超化煤矿22051综采工作面为研究对象，工作面通风系统如图1所示。

图1 22051综采工作面通风系统示意Fig.1 Ventilation system in No.22051 coal face

该工作面为倾斜长壁式一次采全高综采工作面，采用“一进一回”U型通风方式，采用底抽巷穿层钻孔对回采煤层区域进行预抽，开采煤层倾角为3°～5°，煤层厚度为2.67～3.15 m。为了获得准确的22051工作面煤层瓦斯含量，通过钻孔采集煤样进行瓦斯含量测试及瓦斯组分分析，该工作面煤层处于甲烷带，回采工作面采煤工艺参数、煤层赋存参数、煤工业分析结果如下：

采煤工艺综采采出率/%95工作面长度/m195煤层平均厚度/m2.87煤层渗透率/m23.45×10-16煤孔隙率/%4.6煤的密度/(kg·m-3)1.425×103灰分/%9.76水分/%4.62煤层平均倾角/(°)4.5

煤层瓦斯压力与瓦斯组分分析结果见表1，煤对各种瓦斯组分气体的吸附解吸常数测试结果见表2。

表1 煤层瓦斯压力、组分含量Table 1 Gas pressure and component content in coal seam

表2 煤的不同组分气体朗格缪尔吸附常数Table 2 Langmuir adsorption constants of different gas components of coal

2 工作面煤壁多组分瓦斯涌出数值计算模型构建

工作面进回风巷掘进期间、工作面采掘设备布置期间，进风巷与回风巷两帮煤层瓦斯涌出使得进回风巷两帮形成一定宽度的瓦斯排放带，同时采动影响下工作面前方煤体应力重新分布必然导致工作面煤壁前方煤层渗透率发生显著变化，构建回采煤层瓦斯涌出过程简化模型,考虑进回风巷煤壁瓦斯涌出和工作面前方煤层渗透率变化，煤层厚度远小于工作面走向与倾向长度，简化为考虑工作面走向与倾向两方向煤层瓦斯渗流的工作面煤壁瓦斯涌出模型，图2为工作面煤壁瓦斯涌出数值计算物理模型。

图2 22051工作面煤壁瓦斯涌出过程物理模型Fig.2 Physical model of gas emission from coal wall of No.22051 coal face

2.1 工作面煤壁前方煤层多组分瓦斯渗流数学模型

2.1.1 煤层多组分瓦斯赋存模型

煤层中瓦斯赋存形式为吸附瓦斯与游离瓦斯，文献[7]中研究表明采用朗格缪尔等温吸附式描述煤层吸附瓦斯量与瓦斯压力之间关系最为准确，见式(1)。对煤层瓦斯赋存作如下假定：① 煤层瓦斯各组分气体之间吸附解吸过程互不影响，忽略竞争吸附解吸效应；② 不考虑煤层瓦斯吸附解吸过程中热变化，煤层孔隙-裂隙内游离气体遵循道尔顿混合气体分压定律；③ 煤层各处瓦斯吸附解吸过程在瞬间完成。

煤层瓦斯组分比值是不同瓦斯组分的游离量与吸附量之和的比值，且煤层对不同瓦斯组分气体的吸附解吸能力具有明显差异，导致煤层瓦斯组分气体的分压比值与含量比值不相同。基于表2中煤层瓦斯含量、煤层瓦斯组分占比实测数据和表3中煤的不同组分气体朗格缪尔吸附常数测试结果，根据煤的朗格缪尔气体吸附定律、煤层瓦斯质量守恒定律，联立式(1)—式(3)可求解得到煤层瓦斯中甲烷、氮气、二氧化碳各组分气体分压力，按式(4)计算各瓦斯组分分压力在煤层瓦斯总压中占比，表4为煤层瓦斯各气体组分分压计算结果。按式(5)计算煤层孔隙-裂隙内游离瓦斯气体量。

kxpz(z=1,2,3)

(1)

(2)

(3)

(4)

(5)

式中：Qz为单位体积煤瓦斯组分气体z的游离状态与吸附状态含量，m3/m3；az、bz为煤层瓦斯组分z的朗格缪尔等温吸附常数，m3/t，MPa-1；pz为煤层瓦斯组分z的分压力，MPa；ρm为煤的密度，t/m3；W为煤层水分，%；Ad为煤层灰分，%；kx为煤体的孔隙率，%；nz、ng为单位体积煤瓦斯组分z、瓦斯组分g的体积比值；Qg为单位体积煤瓦斯组分气体g的游离状态与吸附状态含量，m3/m3；Q为单位体积煤游离状态与吸附状态瓦斯总含量，m3/m3；mz为瓦斯组分z分压力在煤层瓦斯总压中占比值；p为煤层总瓦斯压力，MPa。

表3 煤层瓦斯各组分分压计算结果Table 3 Calculation results of gas component pressure in coal seam

为了验证煤层瓦斯赋存模型的可靠性，需要对比煤层瓦斯压力直接测试值与煤层多组分瓦斯气体赋存模型计算结果之间的相对误差。为了测得准确的22051工作面煤层瓦斯总压，采用穿层钻孔方式在22051工作面回采区域布置4个瓦斯压力测试孔测试煤层瓦斯总压，煤层瓦斯压力直接测试值与煤层瓦斯总压模型计算结果，相对误差小于5%，验证了该瓦斯含量-压力赋存模型可用于描述22051工作面煤层瓦斯赋存情况。

2.1.2 动态推进过程煤层多组分瓦斯渗流模型

煤层中瓦斯在大裂隙中进行非线性渗流，在小裂隙中进行线性渗流，22051回采工作面煤层裂隙不发育，且未进行煤层增透等措施，假定煤层内部瓦斯流动状态近似为线性层流，符合达西定律[7]，由于该煤层瓦斯由甲烷、氮气、二氧化碳3种气体构成，其中甲烷为主要成分，煤层瓦斯密度按式(6)计算，煤层瓦斯黏度按式(7)计算，利用式(8)计算煤层瓦斯涌出速率。

(6)

1.1×10-8(e1.439ρrm-e-1.111(ρcm)1.858)(Mm)0.5×

(7)

(8)

为了模拟计算在一定推进速度条件下的工作面前方煤层瓦斯渗流过程，需要引入移动坐标系，以回采工作面下端头为坐标原点建立笛卡尔坐标系，该坐标系以工作面推进速度随工作面同步移动(图3)，建立动态回采条件下工作面煤壁前方煤层多组分瓦斯渗流控制方程，见式(9)。

图3 工作面引入移动坐标系模拟工作面开采推进过程Fig.3 Simulation on the mining process of coal face by introducing moving coordinate system

(9)

式中：Vx为工作面推进速度，m/d；x为工作面推进方向上煤层位置坐标值，m；y为工作面倾向方向上煤层位置坐标值，m。

2.2 采动影响下工作面前方煤层渗透率分布模型

采动作用下工作面前方煤体渗透率会发生显著变化[8]，工作面煤壁前方会出现应力三区，卸压区、弹性应力区、原始应力区，应力分布变化必然导致煤体渗透率发生变化[9-12]，文献[13]研究表明工作面煤壁前方卸压区为渗流开放区，弹性应力区为渗流屏蔽区和渗流恢复区，原始应力区对应于原岩渗流区，工作面煤壁前方煤体卸压区渗透率一般增加2个数量级以上，渗流屏蔽区一般为5～8 m，渗流屏蔽区和渗流过渡区一般为5～7 m，如图4所示，忽略采动影响造成的工作面煤壁前方煤层渗透率变化会对计算结果造成显著的误差。

图4 回采工作面煤壁前方煤层渗透率变化模型[12]Fig.4 Permeability change model of coal seam in front of coal wall of coal face[12]

利用工作面煤壁前方煤层渗透率变化模型[12]对22051工作面煤壁前方煤层渗透率进行设置，假定22051工作面煤壁前方煤层渗流开放区范围为0～8 m，渗流屏蔽区和渗流过渡区范围为8～15 m，假定工作面煤壁前方煤层渗流开放区内煤层渗透率随距工作面距离变化满足线性函数，假定渗流屏蔽区和渗流过渡区内煤层渗透率随距工作面距离变化满足抛物线函数，工作面煤壁前方煤层渗透率分布按式(10)计算，将22051工作面原始煤层渗透率代入式(10)可计算得到工作面煤壁前方渗透率分布曲线，如图5所示。

图5 22051回采工作面煤壁前方煤层渗透率分布曲线Fig.5 Permeability distribution curve of coal seam in front of coal wall in No.22051 coal face

(10)

式中：k0为回采工作面煤层原始渗透率，m2；K(x)为工作面推进方向上煤层位置坐标x处煤层渗透率，m2。

3 工作面煤壁多组分瓦斯涌出数值计算模型解算方法

3.1 计算区域网格划分及边界条件设置

图6为回采煤层瓦斯渗流涌出计算模型边界条件设置和计算区域网格划分情况，工作面附近煤层瓦斯压力梯度大，对靠近工作面煤壁区域网格进行加密。空气流过回采工作面期间存在通风阻力损失，进风巷、回采工作面、回风巷巷道空间内存在一定气压差，但巷道内绝对气压值远大于通风阻力值，且煤层倾角较小，故进风巷煤壁、回风巷煤壁、回采工作面煤壁处气压均取工作面中间位置绝对气压值(距进风巷97.5 m处工作面巷道)，表4为回采工作面煤层瓦斯渗流计算区域边界条件设置具体情况。

图6 工作面煤壁前方煤层瓦斯渗流计算区域网格划分Fig.6 Mesh division of gas seepage calculation area in front of coal wall of coal face

表4 工作面前方煤层瓦斯渗流计算区域边界条件设置Table 4 Setting of boundary conditions for calculating gas seepage in coal seam in front of coal wall of coal face

3.2 输运控制方程离散化处理及解算程序

根据计算区域网格划分形式，采用基于有限差分法的C-N离散格式对煤层多组分瓦斯渗流偏微分方程进行离散化处理，该离散格式能够保证计算过程中时刻满足质量守恒定律，解算结果具有物理真实性，回采工作面前方煤层瓦斯渗流计算区域内任一节点瓦斯压力离散模型如图7所示，对每个节点进行编号。在节点瓦斯压力离散模型的基础上，构建每个节点构建其瓦斯压力方程，见式(11)，每个节点的瓦斯压力方程中均包含相邻节点瓦斯压力，最终形成以整个计算区域范围内所有瓦斯节点为未知量的非线性方程组，通过求解该方程组计算得到整个计算区域瓦斯压力场，由于方程组系数矩阵与常数矩阵均包含方程组所求解未知量，采用迭代法求解方程组，计算误差达到一定精度后完成计算。

图7 工作面煤壁前方煤层瓦斯渗流计算区域瓦斯压力节点离散模型Fig.7 Discrete model of gas pressure node in area of gas seepage calculation in coal seam in front of coal wall of coal face

(ai,j+bi,j)Pi,j+ai-1,jPi-1,j+ai+1,jPi+1,j+ai,j-1Pi,j-1+ai,j+1Pi,j+1=0

(11)

(12)

(13)

(14)

(15)

(16)

(17)

式中：Pi,j、Pi-1,j、Pi+1,j、Pi,j-1、Pi,j+1分别为邻近层计算区域内节点(i，j)、节点(i-1，j)、节点(i+1，j)、节点(i，j-1)、节点(i，j+1)位置处瓦斯压力平方，待求解量，MPa2；ai+1,j、ai-1,j、ai,j-1、ai,j+1、ai,j、bi,j均为过程中间值，具体计算公式分别为式(12)—式(17)；Ki,j、Ki-1,j、Ki+1,j、Ki,j-1、Ki,j+1分别为邻近层计算区域内节点(i，j)、节点(i-1，j)、节点(i+1，j)、节点(i，j-1)、节点(i，j+1)位置处煤层渗透率，m2；Δyi,j为节点(i，j)处控制体的宽度，m；Δxi,j为节点(i，j)的控制体的长度，m；yi,j为节点(i，j)的y方向坐标值，m；yi,j+1为节点(i，j+1)的y方向坐标值，m；yi,j-1为节点(i，j-1)的y方向坐标值，m；xi,j为节点(i，j)的x方向坐标值，m；xi+1,j为节点(i+1，j)的x方向坐标值，m；xi-1,j为节点(i-1，j)的x方向坐标值，m。

根据上述离散方程以及求解条件，利用vb语言编制了解算程序，程序结构流程如图8所示。

图8 回采煤层瓦斯渗流解算程序流程Fig.8 Flow chart of calculation program of gas seepage

4 工作面煤壁瓦斯涌出量现场测试

为了对比分析回采工作面煤壁多组分瓦斯气体涌出数值计算模型的现场应用误差分析，需要对回采工作面煤壁瓦斯涌出量进行实测。

4.1 工作面瓦斯涌出量分源测试方法

采用回采工作面瓦斯涌出分源测试方法能够实测工作面巷道内甲烷涌出量与采空区瓦斯涌出量。生产班期间回采工作面瓦斯涌出源主要包括采空区、工作面煤壁、回采落煤，采用源测试方法无法区分计算工作面煤壁甲烷涌出量与回采落煤瓦斯涌出量，检修班期间回采工作面瓦斯涌出源为采空区和工作面煤壁，采用分源测试方法能够实测得到工作面煤壁瓦斯涌出量，因此选择在检修班期间实测回采工作面煤壁瓦斯涌出量。生产班期间割煤产生新暴露煤壁，生产班期间煤壁瓦斯涌出量明显大于检修班期间涌出量，选择夜班生产班结束之后、早班检修班刚开始时刻(具体测试时间为8:00—10:00)进行测试，测试结果介于检修班期间与生产班期间，可近似为工作面煤壁平均瓦斯涌出量。

当前养殖形势和以往有了极大的变化，作为养殖场要顺应时代的变化，通过自身管理水平的提高，做好精细化管理，做好生物安全防控体系，减少鸡群发病并严格按照药物使用规范的要求合理用药，做到健康养殖，这样才能有更加长远的发展。

4.1.1 工作面瓦斯涌出量现场测试

根据分段测定法原理，在对工作面瓦斯分布情况进行测定时，可采取如下5个步骤：

1)将综采工作面沿倾斜方向划分为8个区段，共9个巷道测试断面，在每个测试断面内，沿综采工作面推进方向依次设置6个测点，测点1位置距煤壁最近，测点5位置距采空区最近，测点6位于采空区浅部，由工作面深入采空区0.2～0.4 m，如图9所示。

图9 22051工作面瓦斯风量测点布置Fig.9 Layout of measurement points of gas and air volume in No.22051 coal face

2)对每个巷道测试断面的尺寸进行测定，计算巷道测试断面面积。

3)对每个巷道测试断面内1—6测点均进行瓦斯浓度测定，对每个巷道测试断面内1—5测点平均风速测定。

4)巷道测试断面内测点6处瓦斯浓度近似为该测试断面处采空区漏风流瓦斯浓度。

5)根据瓦斯平衡方程与风量平衡方程，计算每个区段的采空区漏风量、采空区瓦斯涌出量、煤壁瓦斯涌出量。

4.1.2 工作面瓦斯涌出量分源测算方法数学模型

(18)

(19)

(20)

(21)

(22)

由质量守恒定律可得，综采工作面各区段势必遵循瓦斯平衡方程、风量平衡方程，结合式(18)—式(22)得到式(23)，式(23)是以各区段采空区漏风量、采空区瓦斯涌出量、煤壁瓦斯涌出量3个物理量为未知求解量的相互独立方程，将实测数据代入式(23)可求解得到各区段采空区漏风量、采空区瓦斯涌出量、煤壁瓦斯涌出量。

(23)

4.2 工作面煤壁瓦斯涌出量现场实测结果

根据现场实际条件选择分别于2018年3月2日、4月7日、5月11日、6月22日开展工作面瓦斯涌出量现场测试工作，并计算工作面煤壁瓦斯涌出量，同时实测了构建工作面煤壁瓦斯涌出简化模型所需数据，测试数据见表5。

表5 22051工作面煤壁瓦斯涌出量实测结果Table 5 Measured results of gas emission from coal wall of No. 22051 coal face

5 工作面煤壁多组分瓦斯涌出数值计算模型现场应用

为了对比分析回采工作面煤壁多组分瓦斯涌出数值计算模型的现场应用误差，在对回采工作面煤壁瓦斯涌出量进行实测的基础上，采用工作面煤壁多组分瓦斯涌出数值计算模型模拟计算工作面煤壁瓦斯涌出量，对比模拟结果与实测结果，以确定工作面煤壁多组分瓦斯涌出数值计算模型的现场实用性。

5.1 工作面煤壁多组分瓦斯涌出数值计算模型模拟计算

为了设置数值计算模型的边界条件，需测试工作面煤壁瓦斯涌出量现场实测期间工作面中间位置绝对气压，选取在2018年3月2日、4月7日、5月11日、6月22日4个时刻实测了工作面中间位置绝对气压，分别为101.83、99.91、98.29、97.65 kPa。

回采工作面每日的实际推进速度存在差异，工作面煤壁瓦斯涌出取决于煤壁前方煤层瓦斯赋存状况，煤壁前方煤层瓦斯赋存状况与近期一段时间内工作面回采进尺密切相关。第t天工作面瓦斯涌出情况取决于第t天工作面前方煤壁瓦斯分布，第t天工作面前方煤壁瓦斯分布主要取决于第t天工作面推进度、第t-1天工作面前方煤壁瓦斯分布(作为第t天工作面前方煤壁瓦斯分布的初始条件)，而第t-1天工作面前方煤壁瓦斯分布取决于第t-1天工作面推进度、第t-2天工作面前方煤壁瓦斯分布，依次类推，第t天工作面煤壁瓦斯涌出情况与第t天之前一段时期内的工作面推进度密切相关，因此采用一段时期内工作面推进度的平均值作为工作面日推进度进行模拟计算。为了对比分析工作面煤壁瓦斯涌出数值模拟计算结果与现场实测结果，需要以工作面实际日推进度数据为基础，以某一时间周期计算工作面日平均推进速度，见式(24)，第t天工作面平均推进速度为最近一个计算周期内的工作面总推进度除以总时间，平均推进速度计算周期分别取10、15、20、25、30、35、40 d，分别对不同计算周期条件下2018年3月2日、4月7日、5月11日、6月22日这4 d的平均推进速度进行计算，计算结果见表6。

表6 不同计算周期条件下22051工作面日推进速度计算结果Table 6 Calculation result of daily mining speed of No.22051 coal face under different calculation periods

(24)

式中：T为日平均值计算周期，d，取值为30 d；Vt为第t天工作面平均推进速度，m/d；dh为第h天工作面实际推进度，m。

将表1、表7、表8中回采工作面瓦斯基础参数、工作面气压实测值、工作面推进速度实测值代入工作面煤壁多组分瓦斯气体涌出数值计算模型进行数值模拟。图10为工作面煤壁前方煤层瓦斯压力分布数值模拟结果，工作面煤壁瓦斯压力明显大于进回风巷煤壁瓦斯压力。不同工作面参数条件下煤壁瓦斯涌出量数值模拟结果汇总见表7。

图10 工作面煤壁前方煤层瓦斯压力分布模拟结果Fig.10 Simulation results of gas pressure distribution in coal seam in front of coal wall

表7 不同工作面日推进速度计算周期条件下工作面煤壁瓦斯涌出量模拟结果Table 7 Simulated results of gas emission from coal wall of coal face under different calculation periods of daily mining speed of coal face

5.2 现场实测结果与计算模型模拟结果对比分析

将表6中不同回采时刻工作面煤壁瓦斯涌出量实测结果与表10中不同工作面日推进速度计算周期条件下工作面煤壁瓦斯涌出量模拟结果进行对比，工作面日推进速度计算周期为25～35 d时，4次现场测试结果与数值模拟结果之间相对误差均在10%以内(图11)，研究表明22051工作面在当前开采条件和煤层瓦斯赋存渗流条件下，工作面煤壁瓦斯涌出量数值计算模型中的工作面平均日推进速度的最佳计算周期为25～35 d，工作面回采进尺状况对工作面煤壁前方煤层瓦斯渗流运移具有重要影响作用，工作面煤壁瓦斯涌出取决于为25～35 d工作面的回采进尺状况，同时也说明工作面煤壁瓦斯涌出存在时间累加与延迟效应，如果一段时间内工作面回采进尺很小，突然增大推进速度，工作面煤壁瓦斯涌出量也不会突然间增加很大。

图11 不同计算周期条件下工作面煤壁瓦斯涌出量实测与模拟结果相对误差Fig.11 Relative errors between measured and simulated results of gas emission from coal wall of coal face under different calculation periods

6 结 论

1)以22051综采工作面为研究对象，考虑了工作面推进速度和采动影响对工作面煤壁瓦斯涌出的影响，通过引入移动坐标系和采动影响下工作面前方煤层渗透率分布模型，构建了工作面动态推采期间工作面煤壁多组分瓦斯气体涌出数值计算模型，采用有限差分法编制解算程序，可用于模拟计算回采工作面煤壁瓦斯涌出量。

2)以瓦斯平衡方程、风量平衡方程为基础构建工作面瓦斯涌出量分源计算数学模型，实测了工作面煤壁瓦斯涌出量。

3)对比分析了22051工作面煤壁瓦斯涌出量数值模拟结果与实测结果，研究表明，工作面煤壁瓦斯涌出存在时间累加与延迟效应，该工作面煤壁瓦斯涌出数值计算模型中工作面平均日推进速度参数值应该由工作面25～35 d的实际推进度计算得到，该工作面煤壁瓦斯涌出数值计算模型中的工作面平均日推进速度参数取值的最佳计算周期为25～35 d，工作面回采进尺状况对工作面煤壁前方煤层瓦斯渗流运移具有重要影响作用，该工作面在当前开采条件和煤层瓦斯赋存渗流条件下，工作面煤壁瓦斯涌出量取决于25～35 d的工作面回采进尺状况。