杨石刚，蔡炯炜,2，杨 亚，孙文盛，门敬敏

（1. 陆军工程大学爆炸冲击防灾减灾国家重点实验室，江苏 南京 210007；2. 73021 部队，浙江 杭州 310012）

根据中国城市建设统计年鉴[1]的记载，截至2019 年，我国可燃气体管道在城市地下空间的长度达76.7 万公里。然而可燃气体管道时常因阀门失灵、管道破损、操作失误等原因而泄漏，一旦可燃气体在城市地下浅埋管沟内积聚，在遇到点火点时，极易发生严重的可燃气体爆炸事故，同时造成恶劣的社会影响。如2013 年11 月，山东青岛一条输油管道在地下暗渠中发生了严重的可燃气体爆炸事故，造成了62 人死亡、136 人受伤[2]。2014 年7 月，台湾高雄一条输送丙烯液体的地下管道发生了严重的蒸气云爆炸事故，造成了32 人遇难、321 人受伤[3]。因此，对城市地下浅埋管沟内可燃气体爆炸的荷载进行研究，以期减轻爆炸灾害效应，具有重要的经济价值和社会意义。

学者们对长直空间内可燃气体爆炸已开展了不同类型的试验研究。王东武等[4]、司荣军[5]、Ma 等[6]分别在大尺寸的巷道、圆形管道和天然气管道内进行了可燃气体爆炸试验，研究爆炸冲击波在装置内的传播规律，得到了冲击波在受限空间内的超压峰值分布特征和衰减过程。在小尺寸可燃气体爆炸试验方面，Ciccarelli 等[7]在边长7.6 cm、长3.6 m 的管道中进行了不同障碍物拥塞率的多组天然气爆炸试验；Na’inna 等[8]在直径为16.2 cm、长4.5 m 的管道中进行了多组体积分数为10%的甲烷爆炸试验，发现了超压和火焰速度与障碍物尺寸符合平方关系；孙庆文[9]对一段尺寸为0.11 m×0.11 m×6 m、顶部设置了若干个泄爆口的小比例综合管廊进行了试验研究，结果表明：设置泄爆口往往能降低管沟内50%以上的超压和冲量值。另外，Hou 等[10]结合了管沟的实际情况，自主设计了一种Y 形结构的分支管沟，还原了山东青岛“11·22”爆炸事故的破坏场景，分析了管沟内外的超压时程曲线，划分了管沟内超压变化的不同阶段，研究表明管沟盖板会同时受到前驱波和燃烧波的影响。目前，试验研究大多是在长径比较大的小截面管道中进行的。受环境条件、试验操作危险性和外界其他复杂因素的约束，关于大尺寸原型试验研究的报道较少。

数值模拟分析具有经济性和高效性的优点，也是研究可燃气体爆炸的一种重要手段。目前，全球通用的燃气爆炸效应分析数值模拟软件主要有FLACS、AutoReaGas、Fluent、CFX 等，其中FLACS 是一款被众多国际组织和标准认证的气体爆炸分析专业软件。宫广东等[11]应用AutoReaGas 软件定量研究了管道中瓦斯爆炸的特性，分析了不同初始压力、初始温度、障碍物形状、尺寸及空间位置等条件对管道中瓦斯爆炸超压及火焰速度的影响。龚燚[12]利用Fluent 软件对规范中防火分区内一段长为200 m 的综合管廊燃气仓进行了可燃气体爆炸数值模拟分析，结果表明：当泄爆口设置在投料仓和通风口之间时，泄爆效果较好。董浩宇[13]、刘洋等[14]、陈晓坤等[15]、王涛[16]、Hisken 等[17]和Zhang 等[18]利用FLACS 软件分别对综合管廊、密闭方形管道、独头巷道、圆形玻璃管道、带有障碍物的方形管道、隧道进行了可燃气体爆炸数值模拟，研究了爆炸冲击波在不同的长直空间结构内的传播特性。目前的研究多是通过建立长直空间的数值模型分析空间内部的爆炸荷载，而对长直空间内的荷载通过泄爆口向地面传播的规律及特征的研究较少。

本文中，将基于大尺寸管沟可燃气体爆炸试验，利用FLACS 软件对管沟内可燃气体爆炸产生的冲击波通过泄爆口向地面传播的机理、超压时程曲线特征及超压峰值分布规律进行数值模拟研究。

