强制对流换热下竖直冰层融化特性研究

，，，，

(哈尔滨工业大学 建筑学院，黑龙江 哈尔滨 150090)

固液相变问题广泛存在于自然界和工业领域[1-2]。其中，河流冲刷冰面以及相变蓄能材料的凝固[3]和融化现象均属于典型的固液融化问题。该类问题属于第三类边界条件下的对流相变的基础科学问题[4]，对此类问题的研究对工程实际应用有一定的指导意义。

固液相变传热问题即移动边界问题，亦称斯蒂芬问题[6]，以斯蒂芬1891年研究北极厚度而得名。求解这类问题有它固有的困难，这是因为当固相与液相的界面处吸收或释放潜热时，这个界面是移动的[7]。相变过程中的能量变化不仅包含热传导和对流换热，还有相变的潜热变化。由于物理模型的多样性和计算的复杂性，固液相变问题受到越来越广泛的关注。Ho[8]研究了空气外掠半无限大冰层融化换热问题，采用了二次多项式热平衡积分方法对强制对流作用下冰层内的温度场进行了近似求解。Okada[9]采用有限差分法对垂直热壁条件下考虑自然对流作用的冰层融化过程进行了数值模拟，得到了壁面平均Nu和Ra数之间的关系。孟凡康[10]通过热平衡积分法对长圆柱的对流融化过程进行了分阶段求解，并按圆柱内部的温度分布特性和Bi值把传热过程分成了三个阶段。上述研究均是关于封闭腔内的固液相变问题，且在计算分析融化过程把冰层作为整体，并未考虑融冰水液膜的流动和传热特性对冰层融化的影响。

与封闭腔内固液相变问题不同的是，竖直冰层融化过程中的融冰水的流动和传热对冰层的融化产生不可忽略的影响[11]。而融冰水与热空气和冰层的换热属于液膜流动传热问题，更确切的说应该是属于强制对流下的降膜流的流动传热问题。由于表面波的高度很小，仅为1～2 mm[12]，而且通过时间也十分短暂。因此，采用实验方法研究流体动力学特性的局限性很大[13]，仅限于间接地测量壁面剪切应力及膜厚随时间的变化规律。目前主要借助数值模拟方法研究其流动和传热规律，剖析热质传递机理。

冰层的融化传热问题由于非线性的特点，求解十分困难，此外相变层(固体融化)对流等因素的存在使其更加复杂。现有的关于固液相变的研究多数是把融化过程和融化对象整体研究分析，而关于相变液膜的流动特性和传热特性对融化的影响研究较少。本文提出并采用了降膜流模型模拟计算冰层融化过程，首先建立冰层融化模型，然后研究融冰水的流动和传热对冰层融化的影响特性，并分析冰层融化过程中内部温度和融冰水Nu数的变化规律。

1 冰层融化模型建立

为了更好的研究热风融冰条件下融冰水对冰层融化过程的影响，在模拟计算中重点研究融冰水的流动和换热。冰层融化过程中冰外表面产生的融冰水的流动属于自由降膜流，是在表面张力、融冰水粘性力和重力共同作用下形成的。在粘性力作用下流体趋于附着在冰层表面流动；在表面张力作用下，流体趋于平摊形成液膜，以及在不同液相率的情况下形成不同的流型；在重力和粘滞力作用下，向下流动形成表面降膜流动。图1为冰层的物理模型，冰层大小为0.4×0.4 m，厚度为5 cm。

图1 冰层融化模型

利用ANSYS软件对冰层融化过程进行建模计算，采用有限元法求解控制方程组。建立动量方程、连续性方程和能量方程求解热风加热条件下冰层的传热和融化过程。冰层与热风之间的传热方程如公式(1)所示

(1)

式中ρs——冰的密度/kg·m3；

cs——冰的比热容/J·(kg·℃)-1；

Ql——冰的融化潜热/J；

k——冰与空气的对流换热系数/W·(m2·℃)-1，与热空气的风速和温度有关。

冰层与融冰水之间的传热方程如公式(2)所示

(2)

式中Tl——融冰水的温度/℃；

Ts——冰层的温度/℃；

Tm——冰的融化温度/℃；

Hm——冰层的焓值/J；

λs——冰的导热系数/W·(m·℃)-1；

λl——水的导热系数/W·(m·℃)-1；

S——冰层融化的厚度/m。

热空气和融冰水的传热方程如公式(3)所示

(3)

式中ρl——水的密度/kg·m-3；

cl——水的比热容/J·(kg·℃)-1；

t——时间/s；

Ta——热空气的温度/℃；

