黃鑫祥，郭小鋒，齊劍峰

（中原工學院 機電學院， 河南 鄭州 450007）

0 引言

水平軸風力機葉片的氣彈變形包含揮舞、擺振和扭轉(zhuǎn)3 個方面。 常用的風力機動力學分析軟件（如GH-Bladed 和FAST）只考慮了葉片在揮舞和擺振兩個方向的自由度。 然而，隨著風電機組容量的不斷擴大，為降低載荷，大型低風速及海上風電葉片的柔韌度較早期同類型的風電葉片有明顯的增大，如何準確地分析葉片在實際工況下的氣彈扭轉(zhuǎn)變形，進而分析葉片的動態(tài)氣彈扭轉(zhuǎn)變形對風力機氣動功率、載荷及氣彈穩(wěn)定性的影響，是大型風電葉片的氣彈特性評價以及氣彈剪裁設計中必須解決的問題。

郭小鋒[1]運用修正的葉素-動量理論和有限元方法，將風力機的動力學問題轉(zhuǎn)化為靜力學問題， 提出了考慮氣動彈性的葉片性能分析方法，綜合分析了葉片扭轉(zhuǎn)變形對風力機氣動性能的影響。李德源[2]采用“超級單元”模型，應用計算多體動力學理論和風力機氣動模型，建立了受約束的風力機系統(tǒng)氣彈耦合方程。 江術新[3]建立了整機在變槳工況下的氣-彈-控聯(lián)合仿真模型，時域分析了變槳過程中的槳距角變化、功率曲線變化和葉根處受力和力矩的變化。

現(xiàn)有的研究表明，靜氣彈性分析方法沒有反映柔性風電機組的動力學問題本質(zhì)。 風力機的實際工況具有復雜多變性，對于風力機葉片在實際復雜工況下的氣彈扭轉(zhuǎn)變形問題， 及其對風力機氣動性能及氣彈問題性的影響規(guī)律仍須進一步研究。

本文基于葉片及塔筒部件的參數(shù)化表達方法， 采用多體動力學仿真軟件SIMPACK 建立了風力機整機剛?cè)狁詈夏Ｐ停?提出了一種葉片動態(tài)氣彈扭轉(zhuǎn)變形分析的新方法。在幾種典型風況下，對NREL5MW 風力機進行了氣彈動力學分析，采用matlab 軟件對分析結(jié)果進行處理，得到了葉片在典型風況下的動態(tài)氣彈扭轉(zhuǎn)變形數(shù)據(jù)[4]。 在此基礎上，對NREL5MW 風力機在幾種典型工況下的氣動功率、 載荷及氣彈穩(wěn)定性問題進行了分析研究。

1 風力機剛?cè)狁詈夏Ｐ偷臉?gòu)建

1.1 參數(shù)化建模

風力機整機模型是典型的剛?cè)狁詈隙囿w系統(tǒng)。在工作狀態(tài)下，風力機葉片存在明顯的彈性變形，所以將葉片作為柔性體處理。而風力機的其他部件（如機艙、主軸及輪轂等部件）的剛性大，在工作時的彈性變形很小，可以看作剛體。為了縮短計算機的求解時間， 本文將風力機模型進行了合理簡化，只包括葉片、輪轂、塔筒、機艙等主要部件，而不考慮齒輪箱、發(fā)電機等部件。偏航法蘭與塔筒頂端固接，用于連接機艙，雖然機艙中不創(chuàng)建主機架、齒輪箱、發(fā)電機等模型，但是需要把忽略部件的質(zhì)量和轉(zhuǎn)動慣量算入機艙中。 低速軸與機艙為轉(zhuǎn)動連接，輪轂與低速軸進行固接，葉片與輪轂由變槳法蘭連接。

采用EXCEL 表格對葉片的外形和鋪層結(jié)構(gòu)進行參數(shù)化表達，采用matlab 編程調(diào)用PreComp軟件計算得到SIMPACK 葉片建模中所需各個截面的屬性參數(shù)（截面質(zhì)心、彈性中心、剪切中心、扭轉(zhuǎn)剛度、彈性模量等），并自動生成Simpack 軟件Rotorblade Generation 計算模塊的輸入文件[5]，[6]。

