汪 泉，陳 進，郭小峰，李松林，王旭東

（1.重慶大學機械傳動國家重點實驗室，重慶 400030；2.湖北工業(yè)大學機械工程學院，湖北 武漢 430068；3.重慶工商大學機械工程學院，重慶 400067）

0 引 言

隨著風力機風輪直徑的增加，準確的預測其能量輸出及氣動載荷特性變得越來越重要。對于提高風力機風輪功率系數，常見的方法是應用不同的翼型、改變葉片的弦長與扭角分布、進行失速控制及俯仰控制等等。

對于風力機葉片設計研究，Kamoun Badreddinne［1］以風力機輸出功率為優(yōu)化目標，對葉片的升力變化過程進行了研究，并將此與Glauert模型進行了對比，然后對葉片進行了優(yōu)化設計；汕頭大學劉雄［2］等考慮了風場風速的概率分布，以風力機年能量輸出最大為優(yōu)化目標，使用遺傳算法，優(yōu)化得到了1.3MW風力機葉片；華北電力大學韓中合等［3］提出了多次迭代優(yōu)化設定誘導因子初始值的方法，以功率輸出和年發(fā)電量最大為優(yōu)化目標，對1.5MW風力機葉片進行了優(yōu)化設計；重慶大學 Wang Xudong［4］以風力機單位發(fā)電量的成本為目標，對幾種不同裝機容量的風輪進行了優(yōu)化設計，降低了風力機單位發(fā)電量的成本。此外，還有一些學者對此做了很多研究［5－10］，都值得借鑒參考。

然而以上研究大都針對現(xiàn)有葉片進行改進，很少有對具有全新翼型族的葉片進行優(yōu)化設計，而且在對葉片外形的優(yōu)化過程中忽略了減少葉片質量對提高風力機葉片疲勞壽命及降低葉片材料成本這一至關重要的因素。而風力機葉片是由復合材料加工制造而成，葉片密度變化不一，其曲面形狀非常復雜，采用傳統(tǒng)的復合材料層合板理論很難建立葉片的質量計算數學模型，因此，在其他條件一致的情況下（葉片鋪層厚度、密度等），本文提出將葉片的質量計算轉化為葉片曲面面積的計算模型，要使質量最小，意味著對應的曲面面積最小。在整個變槳距風力機葉片設計上，采用作者全新設計的CQU－A翼型族，驗證該翼型族具有較高的氣動性能。基于修正的風力機空氣動力學模型理論，并提出以最大功率系數及最小葉片面積為多目標優(yōu)化模型，以弦長和扭角為設計參數，考慮葉根載荷約束，采用改進的多目標粒子群算法對2MW風力機葉片進行設計與優(yōu)化，并將優(yōu)化結果與初始葉片以及歐盟Tj?re實驗葉片進行對比分析，以驗證該優(yōu)化設計方法的可行性及優(yōu)越性。

1 修正的風輪空氣動力學模型

一維動量理論對風輪進行空氣動力特性分析的前提是假設風輪由無限個葉片組成，而事實上風輪都是由有限葉片來組成的。這樣使得這兩種模型的風力機尾流中的旋渦系不同。普朗特給出了葉尖修正因子F，推導出風輪的推力與扭矩為［11］：

ρ ——空氣密度；

v0——來流風速；

a——軸向誘導因子；

F——普朗特葉尖損失修正因子；

r——葉片展向位置；

ω——風輪轉動角速度；

a′——周向誘導因子；

φ——入流角；

B——葉片個數；

R——風輪半徑；

普朗特葉尖修正模型主要針對軸向質量流量進行了修正，Shen修正模型對風輪葉素中的法向力系數和切向力系數做了進一步的修正［12］，引入修正因子F1并代入經典的葉素理論中得到推力與轉矩為：

g＝exp［－0.125（Bλ－21）］＋0.1

vrel——入流相對風速；

c ——葉片弦長；

λ ——葉尖速比；

Cn——法向力系數；

Ct——切向力系數；

聯(lián)立（1）～（4）可以求出軸向誘導因子a和周向誘導因子a′為：

式中 Y1＝4F sin2φ／（σCnF1）

Y2＝4F sinφcosφ／（σCtF1）

σ＝Bc／（2πr）

