自然循环窄矩形通道内过冷沸腾两相摩擦阻力特性

杨宽，阎昌琪，曹夏昕

（哈尔滨工程大学核安全与仿真技术国防重点学科实验室，黑龙江哈尔滨150001）

引 言

相较于传统的流道，窄矩形通道在流动沸腾条件下单位体积内可以获得更高的传热性能[1]，因而具有窄矩形几何结构特征的流道广泛应用于紧凑式的工业热交换器以及一体化反应堆、实验堆堆芯中[2-3]。采用窄矩形几何结构的板式工业热交换器以及使用板状燃料组件的反应堆堆芯往往具有较小的设备体积，可以适应浮动式核电站以及船用核动力装置的需求。

为了提高船舶核动力装置的安全性，非能动安全技术广泛应用，自然循环仅依靠回路内冷热段流体的密度差即可驱动回路内流体循环流动，可提高核动力装置的固有安全性。同时，这些设备换热面往往具有较大的热负荷，过冷沸腾现象广泛出现在这些设备的流道当中，因此准确计算自然循环条件下窄矩形流道内过冷沸腾两相摩擦阻力对于这些设备的设计和安全评估具有十分重要的意义。

王广飞等[4-6]对多种不同高宽比的窄矩形通道内的常压下空气-水两相流动的摩擦阻力进行了实验研究，对典型的两相摩擦压降计算关系式进行了评价，并根据其实验结果提出了修正公式。秦胜杰等[7]对高压条件下的小宽高比矩形通道内饱和沸腾蒸汽-水两相摩擦阻力进行了实验研究，引入反映小通道对气泡生长限制的无量纲特征参数Nconf，对ChisholmC模型[8]中的参数C值进行修正以进行矩形通道阻力特性预测。刘传成等[9]对常压下窄矩形通道内入口过冷出口达到饱和状态的蒸汽-水两相摩擦阻力特性进行了实验研究，发现Chisholm 系数C与Martinelli 参数X存在指数关系，并且随着全液相Reynolds 数而变化，据此给出了新的分液相摩擦因子的计算方法；在其两相摩擦压降剥离过程中忽略了过冷沸腾段的影响。孙奇等[10]在中压强迫循环工况下开展了ϕ32 mm×3 mm 竖直圆管内的低流速过冷沸腾压降实验，并且同时测量了过冷沸腾空泡率，结合其提出的过冷沸腾真实含气率分布模型[11]剥离出了两相摩擦压降，实验发现低流速过冷沸腾条件下实验段内截面空泡份额不大，总压降中重力压降占总压降份额较大，其值远远大于加速度压降和摩擦压降，同时推荐采用MN 因子计算式或Friedel公式进行摩擦压降计算。颜建国等[12]对高热流条件下内径6 mm 的圆管通道内过冷沸腾流动阻力特性进行了实验研究，并分析了质量流速、热通量、压力、沸腾数、Jakob数等系统参数对流动阻力的影响，提出了一个考虑了管径因素修正项的经验关联式，该关联式的预测误差在±18%范围内；其实验研究中获得的两相流动阻力并未将两相加速压降项和摩擦压降项分离，而是统一处理。赵楠等[13]对常压下宽度分别为3 mm和4 mm的竖直窄缝通道内的过冷沸腾特性进行了可视化实验研究，实验发现随着热通量的增加，进出口压降出现振荡现象，热通量越大振荡幅度越大，气泡的产生、脱离、聚合、湮灭等行为的加剧使得流动的不稳定性产生；窄缝宽度因素对实验段进出口压降有明显的影响，4 mm的振荡幅度要大于3 mm的振荡幅度，其研究中并未对两相摩擦压降进行分离研究。

综上所述，学者们已经对窄矩形通道、常规通道内两相摩擦阻力进行了大量实验研究，实验工质包括空气-水、蒸汽-水以及制冷剂，流动形式包括绝热两相流、过冷流动沸腾和饱和流动沸腾，而对低压低流速自然循环工况下窄矩形通道内的过冷沸腾摩擦阻力特性的研究相对较少。因此，本工作采用去离子水作为实验工质，在低压低流速自然循环工况下开展了单面加热可视化窄矩形通道内的过冷沸腾摩擦阻力特性实验研究，提出了适用于本实验条件的两相摩擦压降的剥离计算方法，对已有的两相摩擦压降计算经验关系式进行适用性评估，并在此基础上提出修正计算模型。

