瓦斯爆炸冲击波在变径管道中的传播特性

Propagation Characteristics of Gas Explosion Shock Waves in Variable Diameter Pipes

摘要: 为了研究管道截面突变对瓦斯爆炸冲击波的传播特性影响，利用ANSYS/LS-DYNA建立变径管道模型，在管道封闭端填充长度0.4 m，9.5%浓度的瓦斯混合气体，对瓦斯爆炸冲击波在变径管道内的传播特性进行了数值模拟。测得了变径管道中心，管道壁面处各测点的速度、超压。结果表明：瓦斯爆炸产生的爆炸波在管道变径区域发生了复杂反射并出现了反流现象，加大了管道变径区域的爆炸强度，使得爆炸波通过变径截面后的短时间内，出现了更高的二次超压峰值，且对管道壁面的冲击也更为严重。因此在对井下巷道进行开拓设计时，为应对可能发生的瓦斯爆炸灾害，应避免井下巷道面积突变或减缓巷道面积的变化程度。

Abstract: In order to study the propagation characteristics of gas explosion shock wave in pipeline section, using ANSYS/LS-DYNA to establish the variable diameter pipe model, the gas mixture gas with a length of 0.4 m and a concentration of 9.5% is filled at the closed end of the pipe, and the propaga-tion characteristics of the gas explosion shock wave in the variable diameter pipe are simulated numerically. The velocity and overpressure of each measuring point at the pipe center and the pipe wall were measured. The results show that: The explosion wave from the gas explosion has a com-plex reflection and a reflux in the pipe diameter area, increased explosive intensity in the pipe di-ameter area. In the short time after the explosion wave passed through the diameter cross-section, there was a higher secondary overpressure peak, and the impact on the wall of the pipe was more serious. Therefore, in the development design of underground roadway, in order to deal with the possible gas explosion disaster, it is necessary to avoid sudden change of roadway area or slow down the variation degree of roadway area.

1. 引言

在我国的煤矿开采进程中，瓦斯爆炸在煤矿事故中占据了很大比重，同时也是煤矿事故中破坏力最大的事故类型。我国学者对于煤矿开采、煤层演化及煤矿安全理论上取得了大量成果 [1] [2] [3] [4]，但由于我国煤矿瓦斯矿井的大量存在，且随着开采深度的进一步加大及高强度机械化采掘和集约化生产，瓦斯涌出量急剧增大 [5]，因此瓦斯爆炸事故时有发生。瓦斯爆炸的破坏效应主要体现在爆炸的传播阶段 [6]，同时煤矿井下巷道结构复杂，巷道截面面积突变处较为常见，发生瓦斯爆炸事故后，爆炸产生的冲击波在巷道面积突变处的传播特性也随之发生巨大变化。因此，通过研究瓦斯爆炸冲击波在变径管道中的传播特性，期望能对井下面积突变巷道、天然气输送管道中变径区域的瓦斯爆炸传播规律提供参考和理论依据，从而有效制定相应防爆抑爆措施，进而减少瓦斯爆炸造成的损害。

在关于变径管道内气体爆炸的研究上，薛景宏发现在温度作用下，变径管道薄弱位置位于管道的变径处周围；变径管道的极限应力和管道温度、温差都有联系 [7]；景国勋将冲击波的基本简短关系式应用到巷道截面积突变面，并推导出了在巷道截面积突变面冲击波波阵面压力的变化式 [8]；印华融研究表明，爆轰波进入突扩管道后会在拐点位置产生涡流并伴随熄爆的现象，且对于不同管径的突扩段，管径越大，形成马赫杆时的位置越靠后且压力越大 [9]。蔺照东发现管道截面突扩明显影响了爆炸冲击波的传输，管道截面积突然变大使气体膨胀导致明显的爆炸冲击波超压衰减 [10]。在关于连通容器中气体爆炸的研究上，尤明伟得出：由于管道的壁面效应产生的拢动与湍流加速火焰传播，使得传爆容器的爆炸压力和强度相较于作为起爆容器时均明显增加 [11]。姚世琪发现：连通系统内的爆炸强度由点火位置、容器及管道的几何尺寸等几个关键因素共同决定了火焰发展以及压力变化，而且各因素之间还会相互影响 [12]。黄佩玉通过改变连通容器中连接管长和传爆方向等条件进行了数值模拟，得出了气体爆炸后容器内的流场变化规律 [13]。

