基于Simtec的近球体瓦斯爆炸数值模拟

苟宝洋，吴 兵，马一飞，张 洋，苏敬亮，贾泽鹏

(中国矿业大学(北京) 应急管理与安全工程学院，北京 100083)

瓦斯爆炸是煤矿发生的主要灾害之一，通常发生于井下受限空间范围内，因此研究受限空间的瓦斯爆炸特性对矿井安全生产具有重要意义。由于爆炸实验具有危险性和重复性差的特点，而数值模拟可弥补上述问题，因此数值模拟目前已经成为研究瓦斯爆炸的重要手段之一。在受限空间的爆炸数值模拟方面，国内许多学者已经做了大量的研究：祝钊[1]通过数值模拟研究了管道内瓦斯爆炸参数的变化规律，并研究了不同因素对爆炸传播的影响；张莉聪等[2]采用19步瓦斯爆炸化学反应机理模拟研究了障碍物形状对瓦斯爆炸传播过程的影响；柳伟等[3]利用LS-DYNA软件研究了5 m长管道内爆炸发生过程流固耦合作用对冲击波的影响；罗振敏等[4]利用FLACS软件对近球形密闭反应容器内的瓦斯爆炸过程进行了模拟，表明采用辐射热换模型模拟的结果与实验结果更为接近；尹彬等[5]利用CHEMKIN软件对采空区煤自燃与瓦斯爆炸的耦合关系进行了模拟研究；王志荣[6]、解北京[7]等利用 Fluent软件分别模拟了连通容器内、分叉管道中的瓦斯爆炸火焰传播规律。然而，目前国内外学者在受限空间内瓦斯爆炸数值模拟研究中很少提及壁面辐射系数对模拟结果的影响。

笔者拟利用 Simtec 软件对近球体密闭容器内甲烷—空气爆炸过程进行数值模拟研究，分析壁面辐射系数对瓦斯爆炸数值模拟精度的影响，以期揭示基于不同甲烷浓度下爆炸过程中的压力、温度等参数变化规律，为受限空间瓦斯爆炸的研究提供理论指导。

1 Simtec软件平台与数值模拟方法

Simtec软件在燃烧和爆炸模拟方面已得到了很好的验证，目前国内有中国科学技术大学、中南大学的学者已在使用该软件[8-9]。Simtec软件采用大涡模型对瓦斯爆炸进行模拟，大涡模型的基本思想是利用滤波函数对动量方程进行过滤，得到不同尺度的涡流加以分类计算[10-13]。其过程遵循如下数学模型：

三大守恒控制方程通用形式如下：

(1)

式中:φ为通用变量，代表u，v，w，T等变量；Γφ为广义扩散系数；Sφ为广义源项。

1.1 大涡模型

用滤波函数对连续不可压缩流的N-S方程进行滤波，得到控制方程如下：

(2)

(3)

1.2 燃烧模型

对于甲烷燃烧爆炸化学反应，其化学反应速率为：

(4)

式中：A为指前因子；T为热学力温度；m为由实验决定的常数；Ea为活化能；R为摩尔气体常数；C、D分别为反应物甲烷、氧气浓度；a、b分别为反应物的反应级数。

反应动力学参数由爆炸实验经验值得到[14-16]。依次给定化学反应动力学参数A、m、Ea、a、b，可得出分步化学反应速率。

2 模拟方案与模型设置

2.1 模拟方案

西安科技大学罗振敏教授在爆炸实验中采用了XKWB-1型20 L近球体密闭气体爆炸特性测试装置，反应容器内部最大直径30 cm，高34 cm，壁面材料为不锈钢。

爆炸前密闭容器内初始温度为14.6～21.0 ℃，初始压力为常压。实验中压力传感器所在位置坐标为(0.295,0.245,0.150)，最终得到了不同甲烷体积分数(6%～14%)下最大爆炸的超压[4]。

笔者在罗振敏爆炸实验的基础上，利用Simtec软件在壁面辐射系数分别为0.9、0.1条件下，模拟了不同甲烷体积分数(6%～14%)下的爆炸过程。模拟过程中的所有参数设定与实验保持一致，后处理时将测点(0.295,0.245,0.150)所在位置的数据提取出来进行分析。

