許波峰，王同光，張震宇，王 瓏

（南京航空航天大學 江蘇省風力機設計高技術研究重點實驗室，江蘇 南京 210016）

0 引 言

風力機葉尖小翼能夠提高風力機的輸出功率，改變風輪的流場特性。風力機葉片旋轉時，由于葉尖壓力面和吸力面的壓力差，氣流會由壓力面繞過尖部端面流入吸力面，形成復雜的三維流動，同時還會產生高強度的葉尖渦。葉尖渦使得葉片效率降低、疲勞載荷增加和葉尖噪聲增大。這種現(xiàn)象在風力機葉片及飛機機翼上普遍存在。早在19世紀末就提出翼尖安裝垂直板可以明顯減輕三維效應從而減小誘導阻力［1］。風力機葉尖加小翼的想法最早由Van Holten（1976）［2］提出。Shimizu Y 等人［3］開發(fā)出一種 Mie型葉尖小翼，成功應用在小型風力機發(fā)電機上，動力輸出增加10%－20%。Johansen J等人［4］研究了分別利用扭轉和后掠的葉尖小翼提高風力機葉片效率的方法。國內內蒙古工業(yè)大學在Mie型小翼的基礎上進行了葉尖小翼的實驗研究［5］和數(shù)值模擬研究［6］。

本文優(yōu)化設計的葉尖小翼為分裂式，以NREL PhaseⅥ葉片為設計對象?？紤]到CFD方法的計算量大而遺傳算法又是要以大量個體計算為基礎，本文采用自由渦尾跡（FVW）方法［7－8］對小翼葉片建模，計算葉片的氣動性能。以風能利用系數(shù)最大和風輪推力系數(shù)最小為目標，采用快速非支配排序遺傳算法［9］（NSGA－Ⅱ），通過不斷遺傳葉尖小翼的形狀參數(shù)，搜索滿足多目標優(yōu)化的Pareto最優(yōu)解集。

1 數(shù)學模型

1．1 設計變量

分裂小翼的俯視圖和左視圖如圖1所示，從前緣至后緣三個小翼順聯(lián)排布，小翼根稍比均為0．75，有相同的后掠角χ，由前往后小翼上反角分別為Ψ1、Ψ2、Ψ3，每個小翼都具有不同的扭角分布。上反角等于0°時，小翼在葉片展向的投影長度為0．59m，保證葉片總長度與NREL PhaseⅥ葉片相同。

圖1 小翼布局示意圖Fig．1 Schematic diagram of winglet

葉片主體部分的翼型為S809，小翼的翼型為NACA0006。定義各參數(shù)的正負和取值范圍：后掠角向后緣方向為正，且0°≤χ≤60°；上反角向吸力面方向為正，且－30°≤Ψ3≤Ψ2≤Ψ1≤60°；扭角定義前緣低頭為正，三個小翼扭角分布分布滿足

1．2 目標函數(shù)

分裂小翼可以減輕三維效應，且能使葉尖渦破碎，從而減小誘導阻力。葉片載荷越大，通過小翼性能提高越明顯。確定小翼優(yōu)化的設計風速為20m／s。設計目標為風能利用系數(shù)最大和風輪推力系數(shù)最小，即

式中，CP為風能利用系數(shù)，CT為風能推力系數(shù)，NG為總的迭代代數(shù)。

2 求解方法

2．1 自由渦尾跡模型

對風力機流場作不可壓和位流假設，氣動模型可以簡化為來流、附著渦線和自由渦面的總和。尾跡建模如圖2所示，葉片主體和小翼均離散成多個葉素段，附著渦線和控制點均在1／4弦線上，從葉素邊界尾緣拖出渦線組成自由渦面。

自由渦尾跡方法允許渦線節(jié)點隨當?shù)亓魉僮杂梢苿?，從而建立渦線的偏微分控制方程。將計算域在時間（Ψ）和空間（ζ）上離散，對控制方程作五點中心差分離散，并取Δψ＝Δζ，則渦線控制方程可以簡化為