对于瓦斯爆炸冲击波在受限空间内的传播特性研究上，国内外学者也取得了许多成就 [14] [15] [16]，同时对瓦斯爆炸在复杂巷道中的传播过程进行了大量工作 [17] [18] [19] [20]。关于瓦斯爆炸冲击波在面积突变巷道传播上的研究也取得了深度的进展。但实验的局限性，使得对获取瓦斯爆炸冲击波的传播过程中的详细信息上有所欠缺，而模拟实验能够较好地还原整个冲击波传播过程，从而得出更为直观的结果。本文采用数值模拟的方法在瓦斯爆炸冲击波在变径管道中传播特性上进行了进一步的分析和完善。由于瓦斯浓度在9.5%时爆炸威力最大，故而，采用显式动力学分析软件ANSYS/LS-DYNA建立变径管道瓦斯爆炸模型，基于多物质ALE算法 [21] 对浓度为9.5%瓦斯爆炸进行仿真分析，期望得出在变径管道内瓦斯爆炸冲击波传播的有益结论，为预防瓦斯爆炸事故的发生、减少瓦斯爆炸事故灾害损失提供理论依据和技术指导。

2. 数学模型

数值实验采用非线性连续介质力学有限元分析方法，根据文献 [22]，所需方程如下：

2.1. 运动方程

位置矢量ξ确定参照坐标系中各点的位置，ξ在空间中的运动速度为υ，质点X在参照坐标系的速度为ω

$\nu =\frac{\partial x\left(\xi ,t\right)}{\partial t}\omega ={\frac{\partial \xi \left(X,t\right)}{\partial t}|}_X$ (1)

式中， $x\left(\xi ,t\right)$ ——参照点ξ在空间中的位置； $\xi \left(X,t\right)$ ——位置矢量；t——所在位置对应的时间。

2.2. 控制方程

2.2.1. 质量守恒方程

取任一连续体作为研究对象，分别用 $\partial {\Omega }_X$， $\partial {\Omega }_x$， $\partial {\Omega }_{\xi }$ 表示连续体的：物质域、空间域、参考域边界，用 ${\rho }_X$， ${\rho }_x$ 和 ${\rho }_{\xi }$ 分别表示连续体中各构形的密度。根据连续介质力学，可得：连续体的质量M在不同构形中可以表示为：

$M={\displaystyle {\int }_{\Delta \xi }{\rho }_{\xi }\text{d}{\nu }_{\xi }}={\displaystyle {\int }_{\Delta x}{\rho }_x\text{d}{\nu }_x}={\displaystyle {\int }_{{\Delta }_X}{\rho }_X\text{d}{\nu }_X}$ (2)

式中， ${\rho }_{\xi }\left(\xi ,t\right)=J1{\rho }_x\left(x,t\right)$； ${\rho }_X=2J{\rho }_X\left(x,t\right)$

2.2.2. 动量守恒方程

由动量守恒定律可知：在t时刻占参考域 ${\Omega }_{\xi }$ 的物体总动量的整体变化率，等于施加在该物体上的所有外力之和，即：

${\frac{\partial }{\partial t}|}_x{\displaystyle {\int }_{\Delta \xi }{\rho }_{\xi }{\nu }_{\xi }\text{d}{\nu }_{\xi }}={\displaystyle {\int }_{\partial \Delta \xi }{t}_i\text{d}{s}_{\xi }}+{\displaystyle {\int }_{\Delta \xi }{\rho }_{\xi }{f}_i\text{d}{\nu }_{\xi }}$ (3)

式中， $t_i$ ——作用在参考域 $\Delta \xi$ 边界 $\partial \Delta \xi$ 上单位表面上的力； $f_i$ ——作用在物体中单位质量体积力。

2.2.3. 能量守恒方程

$E=V_{sij}{\xi }_{ij}-\left(p+q\right)\stackrel{\dot{}}{V}·$ (4)

式中，V——现时构形的相对体积； $\stackrel{\dot{}}{V}·$ ——现时构形的相对体积变形速度；sij，p——偏应力张量和静水压力； ${\xi }_{ij}$ ——应变率张量；q——体积粘性阻力。

3. 数值模拟

3.1. 模型的建立

数值实验所设置的管道相应尺寸和布置如下：总长2.5 m的管道，一端设置为封闭，另一端设置为开口。距离封闭端0.5 m处设置有面积突扩的管道，面积突扩管道长度为0.8 m，截面积突扩管道末端重新与正常尺寸管道相接。其中正常管道的半径设置为0.1 m，变径管道半径为0.2 m。管道壁面厚度设置为0.02 m。在管道封闭端内充填长度为0.4 m，浓度9.5%的瓦斯预混合气体，其余管道内填充为正常空气。用水密薄膜将空气和瓦斯混合气体隔开。

