天线阵面沸腾换热的数值分析及优化

殷翔,钱吉裕,孔祥举,吴华根,曹锋

(1.西安交通大学能源与动力工程学院, 710049, 西安; 2.中国电子科技集团第14研究所, 210039, 南京)

天线阵面沸腾换热的数值分析及优化

殷翔1,钱吉裕2,孔祥举2,吴华根1,曹锋1

(1.西安交通大学能源与动力工程学院, 710049, 西安; 2.中国电子科技集团第14研究所, 210039, 南京)

为研究雷达天线阵面的冷却降温效果,建立了两相流蒸发冷却模型,并着重分析了沸腾换热在局部高热流密度条件下的换热能力,局部热流高达800 kW/m2。采用VOF模型结合用户自定义控制方程,数值计算三维流道内的沸腾传热现象,从气液流动趋势上寻找阻碍换热的因素。结果表明:VOF模型能较好地用于分析气液两相流动中的层状流、泡状流以及气液流动趋势;通过减小弯道处过流面积提高局部流速的方法可以缓解气相在弯道处的汇聚,消除了0.9 K的局部过热;对于单侧高热流蒸发的数值分析,蒸发控制方程的调节系数在500左右为宜。温度计算结果与实验结果趋势一致,总体偏差5K左右,数值分析可以用于研究冷板的沸腾换热。

VOF模型;蒸发冷却;沸腾换热;数值计算

在电子设备领域,电流的发热现象普遍存在且严重危害着设备的正常使用,因而对电子设备的冷却降温研究显得非常重要。尤其是在具有较高热流密度的电子设备上,冷却不足往往导致设备过热损坏,而常规热容冷却方式受冷却极限和运行结构的制约,很难在这一领域满足实际应用[1]。中国科学院电工研究所将蒸发潜热冷却方式应用到大型电机的定子冷却领域,但依旧以实验为主要研究方式[2]。

在低热流密度沸腾换热方面,实验和数值模拟的研究都有所发展[3-4],但数值计算传质传热的精度总体偏低。姜彭采用CFX相变蒸发模型研究水雾冷却的换热性能,数值计算与实验的温差在10 K左右[5]。在高热流领域,数值研究基本限于微细通道内以及微观结构表面的气泡动力行为方面。颛瑞等对高热流微细通道作了数值研究,观察了气泡的产生及湮灭过程,但没有相关实际应用的数值研究[6-7]。国外许多学者采用VOF模型数值研究微观表面的沸腾现象,以及微细通道内的气泡生长[8-9]。总体而言,数值研究蒸发冷却技术在实际三维换热器中的应用相对较少,而在高热流领域,传热传质涉及的问题更为复杂,传热性能以及流动性能较一般流动差异较大。常规流道较微细通道的数值研究条件更为苛刻,而局部具有高热流密度的沸腾传热特性与一般流动的差异性就更大了。

本文考虑雷达天线的运动特性以及局部具有较高热流密度的特征,着重采用数值方法研究冷板换热器的蒸发冷却技术,主要涉及管道内沸腾过程的流动以及传热传质问题。通过对矩形蛇形管流道内的沸腾特性、传热传质、气泡特性以及冷却效果采用FLUENT软件结合用户自定义控制方程进行着重研究分析,同时结合实验验证,为高热流电子设备蒸发冷却换热器的设计提供了指导。

1 沸腾传热控制理论

沸腾传热过程因流动以及传热传质的复杂性,使得其发展处于半经验半理论甚至更多依靠经验的状态,而实验又很难研究流动的内部状态,数值研究相应地弥补了这方面的不足。其中,最复杂部分是在沸腾过程中气相与液相的传热传质问题上。沸腾管内流动状态受管壁面微观结构、气泡形态以及气泡尾流或气相二次流等因素的影响,本文在研究三维流道的换热性能时忽略了这些因素,以保证数值研究在三维实体模型中的可实现性,同时假设气泡为规则球体且气泡直径为定值、气液相界面光滑。本文从微元层面分析某一网格单元体内的传热传质情况,从而模拟计算三维流道内的流动及换热效果。

1.1 气泡动力特性

