基于大涡模拟考虑叶片停机位置大型风力机风振响应分析

柯世堂， 余 玮， 王同光

(1.南京航空航天大学 土木工程系，南京 210016；2.南京航空航天大学 江苏省风力机设计高技术研究重点实验室，南京 210016)

基于大涡模拟考虑叶片停机位置大型风力机风振响应分析

柯世堂1,2， 余 玮1,2， 王同光2

(1.南京航空航天大学 土木工程系，南京 210016；2.南京航空航天大学 江苏省风力机设计高技术研究重点实验室，南京 210016)

为研究停机状态下不同叶片位置对大型风力机塔架-叶片体系风振响应的影响，以某3 MW大型水平轴三叶片风力机为研究对象。基于大涡模拟(LES)方法对叶片八个不同停机位置下的风力机体系流场和气动力性能进行了数值模拟，并通过与国内外实测数据对比验证了该方法的有效性；结合有限元法对考虑叶片不同停机位置的风力机塔架-叶片耦合模型进行了动力特性和风振响应时域分析。主要结论为：不同叶片停机位置对风力机塔架绕流特性和气动力分布影响显著，塔架中上部风振响应受叶片不同停机位置的影响最大，尤其是迎风面和背风面的脉动位移均方差响应较大，相应最大值出现在工况八下塔顶350°位置；最大塔底弯矩出现在工况二的环向330°位置；工况五下三叶片叶尖顺风向位移响应峰值>2.5 m；研究表明在进行大型风力机抗风设计时应考虑不同叶片停机位置的影响。

大型风力机体系；停机状态；大涡模拟；风振响应；气动力分布

大型风力机属于典型的风敏感结构，其塔架-叶片耦合结构具有风场复杂、阻尼小、自振频率低且分布密集等特性，风荷载往往成为其设计的主要控制荷载[1-2]。由于作用在风力机的外部载荷主要是由风况条件决定，高风速条件下风力机将处于停机状态，叶片不同停机位置将引起风力机体系表面流场和气动力分布的差异从而导致风振响应不同。因此，对考虑叶片不同停机位置下大型风力机体系风振响应特性进行研究可为设计人员进行高风速下风力机抗风设计提供参考。

目前，针对大型风力机的风致动力响应分析，王振宇等[3]基于Shinozuka理论，采用谐波叠加法分别模拟良态风和台风时程，对良态风和台风作用下的风力机进行动力时程分析，研究发现台风风向突然偏转90°将会使塔架发生屈服破坏；陈小波等[4]分别分析了屋顶集风效应和离心刚化效应对风力机塔架风振反应的影响，结果表明屋顶集风效应增强了塔架风振反应而离心刚化效应减小了塔架风振反应；KWON等[5]基于风力机塔架-叶片耦合模型提出了一种新的阵风因子计算方法，并通过算例分析表明阵风荷载效应对于大型风力机塔架体系的局部等效静风荷载的分布的影响显著；MURTAGH等[6]采用了时域强迫振动分析方法对风力机耦合体系进行了顺风向响应分析，结果表明不考虑塔架叶片相互影响会导致塔顶的风致响应减小；董辉等[7]基于FLUENT软件对2 MW水平轴风力机在额定工况下遭受暴雨时的风场和风驱雨进行分析，研究了风雨共同作用下风力机风雨荷载的大小和分布规律，结果表明雨滴对塔架的冲击作用不可忽略；金志昊等[8]对风力机叶片采用时域和频域方法进行了顺风向风致振动计算，分析了暴雨停机状态下叶片极限位移最不利情况。已有研究大多是针对叶片固定工况下大型风力机气动力数值模拟和风振响应的分析，鲜有系统考虑叶片不同停机位置下风力机塔架-叶片体系的风振响应特性研究。

