利用叶片元件阻尼方法预测风力涡轮机叶片在后缘迎风状态下的气动弹性不稳定性

《Energy》：Prediction of wind turbine blade aeroelastic instability under trailing edge windward states using blade element damping method

　　李登宁|刘荣祥|郭远林|高传强|张伟伟
西北工业大学航空学院，中国西安710072

摘要
随着风力涡轮机叶片尺寸的增加和结构刚度的降低，气动弹性稳定性和振动已成为叶片设计中的关键技术挑战。特别是，对于尾缘迎风方向的拍动气动弹性不稳定性，缺乏快速准确的预测技术，这成为行业发展的主要障碍。为了解决这个问题，本研究开发了一种用于三维叶片在尾缘迎风条件下的快速气动弹性稳定性预测方法。该方法基于叶片元件理论，将叶片划分为多个参考叶片元件，每个元件被视为一个独立的二维翼型。然后使用降阶模型（ROM）方法计算每个参考翼型的气动弹性阻尼特性。随后，通过整合叶片结构模态信息和对应于参考翼型的模态阻尼系数，采用模态叠加方法来推导整个叶片的气动弹性阻尼。利用这种构建的叶片元件降阶模型（BE-ROM）阻尼方法，成功预测了叶片在177度到183度攻角范围内的气动弹性不稳定区域，并通过风洞测试进行了验证。此外，通过研究几何参数对气动弹性行为的影响，揭示了潜在的不稳定性触发因素。叶片的上部对叶片的整体不稳定性影响最大。

引言
随着风力涡轮机转子直径的增加和结构刚度的降低[1]，气动弹性问题变得越来越突出，成为限制发电效率提升的核心因素之一。在这些问题中，气动弹性稳定性问题尤为关键。Hansen等人[2]对风力涡轮机中的两个主要气动弹性稳定性问题——经典弯曲-扭转耦合颤振和失速颤振——进行了深入全面的分析，阐述了它们的基本原理和特征，并进一步预测随着风力涡轮机尺寸的继续增大，气动弹性问题将变得更加严重。Bossanyi[3]根据风力涡轮机的发展趋势预测，未来的风力涡轮机设计概念将从传统的载荷驱动模式转向稳定性驱动模式。因此，如何准确预测风力涡轮机叶片的气动弹性稳定性已成为风能发电领域需要解决的关键科学和技术挑战。

风力涡轮机气动弹性稳定性的数值预测方法主要分为三种类型：时域仿真方法、频域函数变换判别方法和特征值分析方法。
时域仿真方法通过将流体求解器与结构求解器耦合来实现叶片气动弹性行为的动态仿真。在此过程中，流体求解器通常采用计算流体动力学（CFD）方法，而结构求解器一般基于非线性梁理论，如几何精确梁理论[4]、[5]、[6]。由于具有较高的计算精度和能够可视化叶片不稳定性期间的流体形态变化，该方法在气动弹性研究中得到了广泛应用。例如，Heinz等人[7]和Horcas等人[8]、[9]使用HAWC2和EllipSys3D软件模拟分析了10 MW风力涡轮机叶片在接近90°攻角下的涡激振动。他们研究了攻角、倾斜角和叶尖形状对叶片涡激振动的影响，并比较了不同参数组合下的尾流涡形差异。Dose等人[10]通过将OpenForm与BeamForm耦合实现了叶片的气动弹性耦合仿真。此外，Sessarego等人[11]开发了一种用于三维风力涡轮机叶片时域仿真的粘性-非粘性CFD-CSD（计算结构动力学）方法，验证了该方法在计算效率上的优势。尽管时域仿真方法具有明显优势，但在工程应用中存在显著局限性：它需要大量的计算资源，在使用CFD求解空气动力载荷时对网格质量要求较高，且每次计算通常只能分析一个参数组合。

