平板结冰特性的数值研究及其影响因素

《CMES - Computer Modeling in Engineering and Sciences》：Numerical Study on the Icing Characteristics of Flat Plates and Its Influencing Factors

本研究通过CFD与FENSAP-ICE耦合数值模拟，探究薄板在低温高湿环境中的结冰行为及其气动特性变化。结果表明：冰层厚度随时间线性增长，温度影响结冰形态（融雪冰、硬融雪冰、冰层冰）和分布；风速显著促进冰层厚度与覆盖范围（临界风速低于1m/s时冰 accretion可忽略）；攻角增大使结冰区域向压力面偏移，对称攻角下冰层厚度一致。结冰后薄板前端更平滑，阻力系数略降，升力系数和力矩系数均上升。研究为航空、风电及桥梁等工程结构的冰防护设计提供理论支撑。

　　在寒冷和潮湿的环境中，结构表面如飞机机翼、风力发电机叶片等易发生结冰现象，这对空气动力学性能和运行安全构成了严重威胁。本文研究了在不同环境参数下，薄平板上结冰行为的变化，探索了空气温度、风速和迎角（Angle of Attack, AOA）对结冰形态及空气动力学性能的影响。研究发现，结冰厚度随时间呈线性增长，而低温下会形成雾凇冰（Rime Ice），在结构表面迅速冻结，导致冰层形成；而较高温度下，水滴会在表面形成一层水膜，随后重新冻结形成尖锐的冰锥（Glaze Ice）。在相同条件下，雾凇冰通常会比冰锥形成更厚的冰层。风速的增加显著促进冰层的形成和覆盖，而风速低至1 m/s时，结冰量却相对较少。迎角的变化会导致结冰区域向压力面迁移，相同大小但符号相反的迎角则会形成对称的结冰模式，并产生相同的最大冰层厚度。结冰后，平板的前缘变得更加平滑，轻微降低阻力，同时增加升力和力矩系数。这些发现强调了温度、风速和迎角在决定结冰形态、范围和空气动力学影响中的主导作用，为预测结冰效应和开发防冰策略提供了宝贵的见解。

### 1. 引言

中国地域辽阔，气候条件多样且差异显著。在远北地区，黑龙江省的冬季气温常低于-15°C，而在西部的青藏高原则以显著的日间温度变化和复杂的地形特征著称。这些条件导致了道路结冰，成为重要的交通安全隐患，已被中国气象局列为早期预警系统的重要研究对象。道路结冰会显著降低轮胎摩擦力，影响车辆操控性和制动性能，进而增加事故风险。一个显著的例子是2009年12月28日，杭瑞高速公路鄱阳湖大桥发生的多起碰撞事故，这起事件中低温和浓雾导致了冰层的形成，对生命和财产安全造成了威胁。除了道路基础设施，结冰问题也延伸到其他关键的工程应用中，如航空领域。随着中国航天工业的快速发展，飞行过程中飞机表面的结冰已成为严重的安全问题。结冰会干扰气流的平滑流动，增加阻力，提高失速速度临界值，并降低控制效率，最终威胁飞行安全。在风能领域，极端天气条件，尤其是重庆、云南和湖南等内陆地区，导致了风力发电机叶片上的严重结冰。这不仅减少了发电能力和年均能量产出，还降低了空气动力学性能，甚至在某些情况下导致风力涡轮机停机。芬兰某风力发电场的数据表明，在1996至2002年间，结冰导致了70%的停机时间。因此，结冰问题在多个工程领域都具有极高的重要性。

作为基础空气动力学形状，平板广泛应用于各种工程结构中，如船舶、大型风力发电机叶片、桥梁设计等，结冰行为对这些结构提出了独特的挑战。研究这些挑战对于提高在寒冷和潮湿环境中运行的结构的安全性、优化能效以及推动工程应用至关重要。

### 2. 结冰机制与数值方法

