杨宇鹏，王成龙，张大林，苏光辉，田文喜，秋穗正

(西安交通大学 核科学与技术学院，陕西 西安 710049)

蒸汽发生器是核反应堆动力系统的重要组成部分，是将反应堆一次侧的热量向二次侧传递的关键枢纽。相较于直管式蒸汽发生器，螺旋管式直流蒸汽发生器(HCOTSG)有一些突出优势[1]：1) HCOTSG的特殊结构，使得其相比较直管式蒸汽发生器，在体积相同的情况换热效率大幅提升，同样热负载相同的条件下，其将使得动力系统布置更加紧凑，也更符合核能未来发展的小型化趋势；2) 螺旋管侧与壳侧流体的湍流程度提升，起到强化换热的效果；3) 螺旋管具有自补偿特性，对热膨胀具有良好的调节能力。但相比之下螺旋管式蒸汽发生器存在流动压降增加、设备复杂性增加等缺点。因此，针对HCOTSG的研究中，应综合其传热与流动性能进行全面分析。目前已有较多针对螺旋管内流动过程的实验研究基础，对其内部的流动换热特性以及特殊现象进行了探究[2-5]，对于其壳侧与管侧的耦合流动换热过程通常利用商用软件以及自编程序开展数值模拟计算研究[6-9]。但针对液态金属外掠螺旋管束的流动以及液态金属HCOTSG壳侧与管侧的耦合流动传热过程的研究相对较少。

本文提出利用计算流体力学软件(FLUENT)对液态金属HCOTSG壳侧流体、管侧流体以及管壁固体区域耦合流动传热过程进行数值模拟的方法，对其流动传热特性进行分析。

1 液态金属HCOTSG概述

HCOTSG是一种蒸汽发生器常用形式，其被广泛应用于船用压水堆、小型模块化反应堆[10]，高温气冷堆以及铅铋快堆[11]等。铅基材料(铅、铅铋或铅锂合金等)具有中子性能优良、热工性能优良、安全性好等优点，采用其作为反应堆冷却剂，能使反应堆物理特性和安全运行具有显著优势[12]。在液态金属反应堆中，HCOTSG也得到较为广泛应用。世界第1座液态金属反应堆——法国凤凰堆使用750 MWe螺旋管式蒸汽发生器，日本的280 MWe钠冷快堆(SFR)、徳国的钠冷快堆(SNR-300)以及美国的(ALMR)均采用了HCOTSG。

本文研究对象为铅铋-水HCOTSG，其由分层布置的螺旋管束与同心圆柱套筒组成。圆柱套筒与螺旋管外壁之间的外掠管束区域为一次侧流体域(即壳侧区域)，螺旋管内壁所围区域为二次侧流体域(即管侧区域)，其示意图如图1a所示。螺旋管内为过冷水到过热蒸汽的两相流蒸发过程，壳侧为液态铅铋外掠螺旋管束冷却过程。HCOTSG简图如图1b所示，螺旋管束分三层均匀布置于圆柱套筒中，螺旋方向分别为左旋、右旋、左旋。为保证换热均匀，各层螺旋管长度基本保持一致。各层螺旋管的几何参数列于表1。对此几何结构的液态金属HCOTSG内部的流动传热过程开展了数值模拟计算，并对结果进行了分析。

图1 HCOTSG示意图Fig.1 Schematic diagram of HCOTSG

表1 HCOTSG几何参数Table 1 Gometrical parameter of HCOTSG

2 数值模拟方法

2.1 控制方程与计算方法

1) 控制方程

针对螺旋管内的水蒸发过程，其涉及到气液两相流动以及两相间的传热与界面作用过程。本文选用适合于界面追踪的VOF(volume of fluid)多相流模型对气液两相流过程进行数值模拟，连续性方程、动量方程和能量方程如下。

体积分数方程(连续性方程)，对第q相，有：

(1)

(2)

式中：α为体积分数；S为源项；ρ为密度；v为速度；t为时间。

