环形布膜器垂直管内R113降膜蒸发换热特性数值研究

王为术， 田 苗， 崔 强， 彭 岩， 时小宝

(1. 华北水利水电大学 热能工程研究中心, 郑州 450045; 2. 中信重工机械股份有限公司, 河南洛阳 471003)

环形布膜器垂直管内R113降膜蒸发换热特性数值研究

王为术1， 田 苗1， 崔 强1， 彭 岩2， 时小宝2

(1. 华北水利水电大学 热能工程研究中心, 郑州 450045; 2. 中信重工机械股份有限公司, 河南洛阳 471003)

通过建立垂直管内降膜蒸发物理数学模型，对环形插头型布膜器管内R113的气液两相逆流降膜蒸发换热特性进行二维非稳态数值研究，分析了管内液膜流动分布以及壁面温度和液膜表面温度分布，对比了加热前后液膜厚度的变化.结果表明：随着降膜蒸发过程的进行，液膜下端开始出现液滴飞溅，且不断向上端发展；R113在管内降膜蒸发过程中壁面温度和液膜表面温度沿流动方向逐渐升高，气相温度变化趋势则相反；从壁面到管中心，温度沿径向逐步降低，在近壁面1 mm前后其分布趋势相反；加热后液膜厚度明显减小，且下游液膜厚度变得相对均匀.

降膜蒸发； R113； 壁面温度； 液膜温度； 数值研究

工业余热是指在工业生产过程中可以回收而尚未回收的能量.我国工业余热量大、分布广、品质低、回收利用率低[1-2]，尤其在钢铁、有色金属冶炼、水泥和煤炭等行业更是如此.目前，高品质余热利用技术已逐步成熟，但60～200 ℃的低品质余热回收利用难度较大.因此，低品质余热深度利用、高效余热换热装置的研发是亟待解决的关键技术问题，也是实现节能减排的重要途径.基于有机工质特性研究，采用低沸点工质的微过热度(<3 K)降膜蒸发技术可以实现高回收比的低温余热发电.

降膜蒸发技术传热端差小、效率高，已广泛应用于化工及余热回收领域[3].液膜的均匀状况和管内结垢是影响降膜蒸发换热性能的主要因素，国内外许多学者已对降膜蒸发器的流动特性和传热性能进行了大量实验研究和理论探索.Nusselt[4]对光滑层流及湍流条件下的液膜流动进行了理论研究.王永福[5]提出了二元溶液在垂直圆管内降膜传热传质耦合数学模型及其相应的数值解法.Chun等[6-7]实验研究了竖直管外降膜换热特性，拟合出降膜蒸发换热关联式，且所得公式可用于验证理论模型正确性.Gropp等[8]对R11/R113混合物在竖直管外降膜蒸发传热的研究表明，液膜发生表面蒸发和核态沸腾时，传质阻力对传热系数影响较大.随着计算流体力学(CFD)的迅速发展，CFD数值模拟方法也用于降膜流动与蒸发换热特性的研究.沈涛等[9-10]针对水平管外降膜蒸发，研究了不同管排结构下的成膜特性和传热特性.Ho等[11]对水在波纹板表面的降膜过程进行了二维和三维模拟，发现在不同喷淋密度下液膜的波动范围发生变化，气液界面剪切力对液膜波动有重要影响.刘玉峰等[12]开展了高雷诺数降膜流动特性研究，证实降膜表面波形态与高速摄影拍摄结果符合良好.邱庆刚等[13-14]通过对竖直管内溴化锂溶液进行降膜蒸发传热模拟和试验对比，分析了影响降膜流动的因素，并拟合出适用于层流降膜蒸发的传热系数关联式.

基于高效蒸发换热装置对低品质余热回收利用的迫切需求，笔者以垂直管内降膜蒸发理论为依据，在冷态成膜基础上通过UDF编程，考虑流体物性随温度的变化和气液相变传热传质过程，开展了低沸点有机工质R113在垂直圆管内降膜蒸发流动和传热传质的数值模拟，研究降膜流动液膜分布及温度分布特性，为有机工质低温余热发电关键技术中降膜蒸发器的设计提供指导.