函数变换判别方法通过使风力涡轮机叶片或翼型受到强迫运动来获得结构在特定运动模式下的空气动力响应特性。然后通过分析流体与结构之间的能量转换关系来评估气动弹性稳定性。由于该方法省略了流体与结构之间的耦合，实施起来较为简便，因此在气动弹性研究中得到了广泛应用，例如圆柱体和翼型的涡激振动（VIV）研究[12]、[13]、[14]。例如，Grinderslev等人[15]、[16]通过比较强迫运动方法和完全耦合的流体-结构相互作用（FSI）模拟的结果，证明在特定运行条件下可以忽略结构-流动耦合效应。利用这种方法，他们研究了风力涡轮机叶片在特定模式和频率下的VIV问题。他们还研究了不同来流角度对尾流涡频率的影响，并提供了在不同来流条件下应用CFD方法时优化湍流模型选择的建议。与时域耦合仿真相比，该方法通过排除结构耦合可以减轻计算负担。然而，需要注意的是，对于所有运行条件，不能忽略耦合效应对稳定性分析结果的影响；此外，该方法难以模拟叶片的非线性运动。

特征值分析方法涉及将空气动力模型与运动结构方程耦合，以建立耦合系统的控制方程组。通过分析系统特征值的分布，该方法能够准确评估所有模态阶数的结构阻尼特性，从而深入揭示系统的气动弹性稳定性特性[17]、[18]、[19]。例如，Thirstrup等人[20]基于准稳态空气动力模型推导出了适用于二维翼型失速颤振的气动阻尼公式。他们指出，气动弹性不稳定性更可能发生在边缘方向模式。此外，Hansen[21]将叶片元件动量（BEM）理论与B-L空气动力模型结合，研究了风力涡轮机叶片的弯曲-扭转耦合颤振。该研究揭示了叶片的一阶拍动模式与三阶扭转模式之间的耦合效应可能导致叶片不稳定性。BEM方法结合了叶片元件理论和动量理论，将叶片离散为多个元件单元并评估各自的空气动力贡献。通过迭代求解和积分计算，确定了作用在整个风力叶片上的空气动力。然而，BEM方法本质上忽略了三维流动效应，并缺乏非稳态空气动力数据，如动态失速，这限制了其预测空气动力-结构耦合的能力。此外，Chen[22]通过将结构有限元方法（FEM）与准稳态空气动力建模相结合，建立了叶片空气动力阻尼的计算模型。特征值方法通过矩阵特征值具有较短的计算周期。另一方面，由于该方法从模态角度进行分析，其结果为工程实践中的叶片设计和优化提供了有力指导。

由于尾缘迎风（TEW）条件具有出色的载荷减少能力，因此可以通过调整叶片的尾缘来适应强风天气，从而减少载荷。Cai等人通过风洞测试测量了DU91-W2-180翼型在迎风和尾缘迎风条件下的空气动力系数[23]。他们的结果表明，在尾缘迎风条件下作用在该翼型上的空气动力显著减小，表明具有很大的载荷减少潜力。然而，现有文献指出，风力涡轮机叶片在尾缘迎风条件下运行时可能会经历气动弹性不稳定振动。具体来说，Li等人[24]通过风洞测试发现，DU91-W2-250翼型在接近165度的攻角下表现出拍动方向的气动弹性不稳定性。该实验进一步表明，翼型产生的振动与其结构固有频率一致，且振幅随风速的增加而增大。Lu等人[25]也通过风洞测试确认，IEA 15 MW风力涡轮机的叶片在尾缘迎风条件下表现出类似的拍动方向气动弹性不稳定性。

目前，这种不稳定性的预测主要依赖于两种方法：耦合CFD-CSD（计算流体动力学-计算结构动力学）仿真和风洞测试。然而，前者需要大量的计算资源，而后者面临实验成本高和难以设计弹性缩放模型的挑战。因此，迫切需要一种准确且快速的方法来预测尾缘迎风条件下的气动弹性不稳定性。

气动弹性降阶模型（ROM）方法可以揭示不稳定机制并实现快速的气动阻尼预测。它已被广泛用于研究二维气动弹性问题，如涡激振动、锁相振动、驰振和失速颤振。Li等人[26]使用ROM揭示了圆柱形和方形柱体的气动弹性不稳定机制，指出流动模式与结构模式之间的耦合是这些钝体振动问题的关键。Dou等人[27]使用ROM研究了跨音速预涡流中的频率锁相现象。Lyu等人[28]通过降阶模型分析发现，零频流动模式与结构模式的耦合是NACA0015翼型层流分离颤振不稳定性的触发因素。Luo等人[29]基于ROM开发了一种被动模态控制器，有效抑制了涡激振动问题中的锁相振动现象。最近，Gao等人[30]使用ROM揭示了DU91-W2-250翼型在TEW条件下的气动弹性特性和潜在机制。确定了主要的亚临界流动模式，其特征是一对稳定性裕度较小的特征值。这种流动模式与结构模式之间的耦合在较低的降阶频率下导致了结构模式的不稳定性，从而引发了观察到的拍动方向振动。

