冷表面结霜生长的数值模型

戴明岩,金卫栋

(西安石油大学 石油工程学院,陕西 西安 710065)

自20世纪60年代以来,学者们对霜层生长及其物性变化进行了深入研究,并试图建立数学模型解释这一现象[1]。上海交大赵兰萍等[2]在前期结霜工作基础上,对冷表面上的霜层生长规律、霜层物理性质以及霜层内的热量和质量传递进行了深入的分析。北京理工大学王伟等人[3]建立了强制对流条件下低温表层温度与霜层生长速度、霜层厚度和霜层数量之间的关系,研究了结冰过程中结冰层的密度和导热率对结冰过程的影响。兰州科技学院谢福寿等人[4]在垂直平面上进行了试验,获得了垂直平面上的霜层增长规律,并利用数学方法建立了垂直平面上霜层的热传导方程,对其进行了反复修改,使得该方程能够更为精确地预测霜层的增长,并将其用于实践中。张亮亮等人[5]建立了包含对流、凝结和辐射在内的霜层能量方程式,分析了霜层物理特性、霜层增长速度和管道外传热等的时效性,提出了一套适用于低温环境的霜层增长的数值模拟方法。赵玲倩等人[6]利用可控表面氧化法制备了一种超疏水表面,并将其与裸露的铜质表面进行了比较,对其在天然对流状态下的结霜进行了预测,并对其进行了初步的探讨。

然而,现有的结霜研究主要针对的是机理分析以及实验研究,其对霜层生长特性和霜层物性的认识尚不充分。另外,现有的模式对结霜过程中的结霜换热研究还不够深入,大部分模式也没有反映结霜过程中的结霜换热。针对上述问题,本文拟从传热传质基础理论出发,结合对流换热、相变潜热、辐射换热等因素,采用准稳态的方法,通过物质平衡和能量平衡等手段建立冷表面结霜过程的数学模型,研究在冷表面温度、空气湿度和空气温度不变的情况下,霜层厚度、热阻和密度的变化曲线[7]。

1 结霜模型的建立

霜层表面的霜层形态是决定翅片管与外部空气间的传热性能的重要因素。在冷面上结霜的形成符合物质和能量守恒规律;而在工程实践中,由于受到多种因素的干扰,其传递过程更加复杂。若适当地加以考量,则将提高其数学模型的困难程度。

1.1 结霜的机理

结霜问题的研究已经引发广泛关注,要想抑制结霜就需从结霜机理进行研究[8],现研究表明结霜对温、湿等环境条件有较高的要求,当霜层表层温度大于水汽露点时,仅有显热转移。若霜的表层温度比露点温度低,但比0 ℃高,就会出现水蒸气冷凝现象[9]。另外,若结霜的表面温度在0 ℃以下,水蒸气就有可能结冰。水蒸气在表温及露点温均在0 ℃以下时,可直接变成固态。在两种条件下,冷板面与露点间的温度差都被称为过冷度。

1.2 能量守恒方程

霜层表面的霜层形态是决定翅片管与外部空气间的传热性能的重要因素[9]。在冷面上结霜的形成符合物质和能量守恒规律;而在工程实践中,由于受到多种因素的干扰,其传递过程更加复杂。若加以考量,则将提高其数学模型的困难程度,可以对这些小的变量进行相应的简化或者忽视,并作如下假定:

1)在任何时间点,霜的物理特性参数在霜的厚度方向上是一致的;

2)霜层的增长是一个动态的,但在微元时间内可以视为稳定状态;

3)周围空气的温度和翅片管的表面的温度保持不变化;

4)冷表面上的空气流动看作是自然对流的传热过程;

5)该热传导是一维稳定状态的热传导。

在结霜过程中,其过程中的传热过程由三个方面组成:周围环境温度和霜的表面温度之间的温度差,由其温度差来驱动自然对流换热;周围的空气中的水蒸气在霜中的相变化而产生的蒸发潜热;周围环境和霜层间的热辐射。

q=qs+ql+qr

(1)

式中:q——总热量,W/m2;

qs——显热部分,W/m2;

ql——潜热部分,W/m2;

qr——辐射换热部分,W/m2。

其中,

(2)

因此,能量方程可以描述为:

(3)

式中:Tf——霜层表面温度,K;

δf——霜层表面导热系数,W/(m•K);

hf——霜层表面与周围空气的对流换热系数,W/(m2•K);

Ta——环境空气温度,K;