1 降膜设计要求

降膜蒸发器是基于工质降膜流动的蒸发换热装置，包括布膜器、蒸发元件、气液分离器及排气系统等.降膜流动及换热过程如图1所示，热流体在壁面外侧流动，冷流体经布膜器分配引流在壁面内侧成膜流动，且冷、热流体的流动方向相反.管外热量经管壁导热进入液膜，经液膜导热后，再通过气液界面的对流换热，实现管内气液的传热传质.

降膜蒸发的换热特性优劣主要体现在液体能否均匀分布到每根管上，且能否在管子整个圆周壁面和长度方向上保持连续均匀的液膜分布.如果液体分布不均，势必造成液膜有的厚，有的薄，薄液膜处可能会干涸，甚至烧坏，厚液膜处，热阻较大，可能换热不良.所以液膜过厚或过薄都会导致降膜蒸发器换热特性下降.因此，选择合适的液体布膜器可有效避免换热管换热特性出现较大的波动.不同形式的液体布膜器特点不同，适用的液体种类和工况也不同.对液体分布和布膜装置设计的基本要求是：布膜均匀、操作弹性大、结构简单、制造安装方便以及操作稳定可靠.

竖直降膜蒸发器的布膜器结构形式很多，有溢流型、插头型和喷淋型等.每一种布膜器又可根据不同的操作条件设计成不同结构.图2为环形插头型布膜器的结构示意图.插头型布膜器是在管子顶端放一插头，利用插头与管子内壁间的间隙来促使液体成膜.环形插件式布膜器可以在单根管子内壁形成均匀的液膜，流动阻力小.管径的选用原则是：常压以上操作时采用小管径，一般降膜蒸发器常用管径范围是20～75 mm.布膜装置外缘与加热管内壁间应保留1～1.5 mm的环隙，料液经此环隙均匀分布到各管内壁，此时液膜周向分布均匀，无偏流现象.若环隙过大，则液体呈股状下降，不能均匀润湿管壁；环隙过小时，可能出现干壁现象.

图1 降膜流动及换热过程示意图

2 数值计算方法

2.1计算模型

模拟对象为垂直圆管，采用环形布膜装置实现液体成膜.由于计算模型关于轴线呈中心对称，选取1/2作为计算区域，建立二维物理模型，计算区域如图3所示.插件分内外两侧，外侧与圆管内壁有一定间隙，液体间隙向下流动成膜；内侧为空心柱体，气相由此竖直向上流出，形成与液膜的逆向流动.降膜结构参数见表1.

图3 降膜计算区域结构

对模型进行以下合理简化：(1) 该模拟结构对称，周向各处流动情况相同，故采用径向截面建立二维模型；(2) 主要研究液相降膜换热，略去了环形插件的上端部分；(3) 管子壁厚相对管径较小，且对液膜流动换热影响较小，所以忽略管子壁厚.

表1 物理模型结构参数

2.2网格划分与边界条件

选用高效率的结构化四边形网格离散计算区域，为保证近壁面处计算精度和气液界面的准确捕捉，对近壁处采用均匀边界层加密网格；而在网格精度要求略低的竖直方向网格较为稀疏.计算模型网格划分见图4.

图4 网格划分

利用Fluent 14.0软件求解分析，采用有限体积法对雷诺时均N-S控制方程离散，动量、湍动能和耗散率项均采用二阶迎风格式，应用PISO算法进行求解.R113的物性参数通过查表拟合公式加入到蒸发模型中.收敛残差设置小于10-5.

垂直管内R113降膜蒸发热态的数值研究是基于冷态模拟结果进行的.采用VOF方法捕捉气液两相流自由界面的流动情况，边界条件设置如下：液体入口给定质量入口，液相体积分数为1，入口温度为工质温度；布膜器中心压力入口，液相体积分数为0；圆管出口为压力出口，液相体积分数为0；壁面恒定热流，无滑移；中心线采用对称边界.参考压力为标准大气压.

