复杂荷载环境下海上风力机的建模及动力学特性分析

王青占,赵建中,郭兴明

(上海大学上海市应用数学和力学研究所,上海　200072)

复杂荷载环境下海上风力机的建模及动力学特性分析

王青占,赵建中,郭兴明

(上海大学上海市应用数学和力学研究所,上海200072)

通过分析处于海风、海浪、海流、土壤力等复杂荷载作用下的海上风力机支撑结构,采用赫维赛德阶跃函数和狄拉克δ函数建立了连续统一的、顶端带有集中质量块的悬垂梁风力机动力学模型.基于对控制方程的Galerkin截断,得到离散化的常微分方程组,使用四阶Runge-Kutta方法求解,得到了模型的动力学行为特性云图曲线.通过悬垂梁风力机模型的时程曲线、庞加莱映射对风力机模型进行动力学分析,给出了位移和速度的幅值随激振力频率变化的幅频特性曲线,并研究了垂向激振、自重、变刚度参数对风力机结构振动特性及其稳定性的影响.

海上风力机；变轴力；复杂荷载；变刚度；动力响应

近年来,各国对风能、太阳能等可替代的洁净新能源的开发利用力度加大,尤其是对海上风力机的研究越发引起各国科研工作者的重视[1].以英国、丹麦和德国为代表的欧洲走在了海上风力发电研究的前端,近两年的装机总量不断增加,发电量所占的比重相应提高.图1为海上风力机系统的示意图.海上风力机系统由风力机、支撑结构和地基基础三部分组成,其中风力机由叶片、轮毂和机舱构成；支撑结构包括塔筒和下部结构,下部结构分为固定式和漂浮式两种形式.无论支撑结构采用固定式还是漂浮式,海上风力机与陆地风力机、海上油气平台工程都有较大区别[2],在系统结构、环境条件、荷载特征等方面都具有特殊性.首先,海上风力机与塔架要经受台风的严峻考验,其下部结构还要受到波浪、海流以及与土体之间的耦合作用[3],且风、浪、流是有一定耦合的.国内外的相关学者和研究机构也开始对海上风电工程和其他浮体项目进行研究.南京航空航天大学胡文瑞等[4]进行了大型风力机的空气动力学基础研究,中国科学院力学所周济福等[5]针对海上风电工程和地基的关键力学问题进行了相关研究.国内一些学者通过有限元建模仿真计算对海上风力机进行气动力学、结构静力学和动力学的研究[6-9].国外学者主要对风力机结构进行建模分析[10-11],包括对不同基础、浮动风机以及海床结构的参数化模拟分析.但是考虑复杂海洋环境的综合作用,以及风力机结构刚度变化的报道较少.本工作建立的带有集中质量块的悬垂梁风力机动力学模型较好地考虑了这些因素,研究的结果也可为后续风力机的优化设计和安全评估提供重要参考.

图1　海上风力机系统示意图Fig.1　Sketch of offshore wind turbines

1　海上风力机的结构分析与简化模型

对于海上风力机模型的研究,采用较多的是有限元仿真计算.例如,用有限元进行非线性结构动力学、流固耦合或者土体与结构耦合[12]、叶轮气动力学研究等.王磊等[13]通过计算机仿真,结合水动力学模型和风轮空气动力学模型,建立了“风轮-机舱-塔筒-系泊系统”组成的多柔体系统动力学模型,并对漂浮式风力机系统和近海的定桩式风力机系统进行动力学分析对比,研究结果表明气动载荷与水动力相互耦合对整机结构动力响应及功率波动有着明显影响.空间太阳能电板、飞机机翼、海上风力机塔筒等大尺寸细长物体,在某种程度下都可以简化为梁模型.Piana[14]和Virgin[15]通过建立类似的弹性结构模型研究了上述结构的动态稳定性和振动之间的相互作用,进一步分析了结构承受类似于轴向压力的载荷(尽管可能不会直接导致结构失稳)对结构的固有频率产生的影响,并以此结果为参照对结构施加控制[16].

