程方明，常助川，史 合，高彤彤，罗振敏，葛天姣

（1. 西安科技大学安全科学与工程学院，陕西 西安 710054；2. 西安空间无线电技术研究所，陕西 西安 710100）

管网中的爆炸是发生在燃料储运场所最严重的灾害之一，而障碍物的存在进一步加剧了灾害的严重性[1-2]。在实际的燃料储运场所，管道、罐体内基本上都存在各种障碍物，例如阀门以及管道面积突然缩小的区域等。当火焰在含障碍物的管道内传播时，火焰的形态、火焰传播速度、流场及超压都会受到影响而发生变化，有时因障碍物而从爆燃转为爆轰，造成更严重的后果。因此，研究可燃气体在障碍物管道中的爆炸特性，对减少爆炸事故和降低事故危害具有指导意义[3-4]。

针对障碍物对可燃气体爆炸特性的影响已有大量的研究，学者们在研究中主要关注的是障碍物的形状[5-6]、数量[7]、阻塞率[8-9]、碍物与点火源的距离及相邻障碍间距[10]等边界条件对火焰的湍流加速机理以及爆炸流场特性的影响，但受到现有实验方法和测试手段的限制，很难捕捉到全方位的精细火焰结构变化。随着计算机技术的快速发展，计算流体动力学 (computational fluid dynamics, CFD) 数值模拟等技术被引入来研究火焰[11-17]。雷诺平均模型(Reynolds average Navier-Stokes, RANS)中的k-ε 模型应用最广泛，但该方法模拟所得的火焰形状与实验中的火焰形状不能很好吻合[11]。而大涡模拟(large eddy simulation, LES)将耗散尺度的脉动进行过滤, 只对大尺度脉动进行求解, 其模拟的火焰结构及流场更精细、清晰，特别是在管网中模拟火焰的传播，可以得到较好的火焰形状和流场结构[12]。例如：Wen 等[13]利用大涡模拟，并结合了3 种亚网格燃烧模型对含障碍物的半封闭管道内预混甲烷/空气做了数值模拟研究，发现LES 与Charlette 燃烧模型得到的结果与实验结果最接近，同时障碍物会促使火焰褶皱增强、速度增大。Li 等[14]运用LES 模型对不同障碍物数量的半封闭管道内预混气体爆炸进行了研究，结果表明模拟中的火焰形状与实验有较好的一致性，而且气体爆炸后火焰的压力与速度随着气体浓度与障碍物的增加而增大。Chen 等[15]利用大涡模拟对不同阻塞的障碍物进行了研究，发现火焰在穿过障碍物时被分为4 个阶段，依次是球形火焰、指尖状火焰、喷射火焰和褶皱火焰。王公忠等[4]基于大涡模拟，发现障碍物会使预混火焰失稳，对火焰前锋面有明显的加速作用，而且大涡模拟可以很好地对火焰进行重现。

然而，以往的大部分研究中管道开口和点火位置都是在异端，而且障碍物大多都是单孔的（火焰通过障碍物时路径只有一条），但是实际情况中管道内结构多样，多孔（火焰通过障碍物时路径有多条）障碍物也很常见，且存在点火位置与管道开口同端的火焰回传情况。基于此，本文针对管道开口与点火位置同端的情形，对甲烷爆炸火焰传播受多孔障碍物的影响进行研究：利用高速摄像机对体积分数为10%的甲烷爆炸过程火焰形态进行捕捉，结合大涡模拟对火焰进行重现，通过实验和模拟进行火焰传播的动力学特性分析，揭示多孔障碍物影响预混火焰在管道中的传播规律，以期为该类爆炸事故的防控提供理论参考。

1 实验系统及三维模型

实验及模拟的对象如图1(a) 所示，主要由爆炸管道、高速摄像系统、配气系统、点火系统及数据采集系统等组成。实验管道为圆柱形，高300 mm，内径60 mm，底部设有泄爆口和点火电极，点火电源采用300 VA 电压互感器，电极间距为4 mm，点火放电时间为500 ms，每次试验点火前泄爆口自动打开。高速摄像机采用美国Phantom V411 高速数字摄像机，本次实验选取的拍摄分辨率为1 280×800，拍摄频率为2 000 Hz，每个工况试验至少重复3 次。障碍物布置在距点火位置100 mm 的位置，厚度为15 mm，开设3 个圆形孔洞，每个孔的直径为15 mm，孔与孔的角度为120°。