2.3控制方程

气液两相的质量交换和能量传递对空气流动和气液两相换热有显著影响，尤其是液相蒸发吸收的汽化潜热直接影响着整体的温度分布.通过在Fluent中应用UDF程序添加附加源项，对传热传质过程进行数值模拟.

连续性方程

(1)

式中：ρi为相的密度；u、v、w分别为x、y、z方向的速度分量；Si为源项.

R113工质为牛顿流体，其动量守恒方程如下：

x方向

(2)

y方向

(3)

式中：p为静压；τxyz为应力张量；ρmix为气液两相混合密度；ρmixg、F分别为重力体积力和其他体积力(如两相之间的相互作用力，还可以包括其他模型源场或者用户自定义源项).

能量守恒方程

(4)

式中:E为能量；T为温度；keff为有效热传导系数.

2.4网格无关解分析

为保证数值计算结果的可靠性，需要对网格进行独立性验证，包括近壁面网格尺寸及网格数量等.表2给出了不同边界层尺寸和不同y向节点间距的网格体系.

表2 网格独立性验证

不同算例下的平均液膜厚度如图5所示.从图5可以看出，网格数目由2.1×104增加到5.9×104时，平均液膜厚度变化明显，而在5.9×104～12.3×104内，液膜厚度基本保持不变，表明网格数量大于5.9×104后，继续增加网格数量对平均液膜厚度的影响较小，可忽略不计.综合考虑，选取5.9×104的网格划分方案对垂直管内R113降膜蒸发换热特性进行数值研究.

图5 不同网格尺寸下的液膜厚度

3 计算结果及分析

3.1降膜蒸发管内液膜分布

图6给出了不同时刻降膜蒸发管内液膜流动的分布.由图6可知，在初始时刻，管内R113布膜均匀，覆着率好，液膜较薄且具有波动性，呈现传热效果较好的膜状流动.随着加热的进行，液膜下端开始出现液滴飞溅，这种液滴飞溅不断向上端发展.这是由于一方面下端液膜波动性较大，且随着液膜边流动边吸热，使得下端膜温较高，从而导致R113蒸发传热增强；另一方面，液膜蒸发进入气相，使得上行气体逐渐变成湍流，增强了对液膜的扰动.从图6还可以看出，随着液膜不断受热蒸发，下行液膜逐渐变薄.

3.2典型的降膜蒸发换热特性

液膜吸热方式主要是热传导和对流，辐射换热份额非常小，可以忽略.液膜很薄，仅1 mm左右，所以液膜的蒸发过程不会出现核态沸腾.R113液膜为膜状蒸发，传质过程发生在气液界面.

图7、图8分别给出了喷淋密度Г=1.062 kg/(m·s)、R113入口温度Tin=313 K、热流密度q=5 000 W/m2时，垂直管内流体温度分布云图和壁面温度，以及液膜表面温度沿管长方向的变化曲线.

从图7可以看出，液膜下行过程中，近壁面处液膜温度逐渐升高；随着液膜蒸发和气液换热，上行的气相温度升高幅度较大，升高了13 K.液膜温度沿流动方向逐渐增大，这是因为一方面液膜在近入口段较为平稳，对流换热量相对较小，而进入下面的波动段对流换热增强；另一方面液膜在向下流动时不断受热升温，导致下游温度略高于上游.此外，气相温度沿管长方向与液膜温度趋势相反，这是因为垂直管内气液两相逆流，气相不断受热向上流动所致.

图6 管内R113降膜流动过程

图7 管内流体温度分布图

图8 沿管长方向壁面温度和液膜表面温度的变化

由图8可知，液膜表面温度和壁面温度都是越往下，温度越高.在液膜厚度相对均匀的管子上段，温度值变化幅度较小，而在进入0.4 m之后的液膜波动段，温度有较明显升高，而且温度高低不均，这是因为该段的液膜厚度不均，热阻不同，进而影响温度分布.