當軸向誘導因子a＞0.3時，葉素理論將不再準確可靠，此時就需要對風力機的推力進行修正。那么軸向誘導因子a和周向誘導因子a′為

式中ac＝1／3。

2 葉片設計及優(yōu)化模型的建立

2.1 新型風力機葉片設計

本研究在整個葉片上全部采用作者新設計的CQU－A翼型系列（相對厚度從15%～40%的七種翼型，如圖1、圖2）高性能風力機專用翼型，該翼型系列基于廣義泛函及儒科夫斯基保角變換的理論設計方法。該方法將近似于圓的圖形通過儒科夫斯基保角變換生成翼型，改變其控制方程中的系數就能夠控制翼型型線，通過這種方法，已優(yōu)化出性能較好的翼型［13－14］。

圖1 CQU－A翼型系列Fig.1 CQU－A airfoil series

為了驗證CQU－A翼型系列具有較高的氣動性能，選取葉片主要產生功率區(qū)域的兩種翼型CQUA180及CQU－A210與相當厚度的DU翼型進行對比分析。采用RFOIL軟件計算翼型的氣動性能。其中以自由轉捩工況模擬光滑條件；以固定轉捩工況模擬粗糙條件，轉捩模型采用吸力面（上翼面）處于1%弦長位置固定轉捩，壓力面（下翼面）處于10%弦長位置固定轉捩［15］。圖 3為 CQU－A180翼型與DU96－W－180翼型在相同運行工況下（Re＝3×10 Ma＝0.15）的氣動特性對比圖。從圖中可以看出，無論是光滑條件，還是粗糙條件下，新翼型的升力系數、升阻比均有顯著的提高。圖4為CQU－A210翼型與DU93－W－210翼型在相同工況下的氣動性能對比圖。在翼型失速之前，新翼型的最大升力系數、最大升阻比均有較大提高。表1給出了四種翼型的氣動特性數據對比，結合圖表，表明新設計的翼型無論是在光滑條件還是粗糙條件，均具有高的升力系數和升阻比。

圖2 CQU－A翼型系列Fig.2 CQU－A airfoil series

圖3 CQU－A180翼型與DU96－W－180翼型氣動性能對比Fig.3 The comparison of aerodynamic performance of CQU－A180airfoil and DU96－W－180airfoil

圖4 CQU－A210翼型與DU93－W－210翼型氣動性能對比Fig.4 The comparison of aerodynamic performance of CQU－A180airfoil and DU93－W－210airfoil

表1 翼型氣動性能參數對比Tbale 1 The comparison of aerodynamic performance of airfoils

風力機葉片設計的基本參數如表2所示。其中各個翼型按照對應于2MW葉片展向位置處的氣動條件進行設計，根據額定運行工況條件下，翼型所處展向位置處的弦長及相對來流速度確定雷諾數和馬赫數。翼型沿葉片展向布置及邊界層條件見表3所示，其中相鄰翼型間的過渡翼型，通過相鄰翼型的坐標插值得到。由于翼型沿葉片展向分布時其雷諾數會不一樣，表3中給出了當葉片旋轉速度為20rpm時幾個關鍵葉片段的雷諾數大小，由于馬赫數對葉片翼型的氣動性能影響不大，因此這里取整個葉片的馬赫數為0.15。

表2 葉片設計基本參數Table 2 The basic parameter of the designed blade

表3 翼型沿葉片展向分布及邊界層條件Table 3 The airfoils along the blade span－wise distribution and the boundary layer conditions

2.2 建立多目標優(yōu)化模型

2.2.1 目標函數

對于變槳距風力機而言，風輪轉速可以連續(xù)進行調節(jié)，就有可能使風力機運行在最優(yōu)點Cp，max，因此本文以風力機功率系數最大作為目標優(yōu)化之一：

基于風力機功率系數的定義，我們可以得到以下表達式：

按照修正的葉素動量理論，經過一系列的變換可以得到功率系數另外一種表達式：

風力機葉片是由復合材料加工制造而成，葉片密度變化不一，其曲面形狀非常復雜，采用傳統(tǒng)的復合材料層合板理論很難建立葉片的質量理論計算模型，因此，在其他條件一致的情況下（葉片鋪層厚度及密度），將葉片的質量計算轉化為葉片曲面面積計算模型，要使質量最小，意味著對應的葉片表面面積最小。圖5為葉片面積計算示意圖，葉片曲面面積最小作為目標優(yōu)化之一：

