需要JoVE订阅才能观看此内容。 请登录或开始免费试用

方法文章

水平冷表面上结霜的CFD模拟

583 次观看

⸱

DOI：

10.3791/68133

⸱

2025年9月12日

 ,  ,  ,  ,  , 
1School of Municipal and Environment Engineering, Shenyang Jianzhu University, 2Shenyang Polytechnic College, 3Laboratory of Thermo-fluid Science and Nuclear Engineering, School of Energy and Power Engineering, Northeast Electric Power University

本文内容

摘要

本文提出了一种基于欧拉多相流模型并结合Lee相变方法的数值模型，用于模拟水平低温表面的霜形成过程。该模型通过动态更新霜体积分数以捕捉密度变化，并利用实验数据对霜层厚度、密度及分布进行了验证。

摘要

结霜是一种在制冷、建筑和天然气加工等多个领域中普遍存在的现象。然而，由于其过程的复杂性，建立准确可靠的数值模型仍是一个重大挑战。尽管此前已有多项研究致力于解决这一问题，但现有模型仍存在一定的局限性。本文介绍了一种基于结霜基本机理改进的数值模型。该模型采用欧拉多相流方法，并结合Lee相变模型进行模拟。此外，对最大霜层体积分数的确定方法进行了更新，使模型能够考虑结霜过程中密度的变化。通过与来自不同研究的实验数据在厚度、密度及分布方面的对比，对该模型进行了严格验证。结果表明，霜层厚度的平均绝对相对偏差（MARD）为8.97%，密度的MARD为16.06%。此外，模型预测的霜层形貌与参考文献中报道的实验观测结果高度一致。

引言

结霜是一种在多个不同领域中普遍存在的现象。换热表面霜层的积聚会显著降低传热效率1，阻碍流体流动，并破坏换热器的整体性能2，最终影响其正常运行3。因此，理解结霜的机理与行为对于解决制冷系统中的这一问题至关重要4。近几十年来，已有大量研究致力于探究这些系统中结霜的原因及其特征。

实验研究表明，霜的形成受到多种因素的影响，包括空气温度、湿度以及冷表面的温度5,6,7,8,9,10。大量实验结果表明，较低的进气温度往往导致更厚的霜层6，而较高的湿度水平则有助于形成密度更高的霜层7。Song 等人研究了水平表面上的霜形成过程，发现冷表面的周期性温度变化可引起霜层界面处的融化，从而显著影响霜的形成速率、霜层厚度以及动态霜密度8。其他研究还探讨了霜的形貌与分布特征。Jeong 等人在实验中观察到，霜最初在入口附近形成，导致出现一种称为“霜丘”（frost hill）的现象9。Noorshams 等人研究了水平圆管表面的霜形成，发现圆柱体前表面和后表面的霜层比顶部表面更厚10。此外，已有若干研究11,12,13,14,15,16基于实验观测的霜形成规律，结合理论与经验方法，建立了用于预测一维霜层厚度的模型。Jones 和 Parker 基于分子扩散理论开发了一种霜厚度预测模型11，该模型在3小时内的预测结果与实验数据的偏差保持在30%以下。随着计算技术的发展，越来越多的研究人员转向使用计算流体动力学（Computational Fluid Dynamics, CFD）模拟霜的形成过程。与传统的一维模型相比，CFD 模拟能够显著提供更优的可视化效果，特别是在霜层厚度分布和温度场方面的呈现。Cui 等人基于成核理论开展了霜形成的 CFD 模拟12，其对霜厚度的预测结果与 Lenic 等人提供的实验数据偏差小于13%13。与此同时，针对结构化管内冷凝过程的 CFD 研究表明，诸如凹坑（dimples）14或螺旋节距（helical pitches）15等几何特征能够增强局部传热与传质性能。最近，You 等人16开发了一种基于动态网格的 CFD 模型，将霜层表征为不断生长的多孔介质，并直接纳入水蒸气扩散过程，在保持较低计算成本的同时实现了小于5%的相对偏差。这些研究成果凸显了 CFD 在解决复杂相变现象方面的潜力，为霜形成建模提供了有价值的科学见解。

