强浮力作用下超临界水传热恶化现象数值研究

强浮力作用下超临界水传热恶化现象数值研究

文青龙

1,2

卢冬华

2

杨 珏

2

（1. 上海交通大学核科学与工程学院，上海，200240；2. 中科华核电技术研究院，深圳，518206）

摘 要：本文采用V2F和SST低雷诺数模型对垂直圆管内超临界水的传热恶化现象开展了数值分析，结果表明，这两种模型能较好地预测传热恶化现象发生的起始点且与实验结果的定量差异较小，观察到了超临界流体出现的多次传热恶化现象，SST模型在定性和定量都能较好地预测到这种复杂而有趣的传热恶化现象。对温度场、流场和湍流结构的分析表明，传热恶化的主要原因在于物性和近壁面区域湍流结构的改变，降低了湍动能和能量传输的能力。

关键词：超临界水；传热恶化；湍流模型；数值分析

中图分类号：TK124 文献标识码：A 文章编号：0253-231X（2011）00-0000-00

NUMERICAL INVESTIGATION OF HEAT TRANSFER DETERIORATION IN

SUPERCRITICAL WATER WITH STRONG BUOYANCY EFFECT

WEN Qing-long 1, 2 LU Dong-hua 2 YANG Jue 1

(1. School of Nuclear Science and Engineering, Shanghai Jiao Tong University, Shanghai 200240, China

2. China Nuclear Power Technology Research Institute, Shenzhen 518206, China)

Abstract: In this study, a numerical investigation of heat transfer deterioration (HTD) in supercritical water flowing through vertical tube is performed by using two low-Reynolds number turbulence models, V2F and SST models. These two low-Reynolds models are considerately extentable to reproduce the effect of buoyancy force on heat transfer and show the occurrence of localized HTD. The V2F model and SST model give the best performances in terms of prediction of the onset of HTD. The SST model is able to quantitatively reproduce the two heat transfer deterioration phenomena with low mass flux which have been found in the present study.

Key words: Supercritical Water; Heat transfer deterioration; Turbulence Model; Numerical Analysis

0 前言

超临界水堆(SCWR)因具有设备简化和热效率高等优点而成为世界核能界研究的重点。拟临界点附近急剧的超临界水物性变化使得在特定的热工水力条件下管道内超临界水的流动和湍流结构发生了显著地变化，导致超临界水的传热恶化现象。计算流体力学（CFD）技术由于能弥补高温高压条件下实验技术和费用的限制，从而成为提高对超临界水流动和传热现象认识、发展的重要工具。低雷诺数

湍流模型因其能对近壁面区域的湍流结构进行求解而在超临界流体传热恶化现象研究中也得到广泛地应用。Koshizuka[1]、He[2]、Sharabi[3]、Kim[4]和He[5]等采用低雷诺数湍流模型对圆管内超临界水和CO2等流体的传热恶化进行了数值模拟，分析了强浮力作用和加速作用下的传热恶化现象，然而少有文献从超临界流体的温场、流场和湍流结构等信息出发分析传热恶化的发生机理。

本文将采用V2F和SST两种低雷诺数湍流模型对不同试验工况下超临界水的流动和传热现象开展

收稿日期：2010-XX-XX；修订日期：2010-XX-XX 基金项目：

作者简介：文青龙（1979-），男，四川绵阳人，助理研究员，博士生，主要从事反应堆热工水力研究。

数值研究，从超临界水的温度场和湍流结构的分析出发，进一步探索超临界水传热恶化的发生机理。

本文采用Pis’menny[6]和Shitsman[7]的垂直圆管内超临界水的实验数据进行计算，如表1所示。其中，1#和3#工况均为低流速较高热流密度工况，2#工况实验条件与1#工况相同，只是流动方向相反。

1 计算工况

表1 计算工况 Table 1 Computational cases

Cases 1# 2# 3#

P/MPa 23.5 23.5 23.5

G/kg m-2 s-1)

509 509 430

qw/kW m-2)

390 390 319.87

D/mm 6.28 6.28 8.0

Dir. Up Down Up

Bo*×10-7 5.39~5.75 ---- 9.14~4.57

作者 Pis’menny Pis’menny Shitsman

2 数学模型

假设二维轴对称流动，柱坐标下的质量、动量和能量方程为：

- 质量方程：

