1 引言

冰蓄冷空调系统是在空调负荷很低时制冷储存，而在新空调负荷高峰时化冰取冷，因而可以全部或部分转移制冷设备的运行时间，从而缩小供电网络的峰从电力负荷差。因此，冰蓄冷空调系统在美国、日本、澳大利亚等国家发展非常迅速，已广泛使用十几年之久。

冰蓄冷空调系统由蓄冰设备和空调系统组成。根据ARI《蓄冷设备热性能指南》，蓄冰设备可分为冰盘管、冰球、冰片滑落、冰晶等几种形式[1]。其中应用最广泛的是冰盘管和冰球 [2] 。

近年来许多学者对冰盘管式蓄冰设备进行了研究[3]，生产厂家也开发出许多系列的冰盘管产品。但由于冰盘管管道较细，流动阻力较大，技术要求较高，目前国内还没有一家生产冰盘管的厂家。而冰球式蓄冰设备由于结构简单，已成为蓄冰系统的一个重要发展方向[2]。

由于密封球体内的结冰和融冰过程是一个伴随着相变的导热与自然对流换热的复杂过程，掠过冰球的载冷剂的流动和换热过程也很复杂，因此目前很难看到有关冰球蓄冷和取冷方面的文章和厂家资料。国内较早研究冰球特性的是清华大学的赵庆珠教授，她在1994年发表两篇关于单个冰球蓄冰和融冰过程的文章[4，5]。而实际的冰球蓄冰设备是将许多冰球堆集在一个圆形断面蓄冰罐中，因此了解蓄冰罐的蓄冰和融冰特性对冰球式蓄冰罐的生产、设计和使用具有重要的意义。

2 数学模型的建立和求解

载冷剂在冰球式蓄冰罐中的流动是一个很复杂的过程，涉及载冷剂掠过冰球的换热以及冰球内的要变换热过程。对此作如下简化：①冰球在罐内均匀分布，载冷剂和冰球中介质的热工参数只沿流程发生变化；②冰球的传热性能只与冰球中的冰量有多少有关^[4，5]。根据这两个简化条件即可得到冰球式蓄冰罐的数学模型。

2.1 数学模型

将冰罐中的冰球沿载冷剂流程分成若干排，每排中的冰球个数为N，每个冰球的表面积为A，每排冰球占据的长度为E，冰罐中载冷剂的有效流通面积为F。在每一排中，载冷剂的温度相同，冰球的温度相同。载冷剂侧的能量方程为

　　　　　　　　　　　（1）

　　　　　　　　　　　　　　　　　　（2）

　　　　　　　　　　　　　　　　　　（3）

　　　　　　　　　　　　　　　　　　（4）

式（1）～（4）中，T_f为载冷剂温度，τ为时间，u _f为载冷剂流速，a _f为载冷剂的导温系数，q _ht为单位体积载冷剂与冰球间的换热量，q_fls为单位体积载冷剂通过罐壁的热损失，ρ_f，c _f分别为载冷剂的密度和比热容，K _c为冰球的温度，D为冰罐直径，T_amb为周围环境温度，R_w为冰罐壁热阻，a_w为载冷剂与冰球外表面的换热系数，a_n为冰球中的介质与冰球内表面的换热系数，R₀为冰球壳热阻。

其中a_w由下式确定[6]：

　　　　　　　　　　（5）

　　　　　　　　　　　　　　　　（6）

式（5）（6）中，d_p为冰球直径，λ_f为载冷剂的导热系数，Pr为载冷剂的普朗特数，Re_p为载冷剂的雷诺数，v₀为罐体中载冷剂的表现速度，ν_f为载冷剂的运动粘度。a_n由文献[7]回归得到。冰球侧的能量方程为

　　　　　　　　　　　　（7）

冰球中无冰时： 　　　　　　　　　（8）

冰球中冰水共存时： 　　　　　　　 （9）

冰球中全冰时： 　　　　　　　　　 （10）

初如条件：T_b =T₀ ，T_f =T₀，ε=ε₀　　　　　　　　　(11)

边界条件：T_f | _x=0 =T_fin　　　　　　　　　　　　　 (12)

式（7）～（12）中，Q_ht为单位冰球与载冷剂间的传热量，m为单个冰球中介质的质量，c_w 为水的比热容，ε为冰球中的结冰率，q_lt为冰的常融解潜热，c_i为冰的比热容，T₀为初初始温度，ε₀为初始结冰率，T_fin为冰罐中载冷剂的入口温度。