综上所述，已有大量研究5,6,7,8,9,10,11,12,13探讨了冷表面上的结霜现象，推动了在不同参数条件下对结霜模式的逐步深入理解。尽管已开发出多种数值模型，并纳入了不同维度和机制，但这些模型往往缺乏全面的实验验证。大多数研究6,7,8,9,10,11,12,13主要通过霜层厚度对模型进行验证，限制了其更广泛的适用性16。为克服这些局限性，本文提出了一种融合Lee相变模型与欧拉多相流模型的数值模型，重点聚焦于结霜过程中的核心机制。此外，本文还引入了一种新颖的方法，用于计算霜体积分数的上限，该方法考虑了霜密度随时间变化的特性，从而弥补了以往模型的不足。所提出模型的准确性与可靠性从多个方面进行了评估，包括霜层厚度、密度波动以及在不同实验条件下的结霜形态。这一全面的验证为更准确地预测实际应用中的结霜行为提供了坚实的理论框架。

访问受限。请登录或开始试用以查看此内容。

方案

1. 物理模型与网格

  1. 打开 SpaceClaim，选择 草图 选项卡，然后在创建功能下选择 矩形 选项。
  2. 在 XOY 平面上创建一个二维几何模型，沿 x 轴方向长度为 500 mm，沿 y 轴方向宽度为 15 mm。
  3. 打开 ICEM，进入文件选项卡，选择 几何 选项卡，然后选择 打开几何 并导入二维几何模型。
  4. 选择 二维模型的边界，在部件选项卡下打开 创建部件 功能，并为不同边界分配名称。
  5. 选择 分块选项卡，然后点击 创建分块。在创建分块窗口中，勾选 继承部件名称选项。然后进入创建分块选项卡并选择 初始化分块。
  6. 在初始化分块选项卡中，选择 二维面分块作为类型，并勾选 使用设置初始化。
  7. 在选择选项卡下，方法选项选择 生成面。然后在面选项下，点击 选择面(s)图标，并在图形窗口中选择 二维平面 。
  8. 在面分块选项卡下，方法选项选择 主要映射，自由面网格类型选择 全四边形，自由面网格方法选择 ICEM CFD 四边形。
  9. 在跨曲线合并分块选项卡下，方法选择 全部，忽略尺寸设置为 0.0。
  10. 在预网格参数选项卡中，在网格参数下点击 边参数。在边选项下，点击 选择边(s)图标，并在图形界面中选择 X 方向的边。
  11. 在网格参数下，将节点数设为 500，间距1设为 1e+10，比率1设为 2，间距2设为 1e+10，比率2设为 2，最大间距设为 1e+10。然后勾选 复制参数选项。
  12. 在网格参数选项下，将网格规律设为双几何（BiGeometric），然后点击 应用。在预网格参数选项卡中，在网格参数选项下点击 边参数。
  13. 在网格参数下的边选项中，点击 选择边(s)图标，并在图形界面中选择 y 方向的边。
  14. 在网格参数下，将节点数设为 150，间距1设为 1e+10，比率1设为 2，间距2设为 1e+10，比率2设为 2，最大间距设为 1e+10。然后勾选 复制参数选项。
  15. 在网格参数选项下，将网格规律设为双几何（BiGeometric），然后点击应用。在分块选项下，激活预网格功能，并在提示时点击 是 。
  16. 右键点击 预网格，然后从上下文菜单中选择 转换为非结构化网格 。
  17. 进入输出网格选项卡并点击 求解器设置。在求解器设置中，选择 ANSYS Fluent作为求解器，然后点击 应用。
  18. 在输出选项卡下，点击 写入输入。在新窗口1中，点击 保存；在新窗口2中，点击 是，然后点击 保存；在新窗口3中，点击 打开；在新窗口4中，点击 完成。

