張震宇，王同光，陳 立，許波峰，王 瓏，羅 源

（1．南京航空航天大學(xué) 江蘇省風(fēng)力機(jī)設(shè)計高技術(shù)研究重點(diǎn)實(shí)驗(yàn)室，江蘇 南京 210016；2．中國空氣動力研究與發(fā)展中心 空氣動力學(xué)國家重點(diǎn)實(shí)驗(yàn)室，四川 綿陽 621000）

0 引 言

水平軸風(fēng)力機(jī)的風(fēng)能輸出受到來流風(fēng)況的強(qiáng)烈影響，在偏航以及中等、大風(fēng)速條件下由于葉片背風(fēng)面的邊界層容易發(fā)生流動分離，并會導(dǎo)致扭矩輸出和能量利用效率的明顯降低，以及分離噪聲等一系列影響。風(fēng)力機(jī)葉片的上述特性通過風(fēng)洞實(shí)驗(yàn)和計算流體力學(xué)數(shù)值模擬不斷被證實(shí)［1－2］。

另一方面，水平軸風(fēng)力機(jī)的大展弦比布局一般會產(chǎn)生明顯的葉尖渦結(jié)構(gòu)；在三維旋轉(zhuǎn)主流流動的作用下，螺旋狀發(fā)展的葉尖渦會對葉片產(chǎn)生比較強(qiáng)的誘導(dǎo)影響［3］。

如何準(zhǔn)確捕捉風(fēng)力機(jī)葉片上述流動結(jié)構(gòu)是實(shí)現(xiàn)合理高效氣動設(shè)計的前提。本文采用數(shù)值方法，獲得了以NREL phase VI風(fēng)力機(jī)葉片為原型的流場與葉尖渦結(jié)構(gòu)的準(zhǔn)確模擬結(jié)果，并與縮比模型的PIV實(shí)驗(yàn)結(jié)果形成了良好的對比。在此基礎(chǔ)上，提出了基于多片小翼概念的分裂葉尖形式概念型水平軸風(fēng)力機(jī)葉片設(shè)計方案，并針對其綜合氣動性能利用快速非支配排序遺傳算法（NSGA－II）進(jìn)行了多目標(biāo)優(yōu)化研究并獲得了滿足給定約束的優(yōu)化結(jié)果。

1 數(shù)值模擬

1．1 模型與數(shù)值方法

水平軸風(fēng)力機(jī)葉片流動具有明顯的三維旋轉(zhuǎn)和強(qiáng)烈的湍流特征?；诠こ淘试S的計算能力，除了大渦模擬（LES）與脫落渦模擬（DES）方法外，數(shù)值求解雷諾平均Navier－Stokes（RANS）方程仍然是目前工程上能夠兼顧湍流流動計算量和計算準(zhǔn)確度的一種有效方法。本文采用該方法實(shí)現(xiàn)了風(fēng)力機(jī)葉片基本流場的數(shù)值模擬。

由于暫不考慮側(cè)滑等復(fù)雜風(fēng)況，通過在上述區(qū)域中采用旋轉(zhuǎn)坐標(biāo)系求解三維不可壓RANS方程，原始的非定常葉片繞流可以轉(zhuǎn)化為新坐標(biāo)系下的準(zhǔn)定常流動進(jìn)行求解。其控制方程形式在常規(guī)笛卡爾坐標(biāo)系下的形式基礎(chǔ)上應(yīng)補(bǔ)充附加源項(xiàng)以考慮在非慣性系下的離心力與哥氏力影響。

根據(jù)葉片繞流的周向?qū)ΨQ性，計算區(qū)域采用一半流場并根據(jù)低速不可壓流動特性分別施加速度進(jìn)口、壓力出口，無滑移葉片物面以及周期邊界條件，如圖1所示。

圖1 風(fēng)力機(jī)葉片流場計算區(qū)域Fig．1 Domain for the half flow field of a wind turbine blade

