几何参数对含空气蒸汽冷凝影响的数值分析

肖家禹，孙中宁，李荣绩，李文涛，丁 铭，曹夏昕，边浩志

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

蒸汽冷凝被广泛应用于燃煤电厂[1]、化工[2]、制冷[3-4]和海水淡化[5-6]等工业领域。在核电领域，现役第三代先进核电机组大多采用多道屏障进行防护，并设置非能动安全壳冷却系统(PCCS)以保证事故条件下安全壳的完整性，进而防止放射性物质外泄。系统运行时安全壳内处于大空间自然对流状态，考虑到发生反应堆失水事故(LOCA)或主蒸汽管道破口事故(MSLB)时安全壳内含有大量空气，会大幅抑制蒸汽冷凝传热[7-9]，人们对此进行了大量研究。然而，之前的研究多针对于空气份额、压力及壁面过冷度等热工参数对含空气蒸汽冷凝的影响[9-16]，发现冷凝传热系数与空气份额及壁面过冷度负相关，与压力正相关。而对管径、管长及倾斜角等几何结构参数的分析研究较少。Popiel[17]、Dehbi[18]采用数值模拟方法对管曲率影响进行分析，得到圆管相较平板平均冷凝换热系数的修正系数。Hwang等[19]对影响冷凝现象的项积分得到适用于管外冷凝的传热模型，并对管曲率效应在湍流与层流条件下的影响进行讨论，得出湍流自然对流条件下曲率效应更为显著。Dehbi[13]发现在空气份额较高时曲率对冷凝传热影响较大，局部冷凝传热系数沿冷凝壁面长度方向呈先减小后增大的趋势。Huhtiniemil等[20]进行了倾角对含空气蒸汽冷凝影响的实验研究，并未发现传热能力的增强。Anderson等[21]开展了含不凝性气体蒸汽在平板上的冷凝实验，发现与水平方向倾角小于5°时，冷凝传热系数相较水平板有所增加。Dehbi[22]建立了水平管外冷凝模型，并在不同热工参数和管径条件下进行水平管外含空气蒸汽冷凝的数值模拟研究。Cao等[23]进行了不同管径及倾角条件下的含空气蒸汽冷凝实验，并基于2 276个实验数据点拟合出包含长径比及倾角影响的实验关联式。Bian[24]基于1 500个实验数据点拟合出包含管长及管径影响的实验关联式。值得注意的是，此前对几何参数的研究相对独立，参数选取范围有限，且未进行系统性的评估；理论模型难以反映三维分布；实验测点离散，缺乏对局部现象描述与分析，因此有必要对其进行数值模拟研究。本文采用数值模拟方法在更宽几何参数范围内开展管径、管长及倾斜角对管外含空气蒸汽冷凝传热特性影响的研究，系统性地分析并给出各几何结构参数变化对含空气蒸汽冷凝传热影响的机制，并与实验关联式进行比较。

1 数值计算模型

1.1 基本控制方程

计算过程中流体的流动、传热、传质等物理过程遵循如下控制方程。

连续性方程：

(1)

动量方程：

(2)

能量方程：

(3)

组分方程：

(4)

式中：ρ为密度，kg/m3；w为速度，m/s；Sm为质量源项，kg/(m3·s)；Sρv为动量源项，N/(m3·s)；Sh为能量源项，J/(m3·s)；Sg为组分源项，kg/(m3·s)；fs为表面力，N/m2；fv为体积力，N/m3；E为能量，J；keff为等效导热系数，W/(m·K)；ω为质量份额；D为扩散系数，m2/s；p为压力，Pa；下标g为气体组分。

湍流计算模型选取可实现的k-ε模型，相较于标准k-ε模型，其对旋转流动、强逆压梯度的边界层流动、流动分离和二次流有更好的表现[25]。

1.2 冷凝模型