鉴于此，本文以南京航空航天大学自主研发的某3 MW水平轴三叶片风力机为研究对象，首先基于大涡模(Large Eddy Simulation, LES)拟法对停机状态下风力机体系8种计算工况(由叶片旋转全过程状态下和塔架的相对位置确定)进行了数值模拟，获取了塔架和叶片表面典型测点的气动力时程作为后续风振响应分析的输入参数；然后结合有限元法对考虑叶片不同停机位置的风力机塔架-叶片耦合模型进行了动力特性和风振响应时域分析，最终对比研究了不同停机位置对于大型风力机体系流场特性、气动力分布和风振响应的影响，相关结论对大型风力机抗风的精细化设计具有参考意义。

1 动力特性分析

以某3 MW风力发电机组为研究对象，其中塔架高85 m，顶部直径2.0 m，底部直径2.5 m，塔架为通长变厚度结构，顶壁厚30 mm，底壁厚60 mm；机舱长12.0 m，宽4.0 m，高4.0 m；风轮切入风速3.5 m/s，额定风速12.5 m/s，切出风速25 m/s，风轮倾角5°，风轮转速范围为9～19 r/min，偏航速度为0.5度/秒；各叶片之间成120°夹角，沿周向均匀分布，叶片长度为44.5 m。

根据叶片与塔架的相对位置，并考虑到三叶片体系旋转过程中存在的周期性，以叶片与竖直方向夹角为0°作为初始状态，依次顺时针旋转15°设置8个计算工况。基于ANSYS平台分别建立不同工况下风力机塔架-叶片一体化有限元模型，其中塔架和叶片采用shell63单元，机舱采用beam188梁单元模拟，忽略其内部结构，圆形筏基基础的单元类型为solid65，尺寸为直径24 m，高度2 m，基础底端固结，地基与基础的作用采用combin14弹簧单元模拟。风力机各部件通过多点约束单元耦合形成塔架-叶片一体化有限元模型[9-10]。依据效率与精度均衡的原则，整机模型共划分3 473个单元。不同计算工况下有限元模型见图1。

图1 不同计算工况有限元模型示意图

表1给出了采用分块Lanczos方法计算提取的风力机典型阶的固有频率，可以看出叶片和塔架的耦合效应明显降低了风力机系统的基频，不同计算工况下风力机系统基频均很低，仅为0.209 Hz左右，且各模态之间频率间隔较小。叶片不同停机位置主要影响风力机体系的低阶固有频率，各工况低阶频率出现较小的差别，高阶频率一致。可见停机状态下叶片不同位置对风力机体系的频率和振型影响微弱，低阶振型主要以叶片带动机舱及塔架进行前后挥舞和左右摆动，高阶模态出现塔架及叶片本身的结构变形和失稳形态。

表1 不同计算工况下风力机典型阶固有频率列表

2 数值模拟

2.1 计算域及网格划分

为保证流动能够充分发展，计算流域取12D×5D×5D(流向x×展向y×竖向z，D为风轮直径)。风力机置于距离计算域入口3D，从而保证尾流的充分流动。由于叶片表面扭曲复杂故采用混合网格离散形式，将整个计算域分为内外两个部分：核心区域采用四面体网格，并对风力机周围局部网格进行加密，外围区域采用高质量六面体结构网格，网格总数为795万，计算域及网格具体划分见图2。

2.2 边界条件及计算方法

计算域进口采用速度入口，按照B类地貌设置风速剖面以及湍流强度剖面，其中湍流强度定义式为

(1)

式中：地面粗糙指数为0.15，10 m参考高度处的基本风速为25 m/s。

对应该风力机的切出风速，并通过用户自定义函数实现上述入流边界条件与FLUENT的连接(见图3)；计算域出口设置为压力出口，相对压力为0；计算域

地面以及风力机表面采用无滑移壁面，计算域两侧面和顶面采用对称边界条件。