計算區(qū)域大小沿流向（圖1中－z方向）為（－2R，8R），沿展向?yàn)?R（R為葉片長度）?？偩W(wǎng)格單元數(shù)超過6百萬。為了在葉片下游區(qū)域計算網(wǎng)格局部加密以捕捉尾渦結(jié)構(gòu)，總網(wǎng)格單元數(shù)接近2千2百萬。

上述控制方程的離散形式通過有限體積法進(jìn)行求解。其對流項(xiàng)采用二階迎風(fēng)格式，擴(kuò)散項(xiàng)為中心差分。壓力與速度的耦合采用SIMPLE算法實(shí)現(xiàn)［4］。RANS方程的封閉采用 Menter的k－ωSST模型［5］。

1．2 結(jié)果與討論

利用上述方法對于原始NREL phase VI風(fēng)力機(jī)葉片進(jìn)行了計算，并與NREL盲比實(shí)驗(yàn)以及文獻(xiàn)［2］中S?rensen等人的計算結(jié)果進(jìn)行了對比，見圖2。

圖2 低速軸扭矩結(jié)果對比Fig．2 Comparison of low－speed shaft torque

利用目前方法，尾流加密（22mil）和不加密（6mil）兩套網(wǎng)格算得的低速軸扭矩結(jié)果在主要風(fēng)速范圍內(nèi)與實(shí)驗(yàn)存在良好一致性。除10m／s之外主要計算結(jié)果與實(shí)驗(yàn)相對誤差均在6%之內(nèi)，說明對于該風(fēng)力機(jī)葉片的整體性能把握較好。10m／s風(fēng)速下出現(xiàn)較大偏差的主要原因是此風(fēng)況下葉片吸力面邊界層流動未完全分離，而目前所用RANS算法內(nèi)尚未引入對于流動轉(zhuǎn)捩位置的預(yù)估功能，因此導(dǎo)致該風(fēng)速下對吸力面分離流動的影響把握不準(zhǔn)。

圖3 尾流區(qū)渦量強(qiáng)度分布（周向角60°，風(fēng)速12m／s，尖速比4．91）Fig．3 Vorticity in the wake（azimuth angle of 60°，wind speed of 12m／s，and tip speed ratio of 4．91）

在此基礎(chǔ)上本文利用尾流區(qū)加密網(wǎng)格準(zhǔn)確捕捉了近尾跡區(qū)葉尖渦的發(fā)展過程。圖3中為與風(fēng)洞PIV實(shí)驗(yàn)結(jié)果的對比。該實(shí)驗(yàn)設(shè)置的詳細(xì)介紹見文獻(xiàn)［6］。由圖3可見，在給定風(fēng)速下，計算得到的尾渦渦核位置與渦核強(qiáng)度與實(shí)驗(yàn)非常一致。

2 分裂葉尖布局設(shè)計

2．1 翼尖小翼概念

長期以來關(guān)于固定翼飛行器的研究表明，翼尖小翼能夠有效的起到提高升力，抑制誘導(dǎo)阻力的效果。在此基礎(chǔ)上基于多個小翼的分裂葉尖形式證明在機(jī)翼增升減阻方面具有更為優(yōu)良的效果［7］。小翼設(shè)計在風(fēng)力機(jī)葉片的設(shè)計中在近年來也得到了應(yīng)用［8］。本文提出了基于多片小翼概念的分裂葉尖形式概念型水平軸風(fēng)力機(jī)葉片設(shè)計方案，以研究該布局對于風(fēng)力機(jī)葉片氣動性能的控制效果。其基本形式如圖4所示。

2．2 分裂葉尖數(shù)值模擬

作為初步設(shè)計，在原始葉片85%展長位置改為由三段相互獨(dú)立的小翼組成的分裂葉尖布局。各小翼均無扭轉(zhuǎn)，梢根比均為0．75，后掠角均為45°，從前緣向后緣上反角依次為40°，20°，0°。