Lajimi等[17]用Frobenius幂级数法研究了质量块的质量、旋转惯性以及偏心距对顶端承受集中质量块固支悬垂梁的固有频率的影响,研究结果表明该结构的固有频率随着集中质量块质量的增加而降低.Auciello[18]研究了水对风力机模型固有频率的影响,结果表明水只对较高的固有频率产生显著影响,而对于通常考虑的前三阶固有频率的影响不大.U´sci′lwska等[19]建立了半浸水顶端承受集中质量的悬垂梁模型,并以水面为界限建立了两段梁的运动方程,将集中质量块归类于边界条件,求解六阶奇异矩阵得到了控制其固有频率的特征方程.Andersen等[20]通过弹簧模拟桩体与土壤之间的p-y曲线,研究了海床对风力机支撑机构的影响.

工程设计中海上风力机塔架的固有频率必须要避开风力机转子的自激频率.通过对海上风载以及海浪冲击荷载的观测和统计,找出风力机设计自然频率的安全频域,通过本工作的参数化计算和模拟,可以选择合适的风力机设计参数,提高风力机的效率和使用寿命.通过本工作所考虑的机舱重量和风力机工作动荷载进行调节和振动控制,也可以结合海上环境气候统计[21]分析来调节海上复杂荷载的组合参数,进而更加真实模拟复杂的海洋作用环境.综上所述,本工作所建立的模型和计算方法有很好的适应性,可用于不同海洋环境和风力机参数的普适计算.

2　悬垂梁风力机模型的动力学方程

2.1动力学一般方程的建立

海上风力机支撑结构模型如图2所示.整个舱重和叶轮用一集中质量块M表示,叶片承受的风载等效为作用在轮毂处的集中力f1(x,t),塔体所受的风载为f2(x,t),海浪潮汐的作用为f3(x,t),海流的作用为f4(x,t),海床岩土对塔体的作用为f5(x,t).同时引入随径向和时间变化的轴力,以及由填充钢筋混凝土引起的变刚度高度,依据能量法建立Euler-Bernoulli悬垂梁模型的运动方程[22]：移,M为塔筒顶端(机舱和叶片)质量,ε和σ分别为风力机振动的动质量参数和额定工作频率, l1为塔筒内部混凝土填充高度,m为单位长度空心塔筒的质量,mc为单位长度填充钢筋混凝土的质量,dc为塔筒内径,H(x)表示赫维赛德阶跃函数,δ(x)为狄拉克脉冲函数.

图2　海上风力机支撑结构模型Fig.2　Structure model of offshore wind turbines

2.2控制方程及边界条件的无量纲化处理

为方便计算,对悬垂梁风力机模型进行无量纲化处理.引入无量纲参数,令

式中,Υ为位移标定参量,v(ξ,τ)为横向振动位移w(x,t)的无量纲表示,ξ1为钢筋混凝土填充变量,ξ2为反映海域深浅的量,ξ3为反映海床环境的量,κ为混凝土刚度调节参数,κ1为混凝土与塔筒旋转惯性比参数,η1为塔筒与舱体的质量比(又称固有频率影响因子),η2为旋转惯性参数,κ3为混凝土质量填充参数,α为集中质量块惯性影响参数,τ为无量纲时间尺度,σ为风机激振频率,ζ为结构阻尼比.

风力机模型的无量纲控制方程为

边界条件的无量纲化为

2.3复杂荷载及其无量纲化处理

风力发电机组运行时,其叶片上的风荷载和风力机偏航引起的荷载通过结构和传动机构作用在塔架顶端,所以相关规范规定海上风电机组基础结构设计应考虑风电机组的荷载.这部分荷载包括风轮上的静风压引起的荷载、湍流和尾流引起的荷载、风力发电机偏航引起的荷载和风力发电机组的重力荷载等.动荷载包括风力机运转荷载,以及风、浪、流、冰和地震荷载.为了便于分析计算,分别将风轮荷载简化为轮毂集中力和风力机垂向激振荷载,将塔筒承受的风荷载简化为均匀脉动的荷载,将海面处承受的海浪等冲击荷载用脉冲荷载表示,土壤与塔基之间的作用力[12]用近似简化的p-y曲线表示.具体化简后的荷载表示如下.

海平面处的波浪、海冰船舶冲击荷载(集中)：

