参考GB 50838—2015《城市综合管廊工程技术规范》，构建典型的管廊中掺氢天然气泄漏几何模型（图1）。管廊断面尺寸为2.4 m×1.6 m（高×宽），长度为200 m。管廊内设置一条DN300的燃气管道，管道壁厚4.8 mm ；管道下方每间隔8 m设置300 mm× 300 mm×400 mm（长×宽×高）的混凝土支墩。管廊内设两个通风口，分别负责送风和回风，通风口尺寸为1 000 mm×1 000 mm，距离管廊端部3 000 mm。泄漏口设置在管道中间（距管道两端各100 m）的位置。

利用CFD的Fluent软件模拟管廊泄漏情况时，依据泄漏核心射流区、扩散区等不同的流场特征与主次，将管廊模型划分为5个区域，分别设置不同的网格划分密度。一方面在局部区域加密网格，增强计算的准确性；另一方面控制全局网格数量，节约计算资源。为了验证网格的独立性，设计3套不同网格数量的CFD模型（表1）。

图1 管廊中掺氢天然气泄漏几何模型图（mm） Fig. 1 Geometric model of hydrogen enriched natural gas leakage in utility tunnel (mm)

表1 掺氢天然气泄漏CFD模型网格方案设计表Table 1 Grid scheme of CFD model of hydrogen enriched natural gas leakage

采用3套网格方案在相同的条件下进行迭代计算并收敛，由于该研究更加关注可燃气体体积分数数据，因此重点考察核心区不同网格密度下的掺氢天然气体积分数和射流速度分布（图2）。图2横坐标的“高度”代表射流口到下游的距离，在射流扩散的核心区（即高度0.7~1.2 m的范围内），与中、高密度网格相比，低密度网格模拟结果在射流速度和掺氢天然气体积分数的计算上出现明显的分离，综合考虑计算资源的合理分配与满足工程精度要求，最终选择4×10⁶网格数量进行网格划分。

图2 掺氢天然气泄漏模拟网格独立性检验结果图Fig. 2 Grid independence test results of hydrogen enriched natural gas leakage simulation

考虑掺氢天然气泄漏的可压缩属性及射流后在空间中扩散的湍流特点，采用耦合组分输运的三维可压缩气体动力学湍流模型，选取k-ε realizable两方程模型^[23]、标准壁面函数，组分采用一阶迎风差分模式，能量、密度、湍流等采用二阶迎风格式，设置管道的竖直方向截面为监测面，检测组分体积分数做为收敛判定依据。

依据相关标准与文献设置边界条件如下：①依据GB 50838—2015并由文献[24]中燃气舱室通风速度公式计算得到送风口速度，将其设置为速度入口条件；②回风口连接外部环境，因此设置为大气环境压力的压力出口条件^[24] ；③考虑中压等级的入廊工况，依据GB 50028—2006《城镇燃气设计规范》，设置掺氢天然气管道运行压力为 0.2 MPa、0.4 MPa，泄漏孔口径为6 mm、10 mm^[25]，掺氢比为0、20％^[26] ；④其他壁面设置为无滑移壁面条件；⑤考虑到网格尺寸和库朗数的要求，时间步选取10^-4量级，共计算30 s的泄漏过程。

利用文献[21]中甲烷在管廊中泄漏扩散的实验结果对模型可靠性进行验证。基于文献[21]中的物理实验模型及泄漏孔处的质量流量进行模拟，得到文献中设置的9号、10号采样点处的甲烷体积分数模拟值，并与实验值进行对比（图3）。两者得到的甲烷体积分数变化趋势一致，平均偏差不超过15％，对于工程问题，偏差在可接受的范围之内，验证了所建数值模拟模型的可靠性和准确性。

图3 甲烷泄漏扩散的体积分数模拟值与文献[21]中实验值对比图Fig. 3 Comparison of simulated volume fraction values of methane leakage and diffusion with experimental values in reference [21]

设置掺氢比20％、管道压力0.2 MPa、泄漏孔径10 mm，模拟得到泄漏孔方向分别朝向上、向下、向后及向前工况下，管廊中掺氢天然气泄漏后监测线上可燃气体体积分数随着管廊长度方向的距离和泄漏时间的变化曲线（图4）。

