基于模態(tài)疊加法的大型風力機典型工況動態(tài)特性分析

曹 莉， 孫文磊， 周建星

(新疆大學 機械工程學院，烏魯木齊 830047)

由于全球性能源短缺和環(huán)境污染等問題，風能作為一種清潔的可再生能源，變得越來越重要[1]。風力發(fā)電機組常年工作在變速、變載荷的工作環(huán)境中，在隨機風的作用下會受到頻繁的擾動和動載激勵，對其工作性能和使用壽命有很大的影響[2]。隨著風電機組不斷向大型化、大功率方向發(fā)展，風力發(fā)電機組的動態(tài)特性更加復雜。因此，在風力機的設計過程中，需要考慮風力機在典型工況下的動態(tài)特性。

當前，在風力機動力學分析領域，常用的方法有多體動力學方法MBS (Multibody Systems)、有限元方法FES(Finite Element System)、模態(tài)分析方法及連續(xù)系統(tǒng)COS(Continuous Systems)等[3]，MBS多體動力學方法將實際的機械構件視為剛體來建立系統(tǒng)動力學微分方程(組)，該方法模擬的精度有限[4]；COS連續(xù)系統(tǒng)方法建立的偏微分方程組僅在特殊、簡單的幾何結構及載荷下才可以求解。在風力機整機動力學方面主要采用有限元方法和模態(tài)分析方法[5]。但有限元方法具有較多的自由度，計算和分析成本較高[6]；模態(tài)疊加法通過提取結構的主模態(tài)集(主振型陣)形成坐標變換并實現(xiàn)降價來縮短計算時間，相對于其它方法可使用較少的自由度對風力機動態(tài)行為進行可靠分析，具有計算速度較快、效率高等優(yōu)點。

近年來，國內(nèi)外學者對風力發(fā)電機組的動態(tài)特性問題進行了研究： Asareh等[7]以FAST軟件為平臺，結合氣動力分析了5 MW風力機在不同工況下的動態(tài)響應；劉雄等[8]研究了風力機在湍流風作用下的動態(tài)特性，并揭示了葉片的離心剛化和氣動阻尼對分析結果的影響較大；呂計男等[9]對葉片和塔架進行簡化并建立了風力機整機模型，計算了風力機在啟動工況下的動態(tài)響應及葉片幾何非線性對結果的影響。李明等[10]采用車載法模擬自然風場，對風力機葉片的相似模型進行測試，得到的葉片各截面的載荷數(shù)據(jù)，與BLADED仿真結果較為一致。研究中風力機所處工況多為單一工況，關于風力機在不同工況下的動態(tài)特性研究較少；研究對象多為葉片、塔架、輪轂等單個部件，鮮有文獻對風力機整機的動態(tài)特性研究。

因此，筆者針對復雜工況下的大型風力機，采用了基于模態(tài)疊加法的動態(tài)特性分析方法，建立了風力機動態(tài)特性計算模型，采用該方法對2 MW大型風力機進行分析，同時總結了不同工況下的風力機的動態(tài)特性。

1 風力機動態(tài)特性計算理論

1.1 風力機葉片載荷計算理論

風力機在運行時主要通過葉片捕獲風能，捕獲的風能帶動風輪轉動，再通過傳動系統(tǒng)將風能轉化為機械能，因此葉片是主要的受力部件[11]。葉片所受的力主要包括空氣動力、離心力、和重力。葉片的受力和變形坐標系如圖1所示，其中 為迎風向， 為葉片展長方向， 由右手規(guī)則確定。

圖1 葉片受力和變形坐標系Fig.1 Coordinate system for blade loads and deflections

對葉片的載荷進行綜合可以計算出葉片單位長度所受的載荷[12]：

(1)

式中：FXg，F(xiàn)Ya為葉片單位長度上的氣動載荷；FXg，F(xiàn)Yg，F(xiàn)Zg為葉片單位長度上的重力；FYc，F(xiàn)Zc為葉片單位長度上的離心力。由于筆者所考慮的風力機風輪是靜止的，氣動載荷和離心力可忽略不計，因此迎風向載荷FXa是主要載荷。

1.2 基于模態(tài)疊加法的風力機動態(tài)特性計算理論模型

筆者將有限元法和模態(tài)疊加法這兩種方法進行了結合，其基本思想是以有限元為平臺，反復利用瞬態(tài)分析求解風力機在周期性載荷作用下的穩(wěn)態(tài)響應。瞬態(tài)分析包括完全法、縮減法和模態(tài)疊加法。完全法計算量大且計算時間較長；縮減法不能在時間單元上添加載荷且所有荷載必須加載用戶定義的主自由度上；模態(tài)疊加法相對于其它方法計算速度較快、效率高。因此本文選用模態(tài)疊加法對風力機整機動態(tài)特性進行求解，得到風力機在時變載荷作用下的多自由度運動程[13]：

(2)

[ψ]為系統(tǒng)的振型矩陣，作如下坐標變換：

u=[ψ]η

(3)

式中：η=[u1u2…un]T。

則系統(tǒng)的強迫振動方程為：

