风力机叶片动态气弹变形及其对整机性能的影响

2020-07-23黄鑫祥郭小锋齐剑峰

可再生能源 2020年7期

黄鑫祥，郭小锋，齐剑峰

（中原工学院机电学院，河南郑州 450007）

0 引言

水平轴风力机叶片的气弹变形包含挥舞、摆振和扭转3 个方面。常用的风力机动力学分析软件（如GH-Bladed 和FAST）只考虑了叶片在挥舞和摆振两个方向的自由度。然而，随着风电机组容量的不断扩大，为降低载荷，大型低风速及海上风电叶片的柔韧度较早期同类型的风电叶片有明显的增大，如何准确地分析叶片在实际工况下的气弹扭转变形，进而分析叶片的动态气弹扭转变形对风力机气动功率、载荷及气弹稳定性的影响，是大型风电叶片的气弹特性评价以及气弹剪裁设计中必须解决的问题。

郭小锋[1]运用修正的叶素-动量理论和有限元方法，将风力机的动力学问题转化为静力学问题，提出了考虑气动弹性的叶片性能分析方法，综合分析了叶片扭转变形对风力机气动性能的影响。李德源[2]采用“超级单元”模型，应用计算多体动力学理论和风力机气动模型，建立了受约束的风力机系统气弹耦合方程。江术新[3]建立了整机在变桨工况下的气-弹-控联合仿真模型，时域分析了变桨过程中的桨距角变化、功率曲线变化和叶根处受力和力矩的变化。

现有的研究表明，静气弹性分析方法没有反映柔性风电机组的动力学问题本质。风力机的实际工况具有复杂多变性，对于风力机叶片在实际复杂工况下的气弹扭转变形问题，及其对风力机气动性能及气弹问题性的影响规律仍须进一步研究。

本文基于叶片及塔筒部件的参数化表达方法，采用多体动力学仿真软件SIMPACK 建立了风力机整机刚柔耦合模型，提出了一种叶片动态气弹扭转变形分析的新方法。在几种典型风况下，对NREL5MW 风力机进行了气弹动力学分析，采用matlab 软件对分析结果进行处理，得到了叶片在典型风况下的动态气弹扭转变形数据[4]。在此基础上，对NREL5MW 风力机在几种典型工况下的气动功率、载荷及气弹稳定性问题进行了分析研究。

1 风力机刚柔耦合模型的构建

1.1 参数化建模

风力机整机模型是典型的刚柔耦合多体系统。在工作状态下，风力机叶片存在明显的弹性变形，所以将叶片作为柔性体处理。而风力机的其他部件（如机舱、主轴及轮毂等部件）的刚性大，在工作时的弹性变形很小，可以看作刚体。为了缩短计算机的求解时间，本文将风力机模型进行了合理简化，只包括叶片、轮毂、塔筒、机舱等主要部件，而不考虑齿轮箱、发电机等部件。偏航法兰与塔筒顶端固接，用于连接机舱，虽然机舱中不创建主机架、齿轮箱、发电机等模型，但是需要把忽略部件的质量和转动惯量算入机舱中。低速轴与机舱为转动连接，轮毂与低速轴进行固接，叶片与轮毂由变桨法兰连接。

采用EXCEL 表格对叶片的外形和铺层结构进行参数化表达，采用matlab 编程调用PreComp软件计算得到SIMPACK 叶片建模中所需各个截面的属性参数（截面质心、弹性中心、剪切中心、扭转刚度、弹性模量等），并自动生成Simpack 软件Rotorblade Generation 计算模块的输入文件[5]，[6]。

风力机的控制系统采用基本的变桨变速控制策略，在SIMPACK 中设置风轮转速及桨距角在不同稳态风速下的参数，同时通过设置PID 参数实现对风轮的动态控制。

在SIMPACK 中，刚柔耦合多体系统建模是在浮动参考坐标系法的基础上建立多体动力学方程，并由SIMPACK 求解动力学方程获得风力机的动态响应。在多体系统（MBS）理论中，柔性体看作有限元模型若干个节点的组合，其中每个节点的运动由跟随整体坐标系的整体平动和转动以及在局部坐标系下发生小的线性变形组成。图1 所示为柔性体的一般模型图。

为描述节点的运动，一般采用惯性坐标系下的笛卡尔坐标和表示方位的欧拉角。

图1 柔性体的一般模型Fig.1 General model of a flexible body

1.2 风场模型

