随着能源需求的与日俱增，我国急需进行能源体系转型，并加快发展现代化能源产业，坚持节约能源、保护环境和快速发展一体化。风力发电作为风能利用的主要形式，对于推动能源结构优化、减轻环境压力、应对能源枯竭以及实现可持续发展目标具有深远影响。风力机是实现风力发电的主要设备。风力机通过风力驱动叶片旋转，经过一系列机械装置带动发电机发电，最终将风能转化为电能^［1］。我国的风力发电产业已从初期的探索阶段，迈入商业应用阶段，并在全球市场中占据显著地位。截至2023年，我国风电累计装机已突破4亿kW，同比增长超15%，其中陆上风电装机容量占90%以上，同时风电平均利用率达到97%。我国风电机组占全球装机总量超过 40 %，风力机制造水平已达到国际领先水平。尽管我国风力发电装机容量不断增长，但一些自然环境因素如低温结冰会对风力机的运行状况产生重大影响。在除冰方面现已有较为成熟的方法，包括：热力法除冰，即通过加热使冰融化；机械法除冰，即通过外界机械力使冰破碎。但在结冰方面还有很大的研究空间，因此进行风力机叶片结冰研究是解决风力机结冰及防除冰问题的基础。

1 研究背景 1.1 风力机结冰现象

在欧美等风力发电发展较为成熟的国家和地区，人们发现风力机安装在寒冷潮湿地区时，叶片表面会发生结冰现象^［2］，如图1所示。并且几乎在全球绝大部分区域都可能发生风力机结冰问题。

Laakso 等^［3］对北欧地区风力机在寒冷气候下的运行情况进行了介绍。北欧地区的风电设施主要位于芬兰、瑞典和挪威的山脉及海岸线地区。这些地区拥有优越的自然风力条件，适合安装风力机，但是也存在因极低气温造成的严重的风力机结冰问题。芬兰拉普兰地区风力机结冰最严重的时间是11月至次年2月，每月结冰时间可达200 h。这期间风力机无法正常使用。芬兰几乎所有风电场在冬季都会因为结冰问题停机。挪威沿海地区拥有高风速和较为温和的气候，但当风力机安装在一定海拔高度以上时，整个海岸沿线的风力机都有可能出现结冰现象。瑞典每年都会经历的来自波罗的海的冷风会导致风力机大面积结冰，停机时间长达数周。

战培国^［4］对国外风力机结冰情况进行了分析，通过借鉴国外经验，为我国风力机结冰研究提供了参考。例如：北美地区风电场主要位于寒冷地区；中北部地区明尼苏达州风电场、东西海岸地区佛蒙特州风电场等都偶尔遭遇结冰问题，而靠近北极的偏北部地区风电场更容易受到结冰问题威胁。2021年美国德克萨斯州遭遇大范围极端严寒天气，导致该州超过60%的风力机结冰而无法正常发电，涉及的风力机装机容量达到1 700万kW，使得该州近500万用户停电。

我国风力机主要安装在风资源较为丰富的东北、华北等地区。这些地区冬季温度低、湿度大，极有可能发生结冰问题。而南方地区，由于风力机多安装在高海拔地区，受冬季低温天气影响也容易发生结冰问题。例如2022年受寒潮影响，南方地区风力机普遍结冰，最大影响范围超过40%，从而严重影响电力供应。由此可见，我国大部分地区都可能遭遇风力机结冰问题，且每年导致的发电小时数损失有几百小时，对风力机的安全运行及电力生产造成极大危害。

1.2 结冰对风力机的影响

结冰现象可能会发生在风力机所有迎风表面，例如叶片、整流罩、塔筒、传感器等。叶片是风力机的主要部件，叶片旋转为发电机提供了动力，是风力机能够将风能转化为电能的核心装备。风力机叶片翼型从叶根到叶尖变化较大，同时叶片较长，各处速度差异较大，叶片不同位置结冰后的外形也有较大差别，因此研究叶片结冰在风力机结冰问题中最为关键。