Isv——水蒸气的气固相变焓,kJ/kg;

mv——传质速率,单位为kg/(m2•s);

σ——斯忒藩-玻尔兹曼常数;

ε——辐射常数。

1.3 质量守恒方程

作为一种由冰晶体凝结而形成的孔隙结构,霜层形成与发展的关键是周围大气中水蒸气持续的相变化和水蒸气的浸润。当水蒸气在大气中的自然对流作用下逐渐渗入到霜层中时,霜层的生长质量可以表达为:

mv=hm(ρv,a-ρv,f)

(4)

Mf=mvΔt

(5)

(6)

式中:ρv,f——霜层表面密度,kg/m3;

ρv,a——空气密度,kg/m3;

Mf——单位面积霜层质量,kg/m2;

mv——传质速率,kg/(m2•s);

Δt——单位时间,s。

根据工程热力学原理,可以用下列公式表示远场空气中和霜层表面上的水蒸气密度:

(7)

(8)

式中:φ——环境空气的相对湿度;

Rv——水蒸气的气体状态常数,J/(mol•K);

Psat,a——环境温度水的蒸气饱和压力,kPa;

Psat,f——霜层表面温度水的蒸气饱和压力,kPa。

在环境空气温度下,可用下列公式来求出水蒸气的饱和分压:

(9)

在霜层表面温度下,可用下列公式来求出水蒸气的饱和分压:

(10)

在此,契尔顿-柯尔本相似定律被用于来计算对流传质系数hm,也就是:

(11)

式中:hfs——霜层表面与空气的对流换热系数, W/(m2•K);

cp,a——空气的定压比热容,J/(kg•K);

Le——刘易斯数,取0.854;

α——热扩散系数,m2/s;

D——水蒸气在空气中的扩散系数,m2/s;

hm——传质系数,m/s。

(12)

Ra=Gr•Pr

(13)

(14)

Pr=0.707-3.67×10-4(Ta-273.15)

(15)

(16)

(17)

其中:Ra——空气的瑞利数;

Nu——空气的努塞尔数;

Gr——空气的格拉晓夫数;

Pr——空气的普朗特数;

L——定型尺寸,m;

ν——空气的运动黏度,m2/s。

因为霜面的粗糙度导致了热量传递的加强,因此在凝结后,在霜面上的空气与霜面之间的对流换热系数为:

hfs=1.23ha

(18)

1.4 霜层物性参数

以上的式(1)～(18)构成了能量和质量守恒的方程式,但是这些方程式并不是一个封闭的回路,它必须通过一个辅助方程来解决。本节提出了霜层的热传导率与霜层的表层温度之间的关系。霜层是由空气和冰晶构成的多孔介质,其结构参数直接影响着其热传导特性,关系如下:

1/λf=ζ/λmin+(1-ζ)/λmax

(19)

λmax=(1-ψ)λi+ψλa

(20)

1/λmin=(1-ψ)λi+ψ/λa

(21)

ζ=0.42(0.1+0.995ρf)

(22)

ρf=(1-ψ)ρi+ψρa

(23)

ψ=1-0.710exp[0.228(Tsur-273.15)]

(24)

Tm=(Tw+Tf)/2

(25)

其中:Tm——霜层的定性温度,K;

ψ——孔隙度;

ρi——冰密度,取917 kg/m3;

ρf——霜层密度,kg/m3;

ρa——空气密度,取1.293 kg/m3;

λi——冰的导热系数,W/(m•K);

λf——霜层导热系数,W/(m•K);

λa——空气的导热系数,W/(m•K)。

1.5 结霜发生条件

在对该数学模型进行计算之前,必须对相应的冷板面是否出现结霜情况进行判定[10]。只有在远场大气中的水蒸气含量超过靠近冷板面时,水蒸气含量才会向冷板面上传播,从而导致结霜的形成[10]。如果下式有一个不满足,则结霜不会发生:

(26)

1.6 结霜模型验证

1)输入实验中的参数:冷板面尺寸为0.04 m,风速3 m/s,冷壁面温度tw,空气相对湿度为0.8,环境空气温度为283.65 K;

2)计算霜层物性参数需要使用初始温度,一般可通过公式Tf=(Tw+273.15)/2.0求得,如精度不够则可采用二分法等方法进行重新赋值;

3)确定时间步长Δt;

4)计算初始霜层物性参数得到总结霜能量q1和霜层导热量q2,进行精度比较。如果满足精度要求ε,则输出霜层物性参数;否则使用二分法重新赋值,直至符合精度要求ε;