空气动力学软件AeroDyn 是分析风力机叶片风载的有力工具，为了实现与SIMPACK 软件的联合仿真，需要准备AeroDyn 的输入文件、风场文件、翼型文件、塔筒文件等[7]。本文的风场文件由NREL 提供的IECWIND 软件生成。

AeroDyn 基于BEM 理论建立迭代方程计算风轮气动力和诱导速度，通过迭代求解叶片单元的诱导速度、攻角、推力系数等参数，直到诱导因子和入流角的值收敛到最终值[8]。

AeroDyn 与SIMPACK 动力学联合仿真流程如图2 所示。

图2 动力学仿真流程图Fig.2 Diagram of simulation process

图2 中，AeroDyn 与SIMPACK 两者通过SIMPACK 中241 号力元接口进行数据交换。AeroDyn 同时从输入文件和SIMPACK 求解模块中收集t 时刻的叶片几何形状、工况条件、叶元速度、位置以及流入风等信息，然后利用以上信息计算出叶片各截面的气动数据。 SIMPACK 的求解模块利用241 号力元接口调用气动数据计算得到t+dt 时刻叶片各截面的气弹变形和气动力，计算得到的数据又通过241 号力元接口反馈给AeroDyn 用来计算t+dt 时刻各截面的气动数据，如此往复，直到仿真完成。而控制模块根据输出的反馈力矩和转速对风力机进行变速变桨控制。

为了得到仿真结果的精确性，本文选择SODASRT_2 解算器对模型进行离线积分求解，仿真时长设为120 s，求解步长设为0.01 s。由于SIMPACK 每个时间步长都会调用AeroDyn 计算气动力变化，所以AeroDyn 求解时间步长也设为0.01 s。

为研究阵风、风向角度对叶片扭转变形的影响，选取如表1 所示风场模型。考虑到海上风况的复杂多变性，本文在建立风场模型时将方向变化的极端相干阵风（EDC）和极端工作阵风（EOGR）风况下风速和风向的波动设定为4 个变化周期，并且假定在EDC 和EOGR 风况作用下控制系统故障，研究风力机在极端风况下扭转变形与整机性能的关系。

表1 极端风况Table 1 Extreme wind conditions

1.3 模型验证

通过SIMPACK 中Eigenvalues 模块对风力机（NREL5MW）进行模态分析，计算得到整机的各阶自振频率（表2）。由表2 可知，本文的分析结果同NREL 报告中的数据误差在5%以内[4]，证实了SIMPACK 所建立风力机整机模型的可靠性。

表2 5 MW 风力机整机自振频率对比Table 2 Comparison of natural frequency of 5 MW wind turbine

2 叶片气弹扭转变形分析

风力机在启动阶段的各种变化均不稳定，为更好地研究各风况对叶片气弹变形的影响，本文假设各风况发生变化时风力机处于稳定阶段。假定各风况在风力机启动后40 s 时刻发生变化，并于80 s 之后逐渐恢复为12 m/s 的稳态风。对风力机整机进行了仿真分析，得到4 种工况下叶尖挥舞、摆振的动态响应（图3）。

图3 各风况下的叶尖位移Fig.3 Tip deformation under various wind conditions

为分析叶尖部位翼型的气弹扭转变形数值，本文首先采用SIMPACK 软件得到叶尖某截面各节点的初始坐标以及其相对于初始坐标的变形量，将各节点的初始坐标与对应的位移量相加得到节点在各个时刻的动态坐标值。将各节点的动态坐标值投影到初始时刻的截面上，得到截面扭转变形后的节点坐标值。对于NREL5MW 叶片，距离叶根60 m 处截面的初始形状和在48 s 时的投影形状对比如图4 所示。

图4 在ECD 风况下，0 s 和48 s 时叶尖翼型扭转形变的对比Fig.4 Comparison of tip airfoil torsional deformation at 0 s and 48 s under the ECD wind condition

在风力机的刚柔耦合模型中，假设叶片截面在扭转变形时自身形状不发生变化，叶片截面扭转变形前后两翼型弦长线间的夹角即为该截面在ECD 风况48 s 时刻的气弹扭转变形量。

为分析叶片的动态气弹扭转时域变形趋势，由上述方法得到距离叶根60 m 处截面各个时刻的动态气弹扭转变形曲线（图5）。

图5 各风况下的叶尖扭转变形Fig.5 Tip torsional deformation under various wind conditions

由图5 可知，当风力机运转稳定后，各风况下的叶尖动态扭转角度在2°上下浮动，说明叶片具有的扭转刚度使得叶片在一定风速下的扭转角度会在某一个值上下变动。而叶片各个截面的扭转刚度不同，导致了叶片各截面扭转角度不同，计算叶片实际工作中的扭转角度对叶片的优化设计具有指导作用。

