2019, Vol. 40
2. 华南农业大学 工程基础教学与训练中心，广东 广州 510642;
3. 华南农业大学 电子工程学院，广东 广州 510642;
4. 广东工业大学 机电工程学院，广东 广州 510006
2. Engineering Fundamental Teaching and Training Center, Guangzhou 510642, China;
3. College of Electronic Engineering, South China Agricultural University, Guangzhou 510642, China;
4. College of Mechanical and Electrical Engineering,Guangdong University of Technology, Guangzhou 510006, China
冷藏箱是农副产品贮藏的重要装备。普通贮藏箱自身不带冷源，必须接通电源才能工作[1]，出现断电情况将会造成经济损失，而工业用电又分为峰期和谷期，峰期的电费是谷期电费的2~3倍[2]，因此研究设计一款蓄冷式冷藏箱，既可以防止断电问题造成的农副产品的腐败，也可以避开用电高峰期，为企业节省用电成本，缓解工业用电压力。随着自然资源的日益短缺，既节能环保又保温冷藏两用的蓄冷式冷藏箱成为农产品贮藏的主攻方向[3]。
目前，计算流体动力学(Computational fluid dynamics，CFD)在研究冷藏保鲜环境中得到了广泛的应用[4-8] 。在蓄冷保鲜研究中，谢如鹤等[9]对常见3种不同冷板布置方式下的空载冷藏车厢内温度场进行了10 h的模拟研究，提出对于蓄冷板冷藏车，仅靠自然对流难以形成均匀的温度场及流场，建议通过合理的送风方式，增强箱体内空气对流，从而达到提高冷板冷藏车冷藏效率的目的。詹耀立等[10]对冷板布置在侧部、顶部时车厢内温度的分布和变化规律进行了研究，建议在贮存物和车体底部之间加阻止热传递的隔板或通入循环的冷气。张哲等[11]研究了冷板冷藏车内货物区在不同堆码方式下的温度场及流场，建议在顶部安装适量冷板并加装风机以增强冷空气与货物间的对流换热。需要研究蓄冷冷藏箱在强制对流条件下的情况来丰富冷藏保鲜方面的研究。
蓄冷式冷藏箱内温度的分布决定了农副产品的贮藏特性[12-15]，在降温过程中风机的风速、回风道面积和冷条的初始温度对贮藏室内温度变化速率和温差有着重要影响。本研究以蓄冷式冷藏箱为研究对象，建立三维数值模型，运用软件COMSOL 5.3对模型进行求解，以期获得风机风速、回风道面积、冷条初始温度对箱体内温度场的影响规律，并通过试验验证模型的准确性，研究结果可以为冷藏箱降温参数的设计和优化提供参考。
1 物理模型本研究箱体长×宽×高为2.3 m×2.2 m×1.7 m，如图1所示。保温隔板外皮采用玻璃钢板材料，厚度为2.5 mm，保温隔板中间材质为硬质聚氨酯泡沫板，厚度为0.08 m。冷藏箱分为冷条蓄冷室和货物贮藏室，中间采用硬质聚氨酯泡沫隔开，蓄冷室中有6根冷条，分布在箱体后壁，冷条尺寸为0.18 m×0.10 m×0.12 m，冷条蓄冷室与货物贮藏室通过变频风机进行通风和调节，风机开口的大小为0.22 m×0.22 m，位于箱体隔板的右边，距离箱体底部0.75 m，回风道位置与风机对称，位于隔板的左边。冷条蓄冷后，在风机的作用下，贮藏室的空气进入蓄冷室，与蓄冷室冷条进行热交换后从回风道流入贮藏室，实现贮藏室内的降温。
|
图 1 冷藏箱结构图和实物图 Fig. 1 Structure diagram and photo of refrigerator A、B：平面；1：贮藏室；2：回风道；3：冷条；4：风机 A, B: Plane; 1: Storage room; 2: Return duct; 3: Cold strip; 4: Fan |
为了研究送风系统对冷条贮藏箱内温度场的影响，本研究选取A、B 2个断面作为模拟值的研究平面(图1)，在B平面上布置9个监测点(温湿度传感器)，传感器布置位置如图2所示。运用COMSOL Multiphysics网格剖分功能对结构模型进行网格划分，利用其网格剖分工具和方法，生成四面体网格。同时，将风机和回风道的位置进行局部加密，如图3所示。本研究共涉及3个网格模拟，平均网格划分为366 190个域单元，26 704个边界元，1 975个边单元，网格质量良好。
|
图 2 平面监测点布置图 Fig. 2 Layout of planar monitoring point |
|
图 3 网格模型 Fig. 3 Grid model |
为了便于模拟仿真，对模型进行如下假设[16-19]：箱体密封性良好，不存在漏气、漏冷现象；箱体内气体皆为不可压缩气体，符合Boussinesq假设；箱体内气体为牛顿流体；忽略传感器导热对空气降温过程的影响。
本研究使用k-ε模型研究冷藏过程中流场的分布，k-ε广泛应用于工业领域，相对于更加高端的湍流模型来说，它既有相对较好的鲁棒性，又不需要大量的计算资源[20] 。所以本研究采用k-ε模型以及比较成熟的SIMPLE算法对冷藏箱制冷过程中空气的流场和温度场进行计算。采用有限体积法的控制方程，包括质量、动量和能量守恒方程，得到稳态通用控制方程：
| $ {\rm{div}}\left( {\rho v \varphi } \right) = {\rm{div}}\left( {\varGamma {\rm{grad}}\varphi } \right) + S\text{，} $ | (1) |
式中，φ表示广义变量；Γ表示相对于φ的广义扩散系数；S表示相对于φ的广义源项；ρ表示密度，kg/m3；v表示速度，m/s。
式(1)中各项参数及模型中k-ε的系数详见文献[21] 。
3 边界条件与求解箱体的流体动力来自风机，根据箱体的尺寸条件通过雷诺数求解公式，求得箱体的雷诺数(Re)达到了105，属于高雷诺数湍流模型[21]：
| $ Re{\rm{ = }}UL{\rm{/}}\gamma\text{，} $ | (2) |
式中，U表示速度特征尺度，m/s；L表示长度特征尺度，m；γ运动学黏性系数，m2/s。
3.1 入口条件由于试验箱体采用的是风机驱动循环形式，所以箱体入口条件为速度入口条件。输入湍流参数为湍流强度(I)和水利直径(D)。湍流强度(I)公式[22]为：
| $ I = 0.16{({{\mathop{Re}\nolimits} _{D}})^{ - 1/8}} \times 100 \text{%}\text{，} $ | (3) |
式中，ReD表示以水利直径求出的雷诺数。
通过试验测得风机对应的风速分别为4.5、4.0、3.5、3.0 m/s。
3.2 出口条件出口条件设置为压力出口条件P=0。
3.3 壁面条件影响贮藏箱热负荷的因素主要有：箱体传热的冷消耗、太阳辐射的冷消耗和箱体降温的冷消耗[22]。其中箱体传热和太阳辐射冷消耗设为热通量边界条件，热流密度(q)为：
| $ {q} = h\left( {{\theta_{{\rm{ext}}}} - \theta} \right)\text{，} $ | (4) |
式中，h表示传热系数，取值为0.226，θext表示环境温度 (℃)。
3.4 模型求解模型材料冷条中的蓄冷剂比热容为3.7 kJ/(kg·℃)，密度为1 150 kg·m−3，导热系数为0.37 W/(m·℃)，箱体隔热保温采用聚氨酯材料，密度为45 kg/m3，导热系数0.026 W/(m·℃)，比热容为1.8 kJ/(kg·℃)，箱体内流体设置为空气。结合试验数据，将贮藏室内平均温度初始化为16 ℃，冷条温度初始化为−35 ℃，贮藏室内降温过程是一个瞬态过程，运用comsol瞬态求解器对箱体内部压力、速度和温度的耦合进行求解，添加重力影响因素，设置重力加速度为9.8 m/s2，设置时间步长为1 s，计算60 min内箱体的温度变化。
4 模拟结果与分析 4.1 箱体风速和温度场分布情况对风机风速为4.5 m/s，回风道面积0.014 4 m2，冷条初始温度−35 ℃条件下的贮藏室降温过程进行模拟，当贮藏室温度达到0 ℃，停止运算。获得了降温过程中A、B平面流场分布情况(图4)。图4a、图4b是箱体A截面流场分布图，可以看出在贮藏室内，气流从回风道出来速度较大约为7.0 m/s，顺着气流方向速度逐渐减小，气流达到箱体后壁面受到阻力后向两边流动，形成涡流，贮藏室中间位置速度较小约为0.4 m/s，靠近风机回风口的位置风速逐渐增大，温度场的分布受气流场的影响，可以看出情况基本吻合。图4c、图4d是箱体B截面流场分布图，从图中可以看出正对回风道位置风速较大，温度也较低，其他区域温度分布比较均匀，靠近箱体边缘位置温度较高约为3.5 ℃，和实际情况相符。
|
图 4 冷藏室内流场分布情况 Fig. 4 Distribution of flow field in the refrigerator |
为了验证数值模拟的准确性，对冷藏箱进行降温试验，检测出箱体内温度场的变化情况。试验所用的仪器主要有无纸记录仪(杭州盘古电子公司)和BMWRHTA.S温湿度一体传感器，测量范围为−40~60 ℃，精度为±0.3 ℃。共布置了9路温湿度传感器，将传感器一端布置在监测点，另一端与无纸记录仪连接，设置每隔1 s记录1次数据。试验时间为2017年12月25日，外界平均气温18 ℃，试验开始之前开启压缩机给冷条蓄冷，当冷条蓄冷完成后进行试验，打开箱门，待贮藏室平均温度为(16±1) ℃，关闭箱门，开启风机。为了验证模型的预见性，记录了风机开启后1 h内贮藏室温度的变化情况。重复3次试验取平均值进行分析并与模拟结果进行比较，对比情况如图5所示，其中图5a为箱体贮藏室内平均温度随时间变化情况，图5b为箱体内平均温度为0 ℃时，传感器测量点的温度值与对应模拟值的对比情况。由图5可以看出，试验结果与模拟结果吻合较好，采用绝对误差(E)来评价试验结果和模拟结果的差异，计算方法为：
|
图 5 模拟值和试验值对比 Fig. 5 Comparison of simulated and experimental values |
| $ E = \frac{{\displaystyle\sum\limits_{i = 1}^N {{\rm{|}}{t_i} - {t_{\rm{p}}}{\rm{|}}} }}{N}\text{，} $ | (5) |
式中，ti为试验值，℃；tp为预测值，℃。
通过计算，降温过程中同一时刻平均温度模拟值和对应试验值平均误差为0.68 ℃；同一位置，温度模拟值和对应试验值平均误差为0.26 ℃，出现偏差的主要原因是模型的简化以及传感器的精度。
4.3 模型的实际应用 4.3.1 风机风速对贮藏室温度场的影响风机的风速是影响箱体降温过程的一个重要因素，在回风道面积为0.014 4 m2，冷条初始温度为−35 ℃条件下，依次设置风机风速为3.0、3.5、4.0、4.5 m/s进行数值模拟，当贮藏室平均温度达到0 ℃，停止运算。图6a为贮藏室平均温度变化情况，可以看出风机风速越大降温时间越短，风速为3.0、3.5、4.0和4.5 m/s所对应的降温时间依次为12.0、10.6、9.0和8.0 min，贮藏室内温度变化速率随着时间增加逐渐减小，主要是因为随着贮藏室温度的降低与冷条温度差值逐渐缩小，换热速率降低了；图6b为贮藏室在降温过程中空气温差的变化情况，可以看出风速越大贮藏室内温差越小，在降温前1 min温差迅速增大到8 ℃左右然后达到峰值后再逐渐减小，这是因为冷条温度较低，冷空气进入贮藏室后使温差迅速增大，随着降温时间延长、贮藏室温度与冷条温度差值减小以及贮藏室内空气产生对流，贮藏室内温差逐渐减小。增加风速，使得箱体内部空气的流动加速，温度分布更加均匀。
|
图 6 风机风速、回风道面积和冷条温度对贮藏室温度的影响 Fig. 6 Effects of fan speed、return duct area and cold strip temperature on storage room temperature |
回风道面积是影响流场分布的一个重要因素，在风机风速为4.0 m/s，冷条温度为−35 ℃条件下，依次对回风道面积0.025 6、0.014 4、0.006 4 m2的贮藏室进行数值模拟，当贮藏室平均温度达到0 ℃，停止运算。图6c可以看出，回风道面积越大降温时间越短，回风道面积为0.025 6、0.014 4和0.006 4 m2所对应的降温时间依次为7、9和11 min，这说明增大回风道面积可以缩短降温时间，图6d为贮藏室在降温过程中温差变化情况，可以看出，回风道面积越大贮藏室内的温差越大，这可能是因为增大了回风道面积扩大了单位时间内与冷条的热交换量，从而使温差增大。
4.3.3 冷条初始温度对贮藏室温度场的影响冷条初始温度对降温过程中贮藏室流场分布有着重要的影响，在风机风速为4.5 m/s，回风道面积为0.014 4 m2条件下，依次对冷条初始温度为−35、−30、−25和−20 ℃的贮藏室进行数值模拟，当贮藏室空气平均温度达到0 ℃，停止运算。图6e可以看出冷条初始温度越高降温时间越长，冷条初始温度为−35、−30、−25和−20 ℃所对应的降温时间依次为8、10、14和21 min，降温时间随着冷条初始温度的降低而减少。图6f为贮藏室在降温过程中温差变化情况，可以看出冷条初始温度越低温差越大，这是因为冷条温度越低，贮藏室空气温差越大，换热速率加快，从而导致温差越大。
5 结论通过对蓄冷式冷藏箱降温过程的数值模拟发现：贮藏室温度从16 ℃降低到0 ℃需要8 min，降温较快；正对回风道的区域温度较低，其他区域温度分布比较均匀。
经试验验证，模拟值与试验值随时间变化规律基本吻合，模拟温度和实测温度的绝对误差为0.68 ℃，温度分布平均绝对误差为0.29 ℃，说明所建立的模型比较可靠。
提高风速、增大回风道面积、降低冷条初始温度可以缩短降温的时间，贮藏室内温度变化速率随着降温时间逐渐减小；贮藏室内的温差随风速的增加而减小，随回风道面积的增大、冷条初始温度的降低而增大。采用comsol对蓄冷式冷藏箱降温过程中温度场的数值模拟具有一定的科学性，有助于认清温度场分布规律，为蓄冷式冷藏箱降温过程的参数优化提供参考。
| [1] |
陶文博. 冷板冷藏车的传热模拟及冷板布置方式优化[D]. 广州: 广州大学, 2016.
(  0) |
| [2] |
于娜, 于乐征, 朱景明. 峰谷分时电价下过程工业用户的用电响应模型[J]. 东北电力大学学报, 2014, 34(4): 17-21. DOI:10.3969/j.issn.1005-2992.2014.04.004 ( 0) |
| [3] |
廖明义, 金立芝. 蓄冷板冷藏汽车[J]. 制冷学报, 1987(3): 58-59. ( 0) |
| [4] |
COPERTARO B, PRINCIPI P, FIORETTI R. Thermal performance analysis of PCM in refrigerated container envelopes in the Italian context: Numerical modeling and validation[J]. Appl Therm Eng, 2016, 102: 873-881. DOI:10.1016/j.applthermaleng.2016.04.050 ( 0) |
| [5] |
翁卫兵, 房殿军, 李强, 等. 冷藏运输厢内流场和温度场协同控制[J]. 农业机械学报, 2014, 45(6): 260-265. ( 0) |
| [6] |
刘妍玲, 张岩, 王世清, 等. 果蔬摆放形式对冷库内气流场分布影响的研究[J]. 青岛农业大学学报(自然科学版), 2008(1): 24-27. DOI:10.3969/j.issn.1674-148X.2008.01.007 ( 0) |
| [7] |
张娅妮, 陈洁, 陈蕴光, 等. 机械式冷藏汽车厢体内部气流组织模拟研究[J]. 制冷空调与电力机械, 2007(2): 10-13. DOI:10.3969/j.issn.2095-3429.2007.02.003 ( 0) |
| [8] |
MOUREH J, MENIA N, FLICK D. Numerical and experimental study of airflow in a typical refrigerated truck configuration loaded with pallets[J]. Comput Electron Agr, 2002, 34(1): 25-42. ( 0) |
| [9] |
谢如鹤, 唐海洋, 陶文博, 等. 基于空载温度场模拟与试验的冷藏车冷板布置方式优选[J]. 农业工程学报, 2017, 33(24): 290-298. DOI:10.11975/j.issn.1002-6819.2017.24.038 ( 0) |
| [10] |
詹耀立, 冯国会. 不同冷板布置方式的冷板车内温度场模拟[J]. 沈阳建筑大学学报(自然科学版), 2008(5): 859-862. ( 0) |
| [11] |
张哲, 郭永刚, 田津津, 等. 冷板冷藏汽车箱体内温度场的数值模拟及试验[J]. 农业工程学报, 2013, 29(S1): 18-24. ( 0) |
| [12] |
郭志鹏, 阚安康, 孟闯, 等. 冷藏集装箱内温度场的数值模拟与实验[J]. 上海海事大学学报, 2017, 38(2): 82-87. ( 0) |
| [13] |
何凯强, 庄春龙, 杨小凤, 等. 送风相对湿度对冷藏集装箱内温湿度场的影响[J]. 兵器装备工程学报, 2017(10): 97-103. DOI:10.11809/scbgxb2017.10.020 ( 0) |
| [14] |
刘亚姣, 杨小凤, 庄春龙, 等. 机械式冷藏集装箱内温湿度场的数值模拟及其试验研究[J]. 保鲜与加工, 2016(1): 81-86. ( 0) |
| [15] |
刘杰坤, 陆华忠, 李君, 等. 嵌套真空板冷藏箱体保温性能的建模分析[J]. 西北农林科技大学学报(自然科学版), 2016, 44(8): 226-234. ( 0) |
| [16] |
吕传超, 崔振科, 王乐. 冰箱冷藏室温度场和流场的数值模拟与结构优化[J]. 制冷技术, 2015(5): 56-61. DOI:10.3969/j.issn.2095-4468.2015.05.204 ( 0) |
| [17] |
邓书辉, 施正香, 李保明. 低屋面横向通风牛舍温湿度场CFD模拟[J]. 农业工程学报, 2015, 31(9): 209-214. ( 0) |
| [18] |
李成祥, 陈建东, 李娜, 等. 冷藏箱温度场的数值模拟及其优化[J]. 制冷技术, 2015(3): 39-43. DOI:10.3969/j.issn.2095-4468.2015.03.109 ( 0) |
| [19] |
谢晶, 缪晨, 杜子峥, 等. 冷库空气幕性能数值模拟与参数优化[J]. 农业机械学报, 2014(7): 189-195. ( 0) |
| [20] |
郭嘉明, 吕恩利, 陆华忠, 等. 冷藏运输厢体结构对流场影响的数值模拟[J]. 农业工程学报, 2012, 28(S1): 74-80. ( 0) |
| [21] |
郭嘉明, 方思贞, 曾志雄, 等. 管道式加湿装置湿度场分布的数值模拟及试验验证[J]. 农业工程学报, 2015, 31(16): 57-64. DOI:10.11975/j.issn.1002-6819.2015.16.009 ( 0) |
| [22] |
王广海, 郭嘉明, 吕恩利, 等. 出风道参数对冷藏集装箱温度场的影响[J]. 农业机械学报, 2016(10): 293-301. DOI:10.6041/j.issn.1000-1298.2016.10.037 ( 0) |