1 rurv0rxr- 湍动能方程：

(6)

k1ukrvktrrkxxxtk1rPkGkDrrkr (1)

- 动量方程：

1p2rurvugrxrx1uuv2rererxxrrx1pruvrv2rxrx1rxre其中，剪切产生项为：

(2)

2222uuvvvPkt2xrrrx (7)

重力产生项为：

Gkugx(3)

evvuvr2e22rrxrrku2uvTTcgxT2Tkrxrx3x (8)

其中对于上升流gx-

g，对于下降流gxg，c0.3

其中e为有效粘性系数，定义为：et，t为湍流粘性系数，定义为：

tCfk2湍动能耗散率为：

1urvrrx (4)

t1trxxrrr1Tt(9)

其中C为常数，f为衰减函数。 -

能量方程：

(5)

C1f1PkGkC2f2TtE1ruhrvhrxr1rxThThrrPrTxrPrTr其中，Ttk/

2V2F模型针对v和f还有两个附加方程： -

湍流速度尺度v：

2其中Pr为分子普朗特数，σT为湍流普朗特数，

取常数0.9。

本研究采用两种低雷诺数湍流模型：V2F模型和SST模型，其k–ε方程为：

122uvrvvrrx的加速上升现象，但是这种上升现象十分微弱。无

论是定性上还是定量上， V2F模型和SST模型与实

验结果能很好地吻合。

tv21tv22rkf6vxxrrrk(10)

-

0产生项f方程：

f1frxxrrr

123v2/kC221L2fC1115v/k2TtLkPkL2Tt(11) 模型系数、衰减函数、源项和壁面处理方式参见文献[5]。

计算中，径向和轴向网格结构根据计算工况进行相应的调整，保证近壁面的一个网格的y+<0.2，分析了网格的无关性。对流项采用QUICK离散格式，速度场和压力场采用SIMPLEC算法。超临界水物性根据IAPWS-95程序生成，因管道沿程压降相对系统压力很小，物性只是随流体温度变化。

3 计算结果与讨论

如图1所示为1#工况实验结果和计算值的比较。当壁温低于拟临界温度时，壁温随主流温度升高而平滑增加，但壁温超过拟临界点温度后，出现显著的上升，最大值达到590oC以上，而后又出现下降，这种传热现象为典型的传热恶化现象。V2F和SST两种模型都能在定性上预测到这种传热行为，都能在定量上预测到浮力作用引起的传热恶化起始位置，但是V2F模型低估了传热恶化后期壁面换热能力的恢复，在定量上其壁温的预测值高于实验结果。

图1 1#工况计算与实验的比较结果

Fig.1 Comparison of predictions with experiment for case 1

如图2所示为2#工况的实验结果和计算值的比较。与1#工况相比，下降流的壁温随着主流温度升高而平滑增加，壁温高于拟临界温度时也存在壁温

图2 2#

工况计算与实验的比较结果

Fig.2 Comparison of predictions with experiment for case 2

如图3所示为采用V2F模型时1#和2#工况的实验结果和计算值的比较。在壁温低于拟临界温度的区域，不同流动方向工况的壁温都十分接近，而在壁温高于拟临界温度时，上升流的壁温明显高于下降流工况的壁温。虽然基于平均温度得到的浮力参数Bo*低于发生浮力效应的最低限值(5.6E-7)，但此工况下超临界水的浮力效应却已经十分强烈。为分析这种现象，对5个不同位置的物性和湍流特征进行分析。P1=45D对应正常传热区域，P2=60D和P3=75D分别为传热恶化起始位置和恶化过程区域，P4=89D和P5=94D分别对应传热恶化最严重区域和后期壁面换热能力恢复区域。

图3 V2F模型对1#和2#工况计算与实验的比较结果 Fig.3 Comparison of predictions by V2F model with

experiment for run1 and run2

如图4所示为上升流近壁区域的物性、速度和湍

流特征在不同位置的变化。从图4(a)可以看出在近壁区域出现了温度的最大梯度。因为超临界流体的物性变化剧烈区域为拟临界温度区域(称为大物性变化区)。图中给出了Topc ±5C的区域线，可见在上游P1位置的大物性变化区范围十分狭小，而下游P3和P4位置大物性变化区域的范围明显变大，对应不同位