1 模型验证

1.1 试验

试验系统的主体是管沟装置，总长度为30 m，由10 节长3 m、截面尺寸为0.8 m×0.8 m 的长方体组装而成，如图1(a)所示。试验系统还包含浓度仪、点火系统、燃气装置、气体混合风机、传感器及数据采集装置等[19]，如图1(b)所示。其中管沟端部中间位置预留小孔用以安装点火装置，管沟尾部设置有尺寸为0.8 m×0.8 m 的泄爆口，管沟顶部也设置有尺寸为0.6 m×0.6 m 的泄爆口，泄爆口的数量和尺寸可以根据实际工况灵活调整。

试验以常见的可燃气体甲烷作为研究对象，分别设置了不同的泄爆口方式和甲烷体积分数，具体工况见表1。点火点位于管沟端部中央位置，仅在管沟端部这一节充入甲烷，即可燃气云尺寸为0.8 m×0.8 m×3 m。

表1 管沟可燃气体爆炸工况记录Table 1 Working condition record of combustible gas explosion in pipe trench

1.2 模型验证及网格敏感性分析

基于上述管沟试验装置，利用FLACS 软件建立了相应的数值模型。通过分析以往城市地下管沟燃气爆炸事故可知，爆炸荷载如果通过管沟顶部或尾部的孔口直接传播到地面，通常会造成比较大的危害，因此选取表1 中试验工况1-C 来验证，点火点设置在管沟端部中间位置，管沟内中线上设置10 个压力测点（P1～P10），第1 个测点距端部水平距离为1.1 m，其余相邻测点之间的水平间距均为3 m，如图2所示。

为了验证数值模型的网格敏感性，根据管沟截面尺寸以及泄爆口大小，选取了0.05、0.10 和0.20 m共3 种网格尺寸进行数值模拟，部分测点（P1、P5、P10）的数值模拟结果与试验数据的对比如图3 所示。可以看到，当网格尺寸为0.20 m 时，数值模拟的曲线变化趋势与试验结果基本吻合，可以较好地呈现测点的超压随时间的变化规律，而当网格尺寸为0.10 和0.05 m 时，模拟结果均大于试验数据。这主要是由于预混合燃烧的亚网格模型在网格单元尺寸较小时适用性较差，对燃烧速度的计算往往存在较大偏差，模拟结果会高估爆炸超压。综合考虑计算效率及准确性，选择0.20 m 的网格尺寸来进行数值模拟。

2 地面可燃气体爆炸冲击波传播规律分析

2.1 物理模型

城市管沟通常浅埋在土壤里，由于铺设的长度较长，为了便于检修和后期维护，每隔一定距离都会设置一个检查井。根据国家标准规定，当管沟高度为0.8～1.0 m 时，检查井的最大间距为90 m[20]。因此在实际试验中，建立一段尺寸为180 m×1 m×1 m 的城市浅埋管沟模型，如图4 所示，其中检查井之间的距离为90 m，检查井上方的盖板为无约束泄爆口，为了研究管沟内的爆炸冲击波通过泄爆口向地面传播的分布规律，沿检查井X方向两侧各预留长度为45 m 的管沟，沿管沟Y方向两侧各预留长度为10 m 的地面，管沟顶部与地面平齐、地面上方空气域的高度为20 m（见图5），由于本文的边界条件设为平面波，因此可以避免反射回来的波对结果的影响。

考虑极端的爆炸情况，假设AB两端检查井之间（90 m×1 m×1 m）填充满体积分数为9.5%的甲烷，点火点位于管沟内2 个检查井的中心。环境温度与大气压强分别设置为为20 ℃、100 kPa，边界条件设置为 平面波，湍流模型采用k-ε 模型。为排除其他外界因素对管沟可燃气体爆炸的影响，假设：(1)爆炸产生的能量不会与管沟进行传热而损失；(2)管沟内没有其他障碍物；(3)壁面绝对光滑；(4)在可燃气体爆炸作用下管沟不发生形变。