图1 实验系统及数值模拟网格Fig.1 Experimental system and numerical simulation grid

依据实验装置的实际尺寸进行三维建模，并采用结构化方法进行网格划分，网格数量为1 263 000，如图1(b)所示。

2 计算模型

2.1 大涡模型控制方程

文中对于管道内预混气体燃烧采用大涡模拟，大涡模拟是将全尺度范围上网格尺寸大的湍流运动通过Navier-Stokes (N-S)方程进行直接计算，而对小尺度量的采用亚网格（sub-grid scle, SGS）模型进行建模处理。模拟实验中控制方程被大涡模拟滤波后的三维瞬态方程为[4,16]：

式中：带上标“-”表示该量为LES 滤波的参量，带上标“～”表示该量为质量权重滤波的参量，ρ 为密度，ui和uj为速度分量，t 为时间，p 为压力，hs为显焓，λ 为热导率，T 为温度，R 为混合气体常数，σij为应力张量，τij为亚网格尺度应力。

应力张量σij由分子的黏度μ决定，可表示为：

式中：下标k 为Kolmogorov 尺度。

亚网格尺度应力τij的定义为：

亚网格梯度近似假设为：

式中：μSGS为亚网格黏度，PrSGS为亚网格普朗特数，cp为定压比热。

亚网格模型选用wall-adaptinglocaleddy-viscosity (WALE)，WALE 模型增强了大涡计算中的壁面处理，适合本文中所研究的工况。

2.2 增厚火焰模型

预混火焰的厚度非常小，经常是小于1 mm[17]，而且对网格尺寸要求苛刻。增厚火焰模型（thickened flame model, TFM）适用于较粗的网格来模拟预混火焰的传播[18]，为降低计算成本，本文选用该模型与LES 模型配合来进行计算[15]。增厚火焰模型的基本思想是在不改变层流火焰传播速度的情况下，增加物质扩散率和降低反应速率[4]。层流火焰的速度SL和厚度δL与分子扩散系数D 以及反应速率r 有关[19]：

式中：F 为火焰增厚因子；N 为火焰内指定的网格单元数，本文取4；Δ 为网格尺寸[17]；δL为层流火焰厚度，文中为0.26；Ω 为火焰增厚区域所乘的系数，取值范围为0～1，远离火焰面为0；β 为常数，取值为10，用于控制增厚区与非增厚区之间的过渡区域的厚度；| r¯ |为平均化学反应速率绝对值，max(| r¯ |)为反应区内平均化学反应速率的最大值。

2.3 初始条件及边界条件

实验中的泄爆口在点火前可自动打开，数值模拟中将泄爆口设置为压力出口。实验使用的可燃气体为甲烷，体积分数为10%；在模拟时，为减小计算成本，把甲烷/空气化学反应机理简化为一步反应。管道内的初始温度为300 K，初始压力相对大气压为0 Pa，初始速度为0 m/s。实验中，采用电极进行点火，电极间距为4 mm；模拟中，在管道底端使用直径为4 mm 的高温球形区域进行点火，点火温度设置为1 800 K。使用SIMPLE 算法来求解压力场和速度场的耦合，模拟求解时，每步的求解要求迭代20 步，时间步长设为1×10–5s。

3 结果与讨论

3.1 模拟结果验证

图2 是多孔障碍物条件下，甲烷/空气预混气体火焰传播的高速摄影图像。根据图2 中火焰前锋的位置变化，把火焰传播过程划分为了4 个阶段。阶段Ⅰ，层流快速膨胀阶段（0～21 ms）。该阶段前期，由于未受到管壁的阻碍，火焰前锋光滑且近球形向四周发展。之后火焰前锋横向发展受管壁阻碍，在t =9 ms 从球形转变为指尖形；t = 15 ms 时，火焰前锋裙边接触管壁，受拉伸显著，指尖形明显。阶段Ⅱ，受阻回流阶段（22～27 ms）。随着火焰锋面接近障碍物，障碍物和管顶对锋面的反作用增强，加之底端开口的缘故，在t = 21 ms 后，火焰出现回流现象，高度逐渐降低，光滑锋面开始变形；t = 27 ms 时火焰前锋降到最低点，指尖形消失。阶段Ⅲ，湍流快速发展阶段（28～34 ms）。由于未燃区燃料的持续反应，促使火焰重新向上传播，火焰传播进入该阶段后火焰前锋变得紊乱，褶皱加剧；t = 31 ms 时，障碍物下端火焰出现“空隙”并与泄爆口相通，压力被泄放；同时障碍物上端火焰前锋近似“飘絮状”穿过孔洞（图中障碍物上的黑线为安装螺丝），亮度较弱；之后随着火焰发展，在t = 34 ms 亮度达到最强。阶段Ⅳ，脉动减速阶段（t＞34 ms）。因燃料反应提供的动力和管顶、爆炸产物反作用产生的阻力在不同时刻相互作用，导致火焰前锋出现脉动现象，逐渐减速靠近管顶，在39 ms 时火焰前锋抵达管顶，障碍物下端高温火焰面开始消失，而障碍物上端火焰亮度变弱，爆炸逐渐结束。