风力机叶片结冰会对风力机的运行效率以及设备安全造成极大危害。Fortin等^［5］对风力机叶片因结冰改变外形进行了研究，发现叶片表面结冰容易导致叶片升力下降，阻力上升，升阻比大幅下降，从而严重影响风力机的气动性能；同时，表面结冰改变了叶片原有的粗糙度，使得风力机气动性能恶化；气动性能降低会导致风力机组输出功率减小，甚至停机。在结冰问题严重的地区，每年因结冰造成的电力损失占总电力的比例高达20% ～ 50%^［6］。另外，叶片结冰使叶片质量增加，导致风力机负荷增大，而叶片表面不规则的结冰还会破坏风力机的动平衡，改变叶片固有频率，增加叶片的弯曲力矩，提高叶片的疲劳载荷^［7］。风力机叶片在不停转动过程中，叶尖部分的结冰可能会在离心力、气动力作用下脱落，如图2所示。脱落冰可能会对地面人员及设施造成危害^［8］。

1.3 风力机结冰条件

Frohboese等^［9］介绍了结冰的基本原理。东乔天等^［10］、王聪等^［11］对风力机结冰研究进行了总结。风力机结冰主要分为云内结冰和降雨结冰两类^［12］。云内结冰是指云层中尺寸较小的过冷水滴随气流运动，撞击风力机叶片表面后发生冻结形成结冰，如图3所示。根据风力机外形、运行时风速条件以及结冰气象条件的不同，云内结冰可分为明冰和霜冰。降雨结冰主要是指冻雨等天气条件下，水滴撞击风力机叶片的低温表面而发生结冰。与云内结冰相比，降雨结冰的水滴较大，结冰温度较高。云内结冰的水滴尺寸为μm级，降雨结冰的水滴尺寸通常可达mm级，因此云内结冰的水滴可以悬浮在空气中随气流运动，而降雨结冰的水滴会在空气中沉降。云内结冰形成时温度主要在−20 ～ 0 ℃，甚至低至−40 ℃，而降雨结冰发生时温度通常在−10 ～ 0 ℃之间。

2 研究方法

风力机结冰过程是一个非常复杂的物理过程，涉及流体力学、传热学等多个领域。首先，风力机结冰影响因素较多，包括风力机的几何外形、运行条件和结冰气象条件等。风力机的几何外形包括不同的叶片展长，以及不完全一致的截面翼型。设计条件不同，风力机来流风速、转速、攻角等也不相同，而结冰气象条件涉及过冷水滴和冰晶的含量与直径等。另外，叶片结冰过程是一个复杂的传热传质过程，包括空气、水滴、冰晶、表面液态水、冰层、叶片等固液气三相。结冰是一个非定常过程，会导致叶片表面几何形状发生变化，从而影响周围气流流动和粒子撞击范围，而这又会进一步影响叶片结冰过程。目前研究风力机结冰问题的方法主要包括自然结冰试验、结冰风洞试验和数值模拟方法等^{［13–14］}。

2.1 试验方法

自然结冰试验是指在自然结冰气象条件下对风力机叶片进行结冰试验。试验首先需要选择合适的场地和结冰气象条件，还需确保试验的可靠性和安全性。试验中通过观察和记录风力机在自然气象条件下的结冰情况，可以直接了解风力机表面结冰的分布和厚度等，从而评估风力机的抗冰冻能力和结冰对风力机性能的影响。

结冰风洞试验是通过在结冰风洞中模拟风力机结冰条件的试验方法。与常规风洞不同，结冰风洞中安装的制冷设备和喷雾设备可以用来模拟结冰气象条件所需的低温和水滴。将缩小后的风力机叶片模型安装在结冰风洞中进行试验，从而可较为精确地测量风力机叶片结冰过程中的物理参数，获得叶片结冰分布等数据。我国以及欧美等发达国家和地区都拥有结冰风洞。但受风洞尺寸限制，目前结冰风洞中只能开展风力机翼型的结冰风洞试验。

2.2 数值模拟方法

