降膜蒸发器传热传质与流动过程数值模拟

杨新飞，郭延柱，任丽，王海蕊，于朋玲

(1.山东省环境保护科学研究设计院，山东 济南 250013；2.山东省化工研究院，山东 济南 250014)



降膜蒸发器传热传质与流动过程数值模拟

杨新飞1，郭延柱1，任丽1，王海蕊1，于朋玲2

(1.山东省环境保护科学研究设计院，山东 济南 250013；2.山东省化工研究院，山东 济南 250014)

随着计算机技术的发展，计算流体力学成为深入研究流动传热传质过程的重要方法。文章基于流体动力学基本原理和VOF多相流模型，建立了水平管外降膜流动与相变传热传质过程数学模型，针对降膜蒸发复杂的相变过程，通过用户自定义函数UDF，将编制的计算程序嵌入FLUENT软件相应模块，对降膜蒸发过程进行了模拟研究，阐明了管外液膜、速度、温度分布及其局部传热传质特性沿圆周方向的变化规律。结果表明：在X=0.5～0.8时，液膜厚度较薄，有利于过程传热，圆管底部与顶部局部传热系数较大，并在管顶部出现传热系数最大值；相界面过余温度沿周向距离的增大而增加，圆管下方具有最大的相界面换热温差和最高的蒸发速率；此数值模拟方法和传热传质模型可以较为合理地进行水平管外降膜蒸发过程的理论分析。

降膜; 蒸发; 数值模拟; 传热; 传质

0　引言

在海水淡化过程中，流动介质的工作温度通常在70 ℃以下，开发小温差下的高效蒸发器逐渐成为该领域研究的主流方向[1]。近年来，由于低温多效蒸发系统中的水平管降膜蒸发器在小温差低流量下具有高效的传热传质性能而备受关注[2]，许多学者对降膜蒸发器内传热和流动过程进行了理论分析和实验研究[3-8]。

数值模拟作为传热、传质研究的重要方法，已全面应用到相关领域，作为理论分析、实验研究或工程设计应用的的重要手段，解决了一些实验室或实验装备由于条件受限而难以做到的情景模拟问题。

Jafar等通过CFD软件进行了水平管二维数值模拟并研究了管排数对介质流动形态的影响[9]。随后，针对管间距、管束布置方式，管子直径的模拟研究逐步展开[10-12]，喷淋密度、入口介质温度等对传热的影响也得到了进一步研究[13-14]。为了获得更好的管外成膜特性，诸如椭圆管、蛋形管等异形降膜蒸发管也得到众多学者的关注[15-16]，并且，计算也不再仅仅考虑流动，包括传热过程对膜厚和传热性能的影响也得到了深入研究分析[2.17-18]。然而，液膜蒸发这一相变过程由于涉及传质和传热，增加了数值模拟的复杂程度，因此，传统的降膜蒸发管流动仿真通常不考虑相变蒸发传热过程，从而造成一定程度的计算误差。

文章基于前人的研究，采用ANSYS Fluent6.3.26软件，通过其二次开发接口编制了用户自定义程序udf，将蒸发相变及传热过程考虑进来，分析了降膜蒸发器成膜过程、速度场、温度场及不同喷淋密度下液膜厚度和传热系数的分布规律。模拟结果同已有文献模拟及实验结果进行了分析对比。

1　模型参数与计算方法

1.1物理模型网格及边界条件

文章选取工业常用管径为19 mm的水平圆管作为研究对象。物理模型及其坐标表示如图1所示。选取流动介质为60 ℃的饱和水，外部空间充满此温度下的饱和水蒸气，介质物性参数见表1。

水平管外降膜流动Reynolds数由式(1)表示为

Re=4Γ/μ

(1)

式中：Re为Reynolds数；Γ为单位管长上单侧质量流量，kg/(m·s)；μ为流体黏性系数，pa·s。

文章设定喷淋密度为0.29 kg/(m·s)，此时Re为2469，研究表明，该数下流动处于层流状态[13]。模拟用二维模型，网格及具体边界条件设置如图1所示，网格采用四面体划分方式，壁面进行边界层处理。网格总数量为18973，网格无关性验证结果表明[15]，此时已获得网格无关性解。边界设定速度入口，入口宽度为1 mm，饱和水布液高度为9 mm，管壁为无滑移壁面，计算采用标准壁面函数进行修正，恒壁温Tw为333 K；设定管壁和饱和水接触角为20 °。