圖5 葉片曲面面積計算示意圖Fig.5 The blade surface area calculation schematic

f（x，y）為葉片沿展向r處截面輪廓線長度方程，采用微積分的形式來表達此函數：

xi、yi為三維空間葉片截面數據點坐標，如果數據點足夠多，將會無限接近截面曲線長度。根據葉片空間坐標幾何變換原理，可以將沿展向r處的葉片截面坐標（xi，yi）表示成：

其中x，y為葉片沿展向r處翼型的空間坐標，x0，y0為單位弦長翼型在二維平面的坐標，c為葉片沿展向r處的弦長，β為葉片沿展向r處的幾何扭角。

聯(lián)立方程（13）～（15），采用微積分的方法即可求出整個葉片曲面的面積。

由上面兩個目標函數可知，該葉片優(yōu)化過程實質上是一非線性雙目標函數優(yōu)化數學模型。為求解方便，本文通過設定權值系數將多目標優(yōu)化問題轉化為單目標優(yōu)化函數：

式中，μ1為權值系數，μ1∈［0，1］；為使兩個目標函數為同一個數量級，在這里給定一個比例因子ratio，一般最大功率系數在0.5附近，那么為了使兩個優(yōu)化目標函數重要程度相當，F(xiàn)2（X）／ratio也應該是這一數量級。

2.2.2 設計變量及約束條件

葉片的曲面形狀決定了風力機捕獲風能的效率，而葉片的長度、弦長、扭角及相對厚度直接構成了葉片的形狀。風輪的半徑選取常值為31m，由表3可知，標準翼型沿葉片展向位置是確定的，葉片的厚度分布即可通過插值的方法確定，因此選取葉片的弦長、扭角作為優(yōu)化設計變量。為了避免因翼型族的不同而使優(yōu)化結果對比沒有意義，葉片初始幾何外形參數選擇傳統(tǒng)的歐盟2MW Tj?re型葉片（采用NACA44XX系列翼型）相關參數，該葉片的具體參數（弦長及扭角分布）見參考文獻［4］；同時針對同一功率大小不同翼型族的風輪性能進行研究（假設葉片長度、弦長分布和扭角分布不變）。

針對葉片的弦長，為了保證葉片擁有良好的氣動特性及表面連續(xù)光滑性，本文選取關鍵的8個點（標準翼型沿葉片展向分布的點），采用樣條曲線來對弦長及扭角的展向分布進行控制，這樣設計變量總共就有16個，其取值范圍控制在（表4給出了優(yōu)化設計變量的約束范圍）

表4 優(yōu)化設計變量的約束范圍Table 4 The constrained range for optimized design variables

另外在優(yōu)化設計中，還需對葉片的載荷進行控制約束。葉片的載荷由葉素動量理論計算得出。基于風力機空氣動力學理論知，葉片的法向力系數Cn和切向力系數Ct表示為

式中L為翼型的升力系數，D為翼型的阻力系數。

在風輪工作時，葉根處的彎矩很大，彎矩的增加會直接影響著葉片強度及風力機的疲勞壽命，而葉根處的彎矩主要為揮舞彎矩和擺振彎矩，揮舞彎矩是由推力而產生，擺振彎矩則是由切向力分布產生的在風輪平面內的彎曲力矩，由于擺振彎矩的大小很大程度上受到葉片重力的影響，氣動載荷影響較小，彎矩圖表現(xiàn)出明顯的正弦變化，而揮舞彎矩主要受到氣動載荷的影響［16］，因此，不考慮擺振彎矩約束，對葉根揮舞彎矩進行了約束。

此外，風輪的扭矩由葉片所受的切向力而形成的，扭矩的大小直接決定了風力機的輸出功率，理論上講，輸出功率隨轉矩的增加而增加。但是，轉矩過大會使風力機傳動系統(tǒng)的載荷增加，不利于機械傳動系統(tǒng)的使用壽命。因此，對風輪的扭矩也進行了約束控制

3 優(yōu)化設計結果

