泄爆条件对预混H2/空气燃爆特性影响的数值模拟

周宁，徐莹莹，陈兵，李雪，乔世伟，袁雄军，刘俊，黄维秋，赵会军

（1 常州大学江苏省油气储运技术重点实验室，江苏常州 213164；2 中国安全生产科学研究院，北京 100012）

氢气是一种理想的燃料，但因具有高的燃烧速率，使氢气爆炸过程能够在极短的时间内完成，较其他碳氢燃料具有更严重的爆炸危险性，若在受限密闭空间内不及时泄压，火焰加速燃烧转为爆轰，将对周围造成极大的危害。泄爆能有效降低可燃气体的爆炸强度，减小事故灾害，成为近年来学者们研究的重点[1-2]。郭强等[3]研究了不同泄压面积条件下预混可燃气体爆燃过程中的压力变化，发现小面积泄爆口开启后，压力先下降后上升，且第2峰值较大，大面积泄爆时第2 峰值较小。杜扬等[4]研究了含有双侧分支结构的管道内油气泄压过程中的爆炸超压和火焰演变特性，发现爆炸超压出现了3个典型的超压峰值即p1、p2、pmax，其中p1的形成与管端的密封材料瞬时破裂有关，p2与分支结构泄压有关，而pmax受管内爆炸强度与火焰加速协同效应影响。路长等[5]研究了四氟乙烷（R134a）对甲烷/空气爆炸的影响，发现R134a 的加入会延迟破膜时间，但破膜压力主要受泄爆膜自身影响。李乾等[6]开展不同泄压膜材料、泄压膜层数及泄压口位置对方管内爆炸压力特性的影响，发现随着泄压口位置接近点火端，各泄压膜约束条件下的最大泄爆压力曲线均呈现Z形规律。Wan等[7-8]研究了含有障碍物的通风管道中通风口位置、大小及障碍物对管道内爆炸特性的影响，发现侧面通风口应放置在易燃点附近，并设置在障碍物前面，以发挥理想的安全缓解作用。

目前研究氢气泄爆主要关注泄压面积以及泄爆口位置，对于泄爆口破膜压力的研究大多采取小尺寸实验装置，而大尺寸管道的泄爆口破膜压力的研究较少。因此，本文针对工业生产中的大长径比受限空间泄爆问题进行数值模拟研究，重点分析泄爆口处不同的破膜压力下的火焰传播速度、压力、流场、火焰结构相互耦合作用，揭示H2/空气预混气体燃烧过程中的火焰传播机理，为降低事故灾害提供指导。

1 数值模型及控制方程

1.1 控制方程

预混气体在受限空间内的爆炸过程是一个复杂的燃烧化学反应过程，且具有较强的湍流流动性。大涡模拟（large eddy simulation，LES）能够更清晰地反映出火焰及流场的大涡结构，特别是在受限空间中模拟火焰的传播，可以得到较好的火焰精细结构和火焰内外部流场特性。大涡模拟（LES）较适合燃烧湍流的计算，其模型引入Favre 滤波函数来对Navier-Stokes方程进行滤波处理，分解出只含有大尺度涡的运动方程，而小涡对流场的影响通过在控制方程中引入附加应力项建立亚格子尺度模型（subgrid-scale model）[9]。

大涡模拟控制方程为式(1)～式(3)。

式中，上标横线、波浪线分别为LES、质量权重滤波的参量标注；ρ为密度；p为压力；ui、uj为速度分量；t为时间；μ为分子黏度；τij为亚格子尺度应力；hs为显焓；λ为热导率。

通过求解密度加权平均反应过程变量的输运方程，对火焰锋面传播进行建模，采用Favre 滤波平均的反应进程变量RANS方程为式(4)。

利用预混燃烧模型C方程的Zimont燃烧火焰面亚格子模型来反映预混火焰反应进展变量，描述火焰朝向位置的轨迹，用火焰锋面褶皱和加厚模型计算Zimont湍流火焰速度，如式(6)。

1.2 物理模型和初始条件

