考虑叶片偏航和干扰效应大型风力机体系风振响应与稳定性分析

柯世堂 王晓海

摘 要：叶片偏航和干扰会显著改变大型风力机表面气动力分布模式，进而影响风力机体系的风振响应和稳定性能.以某5 MW大型风力机为研究对象，首先采用大涡模拟（LES）方法进行了最不利叶片位置下考虑6个偏航角（0°、5°、10°、20°、30°和45°）影响的风力机体系流场和气动力模拟，并与规范及国内外实测结果进行对比验证了大涡模拟的有效性.在此基础上，结合有限元方法系统分析了不同偏航角下风力机塔架-叶片耦合模型的动力特性、风振响应和稳定性能.结果表明：不同偏航角下塔架径向位移均值和均方差的最大值均出现在塔架环向0°和180°处，最大塔底弯矩均出现在环向20°处.0°偏航时各叶片顺风向位移响应极值均大于2.7 m，随着偏航角的增大，塔架顶部径向位移、叶片顺風向位移和叶片根部内力的均值及均方差均逐渐减小，而临界风速则呈现先减后增再减小的趋势.综合表明：0°偏航角下风力机体系气动性能和风振响应均最为不利，45°偏航角下风力机体系的稳定性能最为不利.

关键词：风力机；大涡模拟；偏航效应；气动力分布；风振响应；稳定性能

中图分类号：TK83 文献标志码：A

文章编号：1674—2974（2018）07—0061—10

Abstract： The aerodynamic performances of large wind turbine systems are significantly affected by blade yaw and interference， and then the wind-induced response and stability of the wind turbine system are changed. Taking the 5 MW wind turbine as the example，the flow field and aerodynamic forces of the wind turbine systems considering the six yaw angles （0， 5， 10， 20， 30 and 45 degrees） were simulated by large eddy simulation method， and the numerical simulation results were compared with standard curves to verify the validity of the numerical method. On this basis， the dynamic characteristics， wind-induced response and stability of the wind turbine systems under different yaw angles were analyzed by the finite element method. Main conclusions are as follows： The maximum values of the mean value and the mean square deviation of the radial displacement forthe tower under different yaw angles appearat 0 and 180 degrees， the maximum bending moment at the bottom of the tower appears in the circumferential direction of 20 degrees， and the peak value of three blade tip displacement response under 0 degrees yaw is more than 2.7 m. With the increase of the yaw angle， the mean value and mean square deviation of the radial displacement at the top of the tower， the forward displacement of the blade and the internal force of the blade root are gradually reduced， and the critical wind speed decreases first， then increases and decreases again. The results show that the aerodynamic performance and the wind-induced response of large wind turbine systems is the most unfavorable andthe largest under the 0 degree yaw angle， and the stability performance under 45 degrees yaw angle is the most unfavorable.

Key words： wind turbines；large eddy simulation；yaw effect； aerodynamic force distribution；wind-induced response；stability

由于风向连续变化且风力机偏航系统转速较慢，风力机的叶轮转轴无法及时对准来流风向，导致风力机在日常运营时常处于偏航状态[1-2].由于不同偏航角下处于最不利位置[3]的叶片对塔架的遮挡程度不同，风力机体系的表面流场和气动力分布发生明显变化，进而使得风力机体系的风振响应和稳定性发生改变.因此，对不同偏航角下大型风力机体系的风振响应和稳定性进行系统研究具有重要的工程和理论价值.

针对大型风力机的风振响应和稳定性研究，文献[4-5]采用叶素-动量理论和谐波合成法得到了风力机叶片的脉动时程风荷载，并结合有限元方法进行了叶片的动态响应分析；文献[6-7]分别探讨了离心刚化效应和屋顶集风效应对风力机塔架风振响应的影响，结果表明离心刚化效应减小了塔架风振响应而屋顶集风效应增强了塔架风振响应；文献[8]采用计算流体动力学方法对风雨共同作用下风力机的风场和风驱雨过程进行了模拟，研究了风雨共同作用下风力机的受力和变形特性；文献[9]基于HHT方法对风力机启动、紧急停机等非稳定工况进行了时频域特性分析，研究表明塔架和叶片位移响应频率主要集中于一阶振动频率；文献[10]探讨了水平轴风力机叶片气弹稳定的判别数.已有研究主要针对风力机对准来流风向工况下风力机风振响应和稳定性分析，鲜有系统讨论不同偏航角对大型风力机塔架-叶片体系风振响应与稳定性的影响.