2. 结霜形成模拟软件操作

  1. 打开 Ansys Fluent。转到文件选项卡，然后在读取（Read）下选择 Mesh 选项。进入网格缩放（Scale Mesh）界面，将“Mesh Was Created In”设置为 mm 。此处介绍的流程基于表1中案例1所定义的实验条件。
  2. 在求解器设置中，选择压力基（Pressure-Based）类型、绝对速度（Absolute Velocity）公式化方法和瞬态时间（Transient Time）。
  3. 将Y方向的重力加速度设置为-9.81。在模型（Models）中，点击能量（Energy）并启用能量方程（Energy Equation）。
  4. 在模型（Models）中启用黏性（Viscous），并选择k-epsilon (2 eqn)。在k-epsilon模型部分，选择标准（Standard）。
  5. 对于近壁面处理（Near-Wall Treatment），选择标准壁面函数（Standard Wall Functions）。在湍流多相模型（Turbulence Multiphase Model）中选择混合（Mixture）。
  6. 在模型常数（Model Constants）选项中，将Cmu、C1-Epsilon、C2-Epsilon、TKE普朗特数、TDR普朗特数、弥散普朗特数、能量普朗特数、壁面普朗特数和湍流施密特数分别设置为0.09、1.44、1.92、1、1.3、0.75、0.85、0.85和0.7。
  7. 在用户自定义函数（User-Defined Functions）中，将湍流黏度混合物、相1和相2设置为无（none）。在模型（Models）中点击组分（Species）并启用组分输运（Species Transport）。
  8. 在组分输运模型窗口中，勾选扩散能量源（Diffusion Energy Source），并选择相1（phase-1）。
  9. 在模型（Models）中启用多相（Multiphase），并选择欧拉（Eulerian）模型。在多相模型窗口中，选择体积分数参数公式化隐式（Volume Fraction Parameters Formulation Implicit），并将欧拉相数量（Number of Eulerian Phases）设置为2。
  10. 在多相模型窗口中选择相（Phases）选项卡，选择相1 - 主相（phase-1 - Primary Phase），将名称（Name）设为相1（phase-1），并将相材料（Phase Material）设为mixture-template。
  11. 在多相模型窗口中选择相（Phases）选项卡，选择相2 - 次相（phase-2 - Secondary Phase），将名称（Name）设为相2（phase-2），将相材料（Phase Material）设为ice，并启用颗粒（Granular）。
  12. 在相选项卡的相设置窗口中，选择相属性（Phase Property）为颗粒温度模型（Granular Temperature Model），将直径（Diameter）设为0.0001。
  13. 在相选项卡的颗粒属性窗口中，将颗粒黏度设为1e-05，颗粒体积黏度设为30，固体压力设为lun-et-al，颗粒温度设为代数（algebraic），摩擦黏度设为无（none），堆积极限设为用户自定义（user-defined），径向分布设为lun-et-al，弹性模量设为推导（derived）。
  14. 在多相模型窗口的相相互作用（Phases Interaction）选项卡中，选择力（Forces）标签，然后选择相1、相2（phase-1, phase-2），将系数设为wen-yu。
  15. 在多相模型窗口的相相互作用（Phases Interaction）选项卡中，选择力（Forces）标签，然后选择相2、相2（phase-2, phase-2），将恢复系数设为0.9。
  16. 在多相模型窗口的相相互作用（Phases Interaction）选项卡中，选择界面面积（Interfacial Area）选项卡，然后选择ia-symmetric选项。
  17. 在边界条件（Boundary Conditions）中，点击入口（Inlet）并选择相1选项卡（phase-1 tab），将相选项设为相1。
  18. 在相1选项卡的速度入口（Velocity Inlet）窗口中，选择动量（Momentum）选项卡，将速度指定方法设为“大小，垂直于边界（Magnitude, Normal to Boundary）”，参考系设为绝对（Absolute），速度大小设为0.6。
  19. 在相1的速度入口窗口中，选择热（Thermal）标签，将温度设为292.8。在相1标签的速度入口窗口中，选择组分（Species）标签，将h2o设为0.008202。
  20. 在边界条件中，点击入口（Inlet）并选择相2（phase-2）标签，将相选项设为相2。
  21. 在相2标签的速度入口窗口中，选择动量（Momentum）标签，将速度指定方法设为“大小，垂直于边界（Magnitude, Normal to Boundary）”，参考系设为绝对（Absolute），速度大小设为0，颗粒温度设为0.0001。
  22. 在相2标签的速度入口窗口中，选择热（Thermal）选项卡，将温度设为273。在相2标签的速度入口窗口中，选择多相（Multiphase）选项卡，将体积分数设为0。
  23. 在边界条件中，点击出口（Outlet）并选择相1（phase-1）标签，将相选项设为相1。
  24. 在相1选项卡的压力出口（Pressure Outlet）窗口中，选择热（Thermal）选项卡，将回流总温（Backflow Total Temperature）设为300。
  25. 在相1选项卡的压力出口窗口中，选择组分（Species）选项卡，将h2o设为0。在边界条件中，点击出口（Outlet）并选择相2选项卡，将相选项设为相2。
  26. 在相2选项卡的压力出口窗口中，选择热（Thermal）选项卡，将回流总温设为300。
  27. 在相2选项卡的压力出口窗口中，选择多相（Multiphase）选项卡，将回流颗粒温度设为0.0001，将体积分数指定方法设为回流体积分数（Backflow Volume Fraction），并将回流体积分数设为0。
  28. 2.28 在边界条件中，点击壁面（Wall）并选择cold-wall选项卡，将相选项设为混合物（mixture）。
  29. 在cold-wall选项卡的壁面窗口中，选择动量（Momentum）选项卡，将壁面运动选项设为静止壁面（Stationary Wall），将壁面粗糙度模型选项设为标准（Standard），粗糙度高度设为0，粗糙度常数设为0.5。
  30. 在cold-wall选项卡的壁面窗口中，选择热（Thermal）选项卡，选择温度作为热条件（Thermal Conditions），将温度设为252.65，材料设为steel。
  31. 在求解（Solution）中，点击方法（Methods）并打开求解方法窗口。在求解方法窗口中，选择相耦合SIMPLE作为压力-速度耦合方案（Pressure-Velocity Coupling Scheme），选择最小二乘单元基作为梯度空间离散化（Gradient Spatial Discretization），选择二阶作为压力空间离散化（Pressure Spatial Discretization），选择一阶迎风作为密度空间离散化（Density Spatial Discretization），选择一阶迎风作为动量空间离散化（Momentum Spatial Discretization），选择一阶迎风作为体积分数空间离散化（Volume Fraction Spatial Discretization），选择一阶迎风作为湍流动能空间离散化（Turbulent Kinetic Energy Spatial Discretization），选择一阶迎风作为湍流耗散率空间离散化（Turbulent Dissipation Rate Spatial Discretization），选择一阶迎风作为能量空间离散化（Energy Spatial Discretization），选择一阶迎风作为相1 h2o空间离散化（phase-1 h2o Spatial Discretization），选择一阶隐式作为瞬态公式化（Transient Formulation）。
  32. 在求解中，点击控制（Controls）并打开求解控制窗口。将压力欠松弛因子设为0.4，密度欠松弛因子设为1，体积力欠松弛因子设为1，动量欠松弛因子设为0.4，体积分数欠松弛因子设为0.4，颗粒温度欠松弛因子设为0.3，湍流动能欠松弛因子设为0.3，湍流耗散率欠松弛因子设为0.3，湍流黏度欠松弛因子设为0.3，能量欠松弛因子设为0.4，相1 h2o欠松弛因子设为0.4。
  33. 在求解中，点击初始化（Initialization）并打开求解初始化窗口。选择标准初始化（Standard Initialization）作为初始化方法，选择相对于单元区域（Relative to Cell Zone）作为参考系。
  34. 在求解初始化窗口中，将表压（Gauge Pressure）设为0，将湍流动能设为0.00135，将湍流耗散率设为0.001143987，将相1 X速度设为0，相1 Y速度设为0，相1 h2o设为0.008202，相2 X速度设为0，相2 Y速度设为0，相2体积分数设为0，相2颗粒温度设为0.0001，相2温度设为273。然后点击初始化（Initialize）。