数值计算采用3D单精度、分离式求解器，由于风力机所处流场属于非定常且扰流情况复杂，基于大涡模拟法能够对风力机复杂的流场进行更好地模拟[11]。亚格子尺度选用Smagorinsky-Lilly模型，压力速度耦合方程组求解采用SIMPLEC格式，压力项离散采用standard格式，动力离散采用BCD(Bounded Central Differencing)格式，瞬态方程采用二阶隐式，控制方程的计算残差设置为1×10-6，以模型高度和模型高度处的平均风速计算时间步长定为0.001 s。

2.3 结果分析及有效性验证

图4和图5分别给出了不同工况下风力机体系典型截面的涡量分布以及速度流线图。由图4和图5可知，叶片和塔架之间的相互干涉作用将引起风力机体系大范围的涡量增值区域，且叶片不同停机位置引起的涡量分布差异明显，叶片与塔架相对位置越接近时塔架背后形成的尾涡区域越明显；同时叶片与塔架的相对位置不同将引起叶片引起的尾涡不同，塔架的存在破坏了叶片尾涡，而叶片尾涡也将改变塔架涡的形状；来流由叶片流经塔架时，塔架周围形成了明显的速度增大区，并且随着叶片与塔架之间的位置靠近程度增加，塔架周围流体加速越显著。

(a) 网格区域划分

(b) x-y平面网格分布

(c) x-z平面网格分布

图3 平均风速剖面及湍流强度分布Fig.3 Vertical profiles of wind speed and turbulence intensity

风力机体系表面共设置126个测点，其中塔架沿子午向布置8层测点，每层环向12个测点，各叶片沿翼展方向布置10个测点[12-13]，基于以上大涡模拟法获得风力机体系表面风压系数时程曲线，以作为风荷载时程输入参数。限于篇幅，图6仅给出了塔架中上部以及叶片典型测点的压力系数时程曲线，塔架迎风面呈正压分布，塔架两侧及背风面呈现出负压分布且两侧的负压明显大于塔架背风面。由于叶片的阻挡作用使得来流出现撞击、分离、附着以及漩涡脱落等特征湍流的情况，同时塔架与叶片之间的相互干扰作用增大使得塔架表面压力系数的均方差较大。

图4 不同工况风力机塔架典型截面涡量分布图Fig.4 Vorticity distribution of the wind turbine tower under different working conditions

图5 不同工况风力机塔架典型截面速度流线图Fig.5 Velocity contour distribution with streamlines of the wind turbine tower under different working conditions

图6 风力机体系典型测点压力系数时程曲线Fig.6 The history curve of pressure coefficient of wind turbine typical points

图7给出了风力机塔架下部受叶片干扰影响较小的典型截面沿环向分布的平均压力系数和脉动风压系数，并与规范以及国内外实测分布曲线[14-15]进行对比，由图7可知，大涡模拟与规范的平均风压系数沿环向分布规律基本一致，仅在背风区数值略小于规范值；脉动风压分布曲线分布在国内外实测曲线中间，且沿环向分布规律较为接近，考虑到脉动风压分布与实测塔所处的地形、来流湍流和周边干扰密切相关，本文大涡模拟获得的脉动风压分布趋势和数值均在已有实测结果包络中，对比结果也验证了数值模拟的有效性，可用于后续的风振响应时域分析。

(a) 平均风压

(b) 脉动风压

3 风振响应分析

风力机塔架-叶片耦合模型风振响应时域计算基于ANSYS软件平台，采用大涡模拟得到的风压系数作为风荷载时程输入参数，采用隐式方法Newmark和HHT来求解瞬态动力学平衡方程，其中各模态阻尼均为0.02，积分时间步长取为0.5 s，加载时间步数为1 024 步。

3.1 塔架响应