动量方程：

(3)

式中：p为压力；F为源项；μ为黏性系数。

物性参数由各相体积分数决定，如对于密度有：

ρ=∑αqρq

(4)

能量方程：

(5)

(6)

式中：keff为有效热导率；当i取1～3时，ui代表x、y和z方向的速度。

当涉及到固体区域时，其能量方程为：

(7)

式中：ks为导热系数；T为温度。

2) 湍流模型

k-ε模型拥有可靠、收敛性好、内存需求低等优点。在FLUENT中的不同湍流模型中，Realizablek-ε模型对于二次流型是针对低雷诺数以及射流模拟提出的修正模型。考虑到螺旋管内的流动为湍流工况，且伴有二次流的发生，因此选择此模型为数值模拟中采用的湍流模型。考虑到本文的工况，选用的模型为Realizablek-ε湍流模型，其假设流动的情况为完全湍流，从以下的输运方程中可获得湍动能k和湍流耗散率ε：

Gk+Gb-ρε-YM+Sk

(8)

(9)

(10)

式中：Gk为速度梯度引起的湍动能增量；Gb为浮力引起的湍动能增量；Ym为脉动膨胀引起的湍流耗散率增量；C2和C1ε为常数；σk和σε分别为k、ε的湍流普朗特数；Sk和Sε为源项。

2.2 液态金属流动传热模型

相比于常规流体，液态金属的热物理性质存在较大差异，如液态金属的导热能力远大于常规流体。因此，对于液态金属，其热扩散能力远大于动量扩散能力，导致其普朗特数(Pr≈0.01～0.06)与常规流体(Pr≈1)相差巨大。因此，对于液态金属流动的数值模拟，湍流普朗特数(Prt)的影响较为显著。因此在采用CFD方法对液态金属流动传热进行数值模拟时通常需对Prt模型进行修正。本文选择了Cheng等[13]提出的液态金属Prt模型，利用用户自定义文件(UDF)导入FLUENT软件对湍流模型进行修正，公式如下：

(11)

式中：Pe为佩克莱数；A为常数。

2.3 数值模拟方法验证

由于对液态金属HCOTSG的实验数据相对较少，因此对其壳侧液态金属流动、管侧水蒸发两相流动分别进行模型验证，以证明数值模拟方法的正确性。

1) 螺旋管内水蒸发模拟验证

对于管内的水沸腾两相流过程，首先选择Bartolomei等[14]于1967年设计搭建的竖直加热圆管内流动沸腾实验进行了对比验证，圆管长度为2 m，直径为15.4 mm，壁面边界为恒定热流边界，所得数值模拟结果如图2所示。得到的沿轴向高度的空泡份额与实验值误差在25%以内，温度分布结果与实验值误差在3%以内，均符合较好。

在此基础上，将模拟方法应用于螺旋管内水沸腾过程模拟。选取Santini等[15]所做实验的结果与数值模拟结果进行对比。螺旋管内径为12.53 mm，螺旋直径为1 000 mm，节距为800 mm，长为32 mm的螺旋管竖直布置。所得数值模拟结果如图3所示。由不同工况下的换热系数对比可看出，所得结果与实验值符合较好，大部分工况误差在20%以内，因此螺旋管内流动的数值模拟方法正确性得到验证。

图2 Bartolomei实验值与计算值对比Fig.2 Comparison between experimental result of Bartolomei and calculated result

图3 Santini实验值和计算值对比Fig.3 Comparison between experimental result of Santini and calculated result

2) 螺旋管外液态铅铋流动模拟验证

壳侧的铅铋流动过程为液态金属外掠螺旋管束流动，缺乏相关的实验数据。因此，首先选择西安交通大学圆管内的铅铋流动换热实验数据进了对比验证，圆管直径为20 mm，长度为2 m。图4为西安交通大学(XJTU)圆管内铅铋流动实验值和计算值对比。由图4可知，在试验工况下，圆管平均Nu的数值模拟结果与实验值符合较好，大部分工况误差在15%以内。验证了铅铋流动数值模拟方法的准确性。