物理模型见图1，为与实际工业情况接近，拟选用长径比（L/D） 为48、管长为6m、截面为125mm×125mm 的方形管道进行模拟。采用大涡模拟（LES）和Zimont燃烧模型，速度/压力耦合选用PISO 算法，对流项采用二阶迎风格式，扩散项采用中心差分格式，壁面边界条件为无滑移绝热边界条件，时间步长为0.1ms。H2/空气在计算区域内充分混合（当量比为1），形成稳定的常温常压预混气体，在左端进口端面中心位置点火，点火能为1J。在密闭管道中随着燃烧的进行不及时泄压，压力会一直呈上升趋势，因此设置泄爆口位于2m处，任少锋等[11]提出泄压比大于30%时，爆炸压力与火焰速度峰值几乎不受泄压口比率的影响，所以设置泄爆口直径为60mm，泄压比为18%（泄压比是指泄爆口面积与截面面积之比），泄爆条件根据泄爆口的破膜压力改变设置，各工况的泄爆口破膜压力分别为0.1MPa、0.3MPa、0.5MPa、0.7MPa。模拟过程中，先设置泄爆口为管道壁面，当管内压力达到预设的泄爆口破膜压力时，将泄爆口壁面设置为压力出口（pressure-outlet），以达到在实验时泄爆膜会自动打开的效果，从而达到泄压的目的。

图1 物理模型

1.3 网格无关性验证

在大涡模拟（LES）过程中，对于大尺度湍流采用直接求解的方法，而对于尺寸较小的湍流则利用建立的模型求解。其计算特性决定了网格尺寸和时间步长越小，计算结果越精确。但是，由于网格尺寸越小，计算量就越倾向于直接模拟法，因此在网格划分中应尽量平衡计算精度和计算成本。本文在考虑计算成本的基础上，采用4 种网格尺寸（1mm、2mm、4mm、5mm）来分析验证模型的独立性，结果如图2所示。综合考虑精度和计算效费比，最终计算网格为4mm×4mm×4mm，对泄爆口处进行局部加密。

图2 网格无关性验证图

1.4 模型验证

根据文献[12]中有障碍物方形管内的爆炸试验，预混气体燃烧在障碍物的激励作用下能迅速从层流燃烧转变为湍流燃烧，满足数值模拟验证对燃烧过程的要求。管道参数为：截面尺寸为0.1m×0.1m，管道长度为0.5m，堵塞率为0.5，物理模型的建立与实验系统相对应。H2/空气在计算区域内充分混合（当量比为1），设置氢气的摩尔质量20.8g/mol，氢气的摩尔分数占29.6%，氧气占14.8%，氮气占55.6%，氢气的放热量为1.43×108，层流火焰速度设置为2.1m/s[13]。计算结果见图3，大涡模拟（LES）能再现实验趋势。在图3(a)中，数值模拟中火焰结构的演变与实验结果基本一致，且在相近的时刻出现了分叉火焰（红色实线框所示）。在图3(b)和图3(c)中，LES 预测的火焰速度和超压均大于实验结果，火焰速度和超压振荡的计算频率与实验数据误差在7.21%左右。其主要原因是在数值模拟过程中没有考虑管壁的冷却效果。预测的火焰结构、火焰传播速度和超压与实验结果基本吻合，从而验证了大涡模拟（LES）模型的可靠性。以上结果也表明了LES适用于氢气湍流燃烧的计算。

2 破膜压力对管道内燃爆特性的影响

2.1 破膜压力对火焰传播速度的影响

图3 实验与计算结果的对比

图4 火焰传播速度图