1 实验装置和参数

实验装置的结构如图1 所示，其位于一个液压机械摇摆平台上，可进行倾斜以及简谐正弦摇摆运动，本实验中将其调节至水平静止状态。实验装置由主循环回路和辅助冷却回路组成，辅助冷却回路为主回路冷凝器提供二次侧冷却水。

图1 实验装置流程Fig.1 Schematic diagram of experimental apparatus

通过切换主回路阀门开合，主循环回路可以运行在强迫循环和自然循环工况下，将主循环回路调节至自然循环模式运行，去离子水经过预热器加热至指定过冷度后流入窄矩形通道实验本体中被再次加热产生过冷沸腾。

单面加热可视化实验段本体的结构可见文献[14]，窄矩形流道由3 mm 厚度的316 L 不锈钢板和蚀刻在石英玻璃上的一个矩形槽道装配形成，其横截面名义尺寸为71 mm×2 mm，当量直径De为3.9 mm，实验板的加热段长度0.55 m。实验段加热板上沿着流动方向开有四个引压孔，一共采用三个压差传感器测量沿程压降，从下到上三个测压段的长度分别为0.065、0.2、0.215 m。

如图2 所示，在实验段本体支撑架上安装了一个二维导轨，高速摄影仪安装在导轨平台上，可拍摄沿着流动方向上窄矩形通道内不同位置处不同区域范围的气液两相图像。

实验过程中通过电磁流量计测量回路循环流量，两个N型铠装热电偶测量实验段进出口水温，实验过程中使用的流量计和压差传感器的量程等参数如表1所示。

表1 测量仪表参数Table 1 Parameters of measurement instruments

图2 实验段本体结构Fig.2 Schematic diagram of test section

2 实验数据处理

2.1 单相实验验证

为了验证实验测量装置的可靠性，同时为两相摩擦阻力压降的计算提供基础，首先开展了冷态等温强迫循环条件下，在矩形通道内的单相流动阻力特性实验。由于矩形通道内绝热层流流动条件下速度分布存在级数形式解析解，根据Shah 等[15]的研究结果，可在充分发展的层流条件下，高宽比ε在0～1 范围内，将矩形通道内的单相摩擦阻力系数计算关系式简化为与通道高宽比相关的指数形式关系式，其和理论解的最大相对误差在0.05%范围内，该计算关系式为

式中，λ为单相摩阻系数；Re为矩形通道内的Reynolds数；ε为矩形通道截面的高宽比。

如图3 所示，将窄矩形实验段在冷态条件下获得的沿程摩阻系数λ-Re特性曲线在双对数坐标系上绘制，图中λlaminar为采用Shah & London 关系式计算的结果，λBlasius为采用Blasius 关系式的计算结果。Blasius 关系式假设流道壁面到中心的速度为1/7 指数律分布，并且结合部分实验数据拟合得来，最佳适用范围为充分发展湍流区，Re<10000 区域范围。可见，λ-Re曲线出现了典型的层流-湍流转捩现象，层流-过渡区转捩点在Reynolds 数2700 附近，过渡区-湍流转捩点Reynolds 数在3340 附近。Shah &London 关系式和Blaxius 关系式分别能很好地预测层流区和湍流区内的阻力特性，其平均相对误差分别为1.13%和2.04%，可见本实验中的测量装置有较好的精度。对于过渡区内的阻力特性，采用层流区-过渡区转捩点和过渡区-湍流区转捩点两点之间的值进行插值计算，与实验值对比，平均相对误差为2.859%。

2.2 两相摩擦压降计算

在两相实验工况中，入口水温处于过冷状态，实验段同时存在单相段和两相段，为了确定两相段摩擦压降，需要首先确定测压段内两相段长度。实验过程中通过高速摄影仪拍摄沿着流动方向上不同位置处的图像信息，确定过冷沸腾起始点的位置。

图3 窄矩形通道在冷态等温条件下的单相阻力特性曲线Fig.3 Single-phase friction coefficient characteristic in isothermal condition

