1. 引言
随着我国煤矿事故频频发生，煤矿安全问题日益引起人们的关注。瓦斯爆炸事故是煤矿事故中发生频率最高，破坏性最严重的事故。不仅仅造成爆炸地点巷道毁坏，其爆炸威力还可能影响到临近巷道、采区，甚至整个地下巷道系统，巷道修复耗资巨大。
近年来众学者关于受限空间内瓦斯爆炸进行过大量研究，但大多都仅限于单爆源瓦斯爆炸的研究，包括但不限于对瓦斯爆炸的影响因素、传播规律、伤害机理、破坏特性等一系列瓦斯爆炸的物理化学性质进行了研究，例如贾泉升、叶青等 [1] [2] 从各种影响因素出发对矿井瓦斯爆炸进行了研究，发现各种影响因素对于瓦斯爆炸的火焰传播速度和压力均有不同程度的影响并从压力波及湍流等层面解释了以上各因素对瓦斯爆炸产生影响的原因；王海燕 [3] 等基于弱冲击波理论和瓦斯爆炸能量释放，推导出冲击波超压与起爆距离存在非线性反比关系；曲志明、周心权 [4] 等建立了瓦斯爆炸的物理模型和冲击波关系式；徐景德、杜扬、司荣军 [5] [6] [7] 等采用试验方法对巷道中瓦斯爆炸冲击波的传播规律进行了研究；景国勋 [8] [9] 等针对煤尘爆炸传播特性做了实验研究；杨书召、曲志明 [10] [11] 等针对矿井瓦斯爆炸冲击气流的伤害、衰减和破坏做了研究。林柏泉、叶青 [12] 等研究了巷道分岔不同角度对瓦斯爆炸冲击波的影响。然而却鲜有学者对于多爆源的瓦斯爆炸进行过研究，爆炸的介质可区分为固液气三种，目前关于多爆源的研究大多集中于地下爆炸和水下爆炸两种情况，例如张世豪 [13] 模拟了混凝土中多点同步爆炸，分析了多点爆炸的能量聚集现象；发现在炸点附近，多点同步爆炸冲击波峰值压力小于单点冲击波峰值压力，随着传播距离的增加，多点爆炸峰值压强将大于单点峰值压强；孟闻远 [14] 对冰面两点同步爆炸进行研究，发现炸点间距不同，冰面破坏形状不同；胡宏伟 [15] 发现在总药量一定的状况下，炸点数量越多，冲击波和能量的增长幅度越大；韩克华 [16] 对多爆点同步起爆进行研究，发现相对于三点、六点、八点同步起爆，四点同步起爆具有更高的超压，然而在煤矿生产中，如果巷道内多个瓦斯泄露或积聚点在同一或相邻的时间节点被点燃，就将会发生多爆源瓦斯爆炸的情况，其冲击波对煤矿巷道结构和设备的破坏严重影响了煤矿生产和安全，因此，对于巷道内多爆源瓦斯爆炸冲击作用特性的研究具有重要意义。
本研究采用ANSYS/LS-DYNA软件对双爆源瓦斯爆炸条件下的巷道展开模拟研究，对巷道内双爆源瓦斯爆炸冲击作用特性进行了深入探究。通过对模拟结果的分析和与巷道内单爆源瓦斯爆炸冲击作用的对比，以期补足对双爆源条件下的瓦斯爆炸冲击作用特性的研究，并为煤矿安全生产提供有益的参考和指导。
2. 巷道瓦斯爆炸数值模型的建立
2.1. 物理模型的建立与网格划分
本文巷道模型断面以重庆分院大型试验主巷道 [17] 为蓝本等比例构建，为方便计算，根据几何模型的轴对称性，建立10 m长的1/2的两端开口巷道模型，本文中所用材料的模型和参数均采用统一单位制(kg/m/s)，计算模型的长宽高分别为L = 10 m、W = 4 m、H = 8 m，巷道高度为h + r = 2 m，宽度为2 w = 2 m，顶拱半径为r = 1 m，假设巷道内均匀地充满浓度为9.5%的瓦斯预混气体，采用TNT当量法 [18] 来对巷道内瓦斯气体爆炸进行等效代替，通过爆热能、冲量等效关系可以将瓦斯气体当量等效转换为 TNT 当量，经计算10 m长的巷道充满9.5%的瓦斯气体爆炸可等效替代为约30 kg的TNT炸药爆炸，以巷道模型入口底部中点为原点构建空间直角坐标系，假设起爆时间为0 s，巷道内瓦斯气体的两着火点为(0, 4, 3)和(0, 4, 7)，两着火点的点火能量一致，故在这两点填充当量一致的炸药(填充半径一致)，建立起的巷道有限元模型如图1．由于建立的围岩模型是不规则的，难以用映射法来划分网格，所以本文选择扫掠法进行网格划分，网格划分的尺寸设置为5 cm，共划分网格数约133万个，划分网格后的模型如图2，采用Lagrange单元建立围岩模块以实现变形与损伤的可视化，采用Euler单元构建空气与TNT模块以更好地模拟爆炸冲击波的传播。为实现爆炸冲击波对巷道壁面的耦合破坏作用，采用ALE方法建立模型，使固体网格与流体网格重合，并通过关键字*CONSTRAINED_LAGRANGE_IN_SOLID对流固耦合进行控制。
2.2. 边界条件和基本假设
瓦斯爆炸是一种反应极其迅速的化学反应，爆炸过程中产生大量的中间和瞬时产物，在本文中忽略反应的中间过程。为使模拟结果符合巷道真实结果，需在模型四周合理设置边界条件，为降低边界反射效应，在模型除对称面外均设置无反射边界。为防止模型在重力效应下在y轴方向出现不合理位移变形，在模型底部设置全约束刚性边界，在Y-Z面设置对称边界。为了进行合理有效的数值计算，对模型做出以下基本假设：
1) 巷道瓦斯段的瓦斯混合气体处于常温常压状态且混合均匀(9.5% CH4)；
2) 巷道内只存在双爆源瓦斯爆炸两个热源；
3) 瓦斯气体爆炸过程为单步不可逆反应，忽略反应过程中的中间产物；
4) 巷道围岩为均匀连续的；
5) 壁面光滑绝热，不考虑壁热效应。
2.3. 材料参数的选择
2.3.1. 围岩材料参数
围岩采用JOHNSON_HOLMQUIST_CONCRETE (HJC)材料模型描述，该材料模型是在Ottosen模型基础上加入了应变率等参数建立的，能够反映混凝土、岩石等脆性材料在大变形、高应变率和高围压状态下的损伤失效变化。
2.3.2. 空气介质参数
空气介质模型材料选用LS－DYNA中的*MAT_NULL，状态方程选用*EOS_LINER_POLYNOMIAL。
2.3.3. 炸药材料参数
本文采用TNT当量法，将瓦斯爆炸产生的能量转化为一定当量的TNT，炸药采用LS－DYNA中的*INITIAL_VOLUME_FRACTION_GEOMETRY关键字和*ALE_MULTI－MATERIAL_GROUP关键字，以球形填充在空气网格中，每个炸药半径为0.13 m。其本构模型及状态方程分别选用LS－DYNA材料中的*MAT_HIGH_EXPLOSIVE_BURN和*EOS_JWL。
2.4. 观测点的布置
为观测巷道内的压力变化情况，在炸药所处的轴线方向布置1-8号测点，测点的位置为左侧爆源位置向两端每隔0.5 m布置一个，每侧4个共8个，为了观测壁面的动力响应规律，在巷道两爆源间的对称切面上的拱顶、帮部、隅角、底板分别设置1个测点，共4个。具体测点分布位置如图3。
3. 数值模拟结果分析
3.1. 双爆源瓦斯爆炸冲击波超压分析
由前人受限空间单爆源爆炸的研究 [19] 可知，巷道单爆源瓦斯爆炸发生后，冲击波先是以球状的波阵面从爆心向四周扩散开来，在巷道模型的约束下冲击波首先到达巷道拱顶并开始产生反射波，随后再传播至巷道底板和帮部，由于隅角离爆心的距离最远，所以冲击波最后到达隅角并传播至巷道两端。由于反射波的形成和传播，巷道内冲击波的波阵面不再是球形，这种影响随着传播距离的增加更加明显，在爆炸冲击波与巷道围岩的相互作用过程中，由于巷道构造的不规则性，巷道内由围岩产生的反射波的传播方向各不相同，最终反射波与爆炸产生的冲击波汇合，最后以不规则弧面波的形式随着原来的方向沿巷道的纵向传播并逐渐耗散。巷道内双爆源瓦斯爆炸的冲击波超压传播云图如图4所示。从图4可以看出巷道双爆源瓦斯爆炸在冲击波碰撞之前互不干扰独立传播与上述单爆源冲击波的传播过程一致(见图4(a)~(c))，但随着时间的推移，两爆源激发的冲击波将于巷道爆源的对称面相遇，这使得后续双爆源冲击波流场与前述单爆源冲击波流场出现较大差异(见图4(d)~(f))，两爆源所激发的冲击波于对称面相遇发生碰撞进而产生合成波，合成波的传播方向和速度随时间不停改变，经巷道壁面反射作用追赶上最初爆炸产生的冲击波进而发生叠加合成，并且冲击波的碰撞在对称面表现出明显的尖端效应，呈锥形向下扩展，最后由于传播距离的增加，冲击波产生的能量逐渐耗散衰减。
基于流固耦合数值模型，得到双爆源瓦斯爆炸过程中各轴向测点的超压时程曲线如图5所示。从图5可以看出，在双爆源爆炸发生的初期冲击波超压随着时间的推移而各自达到峰值，除8号测点外的其他关于左侧爆源对称的测点(2-7、3-6、4-5)观测到的第一次压力峰值显示相同，并且压力峰值大小随着测点位置离左侧爆源距离的增大而减少，这也验证了前文所述双爆源冲击波碰撞之前互不干扰并随着传播距离的增加而发生能量耗散，而导致巷道双爆源瓦斯爆炸各测点的第一次冲击波压力峰值与近爆源的距离有关的原因是随着传播距离的增加，冲击波波阵面积增大，使得单位面积内冲击波的能量不断减小，所以即使冲击波在传播过程中没有能量损失，能量也会迅速衰减，同时我们可以观察到4-5号这对近爆源测点所观测到的第一次峰值压力为1.743 MPa，与实际巷道瓦斯爆炸能产生的最大爆轰压力相符 [20] ，而8号测点是因为它位于两爆源间的对称轴上，所以在8号测点所观测到的第一次压力峰值为两爆炸冲击波碰撞叠加生成的合成波扩散所形成的压力，通过对比8号测点以及与8号测点关于左侧爆源对称的1号测点的压力峰值可知，1号测点的第一次压力峰值为0.74 MPa而8号测点的压力峰值急剧上涨至7.1 MPa，约是前者的9.6倍，这说明了两爆源所产生的冲击波碰撞叠加前后压力峰值并非简单的叠加而是急剧攀升。
Figure 5. Time history curve of axial shock wave pressure at each measuring point
图5. 轴向各测点冲击波压力时程曲线
巷道内双爆源爆炸所形成的冲击波在两爆源对称面发生碰撞叠加后形成合成波，而合成波会接着向四周发散，这些发散的冲击波会对初始冲击波进行追赶和叠加，从而使各个测点产生第二次压力峰值，我们可以从图中看出，左侧爆源的右侧的各个测点的第二次压力峰值大小为8 > 7 > 6 > 5，这表明各测点第二次压力峰值与测点和两爆源间对称面的距离呈负相关，测点位置离两爆源间对称面越远则可观测到的第二次压力峰值越小，但左侧爆源的左侧测点未发现该规律，这可能是因为合成波的扩散与围岩相互作用更强烈导致能量耗散得更快，这就使得其对于远距离位置影响较小。
3.2. 巷道内双爆源作用下对称切面动力响应
3.2.1. 对称切面的应力变化规律
当巷道内发生瓦斯爆炸时，在冲击波传播及波系演化过程中，由于冲击波与巷道结构相互耦合，这就使得巷道结构产生动力响应而发生破坏，而从上文可知由于冲击波的碰撞叠加位置为两爆源对称面，这就使得此处会发生更为严重的破坏，所以本文关于巷道壁面的响应分析着重以巷道两爆源间的对称切面展开分析，由前人对于巷道内单爆源瓦斯爆炸的研究可知 [21] ，巷道发生单爆源瓦斯爆炸时对称切面各测点出现有效应力峰值的顺序为顶板–帮部–底板–隅角，而有效应力峰值大小的顺序为顶板–隅角–帮部–底板。图6为双爆源作用下巷道围岩对称切面的有效应力响应云图，由图6(a)~(d)可以看出在双爆源瓦斯爆炸作用下，在0.8 ms时巷道对称切面的顶板率先出现高应力区，最大应力达到23.2 MPa，随后在0.88 ms时帮部和底板相继出现高应力区，最大应力21.5 MPa，最后在1.1 ms时隅角处也出现高应力区，最大应力为26.0 Mpa，之后冲击波逐渐衰弱，围岩的有效应力随爆炸时间的增加逐渐减小。
Figure 6. Effective stress cloud map of symmetrical section of surrounding rock
图6. 围岩对称切面有效应力云图
基于流固耦合数值模型，得到双爆源瓦斯爆炸过程中对称切面各测点的超压时程曲线如图7所示。由图7可以看出，对称切面各测点出现有效应力峰值的顺序为顶板–帮部–底板–隅角，与上文有效云图观测结果一致，而有效应力峰值大小顺序为隅角–顶板–帮部–底板。由图7可以看出底板和帮部的有效应力峰值相差不大约为11 MPa，而顶板所观测到的最大有效应力峰值约为17.3 MPa，相较于底板和帮部的有效应力峰值增大了约57.27%，这是由于拱顶的弧形有利于冲击波的反射导致的，而出现最大有效应力峰值的隅角处的有效应力约为25.8 MPa，相较于底板和帮部的有效应力峰值增大了约134.54%，这是由于“犄角结构”同时放大了两侧冲击波在此处的反射效应，导致隅角的有效应力激增，从而成为对称切面所观测到的最大有效应力峰值处，之后各个测点的有效应力随着时间的增加而震荡减小，由此可知单爆源与双爆源对称切面各测点出现有效应力峰值的顺序一致，而在有效应力峰值大小的顺序上，单爆源为顶板处的有效应力峰值最大而双爆源则是隅角处的有效应力峰值最大，所以在巷道发生双爆源瓦斯爆炸的安全设计中虽然与单爆源一样都需要注重顶板和隅角处的防护，但是相较于单爆源瓦斯爆炸更应注重隅角处的防护。
Figure 7. Time history curve of effective stress at each measuring point on the symmetrical section
图7. 对称切面各测点有效应力时程曲线
3.2.2. 对称切面的壁面变形特征
Figure 8. Displacement time history curve of each measuring point on the symmetrical section
图8. 对称切面各测点位移时程曲线
基于流固耦合数值模型，得到双爆源瓦斯爆炸过程中对称切面各测点的位移时程曲线如图8所示。在双爆源瓦斯爆炸的冲击作用下，巷道对称切面各测点都在膨胀变形。各个测点位移震荡增大，最后趋于一个非零的稳定值，这说明壁面结构在双爆源瓦斯爆炸结束后受到了一定损伤。可以看到冲击波到达壁面后顶板的位移峰值最大，达到了1.37 × 10−3 m；其次是帮部，达到了1.31 × 10−3 m；再次是底板，达到了1.24 × 10−3 m；振动幅度最小的位置为隅角处，只有1.2 × 10−3 m。原因是冲击波到达隅角处最晚，当冲击波到达该处时，帮部和底板已经发生膨胀变形而向该处挤压，使该处密度增大，刚性增强，所以难以振动变形，这与由前人对于巷道内单爆源瓦斯爆炸的研究结果 [20] 相一致，所以就巷道损伤而言，双爆源瓦斯爆炸同单爆源瓦斯爆炸一样，需要着重考虑巷道顶板处的防护。
4. 结论
1) 巷道内发生双爆源瓦斯爆炸时，由于两爆源冲击波之间的相互作用使得巷道内的冲击波流场情况与单爆源的情况有着较大差异。
2) 巷道内发生双爆源瓦斯爆炸时，各个测点第一次压力峰值的大小与该处离近爆源的距离呈负相关，第二次压力峰值的大小与该处离两爆源冲击波碰撞面的距离呈负相关，轴向两爆源中心位置第一次压力峰值为7.1 MPa，而其关于左侧爆源对称位置处的第一次压力峰值为0.74 MPa，前者约为后者的9.6倍，这说明两爆源的冲击波碰撞前后压力值峰并非简单地叠加而是急剧攀升。
3) 巷道内发生双爆源瓦斯爆炸时，而有效应力的峰值大小由高到低依次是隅角、顶板、帮部、底板，壁面变形的大小由大到小依次是顶板、帮部、底板、隅角，由此可知在巷道发生双爆源瓦斯爆炸的安全设计中应该注重隅角与顶板处的防护。
基金项目
国家自然科学基金项目“矿井多爆源瓦斯爆炸传播特性及热冲击动力学机制研究”(编号：52174178)、“矿井瓦斯爆炸能量释放转化特性及热冲击能量损耗研究”(编号：52174177)。