風力機的控制系統(tǒng)采用基本的變槳變速控制策略， 在SIMPACK 中設置風輪轉(zhuǎn)速及槳距角在不同穩(wěn)態(tài)風速下的參數(shù)，同時通過設置PID 參數(shù)實現(xiàn)對風輪的動態(tài)控制。

在SIMPACK 中， 剛?cè)狁詈隙囿w系統(tǒng)建模是在浮動參考坐標系法的基礎上建立多體動力學方程， 并由SIMPACK 求解動力學方程獲得風力機的動態(tài)響應。在多體系統(tǒng)（MBS）理論中，柔性體看作有限元模型若干個節(jié)點的組合，其中每個節(jié)點的運動由跟隨整體坐標系的整體平動和轉(zhuǎn)動以及在局部坐標系下發(fā)生小的線性變形組成。 圖1 所示為柔性體的一般模型圖。

為描述節(jié)點的運動，一般采用慣性坐標系下的笛卡爾坐標和表示方位的歐拉角。

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

1.2 風場模型

空氣動力學軟件AeroDyn 是分析風力機葉片風載的有力工具， 為了實現(xiàn)與SIMPACK 軟件的聯(lián)合仿真，需要準備AeroDyn 的輸入文件、風場文件、翼型文件、塔筒文件等[7]。 本文的風場文件由NREL 提供的IECWIND 軟件生成。

AeroDyn 基于BEM 理論建立迭代方程計算風輪氣動力和誘導速度， 通過迭代求解葉片單元的誘導速度、攻角、推力系數(shù)等參數(shù)，直到誘導因子和入流角的值收斂到最終值[8]。

AeroDyn 與SIMPACK 動力學聯(lián)合仿真流程如圖2 所示。

圖2 動力學仿真流程圖Fig.2 Diagram of simulation process

圖2 中，AeroDyn 與SIMPACK 兩 者 通 過SIMPACK 中241 號力元接口進行數(shù)據(jù)交換。AeroDyn 同時從輸入文件和SIMPACK 求解模塊中收集t 時刻的葉片幾何形狀、工況條件、葉元速度、位置以及流入風等信息，然后利用以上信息計算出葉片各截面的氣動數(shù)據(jù)。 SIMPACK 的求解模塊利用241 號力元接口調(diào)用氣動數(shù)據(jù)計算得到t+dt 時刻葉片各截面的氣彈變形和氣動力， 計算得到的數(shù)據(jù)又通過241 號力元接口反饋給AeroDyn 用來計算t+dt 時刻各截面的氣動數(shù)據(jù)，如此往復，直到仿真完成。而控制模塊根據(jù)輸出的反饋力矩和轉(zhuǎn)速對風力機進行變速變槳控制。

為了得到仿真結(jié)果的精確性， 本文選擇SODASRT_2 解算器對模型進行離線積分求解，仿真時長設為120 s， 求解步長設為0.01 s。 由于SIMPACK 每個時間步長都會調(diào)用AeroDyn 計算氣動力變化， 所以AeroDyn 求解時間步長也設為0.01 s。

為研究陣風、 風向角度對葉片扭轉(zhuǎn)變形的影響，選取如表1 所示風場模型?？紤]到海上風況的復雜多變性， 本文在建立風場模型時將方向變化的極端相干陣風（EDC）和極端工作陣風（EOGR）風況下風速和風向的波動設定為4 個變化周期，并且假定在EDC 和EOGR 風況作用下控制系統(tǒng)故障， 研究風力機在極端風況下扭轉(zhuǎn)變形與整機性能的關系。

表1 極端風況Table 1 Extreme wind conditions

1.3 模型驗證

通過SIMPACK 中Eigenvalues 模塊對風力機（NREL5MW）進行模態(tài)分析，計算得到整機的各階自振頻率（表2）。 由表2 可知，本文的分析結(jié)果同NREL 報告中的數(shù)據(jù)誤差在5%以內(nèi)[4]，證實了SIMPACK 所建立風力機整機模型的可靠性。

表2 5 MW 風力機整機自振頻率對比Table 2 Comparison of natural frequency of 5 MW wind turbine

2 葉片氣彈扭轉(zhuǎn)變形分析