基于修正的空氣動力學模型，給定一定的初始條件，設定容許誤差，通過不斷的迭代，由式（5）和（6）可以計算出葉片的軸向誘導因子a和周向誘導因子a′，進而就可以得到葉片上的作用力和功率等。采用改進的多目標粒子群算法，編制程序進行優(yōu)化求解。由于粒子群算法中有些基本參數不定，優(yōu)化時需要根據具體問題不斷調試參數進行迭代，經過多次迭代，作者發(fā)現(xiàn)選取如下參數時該改進粒子群算法優(yōu)化效果及收斂速度良好：慣性權重w取0.9，學習因子C1、C2均取0.5，變量維數為16，種群大小為50，最大迭代次數為200，權值系數μ1為0.25，比例因子ratio為500。葉片具體優(yōu)化設計計算流程圖如圖6所示，明確目標函數、設計變量及約束條件；以葉片弦長和扭角為變量初始化粒子群，采用計算機圖形學中的樣條曲線將粒子群集插值得到弦長、扭角分布曲線；然后基于動量葉素理論及葉片曲面面積積分，計算適應度值；然后根據適應度值更新相應的粒子群參數，判斷是否滿足迭代終止條件，否則繼續(xù)迭代，是則輸出最優(yōu)值。其中翼型的氣動參數采用RFOIL計算，然后采用外插值方法生成翼型±180°攻角范圍的氣動參數。

圖7、圖8給出了優(yōu)化后葉片弦長及扭角分布圖，將該優(yōu)化結果與優(yōu)化前初始葉片進行了比較。由圖7可以看出，新設計的葉片弦長沿展向分布呈現(xiàn)出非線性特性，從整體上看，每個葉素的弦長都要小于初始葉片，在葉片的中部尤為明顯，這是由于將葉片的曲面面積最小作為目標函數之一，因此算法將會朝著有利于減少葉片面積降低葉片質量的趨勢發(fā)展；由圖8可以看出，新設計的葉片扭角分布也呈現(xiàn)出非線性特性，在葉片根部附近扭角要明顯大于初始葉片，而葉尖處相差不太，葉片的扭角增大會減小葉片根部的載荷，進而增加葉片的使用壽命，有利于降低葉片的成本。

圖6 葉片優(yōu)化設計算法流程圖Fig.6 The blade optimized design algorithm flowchart

圖7 弦長分布Fig.7 The distribution of the chord length

圖8 扭角分布Fig.8 The distribution of the twist

圖9為Tj?re葉片、優(yōu)化前葉片和新設計葉片的功率系數隨葉尖速比變化的曲線圖。由圖9可知，初始葉片在設計葉尖速比下（2～7.5）的功率系數要明顯大于歐盟Tj?re實驗葉片，主要原因是由于翼型族的不同而使其氣動性能不同造成的；而新設計的葉片在設計葉尖速比下（6～12）功率系數要明顯大于初始葉片，這主要是由于葉片外形參數不同而引起的。圖10給出了在相同條件下新設計的風輪年均發(fā)電量特性隨風速的變化曲線，相比優(yōu)化前及Tj?re風輪，在額定風速內（V≤12.5m／s），優(yōu)化后的風輪年發(fā)電量分布有了較大的提高，當風速達到額定風速后，通過傳動控制及變槳距調節(jié)系統(tǒng)使功率控制在2MW，這樣年均發(fā)電量及功率將保持不變。

圖9 功率系數特性Fig.9 The characteristics of the power coefficient

圖10 年均發(fā)電量分布Fig.10 The distribution of the annual generating capacity

表5為Tj?re風輪、初始葉片以及新設計的2MW葉片各項性能指標對比。結合上圖及表5可知，初始葉片在最大功率系數、年發(fā)電量方面均優(yōu)于Tj?re型風輪，但是葉根載荷、風輪扭矩均大于Tj?re型風輪；而相比初始2MW葉片，新設計的風輪葉片，最大功率系數為0.50830，提高了4.019%，年發(fā)電量為4.406×106kWh，提高了1.967%，風輪面積為278.906m2（0.2R～1.0R 處），減少了19.005%，這意味著在其他葉片內部結構不變的情況下風輪的質量有了大幅度的減少，提高了葉片的疲勞壽命，降低了葉片的材料成本。由于采用全新的高性能的風力機翼型，較大的氣動力會增加葉根載荷，通過該設計及優(yōu)化方法，風輪葉根的彎矩和扭矩均有所減少，且葉根彎矩減少尤為明顯，表明該設計方法能夠有效地控制葉根載荷。