鉴于此，以南京航空航天大学自主研发的5 MW大型风力机为研究对象，首先基于LES方法对6个偏航角（0°、5°、10°、20°、30°和45°）下叶片处于最不利位置（叶片与塔架重合[3]）时风力机周围风场进行数值模拟，并与规范及国内外实测曲线进行对比验证了数值方法的有效性.然后对不同偏航角下风力机体系进行动力特性分析，并结合完全瞬态法对不同偏航角下大型风力机塔架-叶片耦合体系的风致动态响应和稳定性能进行了对比研究.最后归纳总结了叶片偏航和干扰效应对大型风力机体系气动力、风振响应和稳定性的影响规律.

1 计算方法及工程概况

1.1 大涡模拟理论

由于风力机所处流场属于非定常且扰流情况复杂，基于大涡模拟方法能够对风力机复杂的流场进行更好的模拟[11].大涡模拟法采用滤波函数，将流场中的涡分为大尺度涡和小尺度涡，对大尺度涡进行直接求解，而小尺度涡则采用亚格子模型进行模拟.因此，LES需对N-S方程进行空间过滤，过滤后的连续性方程和N-S方程为：

1.2 完全瞬态法

采用完全瞬态法[12]求解风力机体系瞬态动力学平衡方程，其核心是使用隐式方法Newmark和HHT来直接求解瞬态问题.其中Newmark方法使用有限差分法，在一个时间间隔内，有：

将大涡模拟获得的风压系数作为风荷载时程输入参数，基于完全瞬态法和ANSYS软件平台进行风力机塔架-叶片耦合模型风振响应时域计算，其中各模态阻尼均为2%，积分时间步长取为0.5 s，加载时间步数为2 048步.

2 风场的数值模拟

2.1 工程概况

表1给出了5 MW风力机主要设计参数及模型，其中塔架通长变厚，各叶片沿周向成120°夹角均匀分布.定义风力机叶轮转轴与来流风向的夹角被称为偏航角，设置偏航角为0°、5°、10°、20°、30°和45°六种计算工况[13]，简称为工况一、二、三、四、五和六.大涡模拟中采用固定来流风向、旋转风力机模型的方法实现不同偏航角下大型风力机的风场模拟，如图1所示.

2.2 计算域及网格划分

为保证风力机尾流的充分发展，计算域尺寸取12D × 5D × 5D（流向X × 展向Y × 竖向Z，D为风力机风轮直径），风力机置于距离计算域入口3D处.为了兼顾计算效率与精度，同时考虑到叶片外形复杂，网格划分采用混合网格离散形式，将计算域划分为局部加密区域和外围区域.局部加密区域内含风力机模型，采用非结构化网格进行划分，外围区域形状规整，采用高质量的结构化网格进行划分.表2给出了不同网格方案下网格质量和迎风面压力系数.由表2可知，随着网格总数的增加，网格质量逐渐增加，网格歪斜度和迎风面风压系数呈现逐渐减小的趋势，而1 100万网格数和3 000万网格数的网格质量和计算结果无明显差异.综合计算精度和效率，本文选取1 100万网格总数的方案.计算域及具体的网格划分（30°偏航）如图2所示.

2.3 边界条件及参数设置

计算域进口采用速度入口，按照B类地貌设置风速及湍流强度剖面，其中地面粗糙度指数为0.15，10 m参考高度处的基本风速为25 m/s，对应该风力机的切出风速，并通过用户自定义函数实现上述入流边界条件与FLUENT的连接（见图3）.数值计算采用3D单精度、分离式求解器，计算域边界条件如图3所示.亚格子尺度选用Smagorinsky-Lilly模型，压力速度耦合方程组求解采用SIMPLEC格式，压力项离散采用Standard格式，动力离散采用Bounded Central Differencing格式，瞬态方程采用二阶隐式，控制方程的计算残差设置为10-6.

2.4 有效性验证

根据不同偏航角下游塔架受上游叶片尾迹的干扰程度，将塔架分为未干扰区段（0 ～ 64 m）和显著干扰区段（64 ～ 124 m），定义不同偏航角下塔架正对来流风向的位置为迎风面0°处.图4给出了大型风力机塔架未干扰区段典型截面沿环向分布的平均与脉动压力系数分布曲线，并与规范及国内外实测