風力機在啟動階段的各種變化均不穩(wěn)定，為更好地研究各風況對葉片氣彈變形的影響， 本文假設各風況發(fā)生變化時風力機處于穩(wěn)定階段。 假定各風況在風力機啟動后40 s 時刻發(fā)生變化，并于80 s 之后逐漸恢復為12 m/s 的穩(wěn)態(tài)風。 對風力機整機進行了仿真分析， 得到4 種工況下葉尖揮舞、擺振的動態(tài)響應（圖3）。

圖3 各風況下的葉尖位移Fig.3 Tip deformation under various wind conditions

為分析葉尖部位翼型的氣彈扭轉(zhuǎn)變形數(shù)值，本文首先采用SIMPACK 軟件得到葉尖某截面各節(jié)點的初始坐標以及其相對于初始坐標的變形量， 將各節(jié)點的初始坐標與對應的位移量相加得到節(jié)點在各個時刻的動態(tài)坐標值。 將各節(jié)點的動態(tài)坐標值投影到初始時刻的截面上， 得到截面扭轉(zhuǎn)變形后的節(jié)點坐標值。 對于NREL5MW 葉片，距離葉根60 m 處截面的初始形狀和在48 s 時的投影形狀對比如圖4 所示。

圖4 在ECD 風況下，0 s 和48 s 時葉尖翼型扭轉(zhuǎn)形變的對比Fig.4 Comparison of tip airfoil torsional deformation at 0 s and 48 s under the ECD wind condition

在風力機的剛?cè)狁詈夏Ｐ椭校?假設葉片截面在扭轉(zhuǎn)變形時自身形狀不發(fā)生變化， 葉片截面扭轉(zhuǎn)變形前后兩翼型弦長線間的夾角即為該截面在ECD 風況48 s 時刻的氣彈扭轉(zhuǎn)變形量。

為分析葉片的動態(tài)氣彈扭轉(zhuǎn)時域變形趨勢，由上述方法得到距離葉根60 m 處截面各個時刻的動態(tài)氣彈扭轉(zhuǎn)變形曲線（圖5）。

圖5 各風況下的葉尖扭轉(zhuǎn)變形Fig.5 Tip torsional deformation under various wind conditions

由圖5 可知，當風力機運轉(zhuǎn)穩(wěn)定后，各風況下的葉尖動態(tài)扭轉(zhuǎn)角度在2°上下浮動，說明葉片具有的扭轉(zhuǎn)剛度使得葉片在一定風速下的扭轉(zhuǎn)角度會在某一個值上下變動。 而葉片各個截面的扭轉(zhuǎn)剛度不同，導致了葉片各截面扭轉(zhuǎn)角度不同，計算葉片實際工作中的扭轉(zhuǎn)角度對葉片的優(yōu)化設計具有指導作用。

對比EOGR，EDC 風況作用下葉尖扭轉(zhuǎn)變化曲線可知， 陣風變化對葉尖扭轉(zhuǎn)變形的影響要大于風向角度變化， 其中風向角度變化使得扭轉(zhuǎn)變化頻率增加。 在ECD 風況作用下，葉尖扭轉(zhuǎn)變形量最為明顯，葉尖最大扭轉(zhuǎn)為4.78 °，最小扭轉(zhuǎn)為0.71 °，陣風風速越大，葉尖扭轉(zhuǎn)變形量越大， 而風速和風向角度的共同作用增加了葉片的扭轉(zhuǎn)頻率。

在ECD 風況下，由于風速和風向角度的共同作用造成了葉尖扭轉(zhuǎn)變形的突然增大， 這種短時間內(nèi)葉尖扭轉(zhuǎn)變化的突然加劇容易對葉片造成致命的損壞。 而風向角度多變的風況使得葉尖扭轉(zhuǎn)頻率增大，這種風況下容易造成葉片的疲勞破壞。

3 氣彈扭轉(zhuǎn)變形對整機性能的影響

3.1 氣動性能

根據(jù)BEM 理論計算可得，在11.4 m/s 的穩(wěn)態(tài)風下， 風力機功率為5 MW。 而在實際工作過程中，由于塔筒變形、葉片的氣彈變形等影響，風力機功率并不是一個恒定的值。 不同風況的作用對葉片氣彈變形的影響很大， 導致各風況下的風力機功率存在明顯波動。