图4为各破膜压力条件下火焰传播速度随着点火端距离的变化曲线图。从图4中看出，各工况下火焰传播速度曲线存在3 个波峰、2 个波谷。Guenoche[14]和Gonzalez[15]发现在较长的管道中（长宽比约为20），火焰锋面面积不断发生变化，火焰锋面反转能够反复出现和消失，从而使火焰传播速度也出现相应的加速与减速过程。如图5所示（以密闭管道和破膜压力为0.3MPa 管内火焰结构为例），燃烧初期火焰快速膨胀，火焰传播速度逐渐增至第一个波峰。随后，火焰传播过程中受壁面影响，伴随着火焰结构由指形向“郁金香”转变，火焰开始减速，火焰锋面面积开始减小，火焰传播速度下降至第一个波谷。随着火焰传播，由于“郁金香”形成后火焰前锋逆转而导致火焰锋面面积再次增大，使火焰传播速度增加，直至火焰传播速度到达第二个波峰。此后，由于声波震荡的作用与火焰燃烧的不稳定性，伴随着火焰锋面的反复逆转，火焰传播速度也相应出现反复的加速与减速过程。

此外，由图4 可知，在管道内泄爆口开启后（图中圆圈为泄爆口开启时刻），除破膜压力为0.1MPa 工况火焰传播速度对比密闭管道（closed pipe）增大外，各工况火焰传播速度整体呈下降趋势，且减小了各工况的速度极值。各工况的加速或减速效应主要是由泄爆口开启后气流泄放作用引起的，该分析可由火焰锋面流场图（图6）进行论证。而各工况气流泄放强度的大小主要是泄爆口与火焰峰面的相对位置。当管内泄爆口开启时，火焰峰面在泄爆口的左侧，则泄爆口开启对火焰传播有激励作用，火焰峰面在泄爆口的右侧，则泄爆口开启对火焰传播有抑制作用。对于破膜压力为0.1MPa 的管道，火焰峰面在泄爆口左侧，未燃气体泄放，导致管道内压差增大，且流场未发生反向，同向气流运动激励火焰加速前进，导致管内第二个波峰时的速度峰值达到215m/s，对比密闭管道速度峰值200m/s 反而上升7.5%。而破膜压力为0.3MPa、0.5MPa、0.7MPa 的管道内火焰锋面在泄爆口右侧，燃烧产物从泄爆口泄出，减小管内压差，且管内流场发生反转，未燃气体气流向点火端方向流动，由此说明当泄爆口开启后，气流泄放形成逆向气流运动，抑制火焰传播速度。其中，当管道破膜压力为0.3MPa 时，由于火焰锋面离泄爆口越近，泄流速度较快，火焰传播速度迅速降低，火焰传播速度峰值比密闭管道下降39%，达到了较好的减速效果，说明在密闭管道速度还未达到峰值时提前泄压，能够有效降低火焰传播速度。对于在破膜压力为0.5MPa、0.7MPa 的管道内，火焰锋面距离泄爆口位置越远，泄流泄放作用越弱，火焰传播速度下降幅度越不明显，在第三个波峰时，火焰传播速度有一定程度的下降，在泄爆口开启前，管道内火焰传播速度已达到峰值，泄爆效果较差，泄爆口开启的时间较晚对火焰传播过程有一定的影响。

图5 火焰结构图

图6 泄爆口开启前后时刻的流场图

2.2 破膜压力对管内压力的影响

图7 管道内压力曲线图

图7为管道内压力曲线图。泄爆口开启后，相比于密闭管道，管内压力均呈下降趋势，且破膜压力越小，管内超压峰值越小。在0.21s 时，密闭管道内超压峰值达到1.18MPa。各破膜压力条件下，分别在0.018s、0.046s、0.060s、0.070s 开启泄爆口，出现超压峰值分别为0.35MPa、0.36MPa、0.5MPa、0.7MPa，由此可知，泄爆口开启越晚，管内超压峰值越大。同时，对比密闭管道，超压峰值分别下降了70.3%、69.4%、57.6%、40.7%，并在0.13s 时，各工况条件下开启泄爆口使管内压力能降低到0.07MPa以下，说明在各工况下泄爆口开启均能很好地达到泄爆的效果。泄爆口开启使气流泄放是压力下降的主要原因。如图8所示，对于泄爆口破膜压力为0.1MPa 的管道，压力流场未发生反向，导致泄爆口开启后管内压力是逐渐上升至0.35MPa后才缓慢降压。除破膜压力为0.1MPa工况外，其余各工况泄爆口开启前压力流场都向点火端方向流动[如图8(a)]，泄爆口开启后，从Z轴与X轴切面的压力流场图中可以看出，管内流场均向泄爆口处流动[如图8(b)]。由此可知，泄爆口越早开启，火焰锋面离泄爆口越近，气流泄放强度越大，此时预混气体燃烧还未充分发展，各工况下的峰值超压也就越小，从而说明提前泄放气流，能够有效降低管内的压力，抑制火焰传播。