置的大比热容区域的分布特征图4(b)可十分明显得看出，大物性变化区的中心位置开始远离壁面，而P5位置的大物性变化区域的中心位置却比P4靠近壁面，这种有趣的现象可在后续的分析中得到说明。对于常规流体因近壁面区域的温度变化最为剧烈，其热物性也在该区域变化最剧烈，但超临界流体的物性变化最剧烈的位置由大物体变化区域的位置决定。

浮升力影响可通过密度的变化中得到说明，如图4(c)，在P1和P2位置，壁面区域存在很大的密度变化，该区域存在显著的浮升力，改变了流体速度分布，但其大变化区域仍限制在粘性支层区域(y+<5)，并不可能对湍流特性产生显著影响。在下游P3位置，密度剧烈变化区域拓宽到整个近壁面区域，随着壁面温度上升，促使更下游P4和P5位置近壁区域的流体密度进一步降低。近壁区域低密度流体与中心区域高密度流体之间的密度差产生近壁面区域的浮升力作用，显著改变了该区域的流动和湍流特性。

图4（a）流体温度 Fig.4(a) Fluid temperature

图4(b) 定压比热

Fig.4(b) Specific heat capability

图4(c) 密度 Fig.4(c) Density

图4(c) 流体无量纲速度 Fig.4(d) Velocity

图4(e) 湍流剪应力 Fig.4(e) Turbulent shear stress

图4(f) 湍动能产生项

Fig.4(f) Turbulent shear production

图4(g) 湍动能

Fig.4(g) Turbulent kinetic energy

图4(h) 湍流导热系数

Fig.4(h) Turbulent conductive coefficient

图4 1#工况的近壁面区物性和湍流特征

Fig.4 Fluid properties and flow field near wall region for upward

flow(case 1)

图4(d)给出了速度的分布特征。P1和P2的密度仅局限在粘性支层。其速度表现为典型的强制对流湍流速度分布特征。而在P3位置，浮升力的作用显著增强，近壁面区域的流体速度明显升高，使得整个速度分布变得平坦。而随着浮升力作用的进一步增强，近壁区域的速度近一步增加，P4区域出现典型混合对流的“M”型速度分布特征，而P5 位置因为湍流和壁面换热的恢复，“M”型速度分布的近壁面最大值有所减低。图4(e)、(f)、(g)分布给出了湍流剪应力、湍动能产生项和湍动能分布。从图中可以看出，在P1和P2的位置湍流特征受浮升力的影响局限在壁面的粘性支层区域，湍动能分布受浮力影响较小。而在P3的位置，湍流剪应力，湍动能产生项和湍动能都显著降低，而在P4的位置，浮力影响进一步增大，出现了典型的层流化现象，出现了负值的湍流剪应力。而在P5由于有“M”型的速度分布产生了湍动能和壁面换热能力的恢复。图4(h)给出了不同位置的有效导热系数，可见在随着浮升力的增大，近壁面区域的有效导热系数急剧降低，低有效导热系数的区域也扩大，但在湍流恢复区，有效导热系数得到了增强且低有效导热系数的区域也缩小从而促进壁面换热能力的恢复。

图5(a) 流体温度 Fig.5(a) Fluid temperature

图5(b) 定压比热

Fig.5(b)

Specific heat capability

图5(c) 密度 Fig.5(c)

Density

图5(d) 流体无量纲速度 Fig.5(d)

Velocity

图5(e) 湍流剪应力 Fig.5(e)

Turbulent shear stress

图5(f) 湍动能产生项

Fig.5(f)

Turbulent shear production

图5(g) 湍动能

Fig.5(g) Turbulent kinetic energy

图5(h) 湍流导热系数

Fig.5(h)

Turbulent conductive coefficient

图5 2#工况的近壁面区物性和湍流特征

Fig.5 Fluid properties and flow field near wall region for upward

flow(case 2)

图5给出了下降流的近壁面区域物体和湍流特征。从图5(a)可以看出，在十分靠近壁面的区域存在较为明显的温度梯度，随着壁面温度的升高，在P4和P5 位置流体的大物性变化区域扩大到了粘性支层和过程层区(y+

=30)。在下降流中，近壁面区和中心区域的密度差引起的浮升力方向与流动方向相反，这与密度降低引起的加速效应相反，所以流体速度的改变较小。但密度降低引起的加速效应使得下游(P4,P5)的湍流剪应力和湍动能有较小的较低，一定程度上降低了有效导热系数，从而使得壁面传热有一定的降低。