式中,p5(ξ,τ)=∆5H(ξ3−ξ)tanh(2πξ sin(2ξ))或∆5H(ξ3−ξ)tanh(2πξϕ1(2ξ)),A=1(与循环荷载有关),B与土壤质地(密度)有关,Pu为土体极限抗力.

2.4Galerkin截断选取的模态函数

模态函数为

2.5无量纲控制方程的Galerkin截断

利用选取的模态函数可求得方程(3)的三阶Galerkin离散方程：

式中,φi,j,k(i,j=1,2,3,k=1,2,…,8)与梁的刚度有关,是反映刚度变化的量.

2.6控制方程的矩阵化表示

为了从整体上更加清晰地了解该悬垂梁的计算模型,对控制方程进行如下的矩阵化归一表示

ξ1,ξ2,ξ3为反映风力机海洋环境和风力机刚度控制的设计参量,Λ矩阵清晰反映了风力机所承受的海洋环境荷载.如能获得更加详尽和真实的海洋环境和风力机的资料,理论上可以通过这些参数组合实现对风力机设计的操控和优化,为设计更加合理的风力机提供重要的分析依据.

3　特定条件下求解结果的对比

3.1Galerkin截断误差及有效性验证

为了验证Galerkin截断的收敛性和有效性,取自由振动的均匀悬臂梁进行对比验证.

图3为悬臂梁自由振动控制方程近似解析解与三阶Glerkin截断数值解的时空振动曲线.可以看出,二者整体对比的结果基本一致.图4为悬臂梁自由振动控制方程近似解析解与三阶Galerkin截断数值解的绝对误差和相对误差.可以看出,位移最大的顶端处的相对误差＜3%.因此,该模态下的Galerkin截断完全可以满足工程中的精度要求.

3.2数值结果与有特定解的解析解对比

图3　方程近似解析解与三阶Galerkin截断数值解的时空振动曲线Fig.3　Vibration curves of approximate analytical and numerical solution of 3-order Galerkin truncation

图4　方程近似解析解与三阶Galerkin截断数值解的绝对误差和最大相对误差Fig.4　Absolute and maximum relative errors of approximate analytical and numerical solution of 3-order Galerkin truncation

由此求得无量纲控制方程的解析解为w(ξ,τ)=−cosτ(1+ξ−0.679 7e-0.9934ξ−0.3463e0.9909ξ+e-0.00125ξ(0.0209cos(1.0079ξ)−1.3216sin(1.0079ξ))).此处w(ξ,τ)与v(ξ,τ)类似,均为相同单位标定的无量纲量.

图5　方程数值结果与解析结果的对比Fig.5　Comparisons between numerical and analytical solution for the equation

图5为悬臂梁自由振动控制方程数值结果与解析结果的对比.可以看出,两种方法得出的结果吻合较好.产生误差的原因可能是由于模型顶端集中质量块处边界条件的近似处理,这也从侧面反映了机舱的重量对其振动特性有显著的影响.优化机舱重量或者在机舱中设计减震与控制装置[16],可以提高风力机的使用年限和安全系数.

图6为特殊荷载下悬垂梁风力机模型的幅频特性曲线.根据文献[24]的海上风力机频域设计区域理论,可以看出,在强有力的海上荷载作用下,该参数下的风力机处于Soft-Stiff区域,是比较安全的设计.计算结果表明，风力机在频率fo=0.37 Hz(转化为无量纲参数λo=2.375 6)的惯常荷载作用下是安全的.

图6　特殊荷载下风力机模型的幅频特性曲线Fig.6　Amplitude-frequency curve of the wind turbine model under specific loads

4　数值仿真算例

本工作以文献[24]中给出的海上风力机支撑结构为算例,取结构额定功率为3 MW,塔筒主体材料为Q345钢且塔筒内部填充了钢筋混凝土.

4.1原始物理参数

将以上物理量代入与其对应的无量纲表达式(3)中,并将所得的结果代入式(7)～(9)中,再通过式(10)计算得到控制方程的系数矩阵(刚度矩阵)和荷载矩阵：

4.2动力学数值分析