如图4所示，本实验工况中，核化点首先出现在通道中心处，实验中认为壁面中心区域内首先出现较多的核化点，并且能够持续产生气泡的位置为过冷沸腾起始点，然后记录过冷沸腾起始点位置。图4 中出现的抛物线表示过冷沸腾起始点后壁面上核化点数目的增加和向加热板两侧边缘位置的发展。

实验过程中通过调节系统压力、实验段入口过冷度和实验段加热功率改变工况，图5 为四种不同热工条件下，高速摄影仪拍摄的加热通道中心部分区域的气泡图像。实验中观察到，实验段内的流型主要为泡状流，观察到的气泡大多为在加热壁面上滑移的气泡，其尺寸与通道窄边高度一个数量级。气泡从核化点产生后，当局部过冷度和壁面热通量满足气泡的生存条件，气泡会脱离核化点，沿着流动方向在加热壁面上滑移运动，过程中气泡直径和运动速度不断增大。如图5(a)～(c)所示，随着系统压力减小，入口过冷度的减小和加热功率的增大，通道内的气泡群直径和密度也随之增大，气泡直径和气泡群密度增大到一定程度，聚合过程大量发生，如图5(d)所示，此时大量气泡聚合形成气块，改变通道内的阻力特性和驱动力大小，进而引发流动不稳定，此时测压段内压差产生强烈振荡，因此本实验中获得数据均来自流量较为稳定的过冷沸腾泡状流条件下的实验工况。

图4 过冷沸腾起始点处的气泡图像Fig.4 Bubble images at point of nucleate boiling

图5 不同热工工况下的气泡图像Fig.5 Bubble images in varies set of working condition

实验过程中，在竖直静止条件下，通过压差传感器可获得窄矩形通道测压段内的总压降，而总压降ΔPt由重位压降ΔPg,tp、加速压降ΔPa,tp和摩擦压降ΔPf,tp三部分组成。重位压降和加速压降分别可按式（2）～式（4）计算。

式中，G为通道内的质量流速，kg/(m2·s)；α为局部的空泡份额；ρl和ρg分别为液相和气相密度,kg/m3；xi和xo分别为测压段进出口的质量含气率；αi和αo分别为测压段进出口的空泡份额。

为了剥离出实验段内的两相摩擦压降，需要获得实验段内的空泡份额以及真实质量含气率。由于本实验工况大多处于过冷度沸腾区域，流型为典型的泡状流，此时通道内的质量含气率和空泡份额相对较低，此时认为通道内空泡份额从过冷沸腾起始点开始到净蒸汽产生点之间沿着流动方向近似服从线性分布，该假设和典型加热圆管通道内的流动区域中划分的深度欠热区特性相类似。在圆管内，深度欠热区内的气泡附着在加热壁面上，不能在主流中生存，阻力特性通常作为表面粗糙度的增加处理，而在本实验中观察到，气泡从加热壁面产生后很快就脱离核化点，在加热壁面上滑移过程中不断长大，气泡尺寸也与通道窄边宽度处于一个数量级，气泡行为会显著影响通道的流动特性。

本实验中，采用了孙奇等[16]提出的适用于计算低质量流速条件下的净蒸汽产生点模型来获得净蒸汽产生点位置，该模型通过考虑静蒸汽产生点处气泡的热力学平衡和水力学平衡推导而出。

式中，Xeq,NVG为净蒸汽产生点处的热平衡含气率；q″为壁面热通量，kW/m2;cp为液相的比定压热容，kJ/(kg·k);ilg为汽化潜热,kJ/kg;Prl为液相Prandtl数;σ为表面张力系数,N/m。

在确定了净蒸汽点的位置后，为获得实验段内的空泡份额，还需要获得净蒸汽点处的空泡份额值大小，本实验采用了Rouhani[17]提出的公式计算净蒸汽产生点处的空泡份额，该模型适用于较大压力范围内水介质过冷沸腾条件。

式中，Psys为系统压力，MPa；Ph为通道横截面加热周长，m;Ac为通道横截面面积,m2。

实验段内从过冷沸腾起始点开始各处的空泡份额通过线性插值获得，即通道内从过冷沸腾起始点开始距离x处的空泡份额αx为

摩擦压降的剥离还需要获得通道内的真实质量含气率，而空泡份额通过以上计算已知，所以本实验中该数据通过空泡份额模型反推获得，本实验中采用典型的漂移流模型Chexal-Lellouche 模型[18]计算。根据Coddington 等[19]对空泡份额计算的漂移流模型进行的系统评价研究表明，该模型在非常宽广的工况范围内预测准确度最高，被核反应堆工程中广泛使用的热工水力计算程序RELAP5 和RETRAN所采用。