風力機控制模塊根據(jù)檢測到的實時轉(zhuǎn)速做出變槳響應，使其穩(wěn)定在額定功率附近。 仿真時，在40～80 s，ECD，EOGR 風況的風速和風向角急劇變化，導致風力機變槳速率、葉片翼型攻角變化幅度增大， 葉片氣動載荷增加。 隨著葉片變槳速率增大，葉片在風速上升期卸載越快，對功率的影響也就越明顯。 另外，當控制系統(tǒng)故障時，風力機變槳控制明顯滯后于實際風況， 在風速、 風向變化劇烈，但未到達切出風速的極端風況下，控制系統(tǒng)無法實時調(diào)整到對應風速的槳距角， 從而導致控制系統(tǒng)暫時失效， 造成功率波動異常。 圖6 為在ECD，EOGR 風況下，風力機功率的變化曲線。 由圖6 可知，風力機功率出現(xiàn)了大幅度的波動，峰值功率甚至達到了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 風況，葉片的扭轉(zhuǎn)變形在一定程度上能削減風力機上的極限載荷。

圖7 各風況下輪轂處載荷的變化Fig.7 Load changes at the hub under different wind conditions

3.3 氣彈穩(wěn)定性

在復雜工況下，風速、風向的劇烈變化引起了葉片揮舞、擺振和扭轉(zhuǎn)方向的急劇變化，從而造成葉片氣彈失穩(wěn)。 特別是葉片的扭轉(zhuǎn)變形直接影響攻角變化，從而造成升力系數(shù)的變化，嚴重時會發(fā)生失速現(xiàn)象。 本文通過AeroDyn 對加載風載的葉片進行計算分析， 選取50～60 s 時刻各風況下距離葉根60 m 處截面翼型的氣動數(shù)據(jù)進行處理，得到攻角-升力系數(shù)曲線，并與靜態(tài)翼型攻角-升力系數(shù)曲線進行比較（圖8）。

圖8 各風況下葉尖截面攻角-升力系數(shù)曲線Fig.8 Angle of attack-lift coefficient curve of blade tip section under various wind conditions

由圖8 可知：NWP，EDC 風況下的升力系數(shù)隨攻角呈線性變化，且該時間段內(nèi)EDC 風況下的攻角大于NWP 風況下的攻角，但沒達到翼型的靜態(tài)失速攻角；EOGR 風況下的攻角進一步增大，接近翼型靜態(tài)失速攻角，仍呈線性變化；ECD 風況下的攻角超過靜態(tài)失速攻角并繼續(xù)增大， 升力系數(shù)也隨著攻角的增大總體呈上升趨勢，達到峰值后攻角逐漸減小，升力系數(shù)也隨之減小，曲線開始回歸線性區(qū)，這是典型的動態(tài)失速現(xiàn)象[9]；在ECD 風況作用下，葉尖截面存在明顯的動態(tài)失速現(xiàn)象，而其他風況作用下，葉尖截面均未出現(xiàn)失速現(xiàn)象，所以在風速、風向變化幅度大的風況下，葉片容易出現(xiàn)動態(tài)失速現(xiàn)象，從而造成葉片氣彈失穩(wěn)， 嚴重時葉片會發(fā)生顫振甚至斷裂。

4 結(jié)論

本文基于SIMPACK 軟件建立了剛?cè)狁詈系娘L力機整機模型，對NREL5MW 風力機進行了動力學分析研究，得到以下結(jié)論。

①提出了一種分析葉片動態(tài)氣彈扭轉(zhuǎn)變形的新方法， 該方法可得到葉片在復雜風況下的動態(tài)氣彈扭轉(zhuǎn)變形數(shù)據(jù)。

②極端復雜工況下， 風速和風向的變化會導致葉片產(chǎn)生較大的氣彈變形，進而對風力機氣動、載荷及氣彈穩(wěn)定性產(chǎn)生較大的影響。

③在ECD 工況下， 葉片氣彈扭轉(zhuǎn)變形急劇增大，會導致葉片處于短暫失速狀態(tài)，葉片容易出現(xiàn)失速現(xiàn)象。 風力機功率和載荷波動增大，從而造成葉片氣彈失穩(wěn)，嚴重時葉片會發(fā)生顫振甚至斷裂。