3. 后处理与数据导出配置

  1. 在“结果”下，点击等值线以打开“等值线”窗口。在“等值线”窗口中，启用填充、节点值、边界值、全局范围和自动范围选项。将“等值线类型”选择为相位，并选择体积分数。然后选择相位2作为相位，并点击保存/显示。
  2. 在“计算任务”下，点击解动画以打开“动画定义”窗口。
  3. 在“动画定义”窗口中，将“每……记录一次”设置为1，并选择时间步长。将“存储类型”选择为HSF文件。在“动画对象”选项下选择等值线-1，然后点击确定。
  4. 在“结果”下，选择表面选项，然后点击新建线/阵列面以打开“线/阵列面”窗口。
  5. 在“线/阵列面”窗口中，启用线，将x0 [m]设置为0.21，x1 [m]设置为0.21，y0 [m]设置为0，y1 [m]设置为0.015，然后点击创建。
  6. 在“文件”选项卡下，在“导出”中选择计算过程中选项。然后点击解数据以打开“自动导出”窗口。
  7. 在“自动导出”窗口中，将“文件类型”选择为ASCII，将“位置”选择为单元中心，将“分隔符”选择为空格。将“导出数据间隔”设置为1，并选择时间步长。在“表面”选项下选择线-1。在“物理量”选项下，选择密度（相位-2）和体积分数（相位-2）。然后点击浏览，打开“选择文件”窗口。在“选择文件”窗口中点击确定。在“自动导出”窗口中点击确定。
  8. 在“计算任务”下，点击自动保存（每流体时间）以打开“自动保存”窗口。将“每隔多少秒保存一次数据文件”设置为100，选择流体时间，将“关联案例文件保存类型”选择为仅当修改时保存，然后点击确定。
  9. 在“求解”下，点击运行计算以打开“运行计算”窗口。然后将“时间推进类型”选择为固定，将“时间推进方法”选择为用户指定。
  10. 在“运行计算”窗口中，将“时间步数”设置为7200，将“时间步长”设置为1，将“每时间步最大迭代次数”设置为20，将“报告间隔”设置为1，将“轮廓更新间隔”设置为2。然后点击计算。

