液膜流动有着流量小、换热高的优点,广泛应用于海水淡化[1]、食品加工和核电站安全壳液膜冷却[2]等工业场所,已有许多学者对降膜流动行为进行了研究[3-4]。在降膜流动中,学者们发现液膜受表面张力和固液接触角等因素的影响,当液体流量小时,液膜将收缩成溪流,换热表面出现干斑现象,而干斑会恶化传热[5-6]。为解决此问题,研究人员对波纹板上溪流状态的液膜进行了三维分析,指出适当增加表面横向纹路的倾角可以有效改善沟流,提高填料润湿面积[7];或者在光滑平板表面布置不同形状的微结构,并通过实验研究发现竖向条纹、菱形条纹(条纹方向大致与流动方向相同)较凹凸点、横向细条纹成膜性更好,可以使成膜率从45%~58%提高到74%~80%[8];或者用不同目数的砂粒对铝板表面进行砂磨,发现无论砂粒目数的大小,砂磨过后成膜性均改善非常明显,但随着时间推移成膜性会变差,分析认为可能是铝板表面氧化导致[9]。本文根据已有的研究成果,提出在光滑平板上构造箭型排布的矩形微槽,通过箭型结构的扩张角抑制液膜收缩;探究了箭型排布的矩形微槽平板上液膜流动的特性,为板式蒸发冷凝器的研发提供数据参考。
1 模拟与计算 1.1 物理模型为探究箭型微槽平板上液膜流动及换热特性,本文对平板及箭型微槽平板上的液膜流动分别进行模拟,其物理模型如图 1所示。
本文作如下假设:①流体为连续不可压缩牛顿流体,层流流动;②流体物性参数恒定为常数。
1.2.1 体积分数方程在volume of fluid(VOF)模型中,各相流体共享一个方程组,每一相的体积分数在整个计算区域内被追踪。该方法能准确地计算存在自由液面的多相流动,其控制方程为
$ \frac{{\partial {\alpha _i}}}{{\partial t}} + \mathit{\boldsymbol{u}}\nabla {\alpha _i} = 0 $ | (1) |
$ \mathop \sum \limits_{i = 1}^n {\alpha _i} = 1 $ | (2) |
式中,αi为第i相的体积分数;u为第i相的速度,m/s。
本文只包含气液两相,取n=2。控制方程中的密度、黏度由每个相共同决定,因此,气液两相系统中密度ρ和黏度μ分别表示为
$ \rho = {\alpha _{\rm{L}}}{\rho _{\rm{L}}} + (1 - {\alpha _{\rm{L}}}){\rho _{\rm{G}}} $ | (3) |
$ \mu = {\alpha _{\rm{L}}}{\mu _{\rm{L}}} + (1 - {\alpha _{\rm{L}}}){\mu _{\rm{G}}} $ | (4) |
质量守恒方程为
$ \frac{{\partial \rho }}{{\partial t}} + \nabla \left( {\rho \mathit{\boldsymbol{u}}} \right) = 0 $ | (5) |
式中ρ由式(3)计算得到。
1.2.3 动量方程在整个区域内求解单一动量方程,所得到的速度场被所有相共用。动量方程通过物性参数ρ和u与体积分数相联系
$ \frac{\partial }{{\partial t}}\left( {\rho \mathit{\boldsymbol{u}}} \right) + \nabla \left( {\rho \mathit{\boldsymbol{uu}}} \right) = - \nabla p + \nabla [\mu (\nabla \mathit{\boldsymbol{u}} + \nabla {\mathit{\boldsymbol{u}}^{\rm{T}}})] + \rho g + \mathit{\boldsymbol{F}} $ | (6) |
式中,g为重力加速度,m/s2。
对于存在自由液面的流动,表面张力对流动的影响非常重要,不可忽略。本文采用Brackbill等[10]提出的continue surface force(CSF)模型计算表面张力源项,在动量方程中表面张力源项为
$ \mathit{\boldsymbol{F}} = {\sigma _{ij}}\frac{{\rho {k_i}\nabla {\alpha _i}}}{{({\rho _i} + {\rho _j})/2}} $ | (7) |
式中,σij为气液界面张力系数,N/m;k为界面曲率,由界面处单位法向量n表示
$ k = \nabla \frac{\mathit{\boldsymbol{n}}}{{\left| \mathit{\boldsymbol{n}} \right|}} = \frac{1}{\mathit{\boldsymbol{n}}}\left[ {\left( {\frac{\mathit{\boldsymbol{n}}}{{\left| \mathit{\boldsymbol{n}} \right|}}\nabla } \right)\left| \mathit{\boldsymbol{n}} \right| - \left( {\nabla \mathit{\boldsymbol{n}}} \right)} \right] $ | (8) |
$ \mathit{\boldsymbol{n}} = \nabla {\alpha _i} $ | (9) |
雷诺数定义为
$ Re = \frac{{4\mathit{\Gamma }}}{{{\mu _{\rm{l}}}}} $ | (10) |
式中,Γ为单位宽度质量流量,kg/(m·s);μl为流体动力黏度,N·s/m2。
1.2.5 能量方程能量方程为
$ \frac{\partial }{{\partial t}}\left( {\rho E} \right) + \nabla \left[ {\mathit{\boldsymbol{u}}\left( {\rho E + p} \right)} \right] = \nabla ({K_{{\rm{eff}}}}\nabla T) + {S_{\rm{e}}} $ | (11) |
式中,Keff为有效热传导率,W/(m·K);Se为能量方程源项,W/m2。
1.2.6 比湿面积根据液膜在倾斜板上的铺展效果,定义无量纲参数—比湿面积ap,比湿面积体现了板的润湿程度。
$ {a_{\rm{p}}} = \frac{{{A_{{\rm{wetted}}}}}}{{{A_{{\rm{plate}}}}}} \times 100\% $ | (12) |
式中,Awetted为被润湿的板的面积,m2;Aplate为整个板的面积,m2。
1.3 计算方法与边界条件选择VOF模型,压力项采用body force weighted算法,压力-速度耦合方程的求解采用pressure implicit with splitting of operators (PISO)方法,动量项采用二阶迎风格式,体积分数项采用geo-reconstruct格式,时间步长取10-5 s,平板倾斜角度为45°;计算所用液相介质为饱和水,气相为空气,固液接触角设置为20°,表面张力0.072 N/m。
数值计算的边界设置如下:液相的进出口为速度入口和压力出口,气相的进出口为压力进口和压力出口,气相侧(top)为对称边界条件,平板表面和两侧壁面设置为无滑移壁面条件。
计算参数如表 1所示。
为验证模型的正确性与可靠性,对文献[11]中所给出的实验条件进行计算。文献中的实验条件为:壁面热流密度2 000 W/m2,入口单位宽度流量0.02 kg/(m·s),温度303 K,气相进口速度2 m/s,温度303 K,气液顺流。图 2为沿流动方向平板壁面温度TW,从图中可以看出,模拟数值与实验数据吻合良好。
为验证计算结果对网格的无关性,对光板模型近壁面1 mm高度上划分不同密度的加密网格。图 3是不同网格数量下液体在X=10 mm处的液膜厚度,以液相体积分数等于0.5为气液相界面[12]。从图中可以看出,当网格数量为56×104(即入口1 mm高度上划分16个网格)时,其液膜厚度与网格数量为70×104(入口1 mm高度划分20个网格)时基本吻合。充分考虑计算精度和迭代时间,所有模型均在流体区域上采用数量为56×104的网格密度进行加密处理。
要充分发挥液膜流动强化传热传质的优势,就需要使之润湿传热表面,防止传热表面出现干斑,致使传热恶化。比湿面积ap很直观地表现了液膜对板面的润湿程度。
2.1 箭型角度θ对比湿面积的影响以液相体积分数等于0.5为气液相界面[12],case1、case 2、case 3及光板为研究对象,图 4展示了不同角度的箭型微槽平板和光板在不同雷诺数时,平板表面液膜覆盖的比湿面积。可以看出,液膜在平板上的铺展面积随着雷诺数的增加而增加;带有箭型微槽构型平板上的液膜铺展面积均大于光板表面的液膜面积;在本文研究的雷诺数范围内,60°箭型微槽结构将液膜比湿面积由光板的62%~89%提高到68%~94%,90°箭型微槽结构平板表面上比湿面积为76%~97%,在120°箭型微槽结构平板上,比湿面积为84%~94%。从图中还可以看出,60°箭型结构和90°箭型结构的比湿面积与雷诺数呈线性关系,120°结构与光板表面的比湿面积则随着雷诺数增加先增大后趋于平缓。从图 5可以看出箭型排布的矩形微槽对液膜的铺展有着明显的横向引导效果,从而使得平板表面的比湿面积增大。
图 6展示了case 7、case 8、case 9、case2和光板在不同雷诺数下的比湿面积。在研究范围内,当雷诺数大于318时,所有宽度条件下的微槽平板均比光板有着更大的比湿面积,证明了箭型排布的矩形微槽对降膜流动有着良好的引导作用;同时,当宽度范围在0.3~0.5 mm之间时,平板表面液膜的比湿面积与雷诺数呈线性关系,而当微槽宽度大于0.6 mm时,随着雷诺数的增大,比湿面积先迅速增大后趋于平缓,且当雷诺数足够小(研究中为Re=199)时,其比湿面积较光板模型小。这是因为,当液相流量足够小时,液体从平板自上而下流过微槽,液体没有足够的速度使其越过较宽的微槽。鉴于上述分析,在研究范围内最优的微槽宽度约为0.5 mm,在该宽度下,液膜的比湿面积从光板的62%~89%提高到76%~97%。
以case 4、case 5、case 6、case 2及光板为研究对象,探究了箭型排布的矩形微槽深度对流动在平板表面上液膜的铺展的影响。从图 7可以看出,与微槽宽度对液膜的影响类似,在雷诺数大于318时,case 4、case 5、case 6、case 2液膜的比湿面积均大于光板表面液膜;且在雷诺数足够小(研究中为Re=199)时,微槽深度大于0.4 mm的情况下含矩形微槽的平板上的比湿面积小于光板。探究其原因,作者认为类似于微槽宽度对比湿面积的影响,当雷诺数很小时,液体没有足够的速度使其越过较深的矩形微槽。同样的,微槽深度为0.2 mm和0.3 mm时,比湿面积与雷诺数呈现线性关系;微槽深度为0.4 mm和0.5 mm时,比湿面积随着雷诺数的增大先迅速增大而后趋于平缓。综上所述,在研究范围内最优的微槽深度约为0.3 mm。
平板上箭尾排布的矩形微槽可以引导液膜向外扩展,抑制因表面张力引起的液膜收缩效应。矩形微槽对于液膜流动换热可以起到有益的促进作用,具体表现为以下两个方面:(1)增加了平板上液膜的比湿面积,可以推测,较小的液膜干斑有助于阻止局部换热恶化;(2)增加了液膜与固体的接触面积,相对于平板,凹凸不平的微槽有助于强化对液膜内部流场的扰动,是一种有效的强化换热措施。
尽管本文的研究结果表明箭尾排布的矩形微槽对液膜流动与换热特性是一种有益的促进,但这并不意味着在实际的工业应用过程中应该极力追求最小的液膜流量。在实际工业换热器中,因壁面的污垢及流体物性的影响,过小的液膜流量更容易引起液膜的随机收缩。因此,本文只是为高换热性能工业换热器的研发提供一种技术方案与基础数据支撑,在研发箭尾排布矩形微槽的降膜换热器过程中,仍需要进行大量的实验测试。
[1] |
WUNDER F, ENDERS S, SEMIAT R. Numerical simulation of heat transfer in a horizontal falling film evaporator of multiple-effect distillation[J]. Desalination, 2017, 401: 206-229. DOI:10.1016/j.desal.2016.09.020 |
[2] |
HUANG X G, YANG Y H, HU P. Experimental study of falling film evaporation in large scale rectangular channel[J]. Annals of Nuclear Energy, 2015, 76: 237-242. DOI:10.1016/j.anucene.2014.09.053 |
[3] |
ISHIKAWA H, OOKAWARA S, YOSHIKAWA S. A study of wavy falling film flow on micro-baffled plate[J]. Chemical Engineering Science, 2016, 149: 104-116. DOI:10.1016/j.ces.2016.04.021 |
[4] |
QUAN X Y, GENG Y, YUAN P F, et al. Experiment and simulation of the shrinkage of falling film upon direct contact with vapor[J]. Chemical Engineering Science, 2015, 135: 52-60. DOI:10.1016/j.ces.2015.06.055 |
[5] |
FUJITA T, UEDA T. Heat transfer to falling liquid films and film breakdown-I:sub-cooled liquid films[J]. International Journal of Heat and Mass Transfer, 1978, 21(2): 97-108. DOI:10.1016/0017-9310(78)90212-0 |
[6] |
GANIC E N, ROPPO M N. An experimental study of falling liquid film breakdown on a horizontal cylinder during heat transfer[J]. Journal of Heat Transfer, 1980, 102(2): 342-346. DOI:10.1115/1.3244285 |
[7] |
谷芳.规整填料局部流动与传质的计算流体力学研究[D].天津: 天津大学, 2004. GU F. CFD simulations of the local-flow and mass-transfer in the structured packing[D]. Tianjin: Tianjin University, 2004. (in Chinese) http://cdmd.cnki.com.cn/article/cdmd-10056-2006052134.htm |
[8] |
罗德勤.规整填料表面润湿改进及液膜流动特性研究[D].上海: 华东理工大学, 2010. LUO D Q. Study on wetting improvement and behavior of liquid film on structured packing surface[D]. Shanghai: East China University of Science and Technology, 2010. (in Chinese) http://cdmd.cnki.com.cn/Article/CDMD-10251-1011050567.htm |
[9] |
李迎春.表面改性规整填料塔性能[D].上海: 华东理工大学, 2010. LI Y C. The performance of surface modified structured packing[D]. Shanghai: East China University of Science and Technology, 2010. (in Chinese) http://www.wanfangdata.com.cn/details/detail.do?_type=degree&id=Y1858966 |
[10] |
BRACKBILL J U, KOTHE D B, ZEMACH C. A continuum method for modeling surface tension[J]. Journal of Computational Physics, 1992, 100: 335-354. DOI:10.1016/0021-9991(92)90240-Y |
[11] |
TSAY Y L, LIN T F. Evaporation of a heated falling liquid film into a laminar gas stream[J]. Experimental Thermal and Fluid Science, 1995, 11(1): 61-71. DOI:10.1016/0894-1777(94)00112-L |
[12] |
NAGAYOSHI T, MINATO A, MISAWA M, et al. Simulation of multi-dimensional heterogeneous and intermittent two-phase flow by using an extended two-fluid model[J]. Journal of Nuclear Science and Technology, 2003, 40(10): 827-833. DOI:10.1080/18811248.2003.9715425 |