整个模拟过程使用SOLID164单元，管道采用刚体模型，密度7830 kg/m³，弹性模量201 GPa，泊松比0.25，切向模量为10 GPa。瓦斯混合气体和空气均采用MAT_NULL材料模型，状态方程采用：*EOS_LINEAR_POLYNOMIAL。气体相关参数参见文献 [23]。模拟采用的单位制为统一单位制(kg-m-s)。

由于所研究的变径管道模型为对称模型，且联系到减少计算机计算时间和便于后期观察，故使用ANSYS/LS-DYNA建立1/2的变径管道模型。所建有限元模型如图1所示。

Figure 1. Pipe finite element model

图1. 管道有限元模型

3.2. 网格划分

根据对研究内容的需求，采用映射网格划分方式对模型进行网格划分 [24]，划分网格后，空气单元数量为132,210，瓦斯混合气体单元数量为20,280，管道单元数量为29,172，整个模型共计单元数量为181,662。划分网格前后的有限元模型如图2所示。

Figure 2. Finite element model before and after meshing

图2. 划分网格前后的有限元模型

3.3. 边界条件与初始条件

为了便于研究，对本次数值模型做出以下假设：假设变径管道内除瓦斯爆炸热源外无其它任意热源；假设变径管道壁面绝热，且为光滑壁面；假设变径管道内空气和瓦斯混合气体皆均匀分布；假设模型开口端为无反射边界条件；假设边界的剖面法线方向上的位移为零；忽略瓦斯爆炸化学反应的中间过程。

4. 模拟结果与分析

4.1. 瓦斯爆炸冲击波在整体管道内的传播特性

为了研究变径管道内瓦斯爆炸冲击波的传播特性，在管道内部布置测点如图3所示。

其中G点处于管道面积突扩前位置，H、I、A、C、D、E测点处于管道面积突扩后位置，B、F、J测点处于管道面积突缩后位置。所布置的测点均处于空气域内。

Figure 3. In-pipe measuring points schematic

图3. 管内测点示意

图4为所布置测点随时间变化的超压值，从图中可以看出，最靠近瓦斯混合气体的测点G在瞬间达到了超压峰值，为0.462 MPa。随后爆炸波通过第一个变径截面，依次通过测点H、I，测点H、I处的压力值也迅速达到了各点的超压峰值，对比测点G的超压峰值有递减趋势。随后通过测点A时，测点A的压力值随着爆炸波的到来首先到达了第一个峰值0.205 MPa，随后该点处的超压值在极短时间内有所衰减，之后便以极高的速率达到第二个超压峰值0.604 MPa。爆炸波的传播经过后续测点C时，同样测点C处超压瞬间达到了0.240 MPa，随后该点处的压力值在一段时间内增减后，在0.0024 s处达到了该点的超压峰值0.290 MPa。相较于测点C，测点D、E的压力值曲线图出现相似的情况，超压峰值分别为0.255 MPa和0.433 MPa。爆炸波通过第二个变径截面，依次通过测点B、F、J，在测点B处的超压峰值到达了所有测点的最高值0.693 MPa。测点F超压峰值为0.275 MPa，J点超压峰值为0.488 MPa。

Figure 4. Overpressure time history curve of each measuring point

图4. 各测点超压时程曲线

4.2. 面积突扩

图5为爆炸波通过面积突扩区域后不同时刻的管内压力分布云图，可以间接的看出爆炸波在通过第

(a) 0.0004 s (b) 0.0006 s (c) 0.0007 s (d) 0.0008 s

Figure 5. The pressure distribution after the explosion wave passes through the expanding area

图5. 爆炸波通过突扩区域后压力分布

一个变径截面后向两边扩散，并随着时间的推移与管道壁面进行了接触，随后发生爆炸波的反射。在0.0006 s时刻爆炸波阵面通过测点A，这时候测点A的超压值如上文所述达到第一个峰值，并在0.0008 s时刻，爆炸波与周围管道壁面产生的反射波在测点A处进行冲击。同时由于爆炸波通过变径截面后，与面积突扩处壁面进行反射，产生的反射波往前方进行传播，并且火焰到达后加热气体升压 [25]，几种主要因素的共同作用，使得测点A处在0.0008 s时超压达到了第二个峰值，且大于前一个峰值。该点处的后续压力值虽然在一段时间下降后出现突增，但由于能量的损失和反射波的衰减，超压未高于之前数值。