图4 XJTU圆管内铅铋流动实验值和计算值对比Fig.4 Comparison between experimental result of XJTU and calculated result

在验证液态铅铋模拟方法的基础上，考虑到螺旋管式蒸汽发生器壳侧区域特殊的几何结构，进一步建立更为相近的液态金属横掠管束模型，并利用已有的经验关系式对数值模拟方法进行验证。在建立了液态金属钠横掠管束模型后对其开展数值计算，选择液态金属钠横掠管束的Kalish-Dwyer关系式以及Schad关系式与数值模拟结果进行验证，结果如图5所示，可看出，数值模拟结果与关系式计算结果符合较好，大部分工况误差在15%以内，进一步验证了壳侧液态金属流体域数值模拟方法的准确性。

图5 关系式计算值和模拟值对比Fig.5 Comparison between calculated result by empirical correlation and simultion result

2.4 模型建立与网格无关性分析

1)模型建立

依据表1中的几何参数，利用SOLIDWORKS软件对螺旋管束几何结构进行建模。螺旋管共布置3层，共9根，壳侧区域为套筒与螺旋管束外壁形成的外掠管束区域。模型包含壳侧流体域、管侧流体域以及管壁固体域三部分。在几何模型建立完成后进行网格划分，管侧流体与管壁区域采用六面体网格，壳侧流体域采用四面体网格，部分区域的网格示意图如图6所示。网格划分完成后将建立模型导入FLUENT软件进行数值模拟计算。在进行数值模拟时，需将液态铅铋的物性以及选用的Prt模型利用UDF导入软件并进行相应参数的选取。入口边界采用速度入口边界条件，出口边界采用压力出口边界条件，内外套筒壁面均为绝热边界条件。综合考虑管侧与壳侧流动过程，湍流模型选择Realizablek-ε模型。在监控计算残差的同时对出口温度、含汽率等参数进行监控以保证计算结果的准确性。

图6 部分区域网格示意图Fig.6 Schematic diagram of partial area grid

2) 网格无关性分析

由于本文的数值模拟包含HCOTSG壳侧液态铅铋流动与管侧水蒸发两相流动过程，两侧流动传热过程的差异较大，网格划分方式对数值模拟结果产生的影响也差异较大，因此对两侧的流体流动传热过程分别进行网格无关性分析。对于螺旋管内的两相流过程，本文选取了沿轴向高度不同位置的截面平均压降作为网格无关性的判定依据，所得结果如图7所示。对螺旋管周向节点数量与轴向节点数量分别进行加密处理，发现当周向网格数量为30、轴向网格数量为700时，加密网格对轴向沿程压降的影响较小。因此选择此时的网格划分方式对螺旋管内流体域及管壁固体域进行网格划分。

图7 管侧网格无关性分析Fig.7 Grid independence analysis at tube side

对于HCOTSG壳侧的液态铅铋流动过程，本文选取沿轴向高度的截面平均温度作为网格无关性分析的判定依据，所得结果如图8所示。首先对网格进行加密，发现当网格数量达到580万左右时，继续加密网格对壳侧液态铅铋的沿程温度影响较小。因此选择此时的全局尺寸划分网格。在确定全局网格尺寸的基础上，选取了不同的边界层网格划分方式进行网格划分，当壁面边界y+缩小至34时，继续加密边界层网格轴向温度的影响较小。因此选择边界层网格处理方式对边界层网格进行划分。

3 流动换热特性分析

在方法正确性验证的基础上，选取典型工况的计算结果进行分析，工况参数列于表2。图9示出壳侧区域中截面速度分布。由图9可看出，在不同层螺旋管之间存在明显的流体加速区域，而同层螺旋管的轴向节距之间存在明显的低速区域，这使得液态金属域内的流体搅混明显提升，强化了流体与管壁的传热。且由于各层螺旋管螺旋方向布置方式不同，使得在截面对称区域的速度分布并非对称，因此在后续的敏感性分析中，螺旋方向的布置情况也应作为考虑的敏感性因素。

