一種實用的風力機葉片氣動設計方法研究

昝丙合，楊旭東

（西北工業大學，翼型葉柵國家重點實驗室，陜西 西安 710072）

0 引 言

風能作為環境友好型的可再生能源越來越受到人們的重視，對于風力機葉片的設計，尤其是優化設計顯得越來越重要。目前存在的槳葉設計理論有Betz理論、Glauert模型、Wilson模型和儒科夫斯基模型。其中，Betz理論由于簡化條件過多，因此其設計結果與實際有較大偏差。Glauert優化設計方法考慮了風輪尾流，但忽略了翼型阻力和葉尖損失的影響。Wilson理論對Glauert理論進行了修改，研究了葉尖損失與升阻比對葉片性能的影響。因此，在傳統的風力機槳葉優化設計方法中，最為常用為Glauert優化設計模型以及Wilson修正模型［1］。這些方法都忽略了一些重要的問題，包括風場風況、槳葉各截面之間的約束等。所以采用這些方法設計出來的槳葉是否是最優的還值得商榷［2］。在現代葉片的優化設計過程中都要考慮風場的風速概率分布，通常采用Weibull模型［3－4］，設計目標為在滿足額定功率的要求下，年能量輸出最大［5－7］。本文在采用 PROPID首先對葉片確定一個較優幾何外形的基礎上，采用遺傳算法進一步進行優化。結果發現，與原始葉片相比，風力機的年能量輸出有了較大幅度的提高。

1 PROPID反設計原理

PROPID程序［8－10］是在原性能預估程序PROP［11］基礎上發展而來的。PROPID程序中用于風輪性能預估的翼型數據并不唯一，不同的葉片截面可以選用不同的翼型。而且同一個截面也可以給定不同雷諾數下的翼型數據，根據雷諾數，對升力系數進行線性插值，對阻力系數進行邏輯插值。并通過普朗特損失模型、三維失速延遲模型等對二維翼型數據進行修正，以使得對風輪性能的預估更加準確。此外，PROPID最大的特點是其強大的反設計能力，它允許設計者給定風力機的性能特性，如功率峰值；或者給定葉片的氣動特性，如沿葉片展向的升力系數分布進行多點反設計。PROPID程序主要是通過牛頓迭代來進行反設計。其反設計步驟如圖1所示：

圖1 PROPID反設計步驟Fig.1 Inverse design theory of PROPID program

其中：c0—葉片根部弦長，

?—各個站位相對于根部弦長的差值，

β0.75R—0.75半徑處的扭轉角，

θ—各個站位相對于0.75半徑處扭轉角的差值，

2 遺傳算法

遺傳算法（Genetic Algorithms，GA）［12－13］研究的歷史比較短，20世紀60年代末到70年代初，主要由Michigan大學的John Holland等研究形成了一個比較完整的理論和方法。遺傳算法由于只用編碼和適應度表示問題，并不要求明確的數學方程及導數表達式，因此較之傳統的優化算法通用性更強。

2.1 目標函數

將風輪扭矩M中與葉片幾何外形有關的一切元素統稱為g，則扭矩為：

那么功率為：

其中，ω為風輪的旋轉速度。對于給定風速、轉速與槳距角的情況下，功率Ρ只與葉片的幾何形狀有關，因此，優化過程中的目標函數定義為：

其中pi表示風速為Vi的概率。

2.2 編碼控制

在優化時，選定幾個截面作為控制點，控制點之間通過多項式進行插值。在控制點處給定弦長與扭轉角的上下限。弦長的上限定義為cu（r），下限定義為cl（r）；同樣扭轉角的上限定義為ψu（r），下限定義為ψl（r）。假設種群大小為N，則半徑r處的弦長與扭轉角分別定義為：

其中，bc（r），bψ（r）∈［0，2N－1］。

為了避免在種群中產生使幾何形狀不符合要求的個體，附加約束：

在對功率優化的程序中，如果不附加其他約束的話，弦長總是趨向于上限值。為了避免上述情況的發生，對葉片面積附加約束，在此給定葉片面積的上限，（1＋δA）A0，其中A0為葉片的初始面積。因此有葉片面積的附加約束為：

為了生成光滑的葉片幾何外形，控制點之外的其他半徑處的扭轉角與弦長分布由5階貝塞爾函數差值得到。

3 算例分析

3.1 驗證算例

為了說明PROPID程序對風輪性能預測的準確性，選用美國國家可再生能源實驗室（NREL）開發的Phase VI實驗風力機葉片作為驗證模型，圖2為Phase VI的剖面翼型圖與葉片幾何外形圖，整個葉片均采用S809翼型。圖3為采用PROPID程序對風輪功率的預測結果。

圖2 Phase VI葉片幾何外形Fig.2 Blade geomtry of Phase VI

圖3 RPOPID程序對Phase VI風力機性能預測結果與實驗值的比較Fig.3 Comparison of the performance between the prediction of PROPID and experimental data

發現PROPID程序能夠很好地對風力機性能進行預測。

3.2 采用PROPID進行反設計

本算例對某2葉片，600kW風力機葉片進行反設計，具體的設計參數如表1所示。

