隔离舱室污染气体扩散的大涡模拟研究

刘琦，谢承利，李伟光，刘文川，龙新平

1 中国舰船研究设计中心，湖北 武汉 430064

2 煤矿灾害动力学与控制国家重点实验室，重庆 400044

3 重庆大学 资源与安全学院，重庆 400044

4 武汉大学 水射流理论与新技术湖北省重点实验室，湖北 武汉 430072

0 引 言

船舶舱室是一个相对密闭的室内环境，其空气品质受污染源、颗粒物、温度、湿度等诸多因素影响。日本“钻石公主”号（Diamond Princess）邮轮新冠病毒感染事件使船舶系统设计人员更加清晰地意识到，必须设计可靠的船舶空气调节和通风系统、建立合理的气流组织，以有效保障舱室的空气品质，从而降低病毒通过空气传播导致交叉感染的风险。

为确保船舶传染病隔离病房内含有病毒的污染空气不侵入周围的洁净空间，通常需要设计负压系统以维持病房内部较低的负压环境。美国疾病控制中心认为，控制气流进出房间的最小压差仅为0.25 Pa，但考虑门开闭和热浮升力效应等影响，推荐隔离房间负压大于2.5 Pa[1]。我国《传染病医院建筑设计规范》[2]规定，负压隔离病房与其相邻及相通的缓冲间、走廊压差应保持不小于5 Pa的负压差。此类准则和规范为相应的负压系统设计提供了重要参考，但是无法准确评估病房内压力波动等实际因素对其内部气流组织和污染物分布的影响，同时缺少针对船舶隔离舱室负压系统设计的相关研究。近年来，计算流体动力学（computational fluid dynamics，CFD）技术已在室内污染物扩散分布研究领域得到广泛应用[3-7]。研究人员基于雷诺平均Navier-Stokes方程（RANS）方法并结合动网格技术，对人员走动和门开闭过程的流场特征进行动态模拟，结果表明，室内物体移动诱导的气流运动会对室内流场和污染物的扩散分布产生重要影响[8-10]。由于物体运动与流场特性高度耦合，固体表面诱导的流体运动具有极强的瞬态特性，大涡模拟（large eddy simulation，LES）方法较RANS方法可以更好地捕捉此类流场特性和污染物运动规律[11-12]。

压差是隔离舱室负压系统的重要设计输入，通常可参考借鉴相关标准规范，但目前相关评估指标都是基于静态参数，缺少对隔离舱室门开启的过程中污染物扩散规律的动态评估。试验研究方法具有费用高、周期长、灵活性差等局限性，CFD方法是解决此类问题的有效手段。因此，本文将以某型应急医疗救援船隔离舱室的负压系统设计为研究对象，对其典型隔离舱室开展三维实尺度污染物扩散大涡模拟研究，并对比分析空间离散精度以及时间步长对模拟结果带来的影响；重点对不同压差工况下，模拟分析隔离舱室门动态开启过程中的流场特性以及污染物扩散规律；以期为传染病隔离舱室安全性评估及压差设计提供参考依据。

1 模拟计算方案

1.1 房间三维模型

本文计算采用的三维模型如图1所示，其中隔离舱室的长（L1）、宽（W1）、高（H1）分别为4，3.9和2.6 m，缓冲间长（L2）、宽（W2）、高（H2）分别设置为2.9，3.9和2.6 m，中间舱室门的宽度（WD）为0.9 m、高度（HD）为2.0 m，高度方向即为z轴方向。初始时刻，隔离舱室压力值设定为P1，缓冲室压力设定为P2，则隔离舱室压差为 ΔP=P1-P2。由于本文压差都为负值，下文统一取其绝对值进行表述。舱室门初始时刻处于关闭位置，即舱室门将隔离舱室和缓冲间隔开。计算开始后，舱室门以角速度ω = π/2 rad/s转向缓冲间一侧，在t=1 s时刻停止运动并保持在90°全开位置。

图 1 房间三维模型示意图Fig.1 Schematic diagram of three-dimensional room model

1.2 数学模型

将隔离舱室在压差作用下门开启过程中污染物的输运和扩散过程视为不可压缩流动，且不考虑温度效应以及重力的影响。因此，采用盒式滤波函数对瞬态控制方程进行滤波后，可以得到如下LES控制方程

本文将污染气体视为连续介质，采用CO2作为示踪气体。初始时刻t= 0 s时，仅在隔离舱室存在一定浓度的均匀分布的CO2，而缓冲间内只有洁净空气。在舱室门开启过程中，空气和CO2在房门诱导和房间压差的共同作用下开始进行对流扩散。因此，采用组分输运模型描述CO2的运动如下：

2 结果及分析

2.1 网格尺度及时间步长敏感性分析