(4)

式(4)乘以特征振型[ψ]T可得：

(5)

由于[ψ]具有正交性，則主坐標下的方程為：

(6)

由主振型的正交性可得：

(7)

將式(7)代入式(6)，得到二階微分方程：

(8)

求出每個振型坐標上的位移分量后，采用模態(tài)疊加法可以得到各個節(jié)點自由度上的位移響應：

(9)

2 風力機動態(tài)特性有限元計算方法

2.1 風力機有限元模型

選取陸地某2 MW大型風力機為研究對象，葉輪直徑為40 m，輪轂直徑為2.5 m，塔架高75 m。葉片材料采用玻璃纖維，彈性模量為1.76×1010Pa，泊松比為0.17；塔架為Q345E合金鋼，彈性模量為2.06×1011Pa，泊松比為0.28。葉片結構采用梁單元，塔頂質(zhì)量采用點質(zhì)量單元劃分網(wǎng)格。塔架底部截面采用全約束，約束該截面節(jié)點的所有自由度。風力機有限元模型如圖2所示，單元節(jié)點數(shù)為93 389，單元數(shù)為462 396。

圖2 風力機有限元模型Fig.2 Finite element model of the wind turbine

2.2 風力機動態(tài)響應有限元計算方法

針對上述建立的基于模態(tài)疊加法的風力機動態(tài)特性計算模型，采用有限元法對其動態(tài)特性進行計算，進行求解的具體計算流程，如圖3所示。

圖3 風力機整機動態(tài)特性有限元法計算流程Fig.3 Dynamic characteristics computing process forwind turbines

圖3中的外部程序1是采用葉素動量定理計算葉片上的氣動力[14]，主要是葉片氣動性能計算和速度誘導因子迭代，Wilson等[15]詳細論述了該計算過程，本文不再贅述；主程序主要是通過二次開發(fā)語言APDL將葉片迎風向載荷讀入ANSYS，并將其離散成沖擊載荷，在運用模態(tài)疊加法求解風力機動態(tài)響應；外部程序2主要判斷風力機在求解動響應時是否達到穩(wěn)定狀態(tài)，周期求解動響應所允許的誤差為ε，若滿足式(10)，則表明動響應求解已進入穩(wěn)定狀態(tài)，則可停止計算；通過主程序及外部程序2求得風力機動態(tài)響應后，根據(jù)外部程序3進行擴展可進一步得到風力機的Von Mises動應力。

(10)

2.3 風力機模態(tài)分析

模態(tài)分析是動響應分析的基礎，采用Block Lanczons對風力機進行模態(tài)分析并提取前10階固有頻率，如表1所示。由于篇幅有限。僅列出風力機前五階模態(tài)振型，如圖4。

表1 風力機模態(tài)計算結果

圖4 風力機模態(tài)振型(1～5階)Fig.4 Modal vibration mode (1st to 5st order) of thewind turbine

由圖5可知：一階振型主要時風力機上面兩葉片沿X軸前后彎曲振動；二階振型主要風力機整機沿X軸前后彎曲振動；三階、四階、五階振型主要是葉片一階揮舞振動。對于三階振型，其下葉片保持不動，其余兩葉片沿X軸異向振動；四階振型是下葉片沿X軸反向振動；五階振型為三葉片沿Y軸同向振動。

3 典型工況下風力機動態(tài)特性分析

3.1 穩(wěn)態(tài)風況下風力機動態(tài)特性分析

計算可得風力機輪轂在恒定風速(即額定風速10 m/s)下的動載激勵，如圖5所示。根據(jù)風力機動態(tài)特性計算流程可得到風力機輪轂的位移動響應，如圖6。由圖5、6可知，輪轂載荷的在266.25～274.74 kN之間呈周期性變化，載荷頻率由1/3次諧波、0.9 Hz主頻成分及二倍頻組成；輪轂位移動響應在0.88～0.80 m之間波動，動響應主頻率成分為0.30 Hz，與輪轂1/3次諧波頻率(葉片主頻率)一致。

圖7為風力機整機的Von Mises應力云圖，可以看出塔架底部應力最大，平均應力達152×106Pa，且塔架底部至塔架頂部應力逐漸減??；下葉片的葉跟至葉尖的動載激勵和Von Mises應力如圖8所示，可以看出葉根應力最大，最大值為2.19×107Pa；從葉根至葉片11 m處時Von Mises應力呈拋物線型減小，葉片11 m處至葉尖，Von Mise應力呈減小趨勢但變化不大。

圖5 10 m/s恒定風速下風力機輪轂動載激勵Fig.5 Dynamic loads of hub under a constant wind speedof 10 m/s

圖6 10m/s恒定風速下風力機輪轂位移動響應Fig.6 Displacement dynamic response of hub under a constantwind speed of 10m/s

圖7 10 m/s恒定風速下風力機應力云圖Fig.7 The stress nephogram of the wind turbine under a constantwind speed of 10 m/s

圖8 葉片根部至葉尖動載激勵和Von Mises應力Fig.8 The loads and stress of blade from root to tip