访问受限。请登录或开始试用以查看此内容。

结果

所提出的改进数值模型能够有效捕捉霜层形成的关键特征。该模型基于霜层生长的基本机理，采用欧拉多相流方法，并结合Lee相变模型。该方法使模型能够动态更新最大霜体积分数，从而考虑结霜过程中密度的变化。模拟结果表明，通过将模型对霜层厚度、密度和分布的预测结果与多项研究的实验数据进行对比，模型得到了严格验证。此外，模型预测的霜层形貌与文献报道的实验观察结果高度一致。

图1 展示了模拟物理模型的示意图。计算区域前端包含一个充分发展的入口区域，上壁面设置为绝热边界。左侧边界定义为速度入口，右侧边界作为压力出口。该设置确保了在对流流动条件下霜形成过程具有真实的边界条件。

图2 展示了在3600秒时模拟的霜层体积分数分布。模拟条件总结于表1中，其中壁面温度为252.65 K，空气温度为292.95 K，入口流速为0.6 m/s，相对湿度为58%。最大体积分数达到0.07。霜层分布均匀，中心区域...

访问受限。请登录或开始试用以查看此内容。

讨论

本研究开发了一种数值模型，通过根据时间和运行条件动态调整体积分数的上限，能够模拟低温水平冷表面上的霜层形成过程，从而再现霜密度的变化。尽管本文所展示的验证仅限于传统的冷表面温度条件，但该模型确保了模拟与实验测得的霜层厚度之间的最大偏差不超过-20%（平均绝对相对偏差 MARD = 8.97%），霜密度的最大偏差在29.72%以内（MARD = 16.06%）。此外，模拟得到的霜层分布与实验观测结果高度一致。这些结果对于构建高精度数值霜层模型具有重要意义。

尽管文献中已提出多种数值霜冻模型11,12,13,14,15,16，但本研究通过验证更广泛的参数范围，并实现模拟与实验之间更小的偏差，展现...

访问受限。请登录或开始试用以查看此内容。

披露

作者声明，他们不存在任何已知的相互竞争的经济利益或可能被认为影响本文所报告工作的个人关系。

致谢

本研究由 (XLYC2203184)、(U23A20657) 和 (LJ222410153082) 资助。