从式（3）可以发现，大涡模拟中亚格子涡黏系数与网格尺度相关，较细密的网格尺寸理论上可以求解更多尺度的流体运动，但是其计算量也会相应增加。因此，需要针对网格尺寸对模拟结果的影响进行研究，这里以房间压差ΔP=2.5 Pa的工况为例进行说明。如图2所示，3套网格Mesh 1，Mesh 2和Mesh 3的节点数分别约为148.66万、341.9万和570.9万，并对门框中心线处5个不同时刻的速度分布值进行对比。如图3所示，模拟结果在z轴方向（舱室高度方向）上呈现出相似的分布规律：在门洞上部区域速度值变化较大；在门洞下部区域速度分布较均匀。同时可以发现，基于Mesh 2的模拟结果与Mesh 1的预测结果有较为明显的差异；而随着网格的进一步加密，Mesh 3与Mesh 2的预测结果更加趋于一致。

本文重点研究门开启过程中隔离舱室污染物进入缓冲间的运动规律，因此，可定义污染物总扩散量M作为模拟评估指标，其计算公式为

图 2 模型表面网格节点分布Fig.2 Grids at the room model surface

图 3 不同网格及不同时刻中线处速度分布对比Fig.3 Velocities along the centerline at different time for different meshes

图 4 不同网格污染物总扩散量随时间的变化Fig.4 Time histories of total contaminants dispersion for different meshes

因此，本文将基于Mesh 2针对时间步长对LES模拟的影响开展进一步研究。在LES模拟中，为了确保模拟结果的可靠性以及求解过程的稳定性，时间步长Δt应该满足CFL准则，即

式中：u，v和w分别为计算域中流体速度在x，y，z方向上的分量；Δx, Δy和Δz分别为网格在x，y，z方向的尺度。可根据特征速度和特征长度初步确定LES模拟中的时间步长，因此，本文选取Δt=0.005 s作为参考基准，并设置0.01和0.002 5 s这2种时间步长进行对比分析。

图5所示为采用3种时间步长设置的LES模拟预测的中线处不同时刻速度分布情况。在门洞上部区域，3种时间步长的预测结果一致性较好；但在下部区域，在t=0.9 s和t=1.2 s时刻，0.002 5 s时间步长的预测结果与其他2种时间步长的模拟结果明显偏离。图6所示为采用3种时间步长预测污染物总扩散量随时间的变化规律，可见Δt=0.01 s计算结果出现明显偏差，而Δt=0.005 s和Δt=0.002 5 s这2种时间步长的预测结果能够较好吻合。

由此可见，网格尺度和时间步长都会对LES模拟结果产生一定影响，尽管加密网格以及减小时间步长理论上可以求解更小尺寸的流体流动特征，但也会使得计算量大幅增加，这不利于LES在复杂工程问题上的应用。因此，后续工况的计算基于本节的网格尺寸和时间步长分析结果。

2.2 压差对污染气体扩散的影响分析

选取4种压差工况，即ΔP= 0.25, 2.5, 5.0和10.0 Pa，基于LES对污染物在0～5 s时间段内的扩散规律分别进行模拟计算。计算过程中，分别对图7所示的4个监测点的速度值随时间变化规律进行监测，其中点1距离房间地面高度为0.25 m，点2高度为0.5 m，点3高度为1.5 m，点4高度为1.75 m。

图 5 不同时间步长及不同时刻中线处速度分布对比Fig.5 Velocities along the centerline at different time for different time steps

图 6 不同时间步长污染物总扩散量随时间的变化Fig.6 Time histories of total contaminants dispersion for different time steps

图 7 监测点示意图Fig.7 Illustration of the monitoring points

图8～图11所示分别为4个监测点处在不同压差工况下的速度变化规律。由图可见，4个监测点处的速度分布呈现出相同规律：在ΔP= 0.25 Pa工况下，速度先缓慢增大，在t=1 s时（即门停止运动时），速度迅速降低，随后逐渐降低至最小值并缓慢回升，最终速度值趋于平缓；在ΔP= 2.5 Pa工况，速度值在0～1 s时间段内逐渐降低，随后速度开始回升并趋于平缓；在ΔP= 5.0和10.0 Pa工况，速度值在0～0.5 s之内达到最小值，随后便开始迅速增大。值得注意的是，ΔP= 2.5，5.0和10.0 Pa这3个工况下，在所监测的4个点处都出现了明显的速度峰值，其中ΔP=10.0 Pa时速度峰值出现在t= 1 s附近，随着压差减小，速度峰值也开始减小，且对应的峰值时间开始延后。由此可见，在门开启过程的诱导作用与房间内外不同压差的配合下，流体运动出现了截然不同的规律，这对污染物的扩散也将产生影响。