图2 多孔障碍物管道内火焰传播的高速摄影图像Fig.2 High-speed photographic images of flame propagation in a multi-hole obstacle pipe

为了观察三维预混火焰的结构，选用T=1 800 K 的温度等值面近似作为预混火焰锋面[4]，进行大涡模拟，得到障碍物管道内火焰结构时间序列，如图3 所示。对比实验照片图2，可以看出LES 计算的结果与实验图基本吻合。LES 的模拟中火焰传播也经历了4 个阶段，且实验中每个阶段火焰的特点都被很好地复现。特别是阶段Ⅱ时，火焰出现回流现象，数值模拟可以很好地解释该现象（图4）：因为管顶封闭和孔板的阻碍，引导气流在碰壁后方向改变，向底端泄爆口低压区流动，回流现象显著（图4(b)）；火焰面受到回流以及爆炸产物的阻碍，内部形成湍流涡（图4(b)），最终诱导火焰面内部矢量转向，出现凹陷（图4(c)），火焰逐渐出现回流（图4(d)）。此外，LES 模型计算预混火焰穿过障碍物孔后层流向湍流转变的过程表现也很出色，可以很好地观察到在t = 29 ms 时火焰面穿过孔后呈凹凸不平状，之后随着火焰发展褶皱加剧；t = 36 ms 时火焰面到达管顶，而在障碍物下端因为燃料耗尽，火焰面逐渐破碎消失。

图3 大涡模拟计算多孔障碍物管道内火焰结构时间序列Fig.3 Sequence of three-dimensional flame structure by LES model in a multi-hole obstacle pipe

图4 火焰在第Ⅱ阶段回流时管道内流场变化结构Fig.4 Flow field structure in the duct during the return of the flame in phase Ⅱ

以点火位置为零点，得到火焰前锋随时间变化的数值模拟和实验结果，如图5 所示。图5 中随着爆炸的发展，火焰前锋高度逐渐升高，在接近障碍物后火焰回流导致高度降低。之后因未燃区燃料的继续反应，火焰再次向上传播，在t=40 ms 时高度达到最大值（实验值为0.208 m，计算值为0.211 m，误差小于5%）。图6 给出了火焰前锋速度随时间的变化特性，从图6 中看出，实验和模拟整体上均呈现出加速→减速→二次加速→二次减速4 个过程，这与图2 和图3 火焰传播的4 个阶段特点是完全一致的。

从图6 还可以看出，前4 ms 时，实验值大于模拟值。这是因为实验用高压放电点火，瞬时能量大，初始火核发展速度也大，而模拟采用高温点火，瞬时能量小，因此速度较小。因模拟中管壁假设为绝热，所以4 ms 之后的模拟火焰速度逐渐超过了实验值。在t=22 ms 时，由于火焰回流导致速度出现负值。在t=31 ms 时，火焰穿过孔洞，实测速度骤增至峰值17.8 m/s（模拟值为23.0 m/s）。之后管顶对火焰的约束作用及爆炸产物的吸热和压缩作用，火焰传播速度逐渐降低并呈现出振荡性。

导致误差的原因，除上述之外，可能还包括：LES 模型中的 WALE 亚网格滤波的放大效应；障碍物诱导的湍流-火焰耦合效应；实验采用高速火焰影像测试分析火焰传播速度的精度不高，等。

图5 火焰位置时程曲线Fig.5 Histories of flame front location

3.2 多孔障碍物对火焰结构的影响

图6 火焰速度时程曲线Fig.6 Histories of flame front surface velocity

图7 是火焰前锋和局部流场在障碍物前后的变化图。图中蓝色是火焰前锋；局部放大图中红色区域是已燃区，温度较高；白色区域是未燃区，温度较低；箭头表示流场中的矢量方向。