曲线[14-15]进行对比.

分析可得：1）大涡模拟与规范的平均风压系数沿环向分布规律基本一致，仅在背风区略小于规范值；2）脉动风压曲线包络在国内外实测曲线中，且沿环向分布规律一致，考虑到风力机塔架脉动风压分布与所处的地形、来流湍流和周边干扰密切相关，故认为本文基于大涡模拟得到的风荷载具有一定的有效性.

2.5 结果分析

图5和图6分别给出了塔架未干扰区段和显著干扰区段典型截面速度流线图，对比发现：1）在塔架未干扰区段，不同偏航角下来流均在迎风面0°处发生分流，在塔架两侧出现明显加速，在背风面产生回流及尺度较小的涡旋；2）在塔架显著干扰区段，不同偏航角下叶片的遮挡程度不同，导致塔架绕流出现显著差异.偏航角较小时，叶片的遮挡作用明显，来流在叶片处发生分流并在叶片和塔架之间形成尺度较大的涡旋，塔架两侧流体速度明显减小，塔架迎风面分离点发生偏移，使得尾迹流动产生偏向并出现尺度较大的涡旋；3）隨着偏航角的增大，叶片对塔架的遮挡程度降低，来流在塔架迎风面0°附近发生分离并在背风面形成尺度较小的涡旋，当偏航角达到45°时，塔架绕流特性与未受叶片偏航和干扰时较为接近.

3 风振响应分析

3.1 有限元建模与动力特性分析

基于ANSYS平台[16-17]建立大型风力机塔架-叶片一体化有限元模型，其中塔架和叶片采用SHELL63单元，机舱及其内部结构可作为整体采用梁单元BEAM189模拟，圆形筏基基础的单元类型为SOLID65，基础底端固结，地基与基础的作用采用弹簧单元COMBIN14模拟.通过多点约束单元耦合将各部件连接，形成风力机塔架-叶片一体化有限元模型.依据效率与精度均衡的原则，模型一共划分了4 122个单元.

图7和8分别给出了风力机频率分布曲线和典型模态振型图，由图看出：1）风力机塔架-叶片耦合模型的基频很小，仅为0.197 Hz，第10阶频率为

1.122 Hz，模态之间的间隔非常小；2）第5和10阶为3个叶片做前后挥舞运动，第30和50阶为叶片复杂的前后挥舞和左右摆动，并伴随着塔架的弯曲

变形.通过多阶振型分析表明，风力机塔架-叶片耦合模型低阶振型以叶片变形为主，随着振型频率的增加，叶片挥舞和摆动愈加显著，并伴随着塔架的共同变形.

3.2 塔架响应

径向位移为圆柱壳沿半径方向的位移，可直观反映圆柱壳结构的相对变形，其中负值表示径向朝内，正值表示径向朝外.图9给出了不同工况下塔架径向位移的均值示意图.可以发现：1）不同工况下塔架径向位移均值分布较为一致，不同偏航角主要对塔架中上部的位移产生影响，径向位移随塔架高度的增加逐渐增大，并在塔架迎风面和背风面范围内形成两个极值区域，最大正负位移均出现塔顶0°和180°处；2）随着偏航角的增大，塔顶径向位移均值的最大值逐渐减小，工况一径向位移均值为0.589 m，工况六为0.466 m，塔架位移极值区域呈现逐渐减小的趋势，其中工况一极值区域在塔架80 m以上，工况六在塔架100 m以上.

图10给出了不同工况下塔架径向位移均方差示意图，从图中可以看出徑向位移均方差随塔架高度的增加逐渐增大，并在塔架迎风面和背风面范围

内形成两个极值区域，均方差最大值出现在塔顶0°和180°处；随着偏航角度的增大，径向位移均方差极值区域逐渐减小，塔顶位移均方差最大值呈现逐

渐减小的趋势，工况一均方差最大为0.206 m，工况六均方差最小为0.096 m.

图11给出了不同工况下塔架底部弯矩特征值，由图可以发现不同工况下塔架底部径向和环向弯矩均值分布规律和数值基本一致，且在塔架0°和20°出现弯矩最小值和最大值；径向和环向弯矩均方差分布规律一致，但均方差数值在塔架0°、20°、90°和340°处出现显著差异，其中工况四均方差数值最大.

3.3 叶片响应