由于实际结冰过程较为复杂，开展试验的难度与成本均较高。随着计算机技术的高速发展和计算流体力学（CFD）的进步，利用数值模拟方法对结冰问题开展研究成为可能，进而可以对翼型升阻比等参数进行计算研究^［15］。该方法具有成本低、效率高、可重复性好等优点，迅速成为研究风力机结冰问题的主要方法。

目前，风力机结冰数值模拟方法主要借鉴了飞机结冰数值模拟方法，常用的计算软件有LEWICE和ANSYS软件。结冰数值模拟计算一般分为流场计算、水滴轨迹计算和结冰外形计算等模块，其中：流场计算主要通过模拟风力机周围气流流动情况获得叶片表面压力分布；水滴轨迹计算通过模拟空间中水滴运动过程获得水滴撞击特性，通常采用拉格朗日法和欧拉法描述；结冰外形计算通过建立质量守恒和能量守恒方程模拟表面传热传质，获得表面结冰增长速率和结冰外形^［14］。图4给出了数值模拟得到的风力机结冰外形。结冰数值模拟方法适用于几乎所有的结冰问题研究，但目前该方法对结冰过程中一些物理过程的描述还不够完善，计算结果还需要通过试验结果进行验证。

3 数值模拟研究现状

近年来，国内外针对风力机结冰进行数值模拟的研究对象包括二维风力机翼型和三维风力机叶片，主要内容包括空气−过冷水滴两相流数值计算、结冰生成计算和结冰对风力机气动特性的影响等3个方面。

空气−过冷水滴两相流数值计算是风力机结冰数值模拟的基础。目前，无论是风力机翼型结冰或是三维叶片结冰数值模拟均采用先计算空气流场，再基于流场结果计算水滴轨迹的方法。在风力机二维翼型空气流动和水滴运动研究方面，国外学者主要借助飞机结冰数值模拟软件对风力机翼型结冰进行模拟。美国开发的LEWICE软件^{［16–17］}是目前世界上最知名的飞机结冰数值模拟软件，其通过求解欧拉方程或N−S方程获得翼型的气流场数据，并采用拉格朗日方法进行水滴轨迹计算。其他结冰模拟软件如芬兰的Turbice软件^［18］、英国的TRAJICE2软件^［19］、法国的ONERA软件^［20］等均采用与LEWICE软件相似的计算方法。近年来，国内逐渐开展风力机结冰研究，并取得了一系列研究成果。Li等^［21］、朱程香等^［22］、闫菁菁^［23］、王晓东等^［24］、邓晓湖等^［25］、郑玉巧等^［26］、孔祥逸等^［27］、刘杰等^［28］均对风力机翼型结冰进行了数值模拟，其中空气−过冷水滴两相流数值计算方法与国外LEWICE软件方法类似。对于三维风力机空气流动和水滴运动研究，国外的FENSAP−ICE软件^{［29–31］}通过求解三维N−S方程获得叶片周围旋转流场，再基于欧拉模型求解水滴运动轨迹。国内易贤等^［32］对计算大型风力机水滴收集率的三维数值模拟方法进行了研究，其基本思路为采用多参考坐标系（MRF）方法计算空气流场，并以流场分布的计算结果为基础，求解基于欧拉方法的水滴运动控制方程，进而获得物体表面的水滴收集率。梁健等^［33］、闫菁菁等^［34］的相关研究中采用的方法与易贤等^［32］的类似。