通过四阶Runge-Kutta数值方法求解矩阵微分方程组,可得到稳态时前三阶动力学时程曲线、时空响应曲线(见图7),可以看出u2与u3重合.这些基本曲线反映了梁的动力学特性[25].

图7　稳态下海上风力机顶端处的时程曲线和幅频曲线Fig.7　Time-history and amplitude-frequency curves in stationary state for the top end of offshore wind turbines

从图7(b)可以看出,危险频段为1.8～2.2,折算为实际的频率为0.28～0.34Hz.如果荷载频率在此频段内,说明设计不合理,就需要调节其他影响固有频率的参数.从梁受迫振动的前三阶时程曲线(见图7(a))、相图(见图8)和庞加莱映射(见图9)可以看出,一阶截断存在较大的误差,而二、三阶截断误差要小得多.因此用Galerkin截断来计算梁的动力曲线是收敛的.从图9中可以看出,庞加莱映射是一极限环状的,可见在此参数下梁的运动是准周期的、稳定的响应.

图8　顶端处的稳态前三阶Galerkin截断相图Fig.8　Steady-state first 3-order Galerkin truncation phase diagram at top end

图9　顶端处的稳态前三阶Galerkin截断庞加莱映射Fig.9　Steady-state first 3-order Galerkin truncation Poincar´e map at top end

4.3变刚度参数影响分析

图10为载荷条件相同时无填充和半填充钢筋混凝土塔筒情况下风力机结构的动力特性.

图10　荷载条件相同时无填充和半填充钢筋混凝土塔筒的动力特性Fig.10　Dynamic characteristics of no-filled and half-filled wind turbines under the same load

从图10可以看出,填充混凝土后风力机结构的最大振幅(0.165个标准)要比无填充时(0.465个标准)小得多.因此,塔筒填充钢筋混凝土对风力机结构的振动特性有着显著的影响,通过调节海上风力机塔架设计的固有频率,可以避开海上荷载的常见频段组合.

5　结束语

本工作建立了海上风力机支撑结构的悬垂梁动力学模型.研究结果表明,该模型可以较好地模拟风力机支撑结构的工作现状.通过分析计算复杂载荷下的结构动力特性,研究各种荷载、风力机安装环境参数、混凝土的填充高度、适宜海域的水深和海床条件等对风力机安全性的影响,不仅可以验证风力机在特定海洋环境的安全性能,同时也可以对风力机的参数进行优化设计.相对于有限元计算,这种数值计算模型大大缩减了优化的时间成本；相对于模拟试验,较好地降低了试验成本；而与解析条件下的情况对比,则大大提高了适应求解的范畴.

[1]林鹤云,郭玉敬,孙蓓蓓.海上风电的若干关键技术综述[J].东南大学学报(自然科学版),2011,41(4)：882-888.

[2]尚景宏.海上风力机基础结构设计选型研究[D].哈尔滨：哈尔滨工程大学,2010：11-19.

[3]李益,张宏战,马震岳.基于土-结构相互作用的海上风机模态分析[J].水电与新能源,2013(1)：70-73.

[4]胡文瑞,王同光,李晔.“大型(兆瓦级)风力机的空气动力学问题”专刊前言[J].应用数学和力学, 2013,34(10)：1-2.

[5]周济福,林毅峰.海上风电工程结构与地基的关键力学问题[J].中国科学：物理学力学天文学, 2013,43(12)：1590-1593.

[6]陈前,付世晓,邹早建.海上风力发电机组支撑结构动力特性分析[J].振动与冲击,2012,31(2)：86-90.

[7]顾岳飞.大兆瓦风电机组塔架的有限元分析与优化设计[D].上海：上海交通大学,2012：18-45.

[8]李立岗.风力发电机塔架结构动力特性分析[D].哈尔滨：哈尔滨工业大学,2011：21-29.

[9]严磊.风力发电机支撑体系结构设计研究[D].天津：天津大学,2008：37-41.

[10]LEBLANc C.Design of offshore wind turbine support structures[D].Aalborg：Aalborg University, 2009：16-27.

[11]AGARwAL P.Structural reliability of offshore wind turbines[D].Texas：The University of Texas at Austin,2008.