在已知通道内真实含气率和空泡份额后，可以根据式(3)和式(4)计算出测压段内的加速压降和重位压降，进而可以求得两相摩擦压降的实验值。

3 实验结果和分析

表2给出了竖直两相自然循环实验中的热工参数范围，在该实验范围内通过第2 节给出的计算方法剥离出压降的实验值，并且将实验值和采用典型的均相流和分相流摩擦阻力计算模型的计算结果进行了对比，对其适用性进行了评价。

表2 实验参数范围Table 2 Range of experimental parameters

评价标准为实验值与计算值之间的平均相对误差(MRE)和相对误差在30%以内的数据比例Ω，其中

3.1 计算模型评价

图6给出了采用不同等效黏度计算方法的均相流模型计算结果和实验值的比较，表3 给出了各种均相流模型的平均相对误差MRE 以及相对误差在30%以内的数据比例Ω。由图6 可知，均相流模型计算的结果相较于实验值整体都明显偏小，而且互相之间的差别很小，MRE 均在18%左右。这是由于本实验工况中流型处于泡状流，通道内的真实含气率虽然相对较小，但通道内存在的气泡尺寸与通道窄边宽度相当，其数量密度也相对较大，气泡群对流动的影响较大，而且此时气液两相之间的滑移现象也较为显著，而均相流模型假设中未考虑气液之间相互作用的影响，特别在低质量含气率区域，其等效黏度计算结果未能反映出气相的影响，其计算结果的等效黏度和单向流体的黏度相差不大，所以导致其计算结果整体明显偏小。

表3 均相流模型计算结果和实验值比较Table 3 Comparison of calculation results of homogeneous models with experimental data

图6 均相流模型的计算结果和实验值的比较Fig.6 Comparison of calculation results of homogeneous models with experimental data

图7为分型流模型中采用分液相折算系数方法计算的两相摩擦阻力值和实验值的对比，表4 为这些模型的平均相对误差MRE 以及相对误差在30%以内的数据比例Ω。由图7可知，Sun and Mishiba 模型和Zhang and Mishiba 模型能较好地预测实验段内的两相摩擦阻力值，其MRE 分别为12.21% 和15.38%，且其Ω分别为100%和85.71%，这两个模型中均考虑了在窄矩形通道中较为显著的表面张力的影响，特别是Sun and Mishiba模型，还考虑了分液相或者分气相流态区域的影响，因此其预测结果相对较好，但有轻微高估通道内两相摩擦压降的趋势。ChisholmC模型是基于常规通道内实验数据提出的典型分液相模型，其在窄矩形通道内的预测结果稍差。王广飞等[4]通过其窄矩形通道内空气-水两相流实验数据，基于ChisholmC模型，根据分液相Reynolds 数和分气相Reynolds 数之比30 进行分区拟合修正，其实验工况流动参数范围和本实验不相符，其适用于较大液相和气相质量流量参数范围内计算，故误差较大。Mishiba and Hibiki 模型根据其1～4 mm 小直径圆管内的空气-水两相流数据对ChisholmC模型中的系数C进行了修正，不适用于本实验中的窄矩形通道，而Lee and Lee 模型根据其水平窄矩形通道内空气-水两相流动实验数据，对ChisholmC模型中的系数C也进行了修正，故计算结果和本实验中竖直条件下蒸汽-水的实验数据结果相差较大。

表4 分液相折算系数方法的计算结果和实验值比较Table 4 Comparison of calculation results of ϕ2l based models with experimental data

图8为分型流模型中采用全液相折算系数方法计算的两相摩擦阻力值和实验值的对比，表5 为这些模型的平均相对误差MRE 以及相对误差在30%以内的数据比例Ω。可见，Tran 模型均能较好地预测实验中的两相摩擦阻力值，大部分实验点均落入±30%的误差线内，平均相对误差12.28%；而ChisholmB和Muller-Steinhagen and Heck 模 型 虽 然大部分数据点都落在30%的误差限内且平均相对误差相对较小，但其都明显低估了通道内的摩擦压降。