2.2 模型设置

环境温度为20 ℃，环境压力为1.013×105Pa，求解器选择可压流求解器。湍流模型选用Smagorinsky 大涡模型，Smagorinsky 常数为0.1；采用非等温壁面函数处理壁面传热；燃烧模型为修改的涡耗散模型，辐射模型选择Modak；计算时间设置为200 ms。

密闭近球体容器物理模型如图1所示。

图1 密闭近球体容器数值模拟几何形状

如图1所示，点火源位于容器中心位置，点火方式为电火花点火，点火能量为1 J。利用结构化单元划分网格，x、y、z3个方向上网格数目均为100个，网格总数为1 000 000个，x、z方向单位网格长度为 3.0 mm，y方向单位网格长度为3.4 mm。对于壁面网格分为2层，每一层网格厚度为0.005 mm。在流体初始条件文件中设置近球体内部的甲烷质量分数。对模型做出如下假设：

1)除中心电火花外，容器内没有其他热源。

2)甲烷的初始浓度、温度和压力都均匀分布。

3)密闭空间内气体满足真实气体状态方程。

采用三维瞬态数值方法计算气体爆炸过程，在整个计算区域空间上采用有限体积法离散微分方程，时间上采用显式方法，对流项采用二阶差分格式ALBADA方法离散，扩散项采用中心差分格式离散。

3 模拟结果与分析

3.1 不同壁面辐射系数下的最大爆炸超压与实验结果对比

辐射是影响密闭空间瓦斯爆炸的重要因素，材料表面辐射系数与材料的种类、壁面温度，以及材料表面的粗糙度有关。辐射系数为0.9、0.1时模拟得到不同甲烷浓度(甲烷体积分数)下的爆炸超压并与实验数据对比，结果如图2所示。在甲烷体积分数为6%～12%、辐射系数为0.1时最大相对误差为6.09%，最小相对误差为 0.85%。

图2 壁面辐射系数分别为 0.9、0.1 时模拟爆炸超压与实验结果对比曲线

由图2可知，当壁面辐射系数为0.1时，模拟计算得到的结果与实验数据的偏差最小，在甲烷体积分数为6%～12%时模拟结果与实验数据有较好的吻合关系。表明壁面辐射系数对密闭空间瓦斯爆炸数值模拟精度有重要影响。

对于密闭空间瓦斯爆炸，甲烷浓度与最大爆炸超压呈二次函数关系[17-19]，因此对数值模拟数据进行多项式拟合，函数关系式如下：

y=-0.007 15x2+0.166 07x-0.218 74

对拟合结果进行优度检验，相关系数R2为 0.965 897。拟合二次函数曲线如图3所示。

图3 壁面辐射系数为 0.1时模拟数据多项式拟合曲线

由图3可知，瓦斯爆炸超压最大值在二次函数最高点(11.581 0，0.745 1)处取得。

3.2 壁面辐射系数为0.1时不同甲烷浓度下瓦斯爆炸压力、温度变化规律

壁面辐射系数为0.1条件时，模拟得到了甲烷体积分数为6%～14%下爆炸超压随时间变化曲线如图4所示，温度随时间变化曲线如图5所示。

图4 甲烷体积分数为6%～14%时瓦斯爆炸超压随时间变化曲线

图5 甲烷体积分数为6%～14%时温度随时间变化曲线

由图4、图5可知，最大爆炸超压及温度随时间变化均表现为先迅速上升，随后爆炸超压缓慢下降，温度呈脉冲状震荡缓慢下降趋势。最大爆炸超压在40～75 ms内达到峰值，最大爆炸超压为0.75 MPa，最高温度为2 370 ℃。

分析表明：在 0～40 ms时，爆炸反应处于初始阶段，电火花点燃甲烷气体后，随着化学反应的进行，密闭空间内温度迅速升高，活化能变小，化学反应速率迅速增大，导致爆炸压力迅速升高；在40～75 ms时，化学反应速率达到最大，化学反应充分，爆炸最为剧烈，各种甲烷浓度下甲烷气体均达到了最大爆炸超压；在75 ms以后，化学反应则基本结束，导致密闭空间内爆炸压力缓慢下降。受密闭空间壁面辐射、对流等因素影响，温度表现为震荡缓慢下降[15]。