在土壤及空气域所在的区域设置核心网格，网格尺寸为0.2 m，核心区域外伸展系数为1.2，如图5 所示。在泄爆口上方的空气域中，沿X、Y、Z轴方向上分别设置14 个测点，用于记录超压时程曲线的变化，由于在泄爆口附近的荷载较大、超压变化幅度较大，因此在尾部开孔附近的空气域中测点布置得较密集，在远离开孔的空气域中测点布置得较稀疏，具体设置如表2 和图6 所示。

表2 测点的坐标Table 2 Coordinates of measuring points

表2（续）Table 2 (Continued)

2.2 管沟外冲击波传播过程及机理分析

通过对管沟外部X、Y、Z轴等3 个方向不同测点处超压时程曲线特征的分析，将其划分为3 个阶段：稳定段、超压峰值Δp1段、超压峰值Δp2段，如图7～8 所示。随后研究了超压、燃料、燃烧产物以及速度矢量的动态变化过程，并由此阐述了不同阶段爆炸荷载产生的机理。

阶段Ⅰ—稳定段（t＜200 ms）：管沟内的超压值不断提高，地面上测点的超压值几乎可以忽略且基本不变。图9 为燃料和燃烧产物的体积分数二维分布云图。管沟内的可燃气体遇点火后发生燃烧反应，而此时产生的燃烧产物仅在管沟内部传播，未传播到泄爆口和地面（见图9(b)），同时管沟内的燃料还未通过泄爆口被排放到空气域中（见图9(a)），因此测点未出现超压峰值。

阶段Ⅱ—超压峰值Δp1段（200 ms＜t＜500 ms）：空气域经历第1 个超压峰值Δp1期间。随着燃烧反应的不断加剧，管沟内的火焰不断朝泄爆口方向传播，推动管沟内残余燃料通过泄爆口喷向空气域中（见图9(a)）。当t=300 ms 时，空气域中的残余燃料在流场的作用下（见图10）具有一定的运动速度，波的最大速度为数十米每秒，因此燃料会压缩分界面处的空气，形成一个有限幅度的压力扰动，导致空气域中产生压缩波（见图8）。当空气压缩波的超压作用结束后，测点就受稀疏波的影响，超压值迅速减小乃至低于大气压力。当t=350 ms 时，空气域中的残余燃料不断减少，压缩波继续向周围空气域传播，其强度也逐渐变小直至为零。

阶段Ⅲ—超压峰值Δp2段（t＞500 ms）：管沟外空气域经历第2 个超压峰值Δp2及之后的小幅振荡直至消失。在管沟中由于有泄爆口的存在，当波传播到泄爆口附近时，管沟内部气体与外部大气之间的存在明显的压力差，管沟内部压力能够得到释放，从而会出现负波峰。当t=825 ms 时，火焰沿管沟的方向到达泄爆口附近，同时受管沟尾部壁面约束的作用而向空气域内传播，火焰波阵面产生的燃烧产物具有高温特性（如图11 所示），并与空气域发生热量的交换，使得燃烧产物迅速膨胀，压缩周围空气介质，并推动周围空气向外运动，形成了密度突跃的阵面（如图12 所示），即该强扰动以冲击波形式在空气中传播，因此出现第2 个超压峰值Δp2。此时波的传播速度如图13 所示，在该阶段中波的最大速度能够达到数百米每秒。由于火焰波的强度要远大于阶段Ⅱ中压缩波的能量，同时波的传播速度也远大于阶段Ⅱ中的速度，因此第Ⅲ阶段的最大超压峰值Δp2要远大于第Ⅱ阶段的超压峰值Δp1。随着燃烧产物的膨胀，其压力、密度和运动速度不断下降，随之能量密度就逐渐减小，当燃烧产物内的压力下降到大气压力时，火焰波就脱离燃烧产物向前运动并在尾部形成稀疏区。当t=850 ms 时，火焰波不断衰减并大致呈不规则的半球形在空气域中向外传播，在泄爆口两侧沿Y方向大致呈对称分布，但是在泄爆口两侧沿X轴方向并不呈均匀对称分布（见图14），其可能的原因是受空气流场的影响，当火焰传到泄爆口后，偏向于沿火焰传播的方向分布。同时泄爆口附近还有一部分强度较小的残余火焰会向空气域中传播能量较低的火焰波，导致该阶段会出现多峰值波动状，但是随着振荡的进行，超压峰值在逐渐减小。