图7 分液相折算系数方法的计算结果和实验值的比较Fig.7 Comparison of calculation results ofbased models with experimental data

表5 全液相折算系数方法计算结果和实验值比较Table 5 Comparison of calculation results of ϕ2lo based models with experimental data

3.2 模型的提出与验证

基于本实验中获得的实验数据，基于分相流模型中的分液相折算方法，考虑到影响在窄矩形通道内两相摩擦压降的关键参数，即反映通道内总质量流速影响的全液相Reynolds 数Relo、气液各相阻力占比的Martinelli 参数X和反映窄通道尺寸影响的La的作用，采用以下关系式计算分液相折算系数

图8 全液相折算系数方法的计算结果和实验值的比较Fig.8 Comparison of calculation results ofbased models with experimental data

其中Laplace数La为

图9 显示了拟合关系式和本实验数据的对比，所有数据点均落在正负20%的误差限范围内，平均相对误差MRE为8.54%，与实验数据符合较好。

为验证本文计算方法以及获得的两相摩擦压降关系式的可信度，对过冷沸腾条件下，对窄通道内的阻力特性进行了文献调研，提取Baburajan 等[31]的数据，其实验工况为：系统压力0.12 MPa，质量流速450、600、690 kg/(m2·s)，热通量225～445 kW/m2,入口温度28℃,当量直径5.5 mm,加热段长度0.55 m,与本文实验工况较为接近。采用本文中使用的计算方法对该实验数据进行处理，并且使用本文获得的过冷沸腾两相摩擦阻力关系式(13)计算摩擦压降，进而获得总压降，将其与Baburajan 等[31]的总压降计算关系式结果进行对比，结果如图10 所示，可见所有数据点均在正负20%误差限内，平均相对误差MRE 为6.45%，可见本文方法可以较好地预测过冷沸腾条件下的两相摩擦阻力，具有较宽的适用范围。因此，在过冷沸腾条件下，本文获得的两相摩擦压降计算关系式(13)的适用系统参数范围可拓展为：系统压力0.1～0.3 MPa、质量流速0～700 kg/(m2·s)、热 通 量80～450 kW/m2、入 口 过 冷 度20～70℃。

图9 拟合关系式的计算结果和实验值的比较Fig.9 Comparison of calculation results of new correlation with experimental data

图10 本文计算方法与文献[31]实验数据的对比Fig.10 Comparison of present calculation methods with experimental data from Ref.[31]

4 结 论

(1)本实验通过可视化观察，发现实验工况中两相流型主要为泡状流，通道内存在的气泡尺寸与通道窄边宽度相当，其数量密度也相对较大，大部分在加热壁面沿着流动方向滑移，气泡群行为会对主流流动产生明显的扰动，明显增大了通道内两相摩擦阻力。

(2)基于实验获得的两相摩擦阻力实验数据，对经典的两相摩擦阻力计算模型进行了评价。实验发现，采用均相流模型计算的结果相较于实验值整体都明显偏小，而且互相之间的差别很小，MRE 均在18% 左右。分相流模型中Sun and Mishiba 模型和Tran 模型能较好地预测窄矩形通道实验段内的两相摩擦阻力特性，其MRE 分别为12.21%和12.28%，并且大多数实验点均落入正负30 的误差限内。

(3)针对本文数据，考虑全液相Reynolds 数Relo、Martinelli 参 数X和Laplace 数La等 重 要 参 数 的 影响，给出了分液相折算系数的计算关系式，能较好地预测过冷沸腾下窄通道内的两相摩擦压降，平均相对误差为8.54%，并且通过其他文献中的实验数据对本文中采用的计算方法进行了对比验证，表明修计算关系式具有较高的精度和较宽的适用范围。

符 号 说 明

De——当量直径，m

L——长度，m

ΔPa,tp——两相段加速压降，Pa

ΔPf,tp——两相段摩擦压降，Pa

ΔPg,tp——两相段重位压降，Pa

ΔPt——测压段总压降，Pa

Re——Reynolds数

X——Martinelli参数

Xeq——热平衡含气率

x——质量含气率

ε——矩形通道宽高比

λ——摩擦阻力系数

下角标

cal——计算值

exp——实验值

l,g——分别为液相、气相

i,o——分别为进口、出口

ONB——过冷沸腾起始点

ONVG——净蒸汽产生点起始点