1.1.1 相间传热传质 假设在某一网格区域产生某一球状气泡,气泡的平面截图如图1所示,圆弧为相界面。VOF模型采用线性插值取代实际相界面,当网格步长远小于气泡直径时,这种相界面处理方法产生的数值误差是非常小的。

图1 VOF模型气液界面插值

气相和液相受热流驱使,质量和热量的转换是在球形界面上完成的。初始时刻,流动区域均为温度不大于饱和温度的过冷或饱和液体,受热升温后,贴近壁面部分流体温度首先升高,形成热边界层。当边界层内某一微小区域温度大于对应条件下工质的饱和温度时,可以认为液体“热量多余”,这部分“多余热量”需转换为气相的潜热,从而使得两相状态的温度值在饱和温度附近,保证数值计算结果与实际情况相符。这部分的转换方式是通过源项形式实现的,因而问题的复杂性便落到如何控制气相和液相之间热量和质量的转换上。

基于前文假设,某一光滑相界面在特定条件下的传输质量由Hertz Knudsen方程给出

(1)

式中:M为流体的相对分子质量;R为气体常数;P*为相界面上气相一侧局部压力。

数值计算控制方程一般采用温度值作为控制变量,式中含有压力不利于数值研究。这里的传质问题是发生在饱和状态点附近,Clausius Claperyron方程描述为

(2)

当P*足够靠近饱和状态点时

(3)

此时,三维模型的单位体积质量源项为

(4)

实际流动过程中流道内压降比较小,为了方便数值研究,将流体饱和温度设为定值。令

(5)

在同一状态条件下,γ为定值常数。气泡直径往往难以准确获得,实际操作中发现,针对R134a或R22工质的沸腾传热研究,非恶劣工况下γ的取值在100～500之间,数值计算结果更贴近实验观察。

此时,微小单元内的潜热与显热的转换能量SE=ΔhS,SE为能量控制方程内的能量源项。

1.1.2 气泡动力特性 从传热角度看,过冷流体内的气泡生长过程主要受壁面热流和Manrangoni热流影响[10],流体的表面张力随着温度的升高而降低,在热边界层内会产生表面张力梯度,气液相界面产生Manrangoni驱动力会抑制气泡生长。当壁面过热度进一步提高,壁面热流远高于Manrangoni热流时,气泡急剧生长,流动演化为核态沸腾状态,传热系数和热流密度急剧增大。图2中Fm为Manrangoni驱动力,Fe为表面张力,共同抑制气泡生长;Fb为气泡在热壁面高热流的驱动下产生的气泡膨胀力,促使气泡成长。

图2 气泡受力分布图

此外,气泡在生长过程受到的其他力包括气液相界面压力波动产生的压力、促使气泡脱离壁面的力,以及体积力和实际管道流动中产生的一系列暂不可知的复杂力,对于这些力在本文分析气泡动力特性中均忽略不计。

1.2 数值计算控制方程

为了能更好地数值研究气泡的流动特性,本文采用追踪相界面的方法,即利用VOF模型研究沸腾换热流动。VOF模型处理相界面采用线性插值取代实际界面的方法,这样的处理方式使得数值计算的相界面与实际界面吻合程度较高,但会有微小程度的数值扩散,尤其当网格单元步长与气泡尺寸相当或更大时,这种扩散使得计算结果质量不守恒,产生错误的结果,因而网格单元步长需严格控制。考虑计算机性能限制,一般网格步长取0.5～1 mm为宜,微细通道研究应选得更小。

正是由于这种追踪相界面的处理方法,使得研究对象质量守恒方程是通过求解每一相的体积分数的连续性方程实现的。本文研究对象为液相和气相,其中液相和气相体积分数方程分别为

(6)

(7)

当液相蒸发为气相时,有一部分流体的显热转化为气相的潜热;当气相冷凝为液相时,势必有一部分气相的潜热转化为流体的显热。计算中流体表征温度参数体现为显热,气液焓差表征为潜热。潜热无法直接体现在能量方程中,即而采用源项形式处理显热与潜热之间的转化问题。能量控制方程为

(8)

SE的作用体现为潜热与显热的转化。

VOF模型的动量方程是对混合相求解单一方程实现的,动量的损失体现为压力的损失,本文研究中因压降较小,故忽略动量损失。动量方程为