图1　圆管几何模型、网格及边界条件图

物质密度ρ/(kg·m-3)黏度μ/(kg·m-3·s-1)表面张力σ/(N·m-1)饱和水983.14.699×10-40.06622饱和蒸汽0.13021.0424×10-5-

1.2数学模拟和计算方法

模拟过程为常物性二维不可压缩非稳态层流两相流动，过程具有明显气液相界面，而Volume of Fluid模型(VOF)相对于其他多相流模型更能准确的捕捉两相界面。其引入了流体体积分数α的概念，各项流体体积分数之和等于1，由式(2)表示为

αl+αv=1

(2)

式中:αl、αv分别为饱和水与蒸汽的体积分数。

两相体积分数控制方程可由式(3)表示为

(3)

由Lee提出的相变传质模型[19]可知,蒸发过程由式(4)表示为

(4)

动量方程由式(5)表示为

(5)

其中，ρ、μ分别由式(6)、(7)表示为

ρ=ρlαl+ρvαv

(6)

μ=μlαl+μvαv

(7)

式中:ρv为蒸汽的密度，kg/m3;μl、μv分别为蒸汽和饱和水的黏度，kg/(m3·s)。

(8)

式中：σ为表面张力,N/m;界面曲率κl由式(9)表示为

(9)

(10)

式中：h为单位质量混合相焓值， J/kg；λ为介质导热系数， W/(m·K)；Sh为单位控制体积能量源项，J/(m3·s)。

混合相焓值h由式(11)表示为

(11)

其中，饱和水和蒸汽的焓值hl和hv分别由式(12)、(13)表示为

(12)

式中：Cp,l、Cp,v分别为饱和水与蒸汽的定压热容，J/(kg·k)。

介质导数系数λ由式(13)表示为

λ=λlαl+λvαv

(13)

单位控制体积能量源Sh由式(14)表示为

(14)

式中：hfg为汽化潜热， J/kg。

为方便结果的对比分析，相对坐标X及无量纲过余温度Θ分别由式(15)、(16)表示为

(15)

(16)

式中：θ为柱坐标弧度值；T、Tw及T0分布为介质温度、管壁温度及入流介质温度， K。

为衡量管外液膜传热特性，采用局部传热系数hw由式(17)计算

(17)

式中：qw为局部传热量，J/m2。

设定首相为饱和水，次相为饱和水蒸气，气液界面采用Geo-Reconstruct重构方案。操作压力为0.02 MPa，操作密度为0.1302 kg/m3。控制方程采用有限体积差分法进行离散，动量及能量方程的离散为二阶差分格式，压力离散为Body Force Weighted格式，速度压力耦合方式为PISO。

2　结果及分析

2.1液膜分布

在层流降膜流动过程中，液膜厚度对传热具有显著影响。图2为液膜形成过程及铺展情况，布液高度为9mm时，液膜从触及圆管顶部到覆盖整个圆管表面耗时约为0.06s。图3为不同喷淋密度对液膜厚度的影响，由图知，沿流动方向液膜厚度有逐渐减小的趋势，圆管上侧液膜厚度较下侧大；当X在0.5至0.8之间时，液膜厚度较薄；当X大于0.8时，即液膜在圆管末端由于重力与表面张力的共同作用，液膜厚度增加。随着喷淋密度的增大，管表面液膜厚度增加，喷淋密度过大，液膜厚度太大，不利于传热的进行；喷淋密度同样不易选择太小，较小的喷淋密度易引起下排管束的局部烧干[11]。

图2　液膜形成过程及铺展情况图

图3　不同喷淋密度对液膜厚度的影响

2.2速度分布

膜内流动介质速度分布对液膜导热具有重要影响，较高的膜内速度有利于对流换热的进行。图4为圆管外侧速度分布，取值图右侧0.5m/s矢标。壁面处速度值为0，沿法向迅速增大，具有较大的速度梯度，液膜表面速度值最大，由于重力作用，液膜表面速度沿轴向逐渐增加，最大值出现在管子底部，约为0.75m/s，该位置处液体湍流扰动强烈，传热较强。受表面切应力的影响，界面处饱和蒸汽速度值也较高，然后，沿法向减小至0。

2.3温度分布