对比EOGR，EDC 风况作用下叶尖扭转变化曲线可知，阵风变化对叶尖扭转变形的影响要大于风向角度变化，其中风向角度变化使得扭转变化频率增加。在ECD 风况作用下，叶尖扭转变形量最为明显，叶尖最大扭转为4.78 °，最小扭转为0.71 °，阵风风速越大，叶尖扭转变形量越大，而风速和风向角度的共同作用增加了叶片的扭转频率。

在ECD 风况下，由于风速和风向角度的共同作用造成了叶尖扭转变形的突然增大，这种短时间内叶尖扭转变化的突然加剧容易对叶片造成致命的损坏。而风向角度多变的风况使得叶尖扭转频率增大，这种风况下容易造成叶片的疲劳破坏。

3 气弹扭转变形对整机性能的影响

3.1 气动性能

根据BEM 理论计算可得，在11.4 m/s 的稳态风下，风力机功率为5 MW。而在实际工作过程中，由于塔筒变形、叶片的气弹变形等影响，风力机功率并不是一个恒定的值。不同风况的作用对叶片气弹变形的影响很大，导致各风况下的风力机功率存在明显波动。

风力机控制模块根据检测到的实时转速做出变桨响应，使其稳定在额定功率附近。仿真时，在40～80 s，ECD，EOGR 风况的风速和风向角急剧变化，导致风力机变桨速率、叶片翼型攻角变化幅度增大，叶片气动载荷增加。随着叶片变桨速率增大，叶片在风速上升期卸载越快，对功率的影响也就越明显。另外，当控制系统故障时，风力机变桨控制明显滞后于实际风况，在风速、风向变化剧烈，但未到达切出风速的极端风况下，控制系统无法实时调整到对应风速的桨距角，从而导致控制系统暂时失效，造成功率波动异常。图6 为在ECD，EOGR 风况下，风力机功率的变化曲线。由图6 可知，风力机功率出现了大幅度的波动，峰值功率甚至达到了13 MW。

图6 各风况下的风力机功率变化Fig.6 Variation of wind turbine power under different wind conditions

风力机的功率波动主要受阵风的影响，阵风变化越大，叶片的气弹变形越大，风力机功率变化也就越大，在极端风况下，功率和载荷波动异常，在实际工作中会造成风力机的损坏。

3.2 载荷波动

图7 为各风况下轮毂处载荷的变化曲线。由图7 可知：相比于EDC 风况，ECD 和EOGR 风况有非常明显的载荷波动，最大轮毂推力值分别为1 280 kN 和1 175 kN；ECD 风况作用下的轮毂处最大载荷要小于EOGR 风况，叶片的扭转变形在一定程度上能削减风力机上的极限载荷。

图7 各风况下轮毂处载荷的变化Fig.7 Load changes at the hub under different wind conditions

3.3 气弹稳定性

在复杂工况下，风速、风向的剧烈变化引起了叶片挥舞、摆振和扭转方向的急剧变化，从而造成叶片气弹失稳。特别是叶片的扭转变形直接影响攻角变化，从而造成升力系数的变化，严重时会发生失速现象。本文通过AeroDyn 对加载风载的叶片进行计算分析，选取50～60 s 时刻各风况下距离叶根60 m 处截面翼型的气动数据进行处理，得到攻角-升力系数曲线，并与静态翼型攻角-升力系数曲线进行比较（图8）。

图8 各风况下叶尖截面攻角-升力系数曲线Fig.8 Angle of attack-lift coefficient curve of blade tip section under various wind conditions

由图8 可知：NWP，EDC 风况下的升力系数随攻角呈线性变化，且该时间段内EDC 风况下的攻角大于NWP 风况下的攻角，但没达到翼型的静态失速攻角；EOGR 风况下的攻角进一步增大，接近翼型静态失速攻角，仍呈线性变化；ECD 风况下的攻角超过静态失速攻角并继续增大，升力系数也随着攻角的增大总体呈上升趋势，达到峰值后攻角逐渐减小，升力系数也随之减小，曲线开始回归线性区，这是典型的动态失速现象[9]；在ECD 风况作用下，叶尖截面存在明显的动态失速现象，而其他风况作用下，叶尖截面均未出现失速现象，所以在风速、风向变化幅度大的风况下，叶片容易出现动态失速现象，从而造成叶片气弹失稳，严重时叶片会发生颤振甚至断裂。