核电安全壳内分离式热管蒸发段管束流体流动传热模拟①

， ，，2，

(1.南京工业大学 机械与动力工程学院，江苏 南京 210009；2.华东理工大学 机械与动力工程学院，上海 200237)

安全壳非能动冷却系统是核电非能动安全系统的重要组成部分，在事故发生时为机组提供最终热阱，保障安全壳的完整性。在正常工况下，安全壳内的温度低于50 ℃，压力接近大气压。当事故发生后，反应堆回路中高温、高压水泄漏进入安全壳中，在低压环境中迅速蒸发产生大量水蒸气，促使安全壳中压力、温度迅速增高，严重威胁安全壳的安全性和完整性[1-2]。

文中研究的安全壳非能动冷却系统(图1)，通过大型分离式热管导出安全壳中因事故产生的热量。分离式热管的蒸发段布置在安全壳上部空间，冷凝段布置在安全壳外的大型水池中，蒸发段与冷凝段采用上升管和下降管连接。分离式热管具有无需电力驱动、换热效率高的优点，可以有效地导出安全壳中的多余热量。

图1 基于分离式热管的非能动冷却系统示图

许多学者对传统分离式热管蒸发段进行了实验及模拟研究[3-5]。而大尺寸分离式热管的流体动力学性能和热性能与传统热管有很大区别，不一定适用传统热管实验得到的经验公式。为此，文中采用计算流体力学(CFD)方法模拟大型分离式热管蒸发段管束的基础流动特性，为核电安全壳分离式热管非能动冷却系统的开发提供参考。

1 分离式热管物理模型与数值计算方法

1.1 物理模型及假设

分离式热管蒸发段二维简化模型见图2，管束通过上、下管箱相互连接在一起。图2中下管箱、上管箱直径分别为D1、D2，每根换热管长度为L，换热管内径为d，换热管之间的管间距为t。蒸发段换热管呈倾斜或垂直布置。冷凝段产生的冷凝水从下管箱的入口进入蒸发段，通过壁面加热变成蒸汽从蒸发段的顶部流出到上升管中，最后回到冷凝段，并且一直循环下去。

图2 分离式热管蒸发段二维简化模型

设进口为速度进口，出口为压力出口。利用FIUENT软件对模型内沸腾传热传质流动过程进行数值模拟研究，基于以下假设：①冷凝水进口质量流量恒定，出口压力为定值。②假设流场为湍流，采用可实现k-ε湍流模型。③蒸发段壁温恒定，换热方式为对流换热，不考虑辐射换热及导热。④管束中每个换热管的外部传热系数相同。⑤换热管内气体为理想气体。

1.2 数值方程及计算方法

对管壁采用第三类边界条件，设定管外壁与外部流体的对流传热系数与周围流体温度。压力差值采用body force weight 格式，采用二维非稳态求解器，控制方程的离散采用有限单元体积法，采用SIMPLE速度-压力耦合，动量、能量方程采用二阶迎风算法，湍动能分量、耗散率采用一阶迎风算法，采用流体体积函数(VOF)模型模拟热管蒸发段管束流动沸腾现像[6-7]，采用连续表面力模型(CSF)，工质相变选用蒸发冷凝模型。

热管中的流动沸腾只有气液两相参与，VOF模型的控制方程如下。

连续性方程：

(1)

(2)

动量方程：

(3)

能量方程：

(4)

(5)

(6)

式(1)～式(6)中，SM为质量源项，kg/m3；SF为动量源项，N/m3；SE为能量源项，W/m3；ui、uj分别为xi、xj轴上的速度分量，m/s；δij为克罗内克函数；ρ为流体的密度，ρi为密度分量，kg/m3；E为流体微团的总能，J/kg；τij为黏性应力分量，Pa；qj为j方向的导热热流密度，J/(m2·s)；μ为动力黏度，Pa·s；∂k为各项的体积分数。下标l代表液相，下标v代表气相。

采用Brackbill CSF模型计算表面张力[8]。在CSF模型中，表面张力被视为通过界面的压力跳跃，并且可以表示为动量方程中的1个源项：