[12]LOMBARDI D,BHATTAcHARyA S,WOOD D M.Dynamic soil-structure interaction of monopile supported wind turbines in cohesive soil[J].Soil Dynamics and Earthquake Engineering,2013, 49：165-180.

[13]王磊,何玉林,金鑫,等.漂浮式海上风电机组动力学仿真分析[J].中南大学学报(自然科学版), 2012,43(4)：1309-1314.

[14]PIANA G.Vibrations and stability of axially and transversely loaded structures[D].Turin：Polytechic University Turin,2013.

[15]VIRGIN L N.Vibration of axially-loaded structures[M].Cambridge：Cambridge University Press, 2007.

[16]HE W,GE S S.Vibration control of a nonuniform wind turbine tower via disturbance observer[J].IEEE/ASME Transactions on Mechatronics,2014,20(1)：237-244.

[17]LAjIMI S A M,HEppLER G R.Free vibration and buckling of cantilever beams under linearly varying axial load carrying an eccentric end rigid body[J].Transactions of the Canadian Society for Mechanical Engineering,2013,37(1)：89-109.

[18]AUcIELLO N M.Dynamic analysis of immersed tapered column carrying an eccentric mass at the free[J].Advanced Research in Scientific Areas,2013(1)：412-417.

[19]U´ScI＇LOwSKA A,KO＇LODzIEj J A.Free vibration of immersed column carrying a tip mass[J]. Journal of Sound and Vibration,1998,216(1)：147-157.

[20]ANDERSEN L V,VAHDATIRAD M J,SIcHANI M T,et al.Natural frequencies of wind turbines on monopile foundations in clayey soils—a probabilistic approach[J].Computers and Geotechnics, 2012,43(11)：1-11.

[21]BERAN J,CALvERI L,LANGE B,et al.Offshore wind modelling and forecast[C]//Proceedings of the 6th WRF/15th MM5 Users’Workshop.2005：1-9.

[22]ADHIKARI S,BHATTAcHARyA S.Vibrations of wind-turbines considering soil-structure interaction[J].Wind and Structures,2011,14(2)：85-112.

[23]刘延柱,陈文良,陈立群.振动力学[M].北京：高等教育出版社,1998：130-135.

[24]ADHIKARI S,BHATTAcHARyA S.Dynamic analysis of wind turbine towers on flexible foundations[J].Shock and Vibration,2012,19(1)：37-56.

[25]丁虎,严巧赟,陈立群.轴向加速运动黏弹性梁受迫振动中的混沌动力学[J].物理学报,2013, 62(20)：200502.

本文彩色版可登陆本刊网站查询：http://www.journal.shu.edu.cn

Dynamics modeling and analysis of offshore wind turbines under complicated loads

WANG Qingzhan,ZHAO Jianzhong,GUO Xingming
(Shanghai Institute of Applied Mathematics and Mechanics,Shanghai University, Shanghai 200072,China)

By analyzing the structure and complexity of offshore wind turbines(OWTs), a beam model with free end carrying a concentrated mass body was built to analyze dynamical characteristics of OWTs subject to adverse working environment loads. Factors including vertical excitation,self-weight,rotary inertia and variety of stiffness were considered to study the effect on offshore wind turbines.A general and uniform governing equation was built by taking advantages of the Heaviside step function and Dirac delta function.Based on the Galerkin truncation and the Runge-Kutta time discretization, numerical solutions of the kinematic governing equation were obtained.By comparing with the analytical solution under specified conditions,validity of the method was checked.The time history of the beam’s free end was chosen to represent motion of the beam.Based on the steady time history of the beam,a Poincar´e map was constructed to study its periodic motion.Furthermore,an amplitude-frequency curve was given to find the dangerous frequency range where OWTs exist.

offshore wind turbines；varying axial force；complicated load；varying rigidity；dynamic response

TH 43

A

1007-2861(2016)05-0573-13

10.3969/j.issn.1007-2861.2015.02.016

2015-06-16

国家重点基础研究发展计划(973计划)资助项目(2014CB0462003)

郭兴明(1964—),男,教授,博士生导师,博士,研究方向为连续介质力学与力学中的数学方法. E-mail:xmguo@shu.edu.cn