液膜表面温度能够体现液膜内温度的渗透性，越薄的液膜，温度升高越快，越有利于界面处饱和水的蒸发。图5为不同喷淋密度下周向相界面过余温度Θi的分布。由图知，喷淋密度越小，同位置处液膜表面温度越高，由式(4)，蒸发速率越大。上侧管外液膜蓄积能量，使得温度在圆管底部得到显著提升。液膜内温度的渗透点(即界面液体受壁面加热温度开始升高的点)随喷淋密度的增加逐渐向后推移，喷淋密度为0.49kg/(m·s)时，发生在X=0.8位置处;喷淋密度为0.39kg/(m·s)时，发生在X=0.6位置处;喷淋密度为0.29kg/(m·s)时，发生在X=0.25左右。

图4　圆管外侧速度分布图

图5　不同喷淋密度下周向相界面过于温度分布图

2.4蒸发速率

管外降膜蒸发速率直接关系到蒸发量和降膜器效率，是降膜蒸发过程研究最为关心的问题之一。图6为管外降膜蒸发速率。当X<0.8时，管外的降膜蒸发速率较小;当X>0.8时，蒸发速率迅速升高，最大蒸发速率约为0.3kg/(m3·s)。从图5知，Γ=0.29kg/(m·s)，当X>0.8时，气液截面过余温度值迅速升高，由式(4)知，较大的传热温差可以提高传质速率，因此，降膜管底端蒸发速率始终要明显高于上侧气液截面。

图6　管外降膜蒸发速率图

2.5传热特性

Γ=0.29kg/(m·s)时，管外局部传热系数沿周向变化曲线如图7所示，由图7可知，本次模拟得到平均传热系数为6.58kW/(m2·s)；而圆管底部与顶部局部传热系数较大，传热性能较好，并在管顶部出现传热系数最大值为11.51kW/(m2·s)，这是由于，在管子顶端X介于0～0.3，传热处于非正规状况阶段，局部传热量大，由式(17)知，局部传热系数较大；传热系数沿周向逐渐降低，当X>0.9时，液膜流至圆管低端并开始脱落，由于两侧液膜对冲增大了膜内液体扰动，强化了管壁侧的传热特性，传热量增加，局部传热系数增大。

图7　局部传热系数沿周向变化曲线图

为了验证本研究所用流动与传热传质模型，将现有的相近工况下研究结果同绘于图7中进行比较发现，壁面传热系数沿周向的变化趋势同已有研究结果基本一致，说明本模型及数值方法的合理性和可行性，文章模拟结果在X<0.5时同Ouldhadda等[13]模拟结果和Liu等[22]实验结果相符，但是，当X>0.5时，模拟结果低于Liu等[22]的实验值。这可能与控制相变强度因子r的取值有关，它可以通过相应的实验结果不断修正，进而获取符合降膜流动过程的控制相变强度因子r。该模型可以为后续针对降膜蒸发过程及蒸发量的理论研究提供了一定程度的参考。

3　结论

通过上述研究可知：

(1) 管外液膜厚度沿流动方向逐渐减小，当X大于0.8时，液膜厚度受重力与表面张力的共同作用而增加；

(2) 周向距离X大于0.8时，气液界面过余温度显著升高，相应的蒸发速率明显高于上侧气液截面；

(3) 圆管局部传热系数沿周向逐渐降低，当X>0.9时，由于扰动现象，强化了管壁侧的传热特性，局部传热系数增大，尾部扰动有利于降膜蒸发过程的进行；

(4) 文章采用的数值模拟方法及传热传质模型可行并且较为合理，可以为后续针对降膜蒸发过程及蒸发量的理论研究提供了一定程度的参考，但如何选取合适的控制相变强度因子r尚需进一步的试验及模拟修正。

[1]KronenbergG.,LokiecF..Low-temperaturedistillationprocessesinsingle-anddual-purposeplants[J].Desalination,2001,136(1):189-197.

[2]LuoL.C.,ZhangG.M.,PanJ.H., et al.Flowandheattransfercharacteristicsoffallingwaterfilmonhorizontalcircularandnon-circularcylinders[J].JournalofHydrodynamics,2013(3):404-414.

[3]范延品. 水平管降膜蒸发实验研究[D].大连:大连理工大学,2006.

[4]高国红. 水平管降膜蒸发器传热性能研究[D].天津:天津大学,2013.

[5]RifertV.G.,PodbereznyV.I.,PutilinJ.V., et al.Heattransferinthinfilm-typeevaporatorwithprofiletubes[J].Desalination,1989,74:363-372.