2.3 超压时程曲线特征分析

在三维坐标系中，选取在X、Y、Z轴方向上具有代表性的全过程超压时程曲线，如图15 所示。由图15(a)可以看出，沿X轴方向，各个测点上的第1 个超压峰值Δp1相差较小，第2 个超压峰值Δp2相差较大。在沿X轴方向上，最大超压峰值Δp2位于测点X8，约为13 kPa，且在泄爆口两侧，测点离泄爆口的距离越远，超压峰值Δp1越小。同时，也可以观察到测点X1～X7 一侧的超压峰值小于测点X8～X14 一侧的超压峰值。

由图15(b)可以看出，在沿Y轴方向上，最大峰值超压Δp2位于距离泄爆口最近的测点Y7 和Y8，其值约为7 kPa，并且在泄爆口两侧超压峰值Δp2随着距离的增加呈对称衰减的趋势。

由图15(c) 可以看出，在沿Z轴方向上，第2 个超压峰值Δp2的最大值位于测点Z1，其值约为14.5 kPa，且随着测点高度的上升，火焰波在空气域中逐渐衰减导致超压峰值Δp2也逐渐减小。

2.4 超压峰值分布规律

通过对超压时程曲线的分析，可以看出，城市地下浅埋管沟内可燃气体爆炸冲击波通过泄爆口到达地面后有2 个超压峰值Δp1和Δp2，图16 为这2 个超压峰值在X、Y、Z方向的分布。以泄爆口位置为参照点，管沟内冲击波传播的方向为X的正方向，管沟内冲击波传播的反方向为X的负方向，Z和Y方向则与坐标轴设置保持一致。由图16 可知，与超压峰值Δp2相比，超压峰值Δp1很小，且各个测点之间波动幅度不大。对于超压峰值Δp2，其值随泄爆口距离的增加而衰减。同时，在Y方向上的超压峰值会明显小于X正方向和Z方向上的超压峰值，这是由于爆炸冲击波的传播具有方向性而引起的。

为进一步分析超压峰值Δp2的衰减规律，参照空气域中的超压峰值分布的数据，可以对其进行非线性曲线拟合，得到测点的超压峰值Δp2和泄爆口之间的距离d大致满足指数函数型的关系：

式中：a、b、c为常数。

根据已有测点的数据，通过拟合分别可以得到沿X、Y、Z方向上的Δp2-d关系式：

式中：Δp2为超压峰值，kPa；d为与泄爆口之间的距离，m。该拟合公式适用于甲烷体积分数为9.5%、盖板为无约束的情况。

将式(2)中超压峰值Δp2与距离d的指数函数型方程关系作拟合曲线，如图17 所示。可以得到：曲线的拟合质量较理想，拟合度均在98.8%以上；Δp2拟合曲线在沿Z和X正方向上的超压峰值最大，但同时其衰减速度（斜率）也是最大的，同时Δp2拟合曲线在沿X负方向上的超压峰值最小，其衰减速度（斜率）也是最小的，当测点与泄爆口的距离足够远时，在Z、X、Y方向上的超压峰值大致相等。

3 结 论

利用FLACS 软件建立了城市地下浅埋管沟内燃气爆炸的数值模型，基于管沟可燃气体爆炸的试验数据验证了模型的准确性，随后分析了爆炸冲击波通过泄爆口向地面传播的机理和超压峰值分布规律，具体结论如下。

(1) 爆炸荷载通过泄爆口向地面传播过程可分为3 个阶段：阶段Ⅰ，即稳定段，该阶段测点的超压值很小，基本未受影响；阶段Ⅱ，即超压峰值Δp1段，经历第1 个超压峰值Δp1的期间；阶段Ⅲ，即超压峰值Δp2段，经历第2 个超压峰值Δp2及之后的小幅振荡直至消失。

(2) 阶段Ⅱ的Δp1较小，主要是由残余燃料挤压空气产生压缩波引起的；而阶段Ⅲ的Δp2较Δp1更危险，其主要是由燃料燃烧产生火焰波引起的。

(3) 超压峰值Δp2和与泄爆口之间距离d的关系大致满足 Δp2=a-bcd的关系式，且拟合效果较好。