访问受限。请登录或开始试用以查看此内容。

材料

本文使用的材料清单
姓名公司目录编号评论
FluentANSYS
ICEMDassault system
SpaceClaimANSYS

参考文献

  1. Saygin, A., Basol, A. M., Arik, M. An experimental study on the frost formation over a flat plate: Effect of frosting on heat transfer. Exp Therm Fluid Sci. 144, 110862(2023).
  2. Fang, X., et al. A new frictional pressure drop correlation based on flow patterns for hydrocarbon refrigerants condensation flow. Int J Refrig. 170, 214-223 (2025).
  3. Rong, X., et al. Experimental study on a multi-evaporator mutual defrosting system for air source heat pumps. Appl Energy. 332, 120528(2023).
  4. Jia, Y., Xu, X., Li, Y., Liang, X., Yao, M. Experimental studies on frost and defrost of fine tube bundles under coolant temperature between −20 and −5 °C. Int J Heat Mass Transf. 116, 617-620 (2018).
  5. Song, M., Dang, C. Review on the measurement and calculation of frost characteristics. Int. J. Heat Mass Transf. 124, 586-614 (2018).
  6. Lee, J., Lee, K. -S. The behavior of frost layer growth under conditions favorable for desublimation. Int J Heat Mass Transf. 120, 259-266 (2018).
  7. Lee, Y. B., Ro, S. T. Frost formation on a vertical plate in simultaneously developing flow. Exp Therm Fluid Sci. 26 (8), 939-945 (2002).
  8. Mengjie, S., Shangwen, L., Hosseini, S. H., Xiaoyan, L., Zhihua, W. An experimental study on the effect of horizontal cold plate surface temperature on frosting characteristics under natural convection. Appl Therm Eng. 211, 118416(2022).
  9. Jeong, H., Byun, S., Kim, D. R., Lee, K. S. Frost growth mechanism and its behavior under ultra-low temperature conditions. Int J Heat Mass Transf. 169, 120941(2021).
  10. Barzanoni, Y., Noorshams, O., Basirat Tabrizi, H., Damangir, E. Experimental investigation of frost formation on a horizontal cold cylinder under cross flow. Int J Refrig. 34 (4), 1174-1180 (2011).
  11. Jones, B. W., Parker, J. D. Frost formation with varying environmental parameters. J Heat Transf. 97 (2), 255-259 (1975).
  12. Cui, J., Li, W. Z., Liu, Y., Jiang, Z. Y. A new time- and space-dependent model for predicting frost formation. Appl Therm Eng. 31 (4), 447-457 (2011).
  13. Lenic, K., Trp, A., Frankovic, B. Transient two-dimensional model of frost formation on a fin-and-tube heat exchanger. Int J Heat Mass Transf. 52 (1-2), 22-32 (2009).
  14. Yu, J., Huo, R., Shen, H., Li, X., Zhu, Z. A simulation study on the condensation flow and thermal control characteristics of mixed refrigerant in a dimpled tube. Appl Therm Eng. 231, 120889(2023).
  15. Yu, J., Jiang, Y., Cai, W., Li, X., Zhu, Z. Condensation flow patterns and heat transfer correction for zeotropic hydrocarbon mixtures in a helically coiled tube. Int J Heat Mass Transf. 143, 118500(2019).
  16. You, Y., Wang, S., Lv, W., Chen, Y., Gross, U. A CFD model of frost formation based on dynamic meshes technique via secondary development of ANSYS Fluent. Int J Heat Fluid Flow. 89, 108807(2021).
  17. Cai, W., Fang, X., Li, S., Qiu, G. A modified CFD model for frosting on a horizontal plate. Int J Heat Mass Transf. 229, 125726(2024).
  18. Boyina, K. S., et al. Condensation frosting on meter-scale superhydrophobic and superhydrophilic heat exchangers. Int J Heat Mass Transf. 145, 118694(2019).

访问受限。请登录或开始试用以查看此内容。

重印与许可

标签

Lee
联系我们
向图书馆推荐

隐私

使用条款

政策

版权所有 © 2026 MyJoVE Corporation。保留所有权利。