式中，V∞是來流速度，Ω是風輪旋轉速度，Vind是渦線對當前節(jié)點的誘導速度總和。計算誘導速度運用了Biot－Savart定律，為了避免直線渦元計算出現(xiàn)奇點，且考慮渦的粘性耗散效應，渦核半徑隨時間漸漸增大，變化方式與Lamb－Oseen渦核形式一致，為

式中，r0是初始渦核半徑，αL、δ和ν均為常數(shù)。尾跡迭代過程中，引入松弛因子來修正渦線節(jié)點位置

本文中，上式的松弛因子ω等于0．5。每步迭代后，計算相鄰兩步的尾跡形狀的均方根，當均方根小于10－4表示尾跡收斂。

圖2 尾跡分析示意圖Fig．2 Schematic diagram of wake analyse

利用收斂的尾跡求解葉片附近的流場信息，根據(jù)葉素氣動理論求出葉片的分布載荷以及總體氣動性能。

2．2 NSGA－Ⅱ算法

快速非支配排序遺傳算法（NSGA－Ⅱ）是對NSGA算法進行大幅改進的一種高效分類排序算法，其通過引入Pareto支配關系和個體分布的評價函數(shù)，將復雜的多目標優(yōu)化問題轉化成一類可以比較個體優(yōu)劣的單目標優(yōu)化問題，適用于任意目標數(shù)的優(yōu)化問題。

新算法主要在以下三個方面對NSGA算法加以改進：一、通過使用一種快速的非支配分類排序方法，將目標數(shù)為M，種群大小為N的計算復雜度降低；二、采用擁擠度評價方法取代適應度共享函數(shù)法，改善個體分布的均勻性；三、加入了精英保存策略以改善算法的魯棒性和收斂速度。

NSGA－Ⅱ算法流程見圖3，其中，采用擁擠距離方法作為評估種群分布的主要方法，定義每一層中個體前后相鄰兩點之間在各個目標方向上的距離之和為評價函數(shù)；采用快速非支配排序法根據(jù)個體非劣解的層次關系對種群劃分級別；按照先層級再擁擠距離的順序對個體進行選擇操作，選擇出新的父代種群，再通過遺傳變異得到第i＋1代種群。

圖3 NSGA－Ⅱ算法流程圖Fig．3 Flowchart of NSGA－Ⅱ

3 計算結果與分析

以NREL PhaseⅥ葉片為設計對象，先計算原始葉片的氣動性能，與實驗值做對比，驗證FVW方法的準確性。再使用FVW／NSGA－Ⅱ耦合方法，對小翼形狀進行優(yōu)化。

3．1 氣動模型驗證

圖4為在不同風速下的風輪低速軸扭矩，實心點為實驗結果，空心點為FVW模型計算結果。計算結果與實驗結果在各個風速下吻合度都較好。15m／s和20m／s風速計算與實驗結果存在相反的偏差，這是由于在高風速下葉片大部分剖面都處于失速甚至深失速狀態(tài)下，很難準確模擬翼型的氣動數(shù)據(jù)，而這是準確模擬葉片性能的關鍵。更大風速25m／s計算與實驗結果又很一致可能是一種巧合。從模擬結果可見，F(xiàn)VW模型能夠較準確的模擬葉片的氣動性能，可以為NSGA－Ⅱ方法提供可靠的性能目標。

圖4 低速軸扭矩比較Fig．4 Comparison of low speed shaft torque

3．2 優(yōu)化結果

小翼的優(yōu)化設計采用風能利用系數(shù)CP最大和風輪推力系數(shù)CT最小的兩目標方案，種群大小為36，迭代次數(shù)為300。

表1 計算狀態(tài)Table 1 Computational parameters

圖5給出優(yōu)化得到的分裂小翼葉片的Pareto最優(yōu)解的分布圖。兩目標優(yōu)化的最優(yōu)解并不是唯一解，而是一個解集，且近似分布在一條二次曲線上，推力系數(shù)隨著風能利用系數(shù)的增大而增大。原始葉片（原NREL葉片）用同樣的方法計算CP和CT，計算結果分別為0．02286和0．1635。加小翼的優(yōu)化結果CP最大能達到0．02977，相比原始葉片提高了30%，且CT與原來相當。可見加小翼的葉片氣動性能會有大幅提高。決策者可以根據(jù)意愿和工程需要從Pareto最優(yōu)解中選擇小翼的形狀：若要追求成本最小化，可以選擇A解的小翼形狀；若要追求發(fā)電量最大化，可以選擇C解；若要兩者折中考慮，可以選擇B解。