圖4 帶三個小翼的NREL phase VI葉片模型Fig．4 NREL phase VI model with three winglets

由于采用分裂葉尖形式的氣動布局，多個小翼將原來集中于葉尖附近的環(huán)量進(jìn)行了重新分配。圖5中所示下游近尾跡區(qū)渦核附近的渦量極值隨尾跡周向角的分布情況描述了這種重新分配環(huán)量的效果。分裂葉尖的布局（見圖5中winglet1結(jié)果）使得近尾跡區(qū)渦核附近的渦量極值比原始葉片（見圖5中basic結(jié)果）大為降低，從而有效減弱下游集中渦的誘導(dǎo)作用。

圖5 尾流區(qū)渦量極值沿周向角的分布Fig．5 Azimuthal distribution of vorticity maximum in the wake

另一方面，各小翼的厚度均接近6%，且?guī)в写蠛舐咏?，這種設(shè)計使得在中等和大風(fēng)速下葉片吸力面在前緣附近獲得更大的吸力峰（如圖6中winglet1黑線所示），并體現(xiàn)為相應(yīng)風(fēng)速下軸扭矩輸出的明顯增長，在10m／s－25m／s范圍內(nèi)，增長幅度在7%－40%之間，見圖7中winglet1結(jié)果。

圖6 分裂葉尖布局80%展長位置壓力分布，25m／sFig．6 Pressure distribution at 80%span，25m／s

圖7 分裂葉尖布局低速軸扭矩Fig．7 Low－speed shaft torque of the blade with multiple winglets

3 數(shù)值優(yōu)化

在上述工作基礎(chǔ)上，對葉片關(guān)鍵幾何參數(shù)利用快速非支配排序遺傳算法（NSGA－II）進(jìn)行了多目標(biāo)數(shù)值優(yōu)化。為問題簡化考慮，優(yōu)化風(fēng)速定為20m／s，小翼根稍比均為0．75。

該優(yōu)化問題中需對以下兩目標(biāo)函數(shù)最小化：

其中Cp為風(fēng)能利用系數(shù)，CT為軸向推力系數(shù)。由于CFD方法所需計算量過大，優(yōu)化部分的上述目標(biāo)函數(shù)采用自由渦尾跡方法計算得出，具體算法見文獻(xiàn)［10］。

優(yōu)化變量及約束條件為：

小翼上反角（向吸力面方向?yàn)檎?/p>

小翼扭角：

3．1 遺傳算法

快速非支配排序遺傳算法（Fast and Elitist Nondominated Sorting Genetic Algorithm，NSGAII）是Deb等人于2002年在對NSGA算法進(jìn)行了大幅改進(jìn)的基礎(chǔ)上提出的一種高效分類排序算法［9］。該算法引入了Pareto支配關(guān)系和個體分布評價函數(shù)，從而可以將復(fù)雜多目標(biāo)優(yōu)化問題轉(zhuǎn)化成可比較個體優(yōu)劣的單目標(biāo)優(yōu)化問題，能夠用于進(jìn)行任意目標(biāo)個數(shù)的數(shù)值優(yōu)化工作。該算法根據(jù)Pareto支配關(guān)系對種群內(nèi)個體首先進(jìn)行非支配分層以確定各個個體的層級；再用擁擠距離評價函數(shù)控制各個層級個體分布；然后用個體的層級個體之間距離做選擇操作，并結(jié)合傳統(tǒng)遺傳算法的其它操作過程，最終得到Pareto最優(yōu)解集。

3．2 結(jié)果與討論