含空气蒸汽冷凝过程的求解基于扩散边界层冷凝模型。该模型基于传质理论，能较好地模拟蒸汽在空气中的扩散传质过程。模型通过比较近壁面蒸汽分压下的饱和温度是否低于壁面温度来判断蒸汽是否冷凝，若发生冷凝，则在气液交界面(即近壁面第1层网格)处移除蒸汽所携带的质量、动量、能量。此前研究[26-27]表明，即便在98%蒸汽质量份额的情况下，液膜所产生的热阻也仅约为总热阻的5%，因此模型设置忽略液膜的影响。模拟基于STAR-CCM+软件的液膜模型，等效于在第1层网格处守恒方程中加入如下源项以求解与冷凝相关的流动、传热过程。

质量源项：

Sm=Si=mcond/Δ

(5)

(6)

(7)

式中：mcond为冷凝质量通量，kg/(m2·s)；Δ为近壁面网格单元厚度，m；n为冷凝壁面法向；下标i及0分别表示交界面和标准状态。

动量源项：

Sρv=Smw

(8)

能量源项：

Sh=Smhv

(9)

式中，hv为焓流，J/kg。

2 冷凝模型验证

本文基于COPAIN实验[28]进行冷凝模型的验证。COPAIN装置实验段为高2.0 m、宽0.6 m、深0.5 m的矩形通道，冷凝壁面为高2.0 m、宽0.6 m的侧壁面。由于通道对称，将宽度方向0.3 m处设置为对称平面。几何模型及网格划分如图1所示。

对COPAIN实验中6种典型工况进行模拟验证，模拟结果与实验数据的对比如图2所示，图中h为换热系数。结果表明，模拟值在96%的置信度下分布在±25%误差带内，可认为该模型能较好地模拟含空气蒸汽冷凝的一些传热特性。

图1 COPAIN实验段的几何模型及网格划分Fig.1 Geometric model and meshing of COPAIN experimental section

图2 COPAIN实验验证结果Fig.2 COPAIN experiment verification result

3 模型设置与网格验证

3.1 模型设置

传热管置于一长方体大空间中心，管上方留有一定空间使流动充分发展。冷凝壁面设置为恒温壁面，大空间上方设置为速度进口，下方设置为压力出口，四周设置为绝热壁面。数值计算几何模型与网格划分如图3所示。

研究表明，传热管竖直布置且在主流流速小于1 m/s时，流动处于冷凝换热的自然对流主导区[29]。模拟中入口流速设置为0.1 m/s，总压为0.3 MPa，空气质量份额为0.56，主流温度为对应蒸汽分压下饱和温度114.7 ℃，壁面过冷度为15 ℃。

3.2 网格无关性验证

文献[29]指出边界层层数和壁面Y+值对冷凝换热特性影响大，在Y+值小于5时能较好地对冷凝换热特性进行模拟。对基于5～17层、厚度0.02～0.06 m的边界层网格进行网格无关性验证，并保证壁面Y+值小于5，最终选取15层、厚度0.06 m的边界层网格。

图3 数值计算几何模型与网格划分Fig.3 Geometric model and meshing of numerical calculation

网格无关性验证结果如图4所示。当网格数大于116 000后，冷凝传热系数对网格尺寸的变化不再敏感。

图4 网格无关性验证结果Fig.4 Grid independence verification result

4 计算结果分析

4.1 管径影响

选取管径范围4～60 mm、管长5 m的竖直圆管及高5 m、宽6 m竖直平板进行模拟。图5为管径对平均冷凝换热系数的影响。由图5可知，管径对平均冷凝传热系数有显著影响，管径4 mm时平均冷凝传热系数高达60 mm时的2倍。平均冷凝传热系数随着管径的增大而减小，且减幅迅速减缓，管径大于30 mm后趋于平稳。与Bian[24]和Dehbi[13]所得变化趋势一致，且与实验关联式偏差在±25%以内。Bian-Ding关联式在低壁面过冷区由于该区壁面过冷度指数项绝对值偏大而导致预测值偏高。

图5 管径对平均冷凝换热系数的影响Fig.5 Effect of tube diameter on average condensation heat transfer coefficient