图9 液膜表面平均温度随时间的变化

图9给出了液膜表面平均温度的变化情况.从图9可以看出，随着加热时间的递增，液膜表面平均温度单调升高，在0～0.35 s内升高趋势接近线性，在0.35 s之后温升加速，这是因为随着液膜温度的升高，液膜表面蒸发加快，传质传热量增加，同时随着液膜变薄，导热热阻变小，进一步造成液膜表面升温加快.

3.3径向温度分布及加热前后膜厚对比

图10为距入口分别为0.2 m、0.5 m和0.7 m处混合相沿径向的温度分布.由图10可知，液体整体温度水平沿径向从壁面到管子中心逐步降低，在近壁面1 mm内的液膜区域，距入口越远，温度越高；在距壁面1～19 mm的气相区域则呈现相反的分布趋势.

图10 不同位置沿径向的温度变化图

图11给出了冷态液膜厚度与加热0.45 s后液膜厚度的对比.从图11可以看出，由于蒸发作用，液膜厚度明显减小，整体厚度由冷态的0.7 mm左右减至0.2 mm左右；对比蒸发后管长0.1～<0.4 m和0.4～0.8 m 2段液膜的厚度发现，上游液膜厚度约为下游液膜厚度的2倍，说明下段蒸发换热强度较大.加热前后下游液膜厚度由周期性薄厚不均变为相对均匀，这是因为随着R113的不断蒸发，下游液体流量减小，在液固附着力和表面张力的作用下，液膜厚度就趋于均匀.

图11 加热前后液膜厚度变化

4 结 论

(1) 随着加热时间的增加，液膜下端会出现液滴飞溅，且不断向上端发展，液膜厚度逐渐变薄.

(2) R113在垂直管内气液逆流降膜蒸发过程中，壁面温度和液膜表面温度沿流动方向逐渐升高，气相温度变化趋势则相反.

(3) 流体整体温度沿径向从壁面到管中心逐步降低，且随着距入口距离的增大，近壁面1 mm前后温度呈现相反的分布趋势.

(4) 加热后液膜明显变薄，且加热前后下游液膜厚度由周期性薄厚不均变为相对均匀.

[1] 连红奎, 李艳, 束光阳子, 等. 我国工业余热回收利用技术综述[J].节能技术, 2011, 29(2): 123-133.

LIAN Hongkui, LI Yan, SHU Guangyangzi, et al. An overview of domestic technologies for waste heat utilization[J].EnergyConservationTechnology, 2011, 29(2): 123-133.

[2] 饶文姬, 赵良举, 刘朝, 等. 利用LNG冷能与工业余热的有机朗肯循环研究[J].工程热物理学报, 2014, 35(2): 213-217.

RAO Wenji, ZHAO Liangju, LIU Chao, et al. Research of organic Rankine cycle utilizing LNG cold exergy and waste heat[J].JournalofEngineeringThermophysics, 2014, 35(2): 213-217.

[3] 张猛, 周帼彦, 朱冬生. 降膜蒸发器的研究进展[J].流体机械, 2012, 40(6): 82-86.

ZHANG Meng, ZHOU Guoyan, ZHU Dongsheng. Research progress of falling-film evaporator[J].FluidMachinery, 2012, 40(6): 82-86.

[4] NUSSELT W. Die oberfla chenk odensation des wasserdampfes[J].ZeitschriftdesVereinesDeutscherIngenieure, 1916, 60(27): 541-546.

[5] 王永福. 脂肪酸降膜蒸发传热性能及机理研究[D]. 北京. 清华大学, 2001.

[6] CHUN K R, SEBAN R A. Heat transfer to evaporating liquid films[J].JournalofHeatTransfer, 1971, 93(4): 391-396.

[7] FAGERHOLM N E. International energy engineer[M]. Helsinki, Finland: Helsinki University of Technology, 1985: 3.