優(yōu)化種群大小選取為36，迭代次數(shù)300次。得到的Pareto最優(yōu)解集分布如圖8所示。利用自由渦尾跡方法算得的基本型NREL PHASE VI葉片的風(fēng)能利用系數(shù)Cp與軸向推力系數(shù)CT分別為0．02286和0．1635。比較發(fā)現(xiàn)，目前得到的優(yōu)化解風(fēng)能利用系數(shù)分布均明顯高于該參考值，相對提高27%～31%不等；同時軸向推力系數(shù)CT與NREL PHASE VI原型葉片相比變化不大。

圖8 Pareto最優(yōu)解集分布Fig．8 Pareto solution sets

由NSGA－II算法得到的上述最優(yōu)解集在給定目標(biāo)函數(shù)空間中給出了一組滿足相應(yīng)約束條件的最優(yōu)解，由于該解集任一個最優(yōu)解并不更優(yōu)于其他解，因此在最終的優(yōu)化方案取決于研究者根據(jù)問題性質(zhì)的具體分析與著眼點(diǎn)。由圖8得知，在目前的優(yōu)化范圍內(nèi)，隨著風(fēng)能利用系數(shù)的提高，一般伴隨有軸向推力的增長。選取三個典型狀態(tài)如表1所示，并分別對應(yīng)于軸向推力最小、折中、能量轉(zhuǎn)化率最大三個考慮因素。

4 結(jié) 論

借助于準(zhǔn)確的數(shù)值模擬方法得知，分裂葉尖形式氣動布局能夠有效地分散風(fēng)力機(jī)葉尖渦量的分布，從而減小其誘導(dǎo)影響；通過合理布置小翼，能夠使葉片端部區(qū)域在中等和大風(fēng)速下獲得額外的前緣吸力，實(shí)現(xiàn)更高的氣動效率，并為進(jìn)一步在此風(fēng)速范圍的氣動控制策略提供可行的技術(shù)準(zhǔn)備。

基于NSGA－II遺傳算法的優(yōu)化方法能夠在給定的性能約束范圍內(nèi)給出準(zhǔn)確的全局最優(yōu)解分布。對于上述問題，在不明顯增大軸向推力的前提下，以風(fēng)能利用系數(shù)形式給出的氣動效率得到明顯改善，比原始設(shè)計提高達(dá)30%左右。以上討論表明，NSGA－II遺傳算法與自由渦尾跡方法結(jié)合，是一種行之有效的風(fēng)力機(jī)葉片氣動性能優(yōu)化方法。

［1］ HAND M M，SIMMS D A，F(xiàn)INGERSH L J，et al．Unsteady aerodynamics experiment phase VI：Wind tunnel test configurations and available data campaigns［R］．NREL／TP－500－29955，NREL，2001．

［2］ S?RENSEN N N，MICHELSEN J，SCHRECK S．Navier－Stokes predictions of the NREL phase VI rotor in the NASA Ames 80ft×120ft wind tunnel［J］．Wind Energy，2002，5（2－3）：151－169．

［3］ VERMEER L J，SORENSEN J N，CRESPO A．Wind turbine wake aerodynamics［J］．Progress in Aerospace Sciences，2003，39（6－7）：467－510．

［4］ VERSTEEG H K，MALALASEKERA W．An introduction to computational fluid dynamics－the finite volume method［M］．2nd Ed．，Pearson Education Limited，2007

［5］ MENTER FR．Two－equation eddy－viscosity turbulence models for engineering applications［J］．AIAA Journal，1994，32（8）：1598－1605．

［6］ 肖京平，武杰，陳立，史喆羽，風(fēng)力機(jī)葉尖渦尾跡結(jié)構(gòu)PIV測量實(shí)驗(yàn)研究［J］．應(yīng)用數(shù)學(xué)和力學(xué)，2011，32（6）：683－692．

［7］ 陳明巖，齊孟卜，翼尖帆片的增升減阻研究［J］．航空學(xué)報，1994，15（6）：641－646．

［8］ JOHANSEN J，SORENSEN N N．Aerodynamic investigation of winglets on wind turbine blades using CFD［R］．Ris?－R－1543（EN），2006．

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

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