当泄漏孔朝向上（图4a）时，由于泄漏口位于100 m处的位置，可燃气体泄漏后向管廊顶部射流，在一定范围内，掺氢天然气体积分数分布呈正态对称形式，且随着泄漏时间的增加，掺氢天然气体积分数稳步增大。在泄漏孔的上游约97 m处，掺氢天然气体积分数出现了阶梯式下降，这是由于掺氢天然气泄漏射流后向管廊的上游传递动量，而由国标规定的通风次数和管廊断面面积可知管廊内空气以约0.33 m/s的速度向下游运动，两者形成湍流消耗动量，阻碍了掺氢天然气向上游的扩散；随后，由于掺氢天然气体积分数差异，泄漏掺氢天然气还会继续向上游传递动量，该扩散面在约91 m处稳定，即泄漏扩散向上游传递的极限。掺氢天然气向下游扩散的过程中，约106 m处形成一个突变，这是由于掺氢天然气泄漏在此处的速度方向与管廊内的通风速度同向，但空气的伴流速度不同。泄漏的扩散速度大于下层空气的速度，这种情况下会形成漩涡，在此之后则呈现平缓的发展趋势，随着时间的推进而逐渐增大。

当泄漏孔朝向下（图4b）时，泄漏孔位于管道的正下方，泄漏后射流受到地面的阻挡衰减迅速，虽然此时泄漏量与泄漏孔向上时是相同的，由于速度衰减快，最大甲烷体积分数为3.8％，小于泄漏孔向上时的5.8％，但向下游扩散的速度更快。掺氢天然气体积分数曲线的整体发展趋势与泄漏孔向上时相似，均在泄漏孔100 m处达到最大值，在向上游扩散的过程中，由于与来流通风方向相反，存在大量漩涡使扩散速度与通风速度相抵，在98 m处出现阶梯式下降，同时在95 m处扩散动量完全消耗，是向上游传播的极限位置。在向下游传播的过程中，同样由于射流扩散速度与通风空气存在速度差，同时由于壁面的干扰，出现了大量回流，在向下游扩散过程中出现了第二个峰值，其高度逐渐与第一个峰值接近，是稳态过程后的极限值。

当泄漏孔朝向后（图4c）、向前（图4d）时，由于泄漏后射流速度的方向不是直接通过监测线，而是通过墙壁的反射与空间内的耗散逐渐扩散到监测线位置上，这两种情况下掺氢天然气体积分数曲线的形状更加不规则。但依然可以看出其在泄漏孔坐标处达到最大值，向上游扩散过程中存在阶梯下降、最远扩散点及第二个峰值，同时掺氢天然气体积分数逐渐增加的典型特征。

图4 不同泄漏孔朝向时掺氢比为20%的掺氢天然气泄漏后体积分数变化曲线Fig. 4 Changing curves of volume fraction after leakage of natural gas with 20% hydrogen under different directions of leakage aperture

设置管道压力0.2 MPa、泄漏孔径10 mm、掺氢比为0、泄漏孔朝向分别向上与向下，模拟得到管廊中纯甲烷泄漏后监测线上甲烷体积分数随管廊长度方向的距离和泄漏时间的变化曲线（图5）。与掺氢比20％的掺氢天然气泄漏时的模拟结果（图4a、图4b）进行对比可见：每个对照组中，掺氢比20％的掺氢天然气泄漏后的体积分数均小于纯甲烷的泄漏工况，但曲线的形状及其随时间的发展趋势是相似的，且均具备图4中各个阶段的典型特征。在工程误差范围内可以认为，由于甲烷与氢气在泄漏前已经混合均匀，泄漏后动量大、湍流强度高、扩散发展迅速，两者扩散传质能力的差异不会对甲烷整体分布规律产生明显影响，只是由于甲烷泄漏量的减小，会使空间内分布的甲烷体积分数近似按比例减小。

图5 纯甲烷（掺氢比为0）泄漏后体积分数变化曲线Fig. 5 Changing curve of volume fraction after leakage of pure methane (with hydrogen blending ratio of 0)

设置管道压力0.4 MPa、泄漏孔径10 mm、掺氢比为20％、泄漏孔朝向分别向上与向下，模拟得到管廊中掺氢天然气泄漏后监测线上气体体积分数随着管廊长度方向的距离和泄漏时间的变化曲线（图6）。与管道压力为0.2 MPa时的模拟结果（图4a、图4b）进行对比可见，随着泄漏压力的增大，泄漏口处的流量增大、速度增大，掺氢天然气具有更高的射流动量，从而向周边扩散的结构特征受到小尺度湍流耗散的影响减小，图4中表现的典型特征更加明显，可以清楚地看到其阶梯和峰值的分布演化趋势。