对于含空气蒸汽冷凝过程，传热管壁附近析出聚集所形成的高浓度空气层是传热热阻的主要来源。在此定义高浓度空气层，如式(10)所示：

(10)

式中：ωair,i为气-液交界面处空气质量份额；ωair，b为选取边界层处空气质量份额；ωair，∞为主流空气质量份额；x为选取高浓度空气层差额的比值。

不同管径下空气层厚度及局部冷凝传热系数沿管长方向的变化如图6所示。由图6可知，空气层厚度随着管径的增大而增厚，主要由于管径越大，管壁附近气流对于空气层的扰动越小，更有利于空气层的积累增厚，若管曲率半径无限大，则相当于平板，平板较管对周向扰动有更好的阻碍效果。不同管径下空气层厚度沿管长方向变化规律基本一致，在0～1 m段显著增厚，是由于空气在管壁附近析出和聚集，高浓度空气层形成并迅速发展增厚；在1～5 m段缓慢减薄，主要由于近壁面高浓度空气层的形成，使得近壁区混合气体密度增大，在重力作用下，向下加速流动，近壁区气体与主流气体之间形成流速差，由伯努利方程可知近壁区将形成负压，从而形成主流区向管壁的横向速度(图7)，对近壁面高浓度空气层产生冲刷，利于空气层的减薄。

图6 不同管径下空气层厚度及局部冷凝传热系数沿管长方向的变化Fig.6 Variation of air layer thickness and local condensation heat transfer coefficient along tube length under different tube diameters

由图6还可看出，不同管径下局部冷凝传热系数沿管长方向变化规律基本一致，与空气层厚度呈良好的负相关性。局部冷凝传热系数在0～0.4 m段显著减小，而在0.4～5 m段缓慢增大，与Dehbi[13]所得局部冷凝传热系数沿管长方向变化规律一致。

4.2 管长影响

选取管长范围0.1～7 m、管径38 mm的竖直圆管进行模拟，结果如图8所示，平均冷凝传热系数随管长的增长先迅速减小而后缓慢增大，管长3 m左右达到最小值，在0～1 m范围内管长变化对平均冷凝传热系数有显著影响。由图9可知，主要是由于传热管顶端1 m段内空气层显著增厚，而后则缓慢减薄，局部换热能力先减弱后增强，当减弱和增强达到平衡时，平均冷凝换热系数达到最小值。此前各实验关联式适用范围基本在管长1～3.5 m范围内，而基于前文分析在管长小于1 m时，由于空气层的显著堆叠造成局部传热系数骤降将导致平均换热能力有一定程度的下降，实验关联式并未考虑此段影响，从而得出平均冷凝传热系数与管长呈微弱的正相关性，后续应考虑该段影响。整体上看，模拟结果与实验关联式偏差多在±25%以内。

图7 近壁面加速及横向速度Fig.7 Near-wall acceleration and lateral velocity

由图9可知，短管空气层厚度及局部冷凝传热系数沿管长方向的变化均与长管顶端等长部分高度一致，即沿管长方向下方流动传热并不会对上方产生影响，是因为流动处于自然对流主导区[29]，而流速与重力方向均向下。

4.3 倾斜角影响

极小的横向流速(v=0.1 m/s)可使水平管的换热能力强化50%[22]，后续模拟采用纯自然对流边界条件。

图8 管长对平均冷凝换热系数的影响Fig.8 Effect of tube length on average condensation heat transfer coefficient

选取管长3 m、管径38 mm、倾斜角0°～90°布置的传热管进行模拟。文中倾斜角θ指管与竖直方向的夹角。图10为倾斜角对平均冷凝换热系数的影响。由图10可知，平均冷凝传热系数随着倾斜角的增大逐渐增大，且增幅逐渐减缓，与sinθ呈良好的线性关系，与实验关联式符合良好。