(9)

式中:F为动量源项,体现为表面张力、壁面黏附力对动量的影响。

2 物理模型

实验研究中,采用发热电阻贴合在具有蛇形流道的冷板上来模拟雷达天线阵面的发热及降温过程。如图3a所示,冷板与水平面的夹角记为θ,用以模拟雷达天线在转动过程中重力加速度方向与流动方向的变化关系。本文选取θ=0°,45°,90°,135°,180°共5种工况作针对性研究。图3b给出了10个加热单元的位置编号。

冷板流道区域采用1 mm步长的结构体网格,固体冷板区域采用分块处理的1 mm步长结构体网格,网格总数在250万左右。工质采用R134a,物性参数以NIST查询数据为基础,采用小区间线性插值方式嵌入FLUENT材料库。工作压力为1.68 MPa,饱和温度为333 K,进口过冷度为5K,流量为0.1 m3/h,环境温度为293 K,每个加热单元功率为200 W,热流密度约800 kW/m2。计算模型采用VOF模型,同时考虑Marangoni效应。

(a)冷板结构示意图

(b)加热单元编号

3 结果分析

3.1 整体分析

过冷工质从进口流入冷板换热器,在第一个加热单元附近温度升至饱和温度,温度梯度较大,继续吸热产生相变之后形成两相沸腾流动。在工况1～3下,当壁面受热产生气泡之后,由于相间力不足以克服重力的作用,气相在重力加速度影响下,在贴近发热单元一侧聚集,形成分层流动,使得高热流壁面直接与气相接触,从一定程度上抑制了蒸发传热。图4a为形成层状流动的气相体积分数φ分布图。在工况4、5下,气泡在壁面产生之后,随即上浮至远离加热单元一侧,使得气相能与主流充分进行质、热交换,形成泡状流,充分保证了高热流壁面与冷流体(液相)接触换热,换热性能较好。图4b为形成泡状流的气相体积分数分布图。

图4b取自蛇形管道某一连续弯道的直管段处。本文研究发现,气液两相流体在流经弯道时并不是均匀的层状流,气液相微小单元在相间作用力以及表面张力等综合作用下流经弯道时,气相会在弯道处汇聚。当通过连续弯道时,气相主流区域从上一弯道内壁附近流向下一弯道内壁区域。如图4c所示,阴影部分是气相主流示意图,在连续弯道产生的扰流作用下,气相主流并不是沿某个截面直线流动。图4b为沿管道流动中心的垂直截面图,气相体积分数分布产生了不连续的现象,这种汇聚趋势也正是下文优化分析的着重点。模拟结果同时显示,在过冷度同为5K的情况下,工况4、5的初始沸腾时间要比工况1～3长近4 s,工况4、5的热边界层厚度要明显小于工况1～3。

(a)层状流示意图

(b)泡状流示意图

(c)弯道气相汇聚示意图

为分析不同加热单元附近的换热效果,取工况1、2为研究对象,分析流经每一个加热单元的蒸发量情况,即流经加热单元后气相体积分数的变化值Δφ。图5给出了工质流经每个加热单元时Δφ的变化趋势,沿流动方向整体呈上升趋势,在出口附近有微小差异。工质流经每个加热单元时,Δφ均有所增加,但增加幅度有所不同。

文献[11]指出,气液两相流动的传热性能是核态传热和对流传热的综合作用。气相的导热性能差,阻碍换热,但其流动速度比液相大,对流换热性能好,总体换热性能是这两种因素的综合作用。加热单元10附近对应流道内气相体积分数较大,气膜较厚,气膜直接与高热流壁面接触阻碍了液相的蒸发,同时加热单元10为最后一个弯道,流动下游产生的扰流作用弱于连续弯道流动,故加热单元10的总体蒸发量小于其他地方。蛇形管流道气相体积分数分布图(图4b)显示,在弯道附近容易产生气相的汇聚、漩涡,形成流动死区,从而使得一部分气相直接与高热流壁面长时间接触,产生局部过热现象,换热效果变差,过热程度尤以工况3更为恶劣。

图5 不同加热单元的Δφ趋势图