图6 管道压力为0.4 MPa时掺氢比20%的掺氢天然气泄漏后体积分数变化曲线Fig. 6 Changing curve of volume fraction after leakage of natural gas with 20% hydrogen at pipeline pressure of 0.4 MPa

设置管道压力0.4 MPa、泄漏孔径6 mm、掺氢比为20％、泄漏孔朝向分别向上与向下，模拟得到管廊中掺氢天然气泄漏后监测线上气体体积分数随着管廊长度方向的距离和泄漏时间的变化曲线（图7）。与泄漏孔径10 mm时的模拟结果（图4a、图4b）进行对比可见，6 mm孔径的泄漏流量更小，射流整体的动量低、空间中耗散快，因此向管廊上游扩散的能力远远弱于10 mm孔径泄漏后的扩散。当泄漏孔向上时，泄漏孔径较小时气体仅能扩散至约94 m的位置（图7a），而泄漏孔径较大时气体可以扩散至约91 m的位置（图4a）；泄漏孔向下时，6 mm空间泄漏后气体可以扩散至上游约97 m的位置（图7b），而10 mm孔径泄漏后气体可以扩散至约92 m的位置（图4b）。同时，由于气体泄漏的初始总动量较小，在空间中扩散时容易受到通风伴流的影响，湍流效应加剧，从而使管廊长度方向上甲烷体积分数随时间发展变化的规律性减弱（图7b）。

图7 泄漏孔径为6 mm时掺氢比20%的掺氢天然气泄漏后体积分数变化曲线Fig. 7 Changing curve of volume fraction after leakage of natural gas with 20% hydrogen with leakage aperture of 6 mm

泄漏发生、管廊内的可燃气体检测探头检测到掺氢天然气后，事故通风启动，管廊内的通风频率由6次/h提高至12次/h。泄漏孔下游可燃气体与管廊内空气混合的体积分数仅取决于泄漏量与通风频率的平衡关系，先由小孔泄漏公式计算泄漏量，再由文献[24]中燃气舱室通风速度公式计算通风频率，可构建可燃气体体积分数与流量的平衡公式，用于描述泄漏后混合均匀时的可燃气体体积分数：

式中：C为管廊内可燃气体的体积分数；Q₁为泄漏流量， m³/s；Q₂为通风机的流量，m³/s；φ为孔口系数，对于圆形泄漏孔可以取0.85；d为泄漏孔直径，m；p为管道压力，Pa；ρ为掺氢天然气常压下的密度，kg/m³；R为掺氢天然气气体常数，J/（kg·K）；T为掺氢天然气温度，K；γ为掺氢天然气比热比；n为燃气舱室每小时通风次数；V为管廊体积，m ³ 。

由于掺氢比影响物性参数，可以作为已知量，进一步整理该公式，得到管廊内可燃气体体积分数与泄漏孔直径、管内压力的函数关系：

若掺氢天然气掺氢比为20％，对于不同的泄漏直径、管道压力及管廊内可燃气体体积分数，将式（4）绘制为二维图上的掺氢天然气体积分数等高线（图8）。可见，在现有12次/h的事故通风条件下，不论管道运行压力、泄漏孔直径在何种范围内变化，均可以保证管廊内可燃气体体积分数低于爆炸极限。对于现有的燃气管廊建设而言，设计按国标执行12次/h的通风标准，入廊的管道压力通常为中压0.4 MPa以下，则对于孔径不大于12 mm的小孔泄漏，12次/h的事故通风频率可以满足管廊内可燃气体体积分数不高于爆炸下限的要求。然而，GB 50838—2015规定可以允许1.6 MPa的天然气管道入综合管廊，对于此压力下，发生约7 mm直径的小孔泄漏就会导致管廊内可燃气体体积分数过高，有发生爆炸的危险，此时应当重新审查不同直径小孔泄漏发生的概率和12次/h通风频率的合理性。

图8 12次/h通风条件下掺氢比20%的掺氢天然气泄漏后管廊内可燃气体体积分数等高线图Fig. 8 Contour map of volume fraction in utility tunnel after leakage of natural gas with 20% hydrogen under 12 times/h ventilation