定义与塔架重合的叶片为叶片1，顺时针依次为叶片2和3，表3给出了不同工况下叶片顺风向位移特征值列表，其中峰值因子取2.5.由表可知：1）随着偏航角的增大，叶片1顺风向位移均值和极值均呈现逐渐减小的趋势，工况一位移极值最大为2.79 m，工况六位移最小为1.69 m；2）随着偏航角的增大，叶片2和3的顺风向位移均值、均方差和极值均呈现先减小后增大再减小的趋势，工况一叶片2和3位移极值最大分别为3.51 m和3.73 m；3）工况

一各叶片的顺风向位移极值均大于2.7 m，此时塔架与叶片之间的相互干扰作用最为明显.

表4给出了不同工况下叶片1根部内力特征值，由表可以发现随着偏航角的增大，叶片根部剪力与弯矩的均值逐渐减小，工况一内力均值最大，工况六内力均值最小，工况一和工况四的内力响应均方差较大.

4 风致稳定性分析

4.1 屈曲稳定

图12给出了不同工况下风力机塔架-叶片耦合体系屈曲最大位移和临界风速对比示意图.对比发现不同偏航角对风力机耦合体系屈曲失稳时屈曲最大位移和临界风速影响显著，随着偏航角的增大，屈曲最大位移和临界风速均呈现先减小后增大再逐渐减小的趋势，屈曲最大位移和临界风速最大值分别出现在工况三和工况一. 综合屈曲最大位移和临界风速对风力机体系屈曲稳定性能的影响发现，工况六为最易屈曲失稳工况，工况一屈曲稳定性能最好.

4.2 极限承载力

图13给出了不同工况下风力机位移随风速变化示意图，由图可以发现：1）随着风速的逐级加载，工况一、二和三的最大位移呈现先增大后骤减最后小幅增大的趋势，工况四、五和六呈现先减小后逐渐增大的趋势，工况四、五和六的最大位移显著大于工况一、二和三，工况六位移最大为0.05 m；2）随着风速的增大，工况一、二和三出现“逆向效应（最大位移骤减）”，而工况四、五和六未出现“逆向效应”，偏航角度越大，风力机体系的极限承载力越小.

5 结 论

结合大涡模拟技术和有限元完全瞬态法，对叶片最不利位置下考虑六个偏航角（0°、5°、10°、20°、30°和45°）的大型风力机体系的流场特性和风效应进行了数值模拟，提炼出了叶片偏航和干扰效应对大型风力机体系气动性能、风振响应和稳定性能的影响规律.

结果表明：0°偏航时，叶片和塔架之间的相互干涉作用最为明显，塔架顶部径向位移、叶片顺风向位移和叶片根部内力的均值及极值均为最大.随着偏航角的增大，塔架与叶片的相互干扰作用逐渐减小，塔架顶部径向位移、叶片顺风向位移和叶片根部内力的均值及均方差均呈现逐渐减小的趋势；同时，风力机体系的屈曲最大位移和临界风速均呈现先减后增再减的趋势，风力机体系的极限承载力逐渐降低.综合得到：0°偏航角下风力机体系气动性能和风振响应均最为不利，45°偏航角下风力机体系的稳定性能最为不利.

参考文献

[1] JEONG M S， KIM S W，LEE I， et al. The impact of yaw error on aeroelastic characteristics of a horizontal axis wind turbine blade[J]. Renewable Energy，2013，60（5）： 256—268.

[2] 曹九发，柯世堂，王同光.复杂工况下的大型风力机气动弹性响应和尾迹数值分析研究[J]. 振动与冲击，2016，35（1）：46—52.

CAO J F，KE S T，WANG T G.Numerical analysis for aero-elastic responses and wake of a large scale wind turbine under complicated conditions[J]. Journal of Vibration and Shock， 2016， 35（1）： 46—52. （In Chinese）

[3] KE S T， YU W，WANG T G， et al. Wind loads and load-effects of large scale wind turbine tower with different halt positions of blade[J]. Wind and Structures， 2016，23（6）： 559—575.

[4] 吕计男，刘子强，赵玲，等.大型风力机气动弹性响应计算研究[J].空气动力学学报， 2012，30（1）：125—129.

J N，LIU Z Q，ZHAO L，et al. Large-scale wind turbine aeroelastic responses analysis[J]. Acta Aerodynamica Sinica， 2012，30（1）： 125—129. （In Chinese）

