氢气迁移扩散二维模型分析

2020-06-21侯丽强曹学武

核科学与工程 2020年2期

李颖，侯丽强，杨帆，曹学武，*

(1.上海交通大学机械与动力工程学院，上海 200240;2.中国核动力研究设计院核反应堆系统设计技术重点实验室，四川成都 610041)

自福岛核事故发生以来，核电厂氢气风险一直是国内外研究热点。除堆芯熔化导致的锆水反应会产生大量氢气并引发爆炸风险外，乏燃料水池也存在由冷却系统冷却能力丧失导致乏燃料组件裸露并引起锆水反应的风险，而乏燃料组件熔化后形成的熔融物流入乏燃料水池底部，也可能发生MCCI并产生氢气[1]。除此以外，氢能由于其清洁、储运方便、利用率高等优点，在能源与化工领域具有重要应用前景。而在氢气制取、储运的过程中，储存容器材料长时间与氢气接触会产生氢脆现象，并导致氢气泄漏风险[2,3]。因此，开发氢气迁移扩散模型对于能源安全具有重要意义。

本文在借鉴现有研究的基础上，考虑氢气扩散过程中夹带现象及空气阻力的影响，建立大空间内氢气自源项释放后，迁移上升扩散的二维数学模型，并结合典型实验数据验证了理论模型的可靠性，为氢气风险的评估提供支持。

1 计算模型

1.1 控制方程

在二维圆柱坐标下，忽略温差，单组分氢气自容器下部中心注入充满空气的密闭空间后，在空间内的迁移扩散可以由以下方程描述，定义竖直向上为z方向。

连续性方程：

(1)

z方向动量方程：

(2)

组分方程：

(3)

式中:u、v——分别为轴向速度和径向速度;

χ——氢气的质量份额;

ρ——混合物的密度;

ρa——空气密度;

τd——空气阻力。

研究表明[11]，对于任意沿z方向的高度处，气体速度、密度及浓度组成的函数呈现高斯分布关系，描述如下：

(4)

(5)

(6)

(7)

式中:u0——气体中心速度;

B——描述氢气源运动的特征参数，被定义为射流在某高度平面上的宽度，指同一高度平面上中心点至该平面上速度为中心点速度的1/e的点的径向距离(见图1)；

ρcen——气体中心密度；

ρ0——氢气密度;

χcen——气体中心的氢气质量份额。λ为描述氢气进入空间后沿径向(垂直于z方向)展开的特征参数，一般取1.16[12]。

图1 射流几何示意图

1.2 夹带系数

对连续性方程、动量方程、组分方程在径向(垂直于z方向)进行积分处理：

(8)

(9)

(10)

其中，α为夹带系数，表示氢气源运动至某高度z处，夹带卷吸周围空气的能力。采用Hirst[13]针对动量射流经典理论推导的夹带系数经验关系式描述，如公式(11)所示：

(11)

其中，FrL为中心弗洛德数，定义式为：

(12)

将公式(4)、公式(5)、公式(7)代入公式(10)得：

(13)

1.3 阻力系数

(14)

式中:D——特征尺寸，在喷口处即为喷口直径;

U0——氦气喷入口处速度。

公式(8)、公式(9)可表示为：

(15)

(16)

其中，f为阻力系数，由Blassius方程[14]记为：

f=bRe0-1/4

其中，b取值0.085。

1.4 求解方程

公式(15)和公式(16)即为得到的描述氢气自喷放口处至大空间内迁移扩散特性的控制方程。由于该方程组仍然表现出了非线性，采用差分格式对偏微分方程进行数值求解，得到B与ucen的关系，进而可由公式(14)和公式(7)求得对应的氢气质量份额。

将公式(15)和公式(16)联立变形，得到中心线上以a为基准，高度差为Δz处的ucen和B的计算式：

(17)

(18)

由公式(14)得，

(19)

代入公式(7)得，

(20)

由此，氢气自大空间内中心位置处喷放后，位于中心线上任意高度处的速度、氢气源运动特征参数及浓度均可由公式(17)、公式(18)和公式(20)进行数值求解近似获得，而处于氢气源运动特征参数范围内的偏离中心线位置任意高度处的参量可以通过公式(4)至公式(6)的高斯分布关系获得。

2 模型验证

将耦合Hirst关联式及空气阻力的理论模型计算结果与典型实验数据进行对比分析，以验证模型的可靠性。美国桑迪亚实验室的Schefer等人[15]开展了大量密闭空间内的氢气注射实验研究，获得了速度、浓度在空间内发展变化的数据，以研究动量、浮力、排放等单因素影响机制。计算选用数据如表1所示。

表1 氢气泄漏实验工况

图2展示了采用氢气迁移上升扩散二维模型对氢气迁移上升时气体速度进行预测，并同三种工况下的实验数据比对的结果，其中纵轴表示气体速度，横轴表示以喷口为0 m标高，沿中心线上升的各标高位置与喷口管径的比值。结果反映氢气迁移上升扩散模型能够准确反映Fr数较低的情况下的气体流动，在预测动量主导的气体流动时，存在一定的误差。整体而言并不因夹带系数及空气阻力的引入使气体在轴向方向的运动预测结果产生较大失真。

图2 中心线处气体速度随高度的变化曲线

图3展示了采用氢气迁移上升扩散二维模型对容器中心线处氢气质量份额进行预测，并同三种工况下的氢气质量份额实验数据进行比对后的结果，结果反映与实验结果均存在一定差异，特别是当z/d处于3～50的范围时，即位于容器中部及以下位置时，预测结果失真较大，在同一高度上模型预测结果较实验值明显偏小，这一失真现象在z/d越小时越严重。根据公式(19)可以分析得到这种差异来源于氢气源运动的特征参数B。可以通过对夹带系数及阻力系数修正，从而减缓B值增大，即减缓氢气的稀释过程，使氢气质量份额更接近真实情况。

图3 中心线处氢气质量份额随高度的变化曲线

3 夹带系数及阻力系数修正

由于Hirst关系式中并未对经验系数的适用范围进行明确描述，且Blassius方程多应用于管内流体阻力的计算，针对氢气/空气两组分气体在大空间内由动量或浮力主导的流动迁移工况，通过拟合系数的方法获得优化的夹带系数和阻力系数表达式，以适用于满足大范围Fr数工况下的氢气迁移上升扩散行为。基于工况1(Fr=99，Re=885)实验数据，通过对浮力主导氢气扩散工况进行数据拟合，获得的改进夹带系数关联式及改进阻力系数关联式分别为：

(20)

f=0.04Re0-1/4

(21)

将公式(20)、公式(21)其分别应用于动量主导的氢气扩散工况及处于二者共同作用的工况中，以验证其适用性。如图4至图5所示，与原有模型计算结果相比，采用修正系数的改进模型的预测值与中轴线各标高处的气体速度和氢气质量份额的实验值对比结果得到了明显改善，因此该改进后的氢气迁移上升扩散二维模型可以应用于大空间氢气浓度迁移分布的计算中。