Tb——融冰水靠近热空气的边界层温度/℃；

hl——热空气和融冰水的对流换热系数/(W·m2·℃)-1。

降膜流的流动控制方程如公式(4)所示

(4)

式中u——融冰水流动方向的速度分量/m·s-1；

υ——垂直于融冰水的流动方向的速度分量/m·s-1；

σ——融冰水的表面张力/N；

g——重力加速度/m·s-2；

d——液膜的瞬时厚度/m。

融冰水不仅与热空气有热量交换，其在流动换热过程中还受热空气风压的影响。为了获得融冰水液膜在热风强迫对流作用下的流动特性，在x和y方向对该微元体列力平衡方程，可得

(5)

公式(6)为粘性流体切应力计算公式

τ=2μθ

(6)

式中θ——融冰水的角变形速率。

将边界条件带入式(5)中可得液膜的速度分布

(7)

式中μ——融冰水的动力粘滞系数/kg·(m·s)-1；

τ——融冰水所受切应力/N；

δ——液膜的厚度/m。

由公式(7)可知：液膜内的速度分布与重力和切应力的作用有关。本文所研究的风向垂直与冰层，风压是以正应力对融冰水进行作用的，所以公式(7)是适用的。

根据以往对垂直表面自由降膜流的流动特性的研究表明：

当Re<20时，流动呈平滑的层流；20

根据融冰水的热物性和流动特性，本文研究的融冰水的Re<800，属于层流。冰/水的物性参数如表1所示。

表1 冰/水物性参数值(0 ℃)

在计算过程中，为了方便计算，假设液相流体为不可压缩流体，并忽略液相流动中粘性耗散。

2 模型验证

首先通过相关的实验研究[16]对冰层融化速率和融冰时长进行对比验证。具体的冰层融化速率的实验测试值和模拟计算值如图2所示。图中0时刻为冰层融化的初始时刻。由图中可以看出，冰层开始融化后融化速率快速增加。当7.1 min时冰层的融化速率趋于稳定，表明融冰水的流态变为稳定膜状流，此时冰层融化速率为0.12 kg/min。

融冰时长方面，实验测试的融冰时长为57.2 min，而通过降膜流模型计算的时长为54.4 min，实验值和模拟值的相对误差为4.7%，并且两条曲线除45 min后的走势吻合较好。45 min后实验融冰速率的下降主要是由于风速和风温的不均，而造成的冰层融化不规则，即冰层局部已完全融化。

图2 融冰速率实验测试值和模拟结果对比

图3为冰层融化时长随风温的变化，由图中的模拟结果可知，在风速一定的情况下，融冰时长随着风温的增加而减小，风温30 ℃时融冰时长为130.1 min，风温60 ℃时融冰时长缩短至55.2 min。模拟结果和实验值偏差最大处在风温30 ℃，此时模拟结果为130.1 min，实验值为138.6 min，误差仅为6.5%。

图3 融冰时长随风温的变化(风速15 m/s)

图4为冰层融化时长随风速的变化，融冰时长同样是随着风速的增加而减小，但在5 m/s到8 m/s区间段变化较大，当风速到达11 m/s后，变化趋势渐缓。模拟结果和实验值偏差最大处在风速8 m/s，此时模拟结果为100.2 min，实验值为112.5 min，误差为12.3%低于15%。综上所述采用降膜流模型计算冰层融化和实验结果吻合较好，故该模型在Re<800的范围内适用于冰层融化的计算。

图4 融冰时长随风速的变化(温度50 ℃)

3 冰层融化特性模拟分析

3.1 冰层内温度场

图5所示为融化过程中冰层内部温度变化。冰层融化的初始时刻，与热空气直接接触的冰表面温度迅速升高并首先开始融化。实验测试温度值在0 ℃附近明显高于模拟值是由于随着冰层的融化，在冰层表面的热电偶完全裸露在冰层外部，此时实验值即会随着热气流的加热持续升高，而模拟值是模拟冰层的温度场分布不会高于0 ℃。而中间层和底层的冰层温度也随着时间的推进逐渐升高，当时间到1 761 s时，中间层温度升高到0 ℃。这表明由热空气传递到冰层的热量一部分用于冰层的融化，另一部分通过冰层的导热传递到冰层内部。

图5 融化过程中冰层内部温度变化

由图6可以看出，冰层开始融化后潜热值由0 W/m2迅速升高。到174 s后潜热值逐渐稳定在6 200 W/m2，显热值则由6 400 W/m2快速下降至320 W/m2。

图6 冰层显热量和潜热量随时间的变化