图8给出了不同工况下塔架径向位移响应均值的分布图，其中负值表示径向朝内，正值表示径向朝外。由图8可知，不同工况下塔架径向位移均值分布规律较一致，叶片不同位置对塔架位移主要影响塔架中上部。径向位移随着塔架高度的增加而逐渐增大，最大正负位移均出现在塔顶处350°和170°位置，以工况六位移最大为0.193 m，工况二位移最小为0.133 m；不同计算工况下塔架位移极值在迎风面和背风面的分布范围不同，极值区域随着叶片与塔架位置相对干涉作用增大而增大，工况六极值区域最大在塔架50 m以上，工况二极值区域最小在塔架70 m以上。

图9给出了不同计算工况下塔架位移脉动响应的均方差分布图，由图9可知，① 塔架位移脉动响应的均方差随着塔架高度增加而逐渐增加；② 塔架中上部叶片与塔架之间的相互干涉作用显著，由于叶片和塔架背风面产生涡旋的脱落，同一高度环向分别在塔架迎风面和背风面范围出现两个峰值区域；③ 脉动响应均方差极值区域最小出现在工况三(>80 m)，最大在工况六和工况八(>55 m)；④ 位移响应均方差最大均出现在塔顶350°位置，工况八数值最大为0.292，工况三最小为0.134。

图8 不同计算工况塔架径向位移响应均值Fig.8 The radial displacement of mean response of tower under different calculate working conditions

为对比不同叶片位置引起塔架底部弯矩值的差异，以工况一作为初始状态，将各工况下塔底内力响应值与之作差，图10给出了不同计算工况下塔底内力响应的特征值。由图10可知：① 塔底径向弯矩均值在0°、60°、90°以及330°位置出现差异，其中工况五下较明显，工况二0°位置数值相差最大达到256.93 N·m；② 塔底环向弯矩均值以工况二和工况五的差异较明显，数值相差最大位于工况二330°；③ 塔底弯矩脉动响应均方差以工况六和工况八下0°和330°位置脉动增加更显著。

图9 不同计算工况塔架径向响应均方差Fig.9 RMS of radial displacement response of tower under different calculate working conditions

图10 不同计算工况塔架底部内力响应特征值Fig.10 Characteristic force responses of tower bottom under different calculate working conditions

3.2 叶片响应

表2给出了不同计算工况下叶根弯矩与剪力的特征值，由表可以发现叶片不同停机位置引起了叶根内力特征值的显著差异，剪力与弯矩均值均随着塔架与叶片相对位置的靠近呈现减小的趋势，以工况一内力均值最大，工况五内力均值最小；工况二和工况八内力脉动响应的均方差较大。

定义工况一中与竖直方向夹角为0°的叶片编号为1，顺时针旋转编号叶片2、叶片3。表3给出了不同计算工况下风力机各叶片叶尖顺风向位移响应的峰值，其中峰值因子取为2.5。

表2 不同计算工况叶根内力特征值列表

由表3可知：① 工况五各叶片极大负值均>-2.50 m，工况三各叶片的极值位移较小；② 其余工况下均出现单个叶片的极大负值>-2.00 m，其中最大负值出现在工况八叶片3叶尖处；③ 叶片2不同停机位置下会对塔架产生遮挡效应，以工况五下叶片遮挡效应最显著，极大负值位移达到-2.530 m，此时叶片与塔架之间相互干扰作用最大。因此，叶片不同位置导致叶尖位移出现明显差异，为保证停机状态下叶片与塔架之间的安全性，应考虑最不利工况下叶尖的位移。

表3 不同计算工况叶尖顺风向位移峰值列表

4 结 论

基于大涡模拟和有限元法，对叶片不同停机位置下风力机塔架-叶片体系的流场作用、气动力分布和风振响应特性进行了系统研究，主要结论如下：

(1)停机状态下叶片不同位置对风力机体系的频率和振型影响微弱，低阶振型主要以叶片带动机舱及塔架进行前后挥舞和左右摆动，高阶模态出现塔架及叶片本身的结构变形和失稳形态。

(2)叶片和塔架之间的相互干涉作用将引起风力机体系大范围的涡量增值区域，叶片与塔架相对位置越接近时塔架背后形成的尾涡区域越明显。