3.2 优化分析

本文研究显示:当热边界层内产生的气相能较充分地与主流进行传质传热时,换热性能较好;反之,气相在热边界层内的汇聚部分将抑制换热。为方便研究,截取蛇形管流道某一处弯道作详细分析,重力加速度方向以工况3为准。

(a)截取弯道气相体积分数分布图

(b)截取弯道温度分布图

图6a、6b分别是截取部分的气相体积分数以及温度分布图。通过流动动画显示:图中A区域产生流动漩涡,气相在漩涡中汇聚,使得这部分高热流壁面直接与气相接触,换热效果变差;B区域是在重力加速度作用下形成的分层流动的气相汇聚。A、B区域温度局部最大值为333.9 K,局部过热度为0.9 K。气相微小单元在壁面产生后,在重力、相间力等各种综合力作用下形成漩涡,产生流动死区,不能很好地汇入主流从而使得换热效果变差。若能使得主流对气相的作用力足以克服漩涡产生的综合力作用,那么此处的汇聚情况便可以得到改善。

为此,对弯道作减小过流面积处理,图7a、7b分别是改进模型的数值模拟结果,结果显示图6a中A区域的流动漩涡消失,且A区域温度降至饱和温度,局部过热消失。A区域紧贴着加热单元中心区域,热流密度较高,当这部分区域的壁面与液相充分接触时,则可改善换热效果。

(a)改进模型气相体积分数分布图

(b)改进模型温度分布图

4 实验对比

本文主要从气液相流动特性以及流动的宏观趋势来研究沸腾换热效果,由于流态的实验验证在目前技术条件下还很难进行,因此从加热单元的温度值验证了整个数值模拟系统的可信性。实验冷板水平放置于两相流实验平台,工质由工质泵供给,流经干燥过滤器、流量计,进入实验冷板。冷板进出口均有试液镜观察流动情况,实验观察到冷板出口确为两相状态,且为层状流,与工况1的计算结果一致,随后工质流入冷凝器,最终进入储液罐,完成循环。

实验验证工况下采用R22作为工质,流量为0.618 t/h,进口过冷度约为2 K,加热单元总功率为1 010 W,温度为常温,压力在1 MPa左右。图8为加热单元温度实验值与计算值对比图。加热单元1、2的温度计算值与实验值相差甚微,此处近似为单相管内流动;沸腾现象产生之后,计算值与实验值总体趋势基本一致,但均低于实验值5K左右,在两相蒸发的数值计算中属于可接受的范围,计算结果可以用于冷板沸腾换热研究。数值研究中,加热单元、扩热板以及冷板之间的传热都是理想化的,而实际实验过程中传热存在接触热阻,必然导致加热单元的总体温度略高。由于实验管道中还存在少量不凝结气体以及杂质等,使得换热效果变差,对蒸发换热性能影响比较大,而计算工况都是理想化的纯工质,因此计算工况的冷却效果好于实验工况,加热单元温度值总体偏低。

图8 加热单元温度计算值与实验值对比图

5 结 论

(1)沸腾传热模型中传热传质控制量归结为γ的选择。一般地,对于R134a、R22工质在较高热流密度下,γ选择在500左右较为适宜,可依据网格步长以及热流大小适当地增大或减小γ值。控制方程中,同一状态条件下,影响传质的因素有气泡直径和调节系数β,这些因素均反映在γ的选取上。

(2)VOF模型能较好地模拟计算三维流道内的分层流动和泡状流动,从一定程度上反映两相流动的气液分、聚情况,为改进流道以提高换热系数提供指导。网格步长以及时间步长需慎重选择,否则将会造成质量的数值扩散。步长若与气泡尺寸相当或者更大,则将产生较大的数值扩散,导致计算结果质量严重不守恒,致使计算结果错误。

(3)对于弯道处的气相汇聚问题,通过减小过流面积提高局部流速的方法可以得到缓解,消除漩涡,从而提高换热性能。加热单元温度值的计算与实验结果差异较小,宏观趋势基本一致,数值方法可以用于冷板沸腾换热研究。

[1] 姚涛, 候哲, 顾国彪.蒸发冷却技术应用于大型汽轮发电机的技术可行性 [J].电工技术学报, 2008, 23(2): 1-5.