圖5 Pareto最優(yōu)解分布Fig．5 Distribution of the Pareto－optimal solution

使葉片性能最優(yōu)化的小翼形狀具有一定的規(guī)律性，圖6～圖8給出最優(yōu)解的幾何形狀分布。圖6中顯示，最優(yōu)后掠角基本都在20°附近，當要求高的風能利用系數(shù)時，后掠角會增大至50°左右。從圖7可以看出：三者滿足約束條件Ψ3≤Ψ2≤Ψ1；1號和2號小翼的上反角均在55°至60°之間，3號小翼的上反角則在正負10°范圍內；總有兩個小翼的上反角很接近。三個小翼的扭角約束參數(shù)A和B優(yōu)化結果均比較接近，A在16與20之間，B則在1．5與2之間。取圖5中的B解得到圖8的扭角的徑向分布，小翼根部扭角與葉片主體連接段一致，達到光滑過渡的目的。扭角的徑向變化率約等于6°／m，且變化趨勢與主體部分相反。這些小翼形狀規(guī)律對于今后小翼的設計及改型有很好的指導性作用。

圖6 后掠角分布Fig．6 Distribution of sweepback

圖7 上反角分布Fig．7 Distribution of dihedral angle

圖8 B解的扭角徑向分布Fig．8 Radial distribution of twist angle in solution B

4 結 論

采用FVW／NSGA－Ⅱ耦合方法，對風力機葉尖分裂小翼的幾何形狀優(yōu)化設計進行研究。通過與實驗結果的對比，證明FVW方法計算葉片氣動性能的可靠性。耦合方法優(yōu)化得到一個Pareto最優(yōu)解集，且分布在一條曲線上，在推力系數(shù)相當?shù)那闆r下，風能利用系數(shù)與原始葉片比較能夠提高近30%。小翼的幾何形狀在最優(yōu)解集下分布具有一定的規(guī)律性，這對于今后小翼的設計及改型有很好的指導性作用。

［1］ JUPP J．Wing Aerodynamics and the science of compromise［J］．Aeronautical Journal，2001，105（1053）：633－641．

［2］ VAN HOLTEN．Windmill with aiffuser effect induced by small tip vanes［A］．In：International Symposium on Wind Energy Systems［C］，Camdridge，UK：1976．

［3］ SHIMIZU Y，VAN BUSSEL G J，MATSUMURA S，et al．Studies on horizontal axis wind turbines with tip attachments［A］．EWEC'90Conference proceedings［C］，Madrid，Spain：1990：279－283．

［4］ JOHANSEN J，GAUNAA M，SORENSEN N N．Increased aerodynamic efficiency on wind turbine rotor using winglets［A］．26th AIAA Applied Aerodynamics Conference［C］．Hawaii，USA：2008：12－21．

［5］ 東雪青，汪建文，韓曉亮，等．風力機葉尖加小翼流場的試驗研究［J］．工程熱物理學報，2009，30（10）：1162－1164．

［6］ 賈瑞博，汪建文．小翼對水平風力機流場特性的改變［J］．華北電力大學學報，2010，37（5）：83－87．

［7］ 許波峰，王同光．基于自由渦尾跡法和面元法全耦合風力機氣動特性計算［J］．南京航空航天大學學報，2011，43（5）：592－597．

［8］ WANG T G，WANG L，ZHONG W，et al．Large－scale wind turbine blade design and aerodynamic analysis［J］．Chinese Science Bulletin，2012，57：466－472．

［9］ DEB K，PRATA PA，AGARWAL S，et al．A fast and elitist multi－objective genetic algorithm：NSGA－Ⅱ［J］．IEEE Transactions on Evolutionary Computation，2002，6（2）：182－197．