水滴撞击叶片后冻结成冰过程是风力机结冰研究的核心问题。在风力机二维翼型结冰研究方面，1953年Messinger基于控制体积提出了经典的Messinger结冰模型，并将其用于模拟水滴撞击后的成冰过程，为结冰数值模拟打下了坚实的基础。随后众多学者基于经典Messinger模型提出了各种热力学模型^［35］。LEWICE软件就是利用Messinger结冰模型开展风力机翼型结冰数值模拟研究。Ali等^［36］也在该模型基础上研究风力机翼型结冰。国内各研究团队在研究风力机二维翼型结冰时均采用Messinger结冰模型。Han等^［37］对三维风力机叶片整体进行了冰风洞试验，并利用LEWICE软件进行了数值模拟。结果显示，对于结冰温度较低的霜冰，其模拟结果与冰风洞试验结果吻合良好。当结冰温度升高、液态水含量上升时，水滴撞击叶片后并不会迅速结冰，而是在叶片旋转离心力的作用下，有部分水沿展向流动。此时利用LEWICE软件进行二维翼型结冰数值模拟的结果与试验结果差别较大。目前三维风力机叶片结冰研究方法主要有：第一种是基于动量−叶素理论对叶片结冰进行模拟。将风力机叶片沿展向划分为多个翼型截面，利用二维结冰数值方法计算得到翼型结冰外形，再利用几何方法将各截面连接成三维叶片^［38］。该方法可快速获得结冰后的叶片外形。第二种是在空气−过冷水滴计算结果基础上，建立三维结冰数学模型，对叶片结冰进行直接模拟^［39］。例如FENSAP−ICE软件采用浅水膜模型求解风力机叶片的结冰^［31］。易贤等^［40］采用三维飞机结冰数值模拟方法对叶片结冰进行计算。该模型考虑了过冷水滴的撞击、蒸发、对流换热、冻结以及液态水流动。梁健等^［33］、闫菁菁等^［34］、Hu等^{［41–42］}基于二维Messinger结冰模型建立了三维结冰模型并对叶片结冰过程进行了模拟。Virk等^［43］也采用类似的方法研究叶片结冰过程。

结冰对风力机气动特性的影响一直是风力机设计、制造、安装及运行中极为重要的研究内容。结冰后风力机气动特性研究的对象主要包括二维翼型和三维风力机。在二维翼型结冰后气动特性研究方面，Hochart等^［44］通过试验测试了结冰条件下叶片翼型的气动性能及风力机的一些基本性能。Bragg等^［45］、Mirzaei等^［46］和Yirtici等^［47］分别对结冰后翼型的气动特性进行了研究。国内王文龙等^［48］采用Fluent 软件研究了结冰对NACA 63418翼型气动性能的影响，其中湍流模型采用S−A和LES模型。结果表明，结冰使翼型的失速提前，阻力系数显著增大。朱程香等^［22］对NACA 63618翼型结冰后气动特性进行了研究，发现结冰后翼型升阻比最大降幅达到61%。此外，Li等^［21］、闫菁菁等^［34］、王晓东等^［24］、郑玉巧等^［26］、孔祥逸等^［27］、杜静宇等^［49］也开展了结冰后翼型气动特性研究，得到了相似的规律。在结冰对风力机三维气动特性的影响方面，何玉林等^［50］基于叶素动量理论和广义动态尾流理论研究了结冰对2 MW风力机气动性能、功率以及年发电量的影响。结果表明，结冰使得风力机翼型的升力降低，阻力增加，年发电量仅为无冰情况下的52%。任鹏飞等^［51］通过试验获得叶片三维结冰情况，并利用数值模拟方法对气动特性进行了计算。结果表明，结冰会导致最高达51%的功率损失，同时使得风力机启动风速提高。易贤等^［40］、Wang等^［52］、王强等^［53］基于计算得到的三维结冰外形研究了风力机结冰后的气动特性。

4 结束语

风力发电潜力巨大，但诸如风力机叶片结冰等问题仍需寻找更好的解决方法。数值模拟方法因其低成本、高效率在风力机结冰研究中广泛应用。通过风力机翼型结冰数值模拟，研究翼型结冰外形及其对气动特性的影响，为风力机结冰及防除冰系统设计分析提供理论支持。目前国内外在空气−过冷水滴运动研究中普遍采用松耦合方法，忽略了水滴对空气的作用。在结冰外形计算中认为结冰时间内结冰速率不变，将其再乘以结冰时间来获得结冰外形。同时，研究过程中采用稳定的来流条件，未考虑来流气象参数的空间分布和随时间变化规律，以及风力机叶片变桨等运动形式。这些将导致计算结果与实际结果的偏差。对这些问题进一步开展研究有助于提升解决风力机结冰问题的能力，对提高风能利用效率具有重要意义。