尽管风力涡轮机叶片具有相当大的展向尺寸，但其空气动力特性仍可以通过几个基础翼型（称为叶片元件）来推导，使用叶片元件动量理论和线性插值。换句话说，可以通过二维叠加方法获得三维叶片的复杂空气动力行为，避免了求解复杂的三维问题，显著提高了分析效率。这些方法也广泛应用于具有较大细长比结构的流动诱导振动研究，如细长火箭的分布式空气动力阻尼分析[33]和长柔性竖管的涡激振动预测[34]、[35]。此外，当尾缘迎风时，展向流动较弱，叶片元件之间的空气动力干扰最小，进一步增强了叶片元件叠加方法的适用性。因此，可以使用上述方法计算叶片的空气动力特性。

基于叶片元件理论和翼型气动弹性ROM方法的思想，提出了一种新的叶片气动弹性稳定性预测方法，即叶片元件降阶模型（BE-ROM）方法。该方法采用叶片元件理论假设，将叶片划分为多个元件，并将每个元件视为二维翼型进行分析。构建了叶片元件的ROM，并通过特征值分析在不同流动和结构参数下获得了每个叶片元件的空气动力阻尼。结合叶片结构模态信息，使用模态叠加方法预测三维叶片在尾缘迎风条件下的气动稳定性。为了验证BE-ROM方法的准确性，本文通过三维叶片的风洞测试验证了该方法的预测结果。第2节介绍了该方法理论，包括BE-ROM方法框架和风洞实验设置。第3节介绍了建模过程和元件单元上的阻尼提取。第4节展示了三维叶片的阻尼预测结果，并通过风洞实验进行了验证。最后，第5节讨论了叶片不同展向位置的阻尼贡献。

研究设置和方法
本节首先描述了BE-ROM阻尼预测方法的基础原理，然后介绍了二维翼型的气动弹性ROM方法。

使用ROM提取翼型的阻尼
本节将介绍使用第2.2节介绍的ARX方法为基本叶片元件（即二维翼型）构建气动弹性ROM的过程，并在不同来流条件和结构参数下获得气动弹性阻尼。为了方便说明，本节以DU91-W2-250翼型为例进行详细阐述[31]。
**3D叶片的阻尼预测**
本节将利用第2节提出的BE-ROM方法对叶片的空气弹性阻尼进行预测。首先介绍叶片的空气动力特性和结构参数，然后分析在不同参数条件下该叶片关键部件（包括第3节讨论的DU91-W2-250翼型）的空气弹性阻尼分布图。最后，根据各部分的贡献计算整个叶片的总空气弹性阻尼。

**几何参数和结构参数的影响**
前几节已经实现了在尾缘迎风条件下叶片空气弹性阻尼的预测，并通过风洞试验验证了结果的准确性。本节将探讨几何参数和结构参数（包括展向扭角以及结构质量和刚度分布）的贡献和影响。由于叶片存在展向扭角，不同叶片部位的局部攻角会发生变化。

**结论**
本研究解决了在尾缘迎风条件下风力涡轮机叶片展向方向空气弹性不稳定性的预测问题。提出了一种基于叶片元素理论和空气弹性降阶模型的BE-ROM阻尼方法。利用该方法，在尾缘迎风条件下预测了缩比风力涡轮机叶片的空气弹性阻尼，并通过风洞实验进行了验证。从预测结果来看，所提出的BE-ROM方法具有一定的有效性。

**作者贡献声明**
郭远林：撰写、审稿与编辑、可视化、验证。
高传强：方法论、资金筹集、概念构思。
李登宁：撰写、初步草案、研究、形式分析。
刘荣祥：撰写、审稿与编辑、验证、研究。
张伟伟：监督、资金筹集、概念构思。

**未引用参考文献**[32]

**利益冲突声明**
作者声明他们没有已知的可能会影响本文研究的财务利益或个人关系。

**致谢**
本研究得到了重庆市科技创新重点研发计划（项目编号：CSTB2024TIAD-STX0019）的支持，同时也得到了国家自然科学基金（项目编号：12372271、U23B6009）和中央高校基本科研业务费的支持。