YAO Tao, HOU Zhe, GU Guobiao.Application of evaporative cooling technique to large turbo generator [J].Transactions of China Electrotechnical Society, 2008, 23(2): 1-5.

[2] 熊楠.蒸发冷却技术面对的几个重要问题 [J].电工电能新技术, 2001, 3(1): 36-40.

XIONG Nan.Main tasks for research and practice of evaporation-cooling technology [J].Advanced Technology of Electrical Engineering and Energy, 2001, 3(1): 36-40.

[3] 郭雷.微细通道流动沸腾换热机理及实验研究 [D].济南: 山东大学, 2011.

[4] 曾阳.蒸发冷却换热器的数值模拟与实验研究 [D].长沙: 湖南大学, 2011.

[5] 姜澎, 黄洪雁, 冯国泰, 等.水雾/蒸汽相变冲击冷却的数值模拟 [J].哈尔滨工程大学学报, 2009, 30(10): 1097-1101.

JIANG Peng, HUANG Hongyan, FENG Guotai, et al.Numerically simulating mist/steam impingement cooling with phase change [J].Journal of Harbin Engineering University, 2009, 30(10): 1097-1101.

[6] ZHUAN Rui, WANG Wen.Simulation of subcooled flow boiling in a micro-channel [J].Refrigeration International Journal of Refrigeration, 2010, 34(3): 781-795.

[7] ZHUAN Rui, WANG Wen.Simulation on nucleate boiling in micro-channel [J].Heat Mass Transfer, 2010, 53(1): 502-512.

[8] THOME J R, DUPONT V, JACOBI A M.Heat transfer model for evaporation in micro channels [J].Heat Mass Transfer, 2004, 47(14): 3375-3385.

[9] LEE W, SON G.Numerical simulation of boiling enhancement on a micro-structured surface [J].Heat Mass Transfer, 2011, 38(2): 168-173.

[10]ZHANG N L, CHAO D F.Models for enhanced boiling heat transfer by unusual Marangoni effects under microgravity conditions [J].Heat Mass Transfer, 1999, 26(8): 1081-1090.

[11]RICARDO J, LIMA D S, THOME J R, et al.Two-phase flow boiling in horizontal smooth tubes: new heat transfer results for R134a at three saturation temperatures [J].Applied Thermal Engineering, 2009, 29(7): 1289-1298.

(编辑 荆树蓉)

NumericalInvestigationandOptimizationonBoilingHeatExchangeofAntennaSurface

YIN Xiang1,QIAN Jiyu2,KONG Xiangju2,WU Huanggen1,CAO Fen1

(1.School of Energy and Power Engineering, Xi’an Jiaotong University, Xi’an 710049, China;2.The Fourteenth Institution of China Electronics Technology Group Corporation, Nanjing 210039, China)

The evaporative cooling of two-phase model was established to cool the antenna surface with partial high heat flux.The analysis was focused on the ability of cooling under the local high heat flux of 800 kw/m2.Using VOF model, as well as the governing equation defined by the user, the computation of the boiling heat exchange in the three-dimensional pipe was performed to study the reason of hindering the heat exchange from the flow of the vapor and liquid.Results show that, the model could be well used in the analysis of stratified flow, bubble flow and the flow trend of vapor and liquid in the two-phase flow.Decreasing the flow area of the curved tube and increasing the local velocity could reduce the converging of the vapor.And the local superheat of 0.9 K was dismissed.The coefficient in the evaporation control equation should be increased to about 500 to ensure the accuracy of simulation in the partial high flux model.The overall error of the numerical result was about 5K but the tendency of the temperature was similar with the experiment.So, numerical analysis could be used in the study of the boiling heat exchange for the cold plate.

model of VOF; evaporative cooling; boiling heat exchange; numerical investigation

2014-03-23。

殷翔(1991—),男,硕士生;吴华根(通信作者),男,副教授。

时间:2014-09-01

10.7652/xjtuxb201411011

TM216.1

:A

:0253-987X(2014)11-0064-06

网络出版地址:http:∥www.cnki.net/kcms/detail/61.1069.T.20140901.1009.004.html