表5 新設計葉片與傳統(tǒng)葉片各項性能對比Table 5 The comparison of the performances for the new design blade and the traditional blade

4 結 論

（1）在整個葉片設計上采用作者全新設計的風力機翼型系列CQU－A，驗證了該翼型系列具有較高的氣動性能；基于修正的風力機空氣動力學理論，提出了以最大功率系數及最小葉片面積為多目標優(yōu)化模型，以弦長、扭角為設計變量，建立了新型風力機葉片設計及優(yōu)化數學模型。

（2）應用該設計模型，采用改進的多目標粒子群算法對風力機葉片進行優(yōu)化設計，并將優(yōu)化結果與初始葉片、Tj?re實驗葉片進行對比分析，優(yōu)化結果表明，相比初始葉片及Tj?re葉片，新設計2MW風力機葉片的最大功率系數及年發(fā)電量均有所提高，而且葉片表面面積具有顯著的下降，這意味著在其他葉片內部結構不變的情況下極大的減少了葉片的質量，提高了葉片的疲勞壽命，降低了葉片的材料成本；同時葉根載荷也得到了有效的控制。該新型葉片的研究為設計出高性能輕質量低成本的風力機提供了理論依據。

［1］BADREDDINNE K，ALI H，DAVID A.Optimum project for horizontal axis wind turbines‘OPHWT′［J］.Renewable Energy，2005，30：2019－2043.

［2］劉雄，陳嚴，葉枝全.遺傳算法在風力機風輪葉片優(yōu)化設計中的應用［J］.太陽能學報，2006，27（2）：180－185

［3］韓中合，吳鐵軍.基于遺傳算法的風力機葉片優(yōu)化設計［J］.動力工程，2008，28（6）：955－958.

［4］WANG Xudong，SHEN W Z，ZHU W J，et al.Shape optimization of wind turbine blades［J］.Wind Energy，2009，12：781－803.

［5］JACKSON K J，ZUTECK M D，VAN D C，et al.Innovative design approaches for large wind turbine blades［J］.Wind Energy，2005，8：141－171.

［6］楊從新，宋顯成.一種大型風力機葉片的氣動優(yōu)化設計方法［J］.空氣動力學學報，2011，29（2）：222－229.

［7］王瓏，王同光，羅源.改進的NSGA－Ⅱ算法研究風力機葉片多目標優(yōu)化［J］.應用數學與力學，2011，32（6）：693－700.

［8］RAJAKUMAR S，RAVINDRAN D.Iterative approach foroptimizing coefficient of power，coefficient of lift and drag of wind turbine rotor［J］.Renewable Energy，2012，38：83－93.

［9］KEVIN Maki，RICARDO Sbragio，NICKOLAS Vlahopoulos.System design of a wind turbine using a multileveloptimization approach［J］.Renewable Energy，2012，43：101－110.

［10］RONIT K Singh，RAFIUDDIN Ahmed M.Blade design and performance testing of a small wind turbine rotor for low wind speed applications［J］.Renewable Energy，2013，50：812－819.

［11］ 賀德馨.風工程與工業(yè)空氣動力學［M .北京：國防工業(yè)出版社，2006：87－90.

［12］SHEN W Z，MIKKELSEN R，S?RENSEN J N，et al.Tip loss correction for wind turbines computations［J］.Wind Energy，2005，8：457－475.

［13］陳進，張石強，王旭東，等.基于粗糙度敏感性研究的風力機專用翼型設計［J］.空氣動力學學報，2011，29（2）：142－149.

［14］ SHEN WZ CHEN J CHENG J T et al.Airfoils and methods for designing airfoils［P］.Application No.PCT／EP2010／056810，International patent application.

［15］FUGLSONG P，BAK C.Development of the Ris?wind turbine airfoils［J］.Wind Energy，2004，7：145－162.

［16］肖勁松.風力機空氣動力學［M］.北京：中國電力出版社，2010，87－89.