[5] KONG C，BANG J，SUGIYAMA Y. Structural investigation of composite wind turbine blade considering various load cases and fatigue life[J].Energy，2005，30（11/12）：2101—2114.

[6] 陳小波，李静，陈健云.考虑叶片离心刚化效应的风力机塔架风振反应分析[J].工程力学， 2010，27（1）：240—245. （In Chinese）

CHEN X B，LI J，CHEN J Y.Wind-include response analysis of wind turbine tower considering centrifugal stiffening effect of blade[J].Engineering Mechanics，2010，27（1）：240—245.

[7] 袁行飞， 余俊伟. 屋顶安装型风力机塔架风振反应分析[J]. 浙江大学学报（工学版），2013，47（11）：1911—1916.

YUAN X F，YU J W.Wind-included response analysis of wind turbine tower installed on buiding roof[J]. Journal of Zhejiang University（Engineering Science），2013，47（11）：1911—1916.（In Chinese）

[8] 高乾丰，董辉，邓宗伟，等.大型风力机风雨结构三场耦合分析[J].中南大学学报（自然科学版），2016，47（3）：1011—1016.

GAO Q F，DONG H，DENG Z W，et al.Three-field coupling analysis for large-scale wind turbine with wind-rain-structure[J].Journal of Central South University（Science and Technology），2016，47（3）：1011—1016. （In Chinese）

[9] 杨阳，李春，叶柯华，等.基于HHT方法的非稳定工况风力机结构动态响应时频特性分析[J]. 振动与冲击，2016，35（21）：22—28.

YANG Y，LI C，YE K H，et al. Structural dynamic response characteristics of a wind turbine in time-frequency domain under non-stationary operating conditions based on HHT method[J]. Journal of Vibration and Shock， 2016， 35（21）： 22—28.（In Chinese）

[10] 李仁年，申慧渊.水平轴风力机叶片气动弹性稳定判别数的探索分析[J]. 太阳能学报， 2012， 33（12）： 2152—2158.

LI R N，SHEN H Y.Analysis of discrimination number of aeroelastic stability of horizontal axis wind turbine blade[J]. Acta Energiae Solaris Sinica，2012，33（12）：2152—2158.（In Chinese）

[11] JIMENEZ A， CRESPO A，MIGOYA E. Application of a LES technique to characterize the wake deflection of a wind turbine in yaw[J].Wind Energy，2010，13（6）：559—572.

[12] 柯世堂，曹九发，王珑，等.风力机塔架-叶片耦合模型风致响应时域分析[J].湖南大学学报（自然科学版），2014，41（4）：87—93.

KE S T，CAO J F，WANG L，et al.Time-domain analysis of the wind-included responses of the coupled model of wind turbine tower-blade coupled system[J]. Journal of Hunan University（Natural Sciences），2014，41（4）：87—93.（In Chinese）

[13] 石亚丽，左红梅，杨华，等.偏航角对风力机气动性能的影响[J]. 农业工程学报，2015，31（16）：78—85.

SHI Y L，ZUO H M，YANG H，et al.Aerodynamic performance of wind turbine under different yaw angles[J]. Transactions of the

Chinese Society of Agricultural Engineering， 2015， 31（16）：78—85. （In Chinese）

[14] 李晓娜，陆煜，刘庆宽，等.圆形截面烟囱风致干扰效应试验研究[J]. 工程力学， 2015，32（S1）：159—162.

LI X N，LU Y，LIU Q K，et al.Experimental study on wind-included interference effects of circular section chimneys[J].Engineering Mechanics，2015，32（S1）：159—162.（In Chinese）

[15] NISHIMURA H，TANIIKE Y. Aerodynamic characteristics of fluctuating forces on a circular cylinder[J]. Journal of Wind Engineering & Industrial Aerodynamics，2001，89（1）： 713—723.

[16] KE S T， WANG T G， GE Y J， et al. Wind-induced responses and equivalent static wind loads of tower-blade coupled large wind turbine system[J]. Structural Engineering and Mechanics， 2014，52（3）：485—505.

[17] 鄧露，肖志颖，黄民希，等.考虑流固耦合的近海风机动力响应数值计算[J]. 湖南大学学报（自然科学版）， 2015，42（7）： 1—8.

DENG L，XIAO Z Y，HUANG M X，et al.Numerical simulation of dynamic response for offshore wind turbines including fluid-structure interaction[J]. Journal of Hunan University（Natural Sciences），2015，42（7）：1—8.（In Chinese）