结冰是一种自然界的常见物理现象。在未受扰动的条件下，水滴在-40°C时会自发冻结成冰。当温度在-40°C至0°C之间时，水滴进入一个亚稳态，即过冷水滴。在此状态下，水滴可以在较长时间内保持液态，除非受到外部扰动，如接触冷表面、机械冲击或振动，才会冻结并形成冰层。在自然环境中，低温和高湿度会导致过冷水滴随湿润空气流动并撞击结构表面，从而形成不同程度的结冰。这个过程复杂且依赖于环境因素。整体的结冰过程如图1所示，图中展示了特定条件下冰层形成的示意图。在本研究中，重点研究的是薄平板的结冰现象，这些平板是寒冷和潮湿环境中受到结冰影响的结构的代表性模型。所研究的结冰现象被归类为大气结冰，可以根据两种主要的形成过程进行分类：降雪结冰（降水结冰）和云中结冰（in-cloud icing）。根据国际标准ISO 12494:2017，降雪结冰可以进一步分为冰锥（glaze ice）和湿雪（wet snow），其中冰锥的空气温度在-10°C至0°C之间，而湿雪的空气温度在0°C至3°C之间。相比之下，云中结冰主要受风速和温度影响，进一步分为软雾凇（soft rime）、硬雾凇（hard rime）和冰锥。这些不同的结冰类型如图2所示，展示了环境条件如何决定具体的结冰形式。在本研究中，关注的温度范围低于0°C，特别是云中结冰，尤其是软雾凇、硬雾凇以及冰锥的形成。

为了模拟平板上的结冰行为，本文提出了一个数值框架，包括四个主要步骤，如图3所示。在第一步中，使用前处理软件ICEM CFD生成结构化网格，这些网格作为后续计算的背景网格。第二步使用有限元软件ANSYS FLUENT 2022 R1计算平板周围的流场，应用适当的边界条件和求解器设置。第三步利用FENSAP-ICE中的DROP3D模块，将第二步计算的流场导入，以确定水滴轨迹及其对平板表面的撞击分布。第四步使用“准稳态”方法模拟结冰过程。该方法在小时间间隔内求解稳态控制方程，将前一个时间步（t）的解作为下一个时间步（t + Δt）的初始条件，从而确保更准确地模拟冰层的生长，这与实际结冰现象紧密相关。

### 3. 验证

为验证提出的数值模拟框架，研究结果与Meng等人的实验数据进行了对比。他们的实验主要研究了在不同气象条件下，静态圆柱模型上的雾凇和冰锥结冰现象，实验在定制的、小型的喷雾式结冰风洞中进行。圆柱模型的直径为25.4 mm。在验证中，研究使用了他们的第一组实验条件，即雾凇结冰的条件。关键参数包括风速为60 m/s，迎角为0°，空气温度为-15°C，空气密度为1.225 kg/m3，液态水含量为0.5 g/m3，水滴直径为20 μm，结冰时间为200 s，如表1所示。需要注意的是，在-15°C的温度下，空气密度对结冰模拟的影响可以忽略不计。因此，在本研究中，空气密度在整个数值模拟中被设定为恒定值1.225 g/cm3。计算域和模型尺寸如图4a所示。为了确保高网格质量和计算精度，模型在ICEM CFD中生成了75层网格。第一层网格高度为0.2 mm，网格增长率为1.05。图4a展示了网格设置和计算域尺寸。对比结果如图4b所示，图中x轴表示圆柱表面任意点与负x轴之间的角度，正角度表示顺时针方向。结果表明，在相同的趋势下，数值模拟与实验数据的偏差最大为4.71%，这证实了本文所采用的数值模拟方法在结冰模拟中的可靠性和准确性。

### 4. 数值模拟

#### 4.1 CFD模拟

在验证了所提出的数值框架后，将其应用于薄平板上的结冰模拟。本文中，平板的几何形状定义为宽度B为0.4 m，高度H为0.01 m，宽高比为40。这种宽高比符合工程实践中平板的定义。为了减少前缘的流动分离并确保平板的最优特性，模型的前缘和后缘均被圆滑处理，半径为10 mm。计算域的尺寸为15B×10B，阻塞比为0.25%。由于FENSAP-ICE软件要求使用三维模型，为了减少跨度方向对结果的影响，本文建议将跨度长度扩展为平板宽度。边界条件如图4a所示，左边界和上、下边界被定义为速度入口，右边界为压力出口，前缘和后缘为对称边界，平板部分被建模为固定壁面。