由图7可知，各工况条件下，泄爆口开启附近管内压力出现大幅度波动现象（如图中黑框所示）。在各工况条件下，影响压力波动是由泄爆口附近产生的涡团引起的（如图8 红色圆圈所示），泄爆口开启引起压差使管内流场形成逆向流动，同时在管道内，管道内燃烧产物快速膨胀产生压缩波，压缩波受壁面影响沿管内横向传播，当传到管道末端反射形成反射波，产生逆向流动，加强气体泄放的强度，扰乱管内流场，在泄爆口附近产生强烈的湍流现象，产生的涡团影响局部的热能和物质交换速率，引起压力大幅度波动。

图9为不同泄爆口破膜压力条件下压力上升速率曲线图。爆炸压力上升速率（dp/dt）是指在爆炸过程中测得的爆炸压力随时间变化曲线的最大斜率，是衡量燃烧速率与爆炸强度的标准，在相同压力峰值条件下，爆炸压力上升速率越快，爆炸产生的破坏力越大[16]。从图中发现，在密闭管道内，压力上升速率在0.075s时达到最大。各工况的压力上升速率相比密闭管道的压力上升速率幅度降低，爆炸强度减弱。由前面分析得到泄爆口开启后气流泄放，火焰传播速度较密闭管道呈下降趋势，燃烧速率降低，火焰发展缓慢，压力逐渐下降，导致压力上升速率减小。图9虚线标出了管道内压力上升速率最大和最小数值，在密闭管道内，压力上升速率上升幅度最大，说明爆炸强度最强。相比破膜压力为0.5MPa、0.7MPa 的管道内振幅区间，破膜压力为0.3MPa 的管道压力上升速率上升幅度较小，下降幅度最大，爆炸强度弱，说明提前泄压，能够有效地降低压力上升速率，达到较好的泄爆效果。

图8 泄爆口开启前后时刻的压力流场图

图9 压力上升速率图

3 结论

采用LES模型和Zimont燃烧模型对预混可燃气体进行三维模拟，得到以下结论。

（1）各工况下火焰传播速度曲线存在3 个波峰、2个波谷。由于火焰锋面的反复逆转，火焰锋面面积发生变化，使火焰传播速度存在相应的加速或减速过程。

（2）除0.1MPa 工况外，其余各工况火焰传播速度整体呈下降趋势，有效地降低了各工况的速度极值。泄爆口开启后产生加速或减速效应是由气流泄放引起的。破膜压力为0.1MPa 时，气流泄放出未燃气体增大管内压差，火焰传播速度峰值上升了25%，而在破膜压力为0.3MPa 的管道内，气流泄放作用最强，泄爆口减速效应较好。

（3）泄爆口开启，气流泄放作用占主导地位，导致管内压力较密闭管道整体呈下降趋势，且泄爆口处破膜压力越小，管内超压峰值越小。同时，泄爆口附近产生的涡团使压力在泄爆口开启后出现大幅度震荡。泄爆口开启使管内流线向泄爆口处集聚，气流泄放作用引起减压效应，火焰锋面距离泄爆口位置越近，气流泄放强度越大。

（4）对比密闭管道，各工况的压力上升速率振幅减小，爆炸强度减弱。当管道内压力达到0.3MPa 时开启泄爆口时，压力上升速率的下降幅度最大，比破膜压力为0.1MPa、0.5MPa、0.7MPa时开启泄爆口泄爆效果好。