[6]HuX.,JacobiA.M..Theintertubefallingfilm:Part1—Flowcharacteristics,modetransitions,andhysteresis[J].JournalofHeatTransfer,1996,118(3):616-625.

[7]HuX.,JacobiA.M..Theintertubefallingfilm:Part2—Modeeffectsonsensibleheattransfertoafallingliquidfilm[J].JournalofHeatTransfer,1996,118(3):626-633.

[8]RibatskiG.,JacobiA.M..Falling-filmevaporationonhorizontaltubes—acriticalreview[J].InternationalJournalofRefrigeration,2005,28(5):635-653.

[9]JafarF.,ThorpeG.,TuranO.F..Liquidfilmfallingonhorizontalcircularcylinders[C].16thAustralasianFluidMechanicsConference(AFMC),Sydney.Australia:TheUniversityofQueensland,2007.

[10]何茂刚,王小飞,张颖,等. 水平管降膜蒸发器管外液体流动研究及膜厚的模拟计算[J].热科学与技术,2007(4):319-325.

[11]王小飞,何茂刚,张颖. 水平管降膜蒸发器管外液体流动数值模拟[J].工程热物理学报,2008(8):1347-1350.

[12]邱庆刚,陈金波. 管束排列方式及管间距对水平管外液体成膜情况的影响分析[J].热科学与技术,2011(2):117-122.

[13]OuldhaddaD.,IdrissiA.I.,AsbikM..Heattransferinnon-Newtonianfallingliquidfilmonahorizontalcircularcylinder[J].HeatandMassTransfer,2002,38(7):713-721.

[14]SunF.,XuS.,GaoY..Numericalsimulationofliquidfallingfilmonhorizontalcirculartubes[J].FrontiersofChemicalscienceandEngineering,2012,6(3):322-328.

[15]罗林聪,张冠敏,田茂诚,等. 异形管降膜传热与流动数值模拟研究[J].工程热物理学报,2013(4):710-714.

[16]吕多. 水平椭圆管外液体降膜流动数值模拟[D].辽宁,大连:大连理工大学,2013.

[17]罗林聪, 潘继红, 田茂诚, 等. 管形对水平管降膜圆周膜厚和Nusselt数的影响[J]. 化工学报, 2013(8): 2760-2768.

[18]罗林聪. 水平异形管降膜蒸发流动与传热强化机理及实验研究[D].济南:山东大学,2014.

[19]LeeW.H. .Apressureiterationschemefortwo-phaseflowmodeling[M].Washington:HemispherePublishing,1980.

[20]SandraC.K,,SchepperD.,GeraldineJ.H., et al.Modelingtheevaporationofahydrocarbonfeedstockintheconvectionsectionofasteamcracker[J].Computers&ChemicalEngineering,2009,31(1):122-132.

[21]孙东亮,徐进良,王丽. 求解两相蒸发和冷凝问题的气液相变模型[J].西安交通大学学报,2012(7):7-11.

[22]LiuP.J..TheEvaporatingFallingFilmonHorizontalTubes[D].Masison:UniversityofWisconsin-Madison,1975.

(学科责编：吴芹)

Numerical study on flow, heat and mass transfer characteristics of horizontal-tube falling-film evaporators

Yang Xinfei, Guo yanzhu, Ren li, et al.

(Shandong Academy of Environmental Science, Jinan 250013, China)

A mathematical model with phase transition was set up based on the principles of fluid dynamics and VOF multiphase flow model. The falling film evaporation procedure code was embedded in FLUENT through the user define functions (UDF). The principles in the circumferential direction of wall film thickness, velocity, temperature, mass transfer rate and local heat transfer coefficient was analyzed. Results show that the thickness of falling film was rather thin when the X was about 0.5～0.8, and this is good for the process of the heat transfer. The bottom and the top of the tube have the large value of local heat transfer coefficient, and the biggest value occurs at the top of the tube. In the circumferential direction, the surplus temperature increases with the distance grows, causing the biggest surplus temperature and evaporation rate at the bottom of the tube. The numerical result is in good agreement with the experimental and predicted data in literature, which means that the mathematical model with heat and mass transfer is reasonable and more comprehensive, and can be used in studying the flow, heat and mass transfer characteristics of falling water film outside horizontal tubes.

falling film; evaporation; numerical simulation; heat transfer; mass transfer

2015-07-05

杨新飞(1980-)，男，工程师，硕士，主要从事大气污染防治等方面的研究.E-mail：52789044@qq.com

1673-7644(2016)01-0014-05

X512

A