图12～图15所示为4个不同压差工况下污染物扩散三维视图，以及对应时刻距离地面高度z=1 m平面上的污染物扩散俯视图。污染物浓度采用隔离舱室内初始浓度进行无量纲处理，深蓝色表示初始浓度的污染空气，白色表示洁净空气。如图12所示，ΔP= 0.25 Pa工况下，门开启过程诱导出一团污染物，随着时间推演污染物逐渐向舱室外部移动并逐渐减小，在t= 4.0 s时仍有一部分污染物停留在舱门前部区域，若此时房门关闭，该部分污染物将有滞留在隔离舱室外部的风险。将图13所示ΔP= 2.5 Pa工况的污染物扩散规律与图12对比，可见ΔP=0.25 Pa 时，污染物在舱室外部的分布范围明显减小。在图14和图15所示的ΔP= 5.0 和 10.0 Pa工况下，门诱导出的污染物体积明显减小，且污染物向舱室外部扩散的范围也得到有效控制。本文模拟结果较好地捕捉到如图16所示试验所得的门开启诱导空气流动的典型扩散特征[16]。

图 8 监测点1处不同压差工况下速度随时间变化规律Fig.8 Time histories of velocity under different differential pressures at monitoring point 1

图 10 监测点3处不同压差工况下速度随时间变化规律Fig.10 Time histories of velocity under different differential pressures at monitoring point 3

图 11 监测点4处不同压差工况下速度随时间变化规律Fig.11 Time histories of velocity under different differential pressures at monitoring point 4

图 13 Δ P=2.5Pa时污染物扩散规律Fig.13 Evolution of contaminants dispersion under ΔP=2.5Pa

图 14 Δ P=5.0Pa时污染物扩散规律Fig.14 Evolution of contaminants dispersion under ΔP=5.0Pa

图 15 Δ P=10.0Pa时污染物扩散规律Fig.15 Evolution of contaminants dispersion under ΔP=10.0Pa

图 16 门开启诱导空气流动的试验结果[16]Fig.16 Smoke visualization of door opening induced flow[16]

图 17 不同压差工况下污染物总扩散量随时间的变化Fig.17 Time histories of total contaminants dispersion under different differential pressures

图17所示为基于式（5）计算的不同压差工况下进入隔离舱室外部空间的污染物总扩散量随时间的变化规律。ΔP= 0.25 Pa压差下，进入缓冲间的污染物在t= 3.3 s达到最大值，随后开始缓慢降低，在t= 5 s时仍滞留在外部空间的污染物总量仍处于一个较高水平。ΔP= 2.5 Pa压差下，污染物总扩散量在t= 1.1 s时达到最大值，随后也呈现下降趋势，在t= 5 s时缓冲间内尚有少量污染气体。ΔP= 5.0 和 10.0 Pa工况下，污染物总扩散量都在t= 1 s之内达到最大值，且其峰值明显小于其他2种工况。在ΔP= 5.0 Pa工况下，t= 2 s时污染物总扩散量为0。在t= 3 s时刻，ΔP= 5.0 Pa工况下扩散出的污染物在压差作用下也已经全部返回隔离病房。

3 结 论

本文基于大涡模拟（LES）方法结合示踪气体组分输运方程对隔离舱室门开启过程中污染物的对流扩散规律进行了数值模拟研究。通过对空间离散尺度和时间离散步长开展敏感性分析，并对不同压差工况下隔离舱室门动态开启过程中流场速度分布特性及污染物扩散规律进行对比分析，得到如下结论：

1） 网格尺度和时间步长都会对LES模拟结果产生一定影响，以污染物总扩散量为指标对网格尺寸和时间步长进行敏感性分析，可确保本文LES模拟结果的可重复性。

2） 在门开启过程的诱导与不同压差的协同作用下，流体呈现出不同的运动规律，速度峰值的大小和出现峰值的时间有明显差异，压差越大，速度峰值越大，且出现峰值所需时间越短。

3） 不同压差工况下，污染物扩散速度和扩散范围都呈现出较大差异。当压差 ΔP=0.25Pa时，污染物总扩散量在t=5 s时刻仍维持较高水平；当压差 ΔP=2.5Pa时，污染物总扩散量峰值大幅下降，但污染物在t=5 s时刻仍未完全消失；当压差ΔP=5.0Pa时，污染物总扩散量在t=3 s时即变为0；当压差 ΔP=10.0Pa时，污染物总扩散量在t=2 s时变为0。因此，提高隔离病房负压值一方面可以减小污染物扩散峰值，同时，也可以缩短已扩散污染物在压差作用下返回隔离病房的时间。

4 致 谢

本文的数值计算工作得到了武汉大学超级计算中心的支持和帮助。