(7)

(8)

式(7)～式(8)中，Fvol为动量方程中的1个源项，N/m3；σ为气液相间的表面张力，N/m；k为界面曲率；α为体积分数。

在FLUENT中，通过在连续性方程和能量方程中加入源项来模拟相变。在蒸发和冷凝过程中，用Lee[9-10]提出的方程计算源项。当T>Tsat时，蒸发过程的质量源项和能量源项分别为：

(9)

(10)

(11)

当T>Tsat时，冷凝过程的质量源项和能量源项分别为：

(12)

(13)

(14)

式(9)～式(14)中，T为温度，Tsat为饱和温度，K；ΔH为蒸发焓，J/kg。

1.3 网格独立性验证

为了检查数值计算的准确性，在质量流量90 kg/h、对流传热系数600 W/(m2·K)、饱和温度373 K、外部流体温度395 K的工况下，选取4组不同的网格数对模型进行模拟，计算不同网格数下换热管的努塞尔数Nu，结果见表1。

表1 不同网格数下换热管努塞尔数变化

从表1可以看出，4组网格数下模拟得到的换热管努塞尔数的相对偏差范围在±3%以内。因此，上述数值模型可以模拟单根换热管中流动沸腾过程的一些基本机制。综合考虑模型的可靠性和计算速度，文中采用16 940个网格进行模拟研究。

1.4 模型验证

根据以往研究[11-13]，分离式热管蒸发段换热管中的饱和沸腾主要受核态沸腾影响。Rohsenow根据实验数据整理了核态沸腾的量纲一关系式[14-16]：

(15)

式中，cp为比定压热容，J/(kg·K)；ΔT为壁面过热度，K；r为汽化潜热，J/kg；Pr为普朗特数；s为经验指数，对于水s=1；Cwl为取决于加热表面和液体的经验常数；q为热流密度，J/(m2·s)；g为重力加速度，m/s2；下标l表示饱和液体，下标v表示饱和蒸汽。

为了保证算法的正确性，将Rohsenow根据实验数据整理的核态沸腾的量纲一关系式计算结果与数值模拟结果进行了对比，结果见图3。从图3可以看出，使用本算法计算的结果与实验结果有较高的一致性。

图3 核态沸腾量纲一实验关系式与本文模拟结果比较

2 换热管束流动传热数值模拟结果分析

2.1 管内两相流流型

选取L=2 m、D1=80 mm、D2=130 mm、d=80 mm、t=160 mm、进口管径80 mm、出口管径130 mm的管束模型进行模拟。在入口质量流量90 kg/h、对流传热系数600 W/(m2·K)、饱和温度373 K、外部流体温度395 K工况下，采用垂直布置时分离式热管蒸发段管束内气液两相的分布情况见图4。由图4可以看出，初始状态下管内液态体积分数为1，在前2 s，液态工质通过壁面加热发生相变产生气泡，气泡通过受热膨胀与其它气泡合并，并且不断长大，2 s后气泡基本达到稳定大小。在稳定情况下，从入口到出口，流型依次经历了泡状流、弹状流以及混搅流。蒸发段内的工质除了在进口处存在少量的泡状流外，主要以混搅流的形式存在。换热管中的两相流动行为与文献[17]中管内气液两相流动沸腾传热的过程一致，文中数值模型可以有效模拟管内流体的流动方式和传热性能。

图4 垂直布置时分离式热管蒸发段管束内气液两相分布情况

在入口质量流量90 kg/h、对流传热系数600 W/(m2·K)、饱和温度373 K、外部流体温度395 K工况下，采用90°、70°以及50°倾角布置时分离式热管蒸发段管束中气液两相在同一时刻的分布情况见图5。

从图5可以看出，管束中的流型在入口处为泡状流，随后发展为弹状流和混搅流，且泡状流和弹状流的存在范围较小，混搅流存在范围较大。在垂直管束中，下管箱中工质流型以泡状流为主，上管箱中工质流型以混搅流为主。随着倾斜角度的减小，管束内两相流的流动不对称性变得明显，且在上管壁处出现间歇性的局部干涸。特别是当倾斜角度减小到50°时，在管束上方工质流型呈现为雾状流，出现较大范围的干涸。