#### 4.2 结果与讨论

本节探讨了关键环境因素，包括空气温度、风速和迎角（AOA）对结冰的影响。此外，还分析了不同结冰情景对静态空气动力学系数的影响。

##### 4.2.1 空气温度对结冰的影响

如前所述，结冰类型、风速和温度之间的关系表明，在恒定风速下，温度增加时结冰类型从软雾凇过渡到硬雾凇，最终转变为冰锥。为了更好地理解温度对结冰形态的影响，并区分这三种结冰类型，本文选择了三种温度条件：-15°C、-8°C和-4°C，分别代表软雾凇、硬雾凇和冰锥结冰。风速保持为20 m/s。首先，通过CFD模拟获取平板周围的流场。然后，将这些数据导入FENSAP-ICE，设置液态水含量、水滴直径和结冰时间，以计算结冰。如表3所示，具体条件包括空气温度为-15°C，液态水含量为0.9 g/m3，水滴直径为20 μm，结冰时间为200 s。

图7展示了在三种不同温度条件下，平板上软雾凇、硬雾凇和冰锥的结冰模式。在这些图中，x轴和y轴分别表示平板宽度和高度的位置，如图5a所示。这些坐标系在所有后续图中被一致使用。结果显示，对于所有三种结冰类型，结冰主要集中在平板的前缘，且随时间线性增加。在雾凇条件下，环境温度较低，水滴在撞击平板表面时立即冻结，导致结冰形态的流线型增长。相比之下，冰锥结冰在较高温度下形成，水滴撞击表面后不会立即冻结，而是形成一层薄水膜，随后向前缘流动并重新冻结，形成突出的冰锥。在相同结冰时间下，雾凇结冰的厚度大于冰锥结冰。值得注意的是，硬雾凇在三种类型中形成最厚的冰层，最大厚度出现在前缘，而冰锥的最大厚度出现在冰锥的尖端。这些观察结果证实了温度在通过热传递过程控制结冰形态中的主导作用。

##### 4.2.2 风速和迎角（AOA）对结冰的影响

在现实条件下，风速和迎角（AOA）可能显著变化，因此理解它们对结冰的影响至关重要。为了展示这一点，研究了风速和迎角对雾凇结冰的影响。具体而言，考虑了四种风速（1、5、10、20 m/s）和五个不同的迎角（-10°、-5°、0°、5°、10°），如表5所示。所有其他参数均符合软雾凇结冰条件，即温度为-15°C，液态水含量为0.9 g/m3，水滴直径为20 μm，结冰时间为200 s。

图8展示了在不同风速条件下，平板上结冰模式的变化，迎角分别为0°、5°、-5°、10°和-10°。结果表明，随着风速的增加，水滴更有可能撞击平板的前缘。具体而言，在相同的结冰时间下，较高的风速会导致更多的水滴撞击表面。此外，增加的风速有助于带走冻结过程中释放的潜热，从而加速结冰速率并增加冰层的厚度和覆盖范围。尽管这些变化，风速对冰层整体形状的影响并不显著。例如，当风速从U=1 m/s增加到5 m/s、10 m/s和20 m/s时，在迎角为0°的情况下，平板上的冰层厚度迅速从0.02 mm增加到1.03 mm、2.51 mm和3.51 mm。这表明风速对冰层厚度有显著影响。表6列出了不同风速和迎角组合下的最大冰层厚度。值得注意的是，在正负迎角相同的情况下，平板上的结冰模式在水平轴上对称分布，产生相同的最大冰层厚度值。

##### 4.2.3 软雾凇对空气动力学系数的影响

平板表面结冰会改变其形状，从而影响周围流场并影响空气动力学系数。这些系数对于评估结构的风阻性能至关重要。例如，空气动力学系数在颤振分析中是识别颤振导数的基础。在桥梁风工程领域，这些系数用于计算作用在桥梁结构上的阻力、升力和力矩力，从而进一步分析结构在稳定风荷载下的应力和位移，并评估静态风稳定性。因此，研究结冰前后空气动力学系数的变化对于全面理解结冰对结构整体性能的影响至关重要。