潜热量和显热量的变化可以分为3个阶段。阶段1：该阶段时间较短，且没有发生相变，冰层的热量变化全部来自显热变化；阶段2：此阶段冰层融化潜热量开始增大，显热量降低；阶段3：显热量下降缓慢，融化后期显热量接近为0。

3.2 融冰水对冰层融化的影响

图7为固液相变层的移动距离变化，通过Karapantsios关联式[15]计算值对模拟结果对比分析。边界层的移动距离随着时间的增大而增大，并且其变化趋势随着时间的增大也越来越大，即曲线的斜率(移动距离的斜率即融冰速率)随着时间的增大慢慢变大。模拟结果中当融冰时长为2 000 s时，相变层的移动距离为1.03 cm，此时的瞬时移动速率为0.102 cm/min。

图7 固液相变层随时间的移动

图8为融冰水Nu数随时间的变化，Nu表征融冰水与冰层之间的换热强度。由图中曲线可以看出，融冰的初始时刻，Nu数近似为0，冰层刚开始融化后融冰水的Nu数在冰层融化初期缓慢增大。融冰时间达到185 s时，Nu数升至35并开始急速增大，此阶段融冰水由滴状流逐渐转为稳定的层状流。400 s后，Nu数基本保持在440，此时融冰水的温度边界层和速度边界层都趋于稳定，两侧的温度梯度也趋于稳定。

图8 融冰水Nu数随时间的变化

4 结论

通过采用降膜流模型对竖直冰层融化过程进行模拟计算，分析了冰层融化过程中内部温度变化，计算了不同工况下的融冰速率，研究了融冰水的流动和传热对冰层融化的影响，探讨了冰层融化过程中固液相变层和融冰水Nu数随时间的变化规律。得到的主要结论如下：

(1)采用降膜流模型应用于分析冰层融化过程，模拟结果和实验测量值、Karapantsios关联式结果吻合较好，验证了模型和计算方法的正确性。

(2)融冰水的降膜流动在冰层融化过程中起到重要作用，在冰层融冰水未形成膜状流时冰层融化速率低，在7.1 min后，形成稳定膜状流后冰层融化速率达到0.12 kg/min。

(3)分析了冰层内温度场变化随时间的关系，得到了冰层融化过程中显热量和潜热量随时间变化趋势，其中显热量随时增大而减小，潜热量随时间增大而增大。

[1]曾轶.面向储能应用的纳米复合相变材料固液相变传热特性实验研究[D].杭州:浙江大学,2014.

[2]刘富爽,赵军,胡寿根.格栅型煤粉分配器气固两相流的数值模拟[J].节能技术,2017,35(4):326-329.

[3]马晓军,蒋林滔,黄坤荣,等.强制变相强化传热系统分析[J].节能技术,2006(3):248-249.

[4]Evans JD, King JR. Asymptotic results for the Stefan problems with kinetic under cooling[J].Applied Math,2000,53(3):449-473.

[5]Myers TG, Mitchell SL, Muchatibaya GA. Cubic heat balance integral method for one-dimensional melting of a finite thickness layer[J]. International Journal of Heat and Mass Transfer,2007,50(25):5305-5317.

[6]Ho CD, Yeh HM, Wang WP. Cool thermal discharge obtained with air flowing over melting ice by complete removal of melt[J].International Journal of Heat and Mass Transfer, 2000,27(6):785-794.

[7]OKADA M. Analysis of heat transfer during melting from a vertical wall[J].Heat Mass Transfer,1984,27(11):2057-2066.

[8]孟凡康,于航,褚琦.第三类边界条件下长圆柱对流融化过程的分阶段求解[J].山东农业大学学报,2016(3):441-446.

[9]蒋章焰,宋金田,孔旭静,等.垂直加热管外自由降膜流的破断特性[J].工程热物理学报,1995(2):199-203.

[10]蒋章焰,陶正文,阎维平.垂直自由降膜流表面波的非线性演化[J].华北电力大学学报,1999(1):13-17.

[11]罗勇强,孙博,刘璐,等.降膜流动板结构变化对传热性能的影响[J].化学工程,2014(8):31-34.

[12]叶学民,阎维平,蒋章焰,等.自由降膜表面波流动和传热特性的研究[J].华北电力大学学报,1999(1):7-12.

[13]Karapantsios TD,Paras SV,Karabelas AJ.Statistical characteristics of free falling film at high Reynolds numbers[J].International Journal of Multiphase Flow,1989,15(1):1-21.

[14]Xie T,Dong J,Chen H,Jiang Y,Yao Y.Experimental investigation of deicing characteristics using hot air as heat source[J].Applied Thermal Engineering,2016(107):681-688.