2.2 方程的离散和求解

采用隐式格式对上述方程进行离散，用追赶法求解离散后的方程组。具体计算过程如图1所示。

图１ 框图STAIT程序的计算

根据图1所示的计算框图，笔者编制了计算冰球式蓄冰罐的计算程序STAIT(Software of Time-dependent Analysis for Ice Tube).

3 模型验证

为了验证上述模型是否合理，笔者利用STAIT计算了若干个工况，并与深圳中电大厦门冰蓄冷系统的实测结果进行了比较。

3.1 深圳中电大厦冰蓄冷系统的简介

深圳中电大厦是我国第一家正式运行的冰蓄冷系统，其蓄冰系统部分如图2所示。

图2 深圳中电大厦蓄冰系统示意图

1995年6月，清华大学热能空调教研室赵庆珠教授应邀对该系统进行了实测。此次测试包括单罐的蓄冷取冷过程、双罐串联和并联的蓄冷取冷过程。本文中的实测结果摘自此次实验的测试报告^[8]。

3.2 模拟结果与实测结果的比较

利用STAIT程序，笔者计算了表1所示的8个测试工况下的蓄冷和取冷过程。图3至图10分别给出了这8个工况的计算结果与实测结果的比较。

图3 工况1冰罐入口温度和出口温度随时间的变化	图4 工况2冰罐入口温度和出口温度随时间的变化
图5 工况3冰罐入口温度和出口温度随时间的变化	图6 工况4冰罐入口温度和出口温度随时间的变化
图7 工况5冰罐入口温度和出口温度随时间的变化	图8 工况6冰罐入口温度和出口温度随时间的变化

从图3至图8可以看出，蓄冷开始时，载冷剂的温度下降较快，当降至-4～-6℃后，载冷剂的温度几乎不再下降。这说明冰球结冰前水温下降较快，而一旦开始结冰后，冰球与载冷剂间的换热是非常稳定的。不论是实测还是计算结果，都反映了这一规律。从图中还可以看出，不论是单罐蓄冷还是双罐串联蓄冷，计算结果与实测结果均吻合较好。　

图9 工况7冰罐入口温度和出口温度随时间的变化

图10 工况8冰罐入口温度和出口温度随时间的变化

从图9和图10可以看出，取冷开始的一段时间内，出口水温上升较慢。经过一定时间后，出口水温上升较快，逐渐接近入口水温。计算和实测结果均反映了这一趋势。从图中还可以看出，除取冰开始的段时间里计算和实测有些差别外，二者相互吻合很好。造成上述差别的原因主要在于实际开始取冷时冰罐中的结冰率沿程并不均匀，而模拟计算时却采用了相同的结冰率。

从上述分析可以看出，本文给出的数学模型定性上与实测结果完全一致，定量上与实测结果基本一致，表明本文提出的冰罐模型是合理的。

4 结论

本文针对冰球式蓄冰罐的换热特点，提出了以单个冰球换热特性为基础的冰罐的数学模型。该模型考虑了载冷剂的掠过冰球的换热、冰球内相变换热、冰罐的漏热等因素。并编制了求解冰球式蓄冰罐的换热特性的程序STAIT。利用STAIT程序，计算了8个不同的蓄冷和取冷工况，并将计算结果与深圳中电大厦冰蓄冷空调系统的实测结果进行了对比。结果表明：本文提出的冰球式蓄冰罐的数学模型较好地反映了蓄冰罐的特点，计算结果与实测结果吻合很好，可用于研究冰球式蓄冰罐的热工特性。

5 参考文献

1 AIR. 蓄冷设备热性能指南，1994。

2 李华坚，陈建平，空调系统的蓄冷调荷技术，暖通空调，1992（5）

3 T B Jekel, J W Mitchell and S A Klein. Modeling of Ice-Storage Tanks. ASHRAE Trans. 1993: 1016-1024

4 赵庆副教授，蒋久轶，冰球蓄冷过程的传热研究 第一部分：冰球蓄冰过程的数值模拟。1994年暖通空调年会论文集。

5 赵庆珠，蒋久斩，冰球蓄冷过程的传热研究 第二部分：融冰过程的传热分析，1994年暖通空调年会论文集。

6 W M罗林诺等主编，传热学应用手册（下册），北京：科学出版社，1992：455。

7 蒋久轶，清华大学热能工程系硕士论文，1994。

8 赵庆珠，李吉生等，深圳中电大厦蓄冰空调系统测试报告，1995。