同时对管径面积突扩后管道的壁面处布置测点K如图6，爆炸波通过变径截面后，由于管道空间的突然变大，爆炸波迅速往四周扩散，便形成了如图5所示较为明显的球面波进行传播，同时管中心处的波面在球面波的最前方。在0.0006 s时刻爆炸波通过测点A时，同时波阵面也到达了位于管道壁面的测点K处，两测点的超压时程曲线如图7所示，可以看出两点超压同时达到了第一个峰值，且测点K超压大于测点A超压。这是由于通过突扩截面后，气体流动方向发生复杂改变，在爆炸波往外扩散后诱导产生湍流化程度更高的湍流，加大了爆炸冲击波的强度，从而对测点K产生了更大冲击。图8为爆炸波通过面积突扩截面后在壁面周围的速度矢量图，可以看出明显的反流现象。

Figure 6. Schematic of measuring point

图6. 测点K布置示意

Figure 7. A、K Measuring point overpressure time history curve

图7. A、K测点超压时程曲线

Figure 8. Vector distribution of velocity in suddenly expanded area

图8. 突扩区域速度矢量分布

由于爆炸波波阵面呈球形，在管道壁面发生了反射波往管中心方向运动的斜反射，四周反射波在管中心进行冲击，这也是造成测点A超压达到第二个更大峰值的主要原因。

4.3. 面积突缩

当爆炸波在保持球面波传播一段时间之后，随着传播的继续，球面冲击波波阵面不断增大，球面冲击波和管壁的接触面积越来越大，入射角也逐渐增大。当入射角增加到临界值以上后，开始出现马赫反射。随着冲击波的继续传播，最后逐渐变成均匀的平面正冲击波 [6]，即图9中0.0013 s时刻管内压力分布云图所示。

图9显示在0.0015 s时爆炸波阵面到达第二个变径截面前，与面积突扩管壁面发生接触，随后爆炸波进入变径截面，同时可以看出爆炸波在0.0016 s时刻到达测点B，此时测点B的超压达到上文中所述第一个峰值，随后在0.0019 s时刻达到了所有测点的最高超压值，出现这种情况的原因在于在爆炸波先

(a) 0.0004 s (b) 0.0006 s (c) 0.0007 s (d) 0.0008 s

Figure 9. The pressure distribution after the explosion wave passes through the narrow area

图9. 爆炸波通过突缩截面后压力分布

前在面积突扩管道内传播，由于管道截面面积突然变小，使位于管中心处的爆炸波先行通过变径截面，而两侧爆炸波与面积突缩管道前壁面发生冲击，产生反射波，并在突缩管道壁面处积聚了大量的能量，在短暂时间后，由于后续冲击波的驱使且中间超压小于两侧超压值，使周围反射波往管道突缩截面传播，0.0019 s时刻两侧反射波在管中间进行冲击，同时造成了测点B处的超压达到了第二个更高峰值。

同时火焰气流流经截面面积突缩的管段时，由于边界层较薄，速度梯度值较高，产生很强的湍流脉动。同时由于面积突缩，产生反向流动，导致湍流增大 [26]。同样导致了该处的压力值升高。图10为在面积突缩截面周围测点E，B，F在Z方向(即管道由闭口端向开口端的正方向)上的波速示意图，E，B，F测点位置如图3所示。图11所示为在爆炸波通过突缩截面后的速度矢量图。

Figure 10. Time-history curve of wave velocity in z direction

图10. z方向波速时程曲线

Figure 11. Vector distribution of velocity in suddenly narrowing area

图11. 突缩区域速度矢量分布

综上所述，瓦斯爆炸在变径管道中的传播过程是一个复杂的过程。在变径截面区域不仅产生了反流现象，加强了变径区域湍流程度，同时爆炸波与管道壁面产生了复杂的反射。湍流程度的加强和爆炸波的冲击等因素的共同作用，导致了管道变径区域中心处产生更高的超压峰值。

5. 结论

利用ANSYS/LS-DYNA软件建立三维有限元分析模型，能较好地对瓦斯爆炸冲击波在变径管道中的传播过程进行展现，并能够得到管内所设置测点的波速、超压等信息。通过数值模拟实验对瓦斯爆炸冲击波在变径管道内传播规律的研究，可得出以下结论：

1) 变径管道内发生瓦斯爆炸后，在管道变径区域的反流现象和爆炸波复杂反射的共同作用，导致变径截面区域测点出现了更高的二次超压峰值。

2) 对比规则的等截面直管，变径管道中的瓦斯爆炸传播过程更加复杂，爆炸波波速和峰值超压都要远大于规则管道，对管道壁面的冲击更加严重。

3) 在矿井井下，尤其是易发生瓦斯爆炸的区域，应尽量避免巷道截面面积的突变，或者减小截面面积的变化程度，并利用面积突变巷道瓦斯爆炸传播规律采取对应举措，来降低瓦斯爆炸所带来的损失。

基金项目

国家自然科学基金资助项目(51704111, 51374003)；湖南省教育厅资助科研项目(18A187, 18C0328)。