图7 障碍物前后管内火焰面及流场变化Fig.7 Changes of flame front and flow field before and after the obstacle

图7 中t=28.5 ms 时，火焰恰好穿过孔洞，此时孔中流速较大，孔与孔流出的气体相互扰动作用使流场中间形成不规则涡流；靠近障碍物上端的管壁附近，因回流与来流的“对冲”作用也产生涡流。之后在t=29 ms 时，火焰锋面因涡流的影响向两侧蔓延，呈现出近似“伞”状，之后火焰前锋发展的速度加快并形成卷吸。在31 ms后，火焰携带的低密度高温燃烧产物流入高密度低温未燃气体，密度梯度的产生促使Rayleigh-Taylor 不稳定现象出现，加上涡流之间的挤压，致使火焰扰动增加，大幅度出现卷吸、褶皱，火焰面失稳逐渐转向湍流，该现象在已有的文献中也可以见到[4,20-21]。

图8 所示是陈鹏等[22]对单孔障碍物且管顶开口的管道内，7% 甲烷爆炸火焰传播影响进行的LES 模拟结果。图8(a)显示，穿过障碍物后，火焰面开始褶皱，随着火焰发展褶皱加剧。对应的流场如图8(b)所示。可以看出火焰在前锋面出现湍流涡，而且湍流涡关于火焰面呈现对称现象。将之与本文的火焰面及流场（图7）对比，可以发现本文中的火焰穿过孔后，火焰面褶皱更剧烈，流场内更紊乱，湍流涡更多。分析原因，主要是本研究中甲烷体积分数高（10%），管道内障碍物有三个孔且管顶封闭导致的。

图8 单孔障碍物管道内LES 模拟的火焰和温度场[22]Fig.8 Flame and temperature fields simulated by LES in a single-hole obstacle pipe[22]

3.3 管道内火焰面积的变化

图9 给出不同时刻的火焰面积及其变化率。可以看出：火焰面积随时间先增大后减小，其中t=28.5 ms 是火焰面积增长的一个拐点；在该点后，火焰面积增长率骤升，在32 ms 时达到峰值1.934 m2/s；之后在t=34 ms 时火焰面积也达到峰值0.095 m2；随着燃料减少，38 ms 后，火焰面开始破碎消失，面积开始减小，增长速率开始逐渐稳定，直至燃料耗尽，火焰面最终消失。

图9 火焰面积变化Fig.9 Change of flame area

在图9 中，火焰未通过障碍物前，t=28 ms 时火焰面积达到了最大值0.036 m2，而障碍物下端的管道内表面积为0.034 m2（障碍物上端是0.032 m2），可见褶皱的出现使得火焰面积大于该段管道的内表面积。为了进一步对火焰褶皱情况进行量化对比分析，本文提出了一种衡量火焰褶皱程度的指标，即火焰面褶皱率δ，表达式为：

式中：δ 为褶皱率，Sflame为火焰的面积，Spipe为爆炸管道内表面积。

通过式(11)，可计算得到火焰在未通过障碍物阶段的最大褶皱率为5.56%，通过障碍物后，t=34 ms 时褶皱率达到最大值44.83%，可见火焰穿过障碍物后最大褶皱率增大了39.27%，说明障碍物是促使火焰褶皱程度增大的主要原因。

4 结 论

利用LES 模型对半封闭管道中设置多孔障碍物影响预混火焰传播过程进行了模拟重现，并通过实验验证了模拟结果的可靠性。基于实验和数值模拟结果的分析得到以下结论：

（1）预混火焰在开口与点火位置同端的多孔障碍物管道内传播时，会经历四个阶段，依次是层流快速膨胀阶段、火焰受阻回流阶段、湍流快速发展阶段和火焰脉动减速阶段，火焰传播在这四个阶段中呈现出加速、减速、二次加速、二次减速的波动变化特征；

（2）当火焰传播至障碍物时，由于开口在底端，引导气流受管顶和障碍物影响显著，流动方向被改变，火焰面出现回流凹陷的现象；穿过多孔障碍物后的火焰相互扰动，在涡流作用下呈“伞”状，因气流涡团、Rayleigh-Taylor 不稳定现象的耦合作用，火焰锋面破碎严重，褶皱程度增强，湍流加剧；

（3）障碍物是火焰增速和褶皱程度增强的主要原因。当火焰穿过障碍物后，传播速度与褶皱率明显增大，最大火焰传播速度可增加58.7%，最大褶皱率可增加39.27%。