3.2 不同風況下風力機動態(tài)特性分析

根據(jù)風力機整機動態(tài)特性計算方法計算了風力機在不同風況下的動態(tài)特性，總結不同工況下風力機位移動響應、Von Mises應力最大值及其相應位置，如表2所示。

表2 不同工況下風力機動態(tài)特性

由表2可知，風力機在不同風況下的最大位移動響應位置在都葉尖處；額定風速下風力機的位移響應最大值在所有工況中最小，為2.59 m；正常停機工況下的位移響應最大值為2.89 m，由于停機時產(chǎn)生了強烈的振動使得葉尖位移響應值瞬間減小到了-1.86 m；極端運行陣風作用下的位移動響應最大值在所有工況中最大，為7.19 m，較額定風速下的位移動響應最大值增加可177.22%。風力機在啟動工況下的Von Mises應力最大在在所有工況中最小，為0.44×108Pa；極端運行陣風工況下的Von Mises應力最大值在所有工況中最大，為4.05×108Pa，較額定風速下的Von Mises應力最大值增加了820.45%；且極端運行工況下的Von Mises應力最大值的位置在塔架頂部，其余工況都在塔架底部。在所有工況中，極端運行陣風的位移響應最大值和Von Mises應力最大值最大，在風力機設計中應著重考慮。

3.3 極端運行陣風下風力機動態(tài)特性分析

風力機運行環(huán)境十分復雜，由于丘陵、建筑物和森林等障礙物的阻擋空氣的流動會造成許多不規(guī)則的渦旋，渦旋的流動方向與空氣流動方向一致時，會產(chǎn)生瞬時極大風速，相反，會產(chǎn)生瞬時極小風速。極端運行陣風是指風速突然下降，接著陡然上升，然后再突然下降，最后又上升到初始值的過程。通過風力機動態(tài)特性計算方法可以求得風力機輪轂、葉尖的位移響應(圖9)，輪轂、塔架頂部的Von Mises應力(圖10)。

由圖9可知在穩(wěn)態(tài)風(0～15 s)作用下，輪轂位移響應在0.55～0.80 m之間波動，葉尖位移響應在0.76～2.92 m之間波動。在陣風(15～40 s)作用下，風力機發(fā)生了劇烈的振動，輪轂振幅最大值達2.08 m，較穩(wěn)態(tài)振幅均值增加了208.15%，最小值達在-1.18 m，較穩(wěn)態(tài)振幅均值反向增加了274.18%；葉尖振動幅值最大為7.19 m，較葉尖穩(wěn)態(tài)振動均值增加了290.76%，最小為-4.78 m，較穩(wěn)態(tài)均值反向增加了359.78%。陣風消失后，風力機逐漸回歸到穩(wěn)定狀態(tài)。

圖9 陣風作用下輪轂、葉尖位移動響應Fig.9 Displacement dynamic responses of hub and blade tipunder transient forcing-flurry

由圖10可知在穩(wěn)態(tài)風(0～15 s)作用下，風力機輪轂的Von Mises應力范圍為0.30×106～2.19×106Pa，塔架頂部的Von Mises應力范圍為0.08×108～1.69×108Pa；在陣風(15～40 s)作用下，輪轂的Von Mises應力最大值為5.3×106Pa，相對于其穩(wěn)態(tài)應力均值增加了325.70%；塔架頂部的Von Mises應力最大值為4.05×108Pa，相對于其穩(wěn)態(tài)應力均值增加了357.63%。陣風結束后，輪轂及塔架頂部的Von Mises應力減小至最小值，然后逐漸增加，最終呈周期性波動。

圖10 陣風作用下風力機輪轂塔架頂部的Von Mises應力Fig.10 The stress of hub and the top of tower under transientforcing-flurry

4 結 論

針對典型工況下風力機整機動態(tài)特性問題，提出了基于模態(tài)疊加法的風力機動態(tài)響應計算模型，進而采用有限元方法對其動態(tài)響應進行計算。對2 MW風力機在不同工況下進行動態(tài)特性分析，得到以下結論：

(1)穩(wěn)態(tài)風作用下風力機塔架底部應力最大，平均應力達152×106Pa，且塔架底部至塔架頂部應力逐漸減??；葉片的葉跟處應力最大，最大值達2.19×107Pa；其Von Mises應力從葉根至尖處時呈拋物線型減小。

(2)計算并總結了不同工況下的風力機動態(tài)特性，發(fā)現(xiàn)風力機的動響應最大的位置都在葉尖處，極端運行陣風下的Von Mises應力最大位置為塔架頂部，其余工況都在塔架底部；且極端運行陣風作用下的最大位移響應值、最大Von Mises應力值在所有工況中均為最大，在風力機設計時應著重考慮。

(3)極端運行陣風條件下，葉尖的位移動響應最大值較其穩(wěn)態(tài)響應均值增加了290.76%，最小值較其穩(wěn)態(tài)響應均值反向增加了359.78%；塔架頂部的Von Mises應力最大值較其穩(wěn)態(tài)應力均值增加了357.63%。