对轴向换热特性进行分析，如图11所示。迎流面与侧壁面均产生了不同幅度的传热强化，而背流面则会产生一定程度的传热抑制。除竖直管外，局部冷凝传热系数从迎流面到侧壁面再到背流面均逐渐减小，结合图12可知其形成机制在于：倾斜布置时，管上方气体流经管表面会受到管壁的阻碍，从而在管的迎流面及侧壁面产生横向冲刷造成空气层的减薄，并且会增强壁面附近的扰动，利于传热传质的进行，换热能力显著提高；而背流面下方则会由于传热管的阻挡形成滞流区[30]，利于高浓度空气的聚集，且减小了背流面附近的扰动，导致传热恶化。此外，比较各类传热面换热能力可发现，局部冷凝传热系数均随倾角的增大而逐渐增大，主要是由于倾角的增大，增大了自然对流的冲刷角，更有利于空气层的冲刷脱落。

图9 不同管长下空气层厚度及局部冷凝传热系数沿管长方向的变化Fig.9 Variation of air layer thickness and local condensation heat transfer coefficient along tube length under different tube lengths

图10 倾斜角对平均冷凝换热系数的影响Fig.10 Effect of inclination angle on average condensation heat transfer coefficient

对周向换热特性进行分析，如图11所示。局部冷凝传热系数周向分布极不均匀，周向角α为180°处局部冷凝换热系数最大，向两侧逐渐减小，迎流面局部冷凝传热系数可达背流面3倍以上。对比倾斜管与竖直管可知，仅在周向角α处于0°～40°及320°～360°区域内会产生传热抑制，而其他区域均有不同程度的传热强化，且随着倾角的增大，传热强化能力逐渐增强。图12表明，水平管迎流面空气层较竖直管有明显减薄，并沿周向向两侧逐渐增厚，同时背流面会显著增厚，很好地解释了图11中局部冷凝传热系数沿周向的变化。

图11 不同倾角下局部冷凝传热系数轴向及周向分布Fig.11 Axial and circumferential distribution of local condensation heat transfer coefficient under different inclination angles

图12 冲刷及空气层厚度展示Fig.12 Scouring and air layer thickness exhibition

5 结论

本文采用数值模拟方法研究了管径、管长及倾斜角对管外含空气蒸汽冷凝传热特性的影响。系统性地分析并给出几何参数变化对含空气蒸汽冷凝传热特性影响的原因，得到如下结论。

1) 管径、管长及倾斜角对管外含空气蒸汽冷凝传热特性均有显著影响，减小管径和管长及增大倾斜角均有利于提高平均冷凝传热系数。但当管径超过30 mm、管长超过1 m及倾斜角超过60°时，平均冷凝传热系数对几何参数的变化不再敏感。

2) 平均冷凝传热系数随管径的增大而减小，且减幅迅速减缓，4 mm时平均冷凝换热系数可达60 mm时的2倍。与Bian和Dehbi所得趋势相一致。结果表明，冷凝传热系数与高浓度空气层厚度呈良好的负相关性，周向对空气层的扰动是管径影响的主要原因。局部冷凝传热系数沿管长方向先迅速减小、后缓慢增大，与Dehbi所得趋势相吻合，且变化趋势基本不受管曲率的影响。

3) 平均冷凝传热系数随管长的增长先迅速减小、后缓慢增大，管长3 m左右达到最小值。原因在于空气层沿管长方向先迅速增厚后缓慢减薄所引起局部冷凝传热系数的相反变化。已有实验关联式的适用范围往往不包含管长0～1 m段，导致其呈简单的正相关性，后续应考虑该段的影响。

4) 平均冷凝传热系数随倾斜角的增大而增大，且增幅逐渐减缓，与sinθ呈良好的线性关系。迎流面与侧壁面均产生了不同幅度的传热强化，而背流面则会产生一定程度的传热抑制。周向角α=180°处局部冷凝传热系数最大(即迎流面中心线处)，向两侧逐渐减小，周向角α在40°～320°范围内均会产生不同幅度的传热强化。传热强化主要源于迎流面及侧壁面的自然对流冲刷及脱落所引起的空气层减薄，而传热抑制则主要源于背流面滞流区形成导致的空气层积累增厚。