甲烷体积分数为10%时，不同时刻的爆炸压力、温度云图如图6～7所示。

图6 甲烷体积分数为10%时不同时刻瓦斯爆炸压力云图

图7 甲烷体积分数为10%时不同时刻爆炸温度云图

3.3 壁面辐射系数为0.1时不同甲烷浓度下最大压力上升速率、爆炸强度指数变化规律

在Origin中将图4中不同甲烷浓度下爆炸超压曲线对时间进行求导，得到甲烷体积分数6%～14%条件下的dp/dt，求解每一浓度下dp/dt的最大值，即可得到甲烷体积分数为6%～14%条件下的最大压力上升速率(dp/dt)max。

最大压力上升速率随甲烷体积分数的变化曲线如图8 所示。

图8 最大压力上升速率随甲烷体积分数变化曲线

由图8可以看出，当甲烷体积分数为10%时，(dp/dt)max达到最大值32.36 MPa/s。在峰值左侧，(dp/dt)max逐渐上升；在峰值右侧，当甲烷体积分数为12%时出现台阶式上升，随后逐渐下降。

结合图4可知，最大爆炸超压pmax与最大压力上升速率(dp/dt)max是在峰值左侧正相关的，都是呈现逐渐增高的趋势；在峰值右侧是负相关的，都是总体呈现逐渐下降的趋势[20]。

由于最大压力上升速率和甲烷体积分数呈二次函数关系[17]，对图8进行多项式拟合，结果见图9。

图9 最大压力上升速率随甲烷体积分数变化拟合曲线

函数关系式如下：

y=-0.446 9x2+9.625 8x-20.373 3

对图9中拟合结果进行优度检验，相关系数R2为0.914 96，表明曲线拟合较好。

除了温度、压力,以及最大压力上升速率等参数外，爆炸强度指数也是描述密闭空间爆炸特性的一个重要参数。爆炸强度指数定义为最大压力上升速率和装置体积三次方根的乘积，其计算公式如下：

(5)

式中：KG为爆炸强度指数，MPa·m/s；V为密闭容器体积，m3。

甲烷体积分数为6%～14%时，最大压力上升速率及爆炸强度指数计算结果如表1所示。

表1 甲烷体积分数为6%～14%时最大压力上升速率及爆炸强度指数计算结果

由表1可知：最大压力上升速率及爆炸强度指数随着甲烷体积分数增大先增加，随后逐渐降低。在甲烷体积分数为10%时，爆炸强度指数达到最大值87.84 MPa·m/s。结合图8可知，甲烷体积分数为10%是最佳甲烷爆炸浓度。

4 结论

1)利用Simtec软件模拟壁面辐射系数分别为0.9、0.1时，密闭近球体容器中基于不同甲烷体积分数(6%～14%)的预混气体爆炸过程。在辐射系数为0.1时，模拟计算得到的结果与实验数据的偏差最小；在甲烷体积分数为6%～12%时模拟结果与实验结果较吻合。表明壁面辐射系数对密闭空间瓦斯爆炸数值模拟精度有重要影响。

2)甲烷爆炸压力、温度随时间变化均表现为先迅速上升，在40～75 ms时压力与温度均达到最大值。最大爆炸超压为0.75 MPa，最高温度为2 370 ℃；在 75 ms 以后，化学反应基本结束，在密闭空间内爆炸压力缓慢下降，受密闭空间壁面辐射、对流等因素影响，温度表现为震荡缓慢下降。

3)最大爆炸超压与最大压力上升速率均随着甲烷体积分数增加先增大到峰值后逐渐减小，最大爆炸超压与最大压力上升速率二者在峰值左侧成正相关，都是呈现逐渐增加的趋势；在峰值右侧成负相关，都是总体呈现逐渐下降的趋势。最大压力上升速率和甲烷体积分数的二次函数关系式为：y=-0.446 9x2+9.625 8x-20.373 3；在甲烷体积分数为10%时，爆炸强度指数达到最大值87.84 MPa·m/s，此时的甲烷体积分数是最佳甲烷爆炸浓度。