图8 壳侧网格无关性分析Fig.8 Grid independence analysis at shell side

表2 工况参数Table 2 Working condition

图9 壳侧区域中截面速度分布Fig.9 Velocity distribution of middle section in shell side area

在壳侧液态铅铋域内，其沿高度方向的截面平均温度与压降分布如图10所示。进出口的直管段区域，管束对壳侧流体扰动较小，因此换热器内流体的换热效果明显弱于螺旋管束区域。相应地，螺旋管束区域的压降也明显高于直管区域。在利用HCOTSG的优点时，应同时关注其特殊结构带来的压降损失。

对于管内的水蒸发两相流动，各层内的螺旋管沿程温度与热流密度分布如图11所示。可看出，各层螺旋管在同一轴向高度下的温度与壁面热流密度分布一致，因此在此几何布置方式下，螺旋管束传热均匀性较好。在入口段与蒸发段，较大的温差以及水蒸发的汽化潜热使得壁面的热流密度维持较高的水平。而在水完全蒸发后，螺旋管内的换热过程为水蒸气的过热过程，此区域内热流密度处于较低水平。因此，为更好发挥蒸汽发生器的换热能力，针对此结构的HCOTSG，应在分析工况参数影响的同时选择出合适的工况范围。

图10 壳侧区域沿程温度和压降分布Fig.10 Distribution of temperature and pressure drop along shell side area

图11 管侧温度与热流密度分布Fig.11 Distribution of temperature and wall heat flux along tube side

图12 螺旋管内温度与压降沿程分布Fig.12 Temperature and pressure drop distribution in helical coil tube

螺旋管内温度与压降沿程分布如图12所示，由于入口段和出口段与接管处存在较大的弯曲角度，因此在这两处压降损失较大，螺旋区域内压降变化较为均匀。沿程温度在入口段快速升高，上升趋势逐渐平缓。螺旋管内部流体除了受到黏性力作用，还受到惯性力作用。最内层螺旋管内的沿程截面液相分布变化如图13所示。在螺旋管内，液态水所受离心力较大，因此其位于螺旋管外侧；而在部分区域，在螺旋管内会出现特殊的液膜倒置现象，即液态水集中在螺旋管内侧的现象，该现象与已有的螺旋管内两相流动实验情况一致。

图13 管侧区域截面液相体积份额Fig.13 Liquid volume fraction in tube side area

4 总结

本文提出了利用计算流体力学对液态金属HCOTSG壳侧流体、管侧流体耦合流动传热过程进行数值模拟的方法，分析了其流动传热特性，得出以下结论。

1) 本文提出的数值模拟方法可对液态金属螺旋管式蒸汽发生器内部的流动换热情况进行模拟，并将所得计算结果与实验结果和经验关系式结果对比：对于管内的水蒸发两相流模拟，通过与Bartolomei的直管沸腾实验以及Santini的螺旋管内沸腾实验进行对比，误差均在25%以内；对于壳侧液态金属模拟，通过与西安交通大学圆管铅铋流动换热实验以及Kalish-Dwyer关系式和Schad关系式进行对比，误差均在15%以内，数值模拟方法的正确性得到验证。

2) 在本文提出的螺旋管式直流蒸汽发生器几何结构下，针对典型工况的流动换热情况进行了分析。壳侧铅铋域存在明显的加速区域，且由于各层螺旋管的布置方式存在差异，因此在各层间的流体速度也存在差异，增加了流体域内的搅混。对于螺旋管内部的两相流动过程，过冷水被加热最终变为过热蒸汽，在过热蒸汽区的热流密度明显下降，因此应针对不同几何结构，确定合适的工况范围。同时发现管内存在液膜倒置现象，与已有实验所得结果一致。