图5 不同倾角布置时分离式热管蒸发段管束内气液两相在同一时刻的分布情况

改变管束管间距，选取L=2 m，D1=80 mm，D2=130 mm，d=80 mm，进口管径80 mm，出口管径130 mm，管间距分别为140 mm、160 mm、180 mm的管束模型进行模拟。在质量流量90 kg/h、对流传热系数600 W/(m2·K)、饱和温度373 K、外部流体温度395 K工况下，不同管间距管束中气液两相的分布情况见图6。从图6中可以看出，对于不同管间距的管束，管内流型都是以混搅流为主。

图6 不同管束管间距下垂直布置时分离式热管蒸发段管束内气液两相分布情况

2.2 管内流场

倾角分别为90°、70°、50°布置的管束模型速度分布云图见图7，管间距分别为140 mm、160 mm、180 mm管束模型的速度分布云图见图8。从图7和图8的模拟结果得出，90°、70°以及50°倾角布置时管束中工质的平均速度分别为1.46 m/s、1.52 m/s、1.28 m/s，管间距分别为140 mm、160 mm、180 mm的管束模型平均速度分别为1.89 m/s、1.46 m/s、1.35 m/s。可以看出，在不同的倾斜角度以及不同的管间距下，管束中工质流速相差不大。管内流场的速度分布较为均匀，管子中心的流速要大于管壁附近的流速。在加热管中，蒸汽带动液态水向上流动，部分液态水回流并冲刷壁面。对比气液相分布图可以发现，上升的气泡速度要大于液态水的流速。

图7 不同倾角布置时分离式热管蒸发段管束内速度分布云图

图8 不同管间距下垂直布置时分离式热管蒸发段管束内速度分布云图

2.3 管内换热特性

事故发生时安全壳中温度会瞬间超过100 ℃，需要通过分离式热管散热来保证安全壳中的温度在设计值(148.89 ℃)以下[1]。为考察蒸发段管束的换热特性，采用努塞尔数Nu评估不同工作条件下的换热强度。努塞尔数是表示对流换热强烈程度的一个准数，又表示流体层流底层的导热阻力与对流传热阻力的比。

不同倾角布置时分离式热管蒸发段管束努塞尔数随时间变化情况见图9。在热管启动后努塞尔数快速增加又短暂下降，0.5 s后继续增加，1 s后基本达到稳定值。取管束在1～3 s的努塞尔数平均值，倾角为50°、70°、90°时对应的努塞尔数分别为359.67、378.91、396.16。倾角为50°～90°时，管束努塞尔数整体上随倾角的减小而减小。

图9 不同倾角布置时分离式热管蒸发段管束努塞尔数随时间变化情况

不同管间距下分离式热管蒸发段管束努塞尔数随时间变化情况见图10。取管束在1～3 s的努塞尔数的平均值，管间距为140 mm、160 mm、180 mm时对应的努塞尔数分别为383.61、396.16、411.01。可以看到，随着管间距的增大，管束努塞尔数也随之增大，但总体相差不大。

图10 不同管间距下分离式热管蒸发段管束努塞尔数随时间变化情况

3 结语

提出了一种基于VOF的安全壳中分离式热管蒸发段管束沸腾流动模型，该模型可以有效模拟流型的演化，获得良好的流体动力学和热力学模拟结果。数值模拟结果：①大型分离式热管蒸发段管束中的工质流型经历了泡状流、弹状流以及混搅流。50°倾角布置的管束上端出现了雾状流。②在换热管中，蒸汽带动液态水向上流动，部分液态水回流并冲刷壁面。管内流场的速度分布较为均匀，换热管中心的流速要大于管壁附近的流速。③在相同的换热条件下，倾角为50°～90°时，随着倾角的增大，努塞尔数整体上随着倾角的增大而增大，倾角越大对应的换热效果越好。在管间距取1.75～2.25倍的管径时，随着管间距的增大，努塞尔数随之增大，换热效果也越好，但总体相差不大。