(3)叶片不同停机位置显著影响塔架中上部的位移响应，其中工况六的位移响应均值最大，各工况下位移最大值均出现在塔顶170°和350°处，脉动响应均方差在迎风面和背风面出现峰值区域，最大值出现在工况八塔顶350°处。

(4)塔底弯矩响应特征值在环向0°和330°处差异显著，其中均值以工况二和工况五最明显，脉动响应均方差以工况六和工况八差异最显著。

(5)叶根剪力与弯矩均值随着塔架与叶片相对位置的靠近呈现减小的趋势，工况五中三叶片叶尖位移峰值均>-2.50 m，叶尖最大位移峰值出现在工况八。

[ 1 ] 柯世堂，王同光. 偏航状态下风力机塔架-叶片耦合结构气弹响应分析[J]. 振动与冲击, 2015, 34(18): 33-38. KE Shitang， WANG Tongguang. Aero-elastic responses analysis of wind turbine tower-blade coupled model in yaw condition[J]. Journal of Vibration and Shock, 2015, 34(18): 33-38.

[ 2 ] KE S T, GE Y J, WANG T G, et al. Wind field simulation and wind-induced responses of large wind turbine tower-blade coupled structure[J]. The Structural Design of Tall and Special Buildings, 2015, 24(8): 571-590.

[ 3 ] 王振宇, 张彪, 赵艳, 等. 台风作用下风力机塔架振动响应研究[J]. 太阳能学报, 2013, 34(8): 1434-1442. WANG Zhenyu, ZHANG Biao, ZHAO Yan， et al. Dynamic response of wind tuebine under typhoon[J]. Acta Energiae Solaris Sinica， 2013, 34(8): 1434-1442.

[ 4 ] 陈小波, 李静, 陈健云. 考虑叶片离心刚化效应的风力机塔架风振反应分析[J]. 工程力学, 2010, 27(1): 240-245. CHEN Xiaobo, LI Jing, CHEN Jianyu. Wind-include response analysis of wind turbine tower considering centrifugal stiffening effect of blade[J]. Engineering Mechanics， 2010, 27(1): 240-245.

[ 5 ] KWON D K, KAREEM A, BUTLER K. Gust-front loading effects on wind turbine tower systems[J]. Journal of Wind Engineering & Industrial Aerodynamics, 2012, 104(3):109-115.

[ 6 ] MURTAGH P J, BASU B, BRODERICK B M. Along-wind response of a wind turbine tower with blade coupling subjected to rotationally sampled wind loading[J]. Engineering Structures, 2005, 27(8):1209-1219.

[ 7 ] 董辉, 高乾丰, 邓宗伟, 等. 大型风力机风雨荷载特性数值研究[J]. 振动与冲击, 2015, 34(15): 17-22. DONG Hui, GAO Qianfeng, DENG Zongwei, et al. Numerical simulations for wind and rain loads of large-scale wind turbines[J]. Journal of Vibration and Shock, 2015, 34(15): 17-22.

[ 8 ] 金志昊, 范宣华, 苏先樾,等. 风力机叶片顺风向风致振动研究[J]. 南京航空航天大学学报, 2011, 43(5): 677-681. JIN Zhihao, FAN Xuanhua, SU Xianyue, et al. Research on blade wind-included vibration in wind direction[J]. Journal of Nanjing University of Aeronautics & Astronautics， 2011, 43(5): 677-681.

[ 9 ] 柯世堂, 曹九发, 王珑, 等. 风力机塔架-叶片耦合模型风致响应时域分析[J]. 湖南大学学报(自然科学版), 2014, 41(4): 87-93.KE Shitang, CAO Jiufa, WANG Long， 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 Science)， 2014, 41(4): 87-93.