5)计算整个时间长度上霜层的各个物性参数;

6)输出结果,计算结束。

在计算中,对结霜模型使用固定的冷表面温度,分别计算冷板温度在-15.5 ℃和-19.7 ℃的冷表面温度下,1 h内霜层厚度,并将冷表面温度对动态结霜过程霜层特征影响的实验数据[3]与本研究中的结霜模型进行了比较。

如图1所示,随着时间,霜层厚度缓慢增加。此外,冷表面温度降低,结霜厚度越大。从图中可以看出,这个模型的预测结果与实验数据的偏差很小,所以这个模型被视为是可信的。

图1 模拟霜层厚度与实验霜层厚度对比

2 分析与讨论

2.1 冷板温度对霜层特性的影响

根据模拟条件,选取150～273 K的冷表面温度,进行计算。结果显示,冷表面温度对霜的形成有很大影响,随着温度的增大,霜层厚度逐渐减小,273 K的低温下几乎没有结霜。同时表明,结霜时间的延长会导致霜层厚度和热阻的增加,且两者的增长趋势一致。在8 h内,霜层厚度最大值为73.1 mm,热阻最大值为4.2 (m2•K)/W。这些结果表明,结霜时间和冷表面温度是影响霜层厚度和热阻的重要因素。随着冷表面温度的增加,霜层密度呈现出增加的趋势,并且在温度较低的时候,其生长速率变得缓慢,以霜层厚度的生长为主。而在接近238.1 K时,由于内部结霜的致密,结霜的厚度增加缓慢。同时,冰霜的浓度也在不断地增加,在8 h时达到了最高值。当冷表面温度增加时,霜层厚度和热阻减少,而霜层密度增加。随着时间的推移,结霜的厚度增加,其热阻和致密程度也随之增加。

2.2 环境湿度对霜层特性的影响

根据模拟条件,选取0.15～0.95的空气湿度,将空气温度和冷面温度分别设置为293 K和233 K,进行计算。空气湿度升高会导致霜层厚度增加,但增长速率逐渐减缓。同时,随着时间的推移,结霜量不断增加。与1 h时的5.8 mm相比,8 h时霜层厚度更大,达到了13.2 mm。结果表明,随着空气湿度的增加,霜层热阻开始上升,后来又下降,在湿度为0.42～0.52时达到最高点。霜层厚度增长速率较慢,说明传质速率变缓。在这个湿度范围内,霜层热阻达到最大值,并且会随着时间的推移而增加。结果表明,随着空气湿度的增加,霜层密度也呈增加趋势,并且根据实验数据可知,8 h时,霜层密度增长得最快,其值可达到105.23 kg/m3。从总体来看,空气湿度对霜层的厚度、密度影响较大,但对霜层热阻影响较小,其值在0.097～0.18 (m2•K)/W进行变化。因此,为了减少结霜,必须控制空气湿度不要在0.4～0.5。

2.3 环境温度对霜层特性的影响

根据模拟条件,选取233 K的冷表面温度和0.75的空气湿度条件,对不同空气温度下霜层物性进行了计算。结果表明,存在一个空气温度,使得结霜厚度达到最大值。随着时间的延长,霜层厚度逐渐增加,在8 h时达到最大值22.8 mm。这说明时间对霜层的厚度有着显著的影响。同时,在不同时间内,空气温度对霜层厚度的影响也不同,这些结果对于设计和优化制冷设备以及控制结霜过程具有重要意义。霜层热阻与霜层厚度的变化趋势基本相同,这表明存在一个空气温度,使得霜层热阻达到最大值。在温度范围为253～263 K时,霜层的热阻达到最大值,最大值为1.1 (m2•K)/W,同时,随着时间的延长,霜层热阻也会逐渐增大。

可以发现在253～263 K范围内,霜层的厚度和热阻都达到最大值。根据数据,当空气温度升高到263 K以上时,霜层密度的增长速度变得更加迅猛,8 h内,密度达到了199.3 kg/m3。当气温较低时,霜层增长的重点在于厚度,因此霜层较薄,热阻较小;而当气温升高时,霜层的消融和增长速度加快,导致其致密性和导热性能增强,热阻也随之减小。因此,霜层厚度和热阻都在气温为253～263 K时达到最大值,这是由于此时霜层的增长速率和致密性、导热性能达到了平衡的结果。