[8] GROPP U, SCHLÜNDER E U. The influence of liquid-side mass transfer on heat transfer and selectivity during surface and nucleate boiling of liquid mixtures in a falling film[J].ChemicalEngineeringandProcessing:ProcessIntensification, 1986, 20(2): 103-114.

[9] 沈涛. 水平管外降膜蒸发流动与传热特性的数值模拟[D]. 大连: 大连理工大学, 2013.

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

QIU Qinggang, CHEN Jinbo. Influential analysis of tube arrangement and tube spacing on film distributions[J].JournalofThermalScienceandTechnology, 2011, 10(2): 117-122.

[11] HO C D, CHANG H, CHEN H J, et al. CFD simulation of the two-phase flow for a falling film microreactor[J].InternationalJournalofHeatandMassTransfer, 2011, 54(15/16): 3740-3748.

[12] 刘玉峰, 童一峻, 任海刚, 等. 高雷诺数竖壁降膜流动特性的数值研究[J].红外技术, 2010, 32(10): 567-571.

LIU Yufeng, TONG Yijun, REN Haigang, et al. Numerical investigation on flow characteristics of falling water film down a vertical plate at high Reynolds number[J].InfraredTechnology, 2010, 32(10): 567-571.

[13] 邱庆刚, 陈金波. 水平管降膜蒸发器管外液膜的数值模拟[J].动力工程学报, 2011, 31(5): 357-361, 374.

QIU Qinggang, CHEN Jinbo. Numerical simulation of film formation on horizontal-tube falling film evaporators[J].JournalofChineseSocietyofPowerEngineering, 2011, 31(5): 357-361, 374.

[14] 杨颖, 王洋, 石程名, 等. 竖管内溴化锂溶液降膜蒸发数值研究[J].太阳能学报, 2012, 33(4): 658-662.

YANG Ying, WANG Yang, SHI Chengming, et al. Numerical research of falling film evaporation of lithium bromide solution in vertical tube[J].ActaEnergiaeSolarisSinica, 2012, 33(4): 658-662.

Numerical Investigation on Heat Transfer Characteristics of R113 Falling Film Evaporation in Vertical Tube of an Annular Membrane Evaporator

WANGWeishu1,TIANMiao1,CUIQiang1,PENGYan2,SHIXiaobao2

(1. Thermal Engineering Research Center, North China University of Water Resources and Electric Power, Zhengzhou 450045, China; 2. CITIC Heavy Industries Co., Ltd.,Luoyang 471003, Henan Province, China)

By establishing physical and mathematic model for R113 falling film evaporation in vertical tube of an annular membrane evaporator, two-dimensional unsteady numerical simulation was conducted on heat transfer characteristics of the gas-liquid countercurrent falling film evaporation, so as to analyze the distribution of liquid film flow, wall temperature and film surface temperature, and to compare the thickness of liquid film before and after heating. Results show that in the process of falling film evaporation, the splash of droplets starts from the lower end of liquid film, which propagates gradually toward the upper end; for R113, both the wall and film temperature rise along the flow direction during falling film evaporation, while the gas temperature falls; the temperature reduces gradually from wall to center along the radial direction, and the temperature distribution tends to be opposite at the position 1 mm away from the wall; the thickness of liquid film reduces greatly after heating, when the thickness of downstream film becomes relatively uniform.

falling film evaporation; R113; wall temperature; film surface temperature; numerical simulation

2016-09-08

：2016-11-23

国家重点研发计划资助项目(2016YFB0601404)；河南省科技创新人才计划资助项目(154100510011)；河南省高校科技创新团队支持计划资助项目(16IRTSTHN017)

王为术(1972-)，男，重庆开县人，教授，博士，主要从事多相流动与传热方面的研究. 田 苗(通信作者)，女，硕士研究生，电话(Tel.):13083692587; E-mail：799209047@qq.com.

1674-7607(2017)09-0710-06

：TK121

：A

：470.10