图6 3＃工况的计算和实验结果比较

Fig.6 Comparison of predictions with experiment for case 3

图6给出了3#工况计算和实验结果比较。在3#

实验工况中发现在壁面温度高于拟临界点温度，而主流温度低于拟临界点温度的范围内出现两次换热恶化现象，其中第一次传热恶化的壁面温升明显高于第二次恶化的壁面温升。所有的低雷诺数湍流模型都能较好定性预测到这种现象。但除了SST 模型外，其他模型在预测第一次传热恶化的换热能力恢复和第二次传热恶化过程中都明显高估了壁面温度。而SST 模型在定量上能较好地预测到这种有趣的传热恶化现象。采用SST 模型的计算结果，对四个不同位置分析二次传热恶化现象，P1、P2为第一次传热

恶化位置，P3、P4 为第二次传热恶化位置。

图7(a) 定压比热

Fig.7 (a) Specific heat capability

图7(b) 速度

Fig. 7(b) Velocity profile

图7(c) 湍流剪应力 Fig.7(c) Turbulent shear stress

图7(d) 湍动能

Fig. 7(d) Turbulent kinetic energy 图7 3#工况的近壁面区物性和湍流特征

Fig.7 Fluid properties and flow field near wall region for case 3

图7为3#

工况不同位置的物性和湍流特征。两次恶化的过程与1#工况的是类似的。但是第一次传热恶化后期的湍流和壁面换热能力的恢复十分充分，近壁面区域的温度显著降低，大物体变化区域再次回

到十分靠近壁面的位置，随着主流温度的升高，近壁面流体温度再次上升，大物性变化区域再次扩大，并离开壁面，在近壁面区域再次形成较大范围的低密度区域，与主流的高密度区域形成密度差，产生明显的浮升力作用，出现了第二次传热恶化现象。但第二次传热恶化中主流温度和近壁面区域的密度差比第一次传热恶化时低，浮升力的作用相对较弱，所以第二次恶化现象的壁面上升比第一次要低。

4 结论

本文采用V2F和SST两种低雷诺数湍流模型对超临界流体的强制高热流密度低流率工况进行了分析计算。基于数值结果对浮力作用下的传热恶化现象进行了分析，得到了以下结论：

1） 低雷诺数湍流模型能对物性变化不剧烈的强制对流换热进行较好的预测； 2） V2F模型和SST模型能较好地预测到强浮力作用下的超临界流体传热恶化现象起始位置，但在所有模型都高估了浮力作用引起的壁面传热能力恶化，并低估了传热恶化后期的传热恢复。

参考文献

[1] Koshizuka, S., Takano, N., Oka, Y., 1995. Numerical

analysis of deterioration phenomena in heat transfer to supercritical water[J]. Int. J. Heat Transfer 38, 3077–3084. [2] He, S., Kim, W.S., Jiang, P.X., Jackson, J. D.,

2004.Simulation of mixed convection heat transfer to carbon dioxide at supercritical pressure[J]. J. Mechanical Engineering Science, 218, 1281-1296

[3] Sharabi, M., Ambrosini, W., He, S., Jackson, J.D., 2008.

Prediction of turbulent convective heat transfer to a fluid at supercritical pressure in square and triangular channels[J]. Annals of Nuclear Energy 35, 993-1005. [4] Kim, J.K., Jeon, H.K., Yoo, J.Y., Lee, J.S., 2005.

Experimental Study on Heat transfer Characteristics of Turbulent Supercritical Flow in Vertical Circular/Non-Circular Tubes[C]. In: Proc. of the 11th NURETH-11, Oct. 2-6. Avignon, France.

[5] He, S., Kim, W.S., Jackson, J.D., 2008. A computational

study of convective heat transfer to carbon dioxide at a pressure just above the critical value[J]. Applied Thermal Engineering 28, 1662-1675

[6] Pis’menny, E.N. A.E., Pioro, I.L., 2006. Heat Transfer to

Supercritical Water in Gaseous State or Affected by Mixed Convection in Vertical Tubes[C]. In: Proc. of the ICONE14, July 17-20, Miami, USA.

[7] Shitsman, M.E., 1963. Impairment of the heat transmission

at supercritical pressures. High Temperatures 1, 237-244.