[10] 柯世堂, 王同光, 陈少林, 等. 大型风力机全机风振响应和等效静力风荷载[J]. 浙江大学学报(工学版), 2014, 48(4): 686-692. KE Shitang, WANG Tongguang,CHEN Shaolin, et al. Wind-included responses and equivalent static wind load of large wind turbine system[J]. Journal of Zhejiang University(Engineering Science), 2014, 48(4): 686-692.

[11] LUO Kun, ZHANG Sanxia , GAO Zhiying, et al. Large-eddy simulation and wind-tunnel measurement of aerodynamics and aeroacoustics of a horizontal-axis wind turbine[J]. Renewable Energy, 2014, 77: 351-362.

[12] 彭文春, 邓宗伟, 高乾丰, 等. 风机塔筒流固耦合分析与受力监测研究[J]. 工程力学, 2015, 32(7): 136-142. PENG Wenchun, DENG Zongwei, GAO Qianfeng, et al. Fluid-solid interaction analysis and stess monitoring research of wind turbine tower[J]. Engineering Mechanics, 2015, 32(7): 136-142.

[13] 黄帅, 宋波, 贺文山, 等. 考虑脉动风影响的风电塔风致动力响应数值分析与现场监测比较[J]. 土木工程学报, 2012(增刊1): 102-106. HUANG Shuai, SONG Bo, HE Wenshan, et al. Numerical analysis of dynamic response and in-site monitoring of wind power tower considering the influence of the fluctuating wind load[J]. China Civil Engineering Journal, 2012(Sup1): 102-106.

[14] 李晓娜,陆煜,刘庆宽,等. 圆形截面烟囱风致干扰效应试验研究[J]. 工程力学, 2015(增刊1): 159-162. LI Xiaona, LU Yu, LIU Qingkuan, et al. Experimental study on wind-included interference effects of circular section chimneys[J]. Engineering Mechanics, 2015(Sup1): 159-162.

[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.

Wind-included vibration response analysis for a large wind turbine blade-tower systembased on large eddy simulation

KE Shitang1，2, YU Wei1，2, WANG Tongguang2

(1.Department of Civil Engineering, Nanjing University of Aeronautics and Astronautics, Nanjing 210016, China；2.Key Laboratory of Jiangsu Province High Tech Design of Wind Turbine, Nanjing University ofAeronautics and Astronautics, Nanjing 210016, China)

In order to study the effects of different blade positions on wind-included vibration responses of a large wind turbine tower-blade system under a stopped status, a 3 MW wind turbine with 3 blades developed by Nanjing University of Aeronautics and Astronautics was studied. First of all, the effects of eight different blade positions on the flow field and aerodynamic performances of the wind turbine system under stopped status were simulated with the large eddy simulation(LES). In order to verify the effectiveness of the LES method, the simulation results were compared with the actual test data of home and abroad. Then the finite element method was used to analyge dynamic characteristics and wind indaced vibration responses of the wind turbine tower-blade coupled model with different blade positions under a stopped status. The results showed that different blade positions have obvious effects on the flow characteristis of the wind turbine tower frame and the aerodynamic force distribution;the effects of different blade positions on the wind-induced vibration responses of the upper-middle part of the tower are the maximum;the mean square deviation between fluctuating displacements of windward side and leeward side is larger,its maximum value occurs at the tower top 350°position under the operation condition of 8;the maximum bending moment of the towerbottom appears at the annular 330°postion under the operation condition of 2;the peak displacements of 3 blade tips along the wind reach more than 2.5 m under the operation condition of 5.The study indicated that the effects of different blade positions under a stopped status should be considered for the anti-wind design of large wind turbine systems.

large wind turbine system；stopped status；large eddy simulation；wind-included vibration response；aerodynamic force distribution

国家973计划项目(2014CB046200)；中国博士后科学基金(2015T80551)

2015-10-28 修改稿收到日期：2016-02-25

柯世堂 男，博士，副教授，1982年生

TK83

A

10.13465/j.cnki.jvs.2017.07.014