表1 設計風力機參數Table 1 Parameters of wind turbine

對于變速、變槳距的風力機而言，可以通過改變轉速與槳距角來控制葉尖速比的大小，使其始終在功率系數最大的條件下進行工作。因此，設計目標就是使風力機的功率系數峰值最大。在本算例中，通過給定沿葉片展向的升力系數分布和周向誘導因子分布，來使得風輪功率系數峰值最大。通過平移和扭轉的過程來確定最佳的分布曲線。圖4為平移過程對風輪功率峰值的影響，圖5為扭轉過程對風輪功率峰值的影響。由于軸向誘導系數的平移與扭轉擾動過程與升力系數相同，因此只畫出升力系數的平移與扭轉擾動關系圖。

最終確定軸向誘導因子和升力系數的分布如表2所示。

圖4 平移過程對功率系數峰值的影響Fig.4 Variation of CPmaxbecause of shifting of target CL／axial induction distribution

圖5 扭轉擾動下對功率系數峰值的影響Fig.5 Variation of CPmaxbecause of rotating of target CL／axial induction distribution

表2 軸向誘導因子與升力系數分布Table 2 Optimized CL／axial induction distribution

3.3 運用遺傳算法進行優化設計

以PROPID反設計得到的葉片為基礎進行優化，優化后得到的葉片幾何如圖6所示，從圖中可以看到，在葉片根部弦長變化梯度緩和，沿葉片展向，弦長有不同程度的減小，這使得葉片的結構安全性更高。當地扭轉角的分布變化平緩，有利于葉片的制造。

圖7為優化前后的功率曲線對比，可以看出在風速達到12m／s之前，優化后的風輪功率優于原始風輪。

圖8為原始葉片、PROPID反設計葉片以及優化后葉片的年能量輸出對比。由于當風速達到12m／s時，風輪達到了額定功率600kW，因此在計算年能量輸出時，只計算風速在12m／s之前的情況。在計算能量輸出時風速分布采用Weibull分布，尺度參數選為8.6，形狀參數選為2.66。原始風輪年能量輸出為1063169kWh／年，PROPID反設計得到的風輪的年能量輸出為1093854kWh／年，優化后風輪年能量輸出為1126017kWh／年。比原始葉片提高了5.9114%。

圖6 優化前后幾何形狀對比Fig.6 Comparison of blade parameters between the optimized blade and initial blade

圖7 優化前后功率對比Fig.7 Comparison of power between the optimized blade and initial blade

圖8 優化前后年能量輸出對比Fig.8 Comparison of the Annual Energy Production

4 結 論

本文提出了水平軸風力機葉片的優化設計模型，在PROPID進行反設計的基礎上，采用該模型進行優化設計。以某600kW風力機為設計目標，結果發現，結合PROPID程序與該優化設計模型，能夠很好地提高風力機的性能，年能量輸出有了較大幅度的提高。采用5階貝塞爾函數對控制點之間的幾何進行插值獲得光滑的葉片幾何外形。因此，本文的優化設計方法能夠高效地獲得實用的高性能的風力機葉片。

［1］SPERA D A.Wind turbine technology［M］.New York：ASME Press，1994.

［2］MORGAN C A，GARRAD D.The design of optimum rotors for horizontal axis wind turbines［M］.London：Garrad Hassan and Partners，1998.

［3］DET NORSKE VERITAS，RISφNATIONAL LABORATORY.Guidelines for design of wind turbines［M］.Det Norske Veritas and RisφNational Labaratory，2002.

［4］SPERA D A.Wind turbine technology［M］.New York：ASME Press，1994.

［5］MORGAN C A，GARRAD D.The design of optimum rotors for horizontal axis wind turbines［M］.London：Garrad Hassan and Partners，1998.

［6］GIGUERE P，SELIG M S.Blade design trade－offs using low－lift airfoils for stall－regulated HAWTs［A］.Proceedings of ASME／AIAA Wind Energy Symposium［C］.Reno，Nevada Jan，1999.

［7］BENINI E，TOFFOLO A.Optimal design of horizontal－axis wind turbines using blade－element theory and evolutionary computation［J］.ASME Journal of Solar Energy Engineering，2002，124：357－363.

［8］SELIG M S，TANGLER J L.Development and application of a multipoint inverse design method for horizontal axis wind turbines［J］.Wind Engineering，1995，19（2）：91－105.

［9］GIGUERE P，SELIG M S.Aerodynamic blade design methods for horizontal axis wind turbines［A］.13thCanadian Wind Energy Association Conference and Exhibiton［C］.Quebec City，Canada，1997.

［10］SELIG M S.PROPID user manual（beta version 3.0）［M］.1998.

［11］HIBBS B，RADKEY R L.Small wind energy conversion systems（SWECS）rotor performance model comparison study［M］.Aerovironment，Inc.，Prepared for Rockwell International Corporation，1981.

［12］GOLDBERG D.Genetic algorithms in search，optimization and machine learning［M］.Addison－Wesley，1989.

［13］DAVIS L.Handbook of genetic algorithms［M］.Van Nostrand Reinhold，1991.