图10展示了平板在形成软雾凇结冰（参考表3中的案例1）前后的空气动力学系数。结果显示，在所选的软雾凇结冰条件下，平板前缘的形状变化不大，导致结冰前后的空气动力学系数相对相似。对于阻力系数，结冰后观察到轻微的降低。这归因于前缘的平滑和流线型变化，如图7a所示的结冰分布模式。值得注意的是，在迎角为-5°时观察到负的阻力系数。这种异常是由于在该迎角下，较高的风速导致前缘形成小的冰锥，改变了流场，从而产生负的阻力系数。相反，升力系数和力矩系数都显示出轻微的增加，且在迎角增加时，这两个系数都会增加。

### 5. 结论

本研究开发并应用了一个集成的CFD-FENSAP-ICE数值框架，用于模拟平板上的结冰现象，并量化环境参数的影响，评估空气动力学影响。研究主要关注在不同温度下结冰模式和厚度的变化，分析风速和迎角对软雾凇结冰的影响，以及研究雾凇结冰对结构性能的空气动力学后果。主要结论如下：

1. 平板上的结冰厚度随时间呈线性增长。在雾凇结冰情况下，相对较低的环境温度导致水滴在撞击时立即冻结，形成流线型的结冰。相比之下，冰锥结冰在较高温度下形成，水滴在表面形成一层薄水膜，随后向前缘流动并重新冻结，形成突出的冰锥。在相同条件下，雾凇结冰通常会比冰锥形成更厚的冰层，最大厚度出现在前缘，而冰锥的最大厚度出现在冰锥尖端。这些观察结果证实了温度在通过热传递过程控制结冰形态中的主导作用。

2. 风速的增加导致更多的水滴撞击表面，加速结冰速率，增加冰层的厚度和覆盖范围。在非常低的风速（如1 m/s）下，由于水滴撞击量较少，结冰量在测试条件下几乎可以忽略不计，这表明存在一个显著结冰的风速阈值。

3. 迎角的增加扩大了压力面的撞击区域，而减少了吸力面的撞击区域，使得结冰区域向压力面迁移。对于相同大小但符号相反的迎角，结冰模式在水平轴上对称分布，产生相同的最大冰层厚度值。

4. 结冰改变了前缘的几何形态，使其更加平滑和流线型。这略微降低了阻力系数，与结冰前相比。在-5°迎角和较高风速的情况下，局部冰锥的形成改变了流场，导致阻力系数为负。升力和力矩系数则显示出轻微的增加，且在迎角增加时，这两个系数都会增加。

本研究探讨了平板上的结冰现象及其对结构性能的影响，为后续研究提供了基础见解和可行性证据。其局限性在于专注于特定计算条件下的基本结冰情景。未来的工作将扩展该框架，应用于相关的工程系统，包括飞机（机翼/机翼）、风力涡轮机（旋转叶片）和桥梁（缆索/箱型梁）等代表性几何形状。此外，将纳入完全时间依赖的结冰模拟，以捕捉非稳态结冰和主动/运行中的除冰情景，通过动态网格更新来解析结冰生长和脱落。还将评估其他湍流模型，并进行针对性的敏感性和不确定性分析。最后，将通过在符合相似性定律的条件下校准和验证预测结果，以提高框架的通用性和实际应用价值。

### 资助声明

本研究得到了国家自然科学基金（52278532）和四川省科技计划（2024NSFSC0153）的支持，特此致谢。

### 作者贡献

Jin Zhu：概念化、监督、资助获取、审阅与编辑。Yanxin Xu：形式分析、调查、方法论、初稿撰写。所有作者均审阅了研究结果，并批准了手稿的最终版本。

### 数据与材料的可获取性

支持本研究的所有数据、模型或代码均可在合理请求下获得。

### 伦理审批

本研究不涉及伦理审批。

### 兴趣冲突

作者声明在本研究中没有需要报告的利益冲突。