覆冰對風力機專用翼型氣動性能影響的數值研究

張聘亭, 楊 濤, 盧緒祥

(1.華中科技大學能源與動力工程學院430074;2.長沙理工大學能源與動力工程學院410114)

風能是一種清潔可再生能源,面對越來越嚴重的能源危機和環境問題,利用綠色風能是普遍的趨勢和潮流.截至2010年年底,中國風力發電的裝機容量已經達到了44 730 MW,占全球裝機容量的22.4%[1].我國風能資源主要集中在北部(東北、華北、西北)地區、東南沿海及其島嶼以及內陸個別地區,不同地域的環境差別很大,大部分地區氣候變化急劇,溫差幅度較大,都有不同程度的冰凍.處于臨界狀態的雨、雪、霧或露遇到低溫的設備和金屬結構表面時會在設備以及結構的表面結冰,因此運行在寒冷氣候環境中的風力機經常面臨嚴重的結冰問題[2-3].風力機葉片的覆冰問題在20世紀90年代就引起了國外學者的關注:Bose[4-5]研究了小型水平軸風力機覆冰的過程和形態,并觀測了覆冰過程中風力機輸出功率的變化,對壓力面和吸力面的覆冰狀態做了分析;隨著CFD技術的發展,國內外的學者開始利用CFD技術對風力機葉片的覆冰問題進行研究,Fortin等[3]和Imamura等[6]利用CFD技術對風力機旋轉過程中的覆冰過程進行了模擬,Fortin等[3]還對風力機葉片融冰過程進行了模擬,并對覆冰后風力機葉片的性能進行了分析;李聲茂等[7]對寒冷氣候條件下直線翼垂直軸風力機葉片表面結冰進行了數值模擬.

風力機葉片是風力機獲取風能的關鍵部位,其氣動性能直接關系到風力機的運行效率,甚至關系到風力機運行的安全性.風力機翼型的升阻力曲線直接關系到風力機的性能,其斜率也影響著風力機運行的穩定性[8],覆冰所導致的翼型失速被推遲也會致使風力機過載[2].因此,研究風力機翼型在覆冰情況下的氣動性能改變,對風力機覆冰狀態下運行決策起著重要的作用.

筆者采用通用CFD軟件Fluent 6.3,利用S-A湍流模型,對風力機專用翼型S809不同覆冰狀態下的靜態失速特性進行了研究,得到了靜態失速特性下的升阻力曲線,并利用 Fluent的自定義函數UDF以及動網格對覆冰翼型的動態失速特性進行了研究.通過和非覆冰翼型的比較,得到了覆冰對翼型性能和動態失速特性的影響.

1 數值模型

1.1. S-A湍流模型

Spalart-Allmaras湍流模型由Spalart和 Allmaras于1992年提出[9].該模型是一種相對簡單的通過求解工作變量ν～的輸運方程得到湍流運動黏度νt的單方程湍流模型,其具體形式如下:

S-A模型在計算受制于壓力梯度的邊界層流動中取得了很好的結果,因此,S-A模型被普遍應用于透平機械的數值計算中.文獻[10]～文獻[12]分別利用S-A模型對翼型動態失速特性進行了仿真,得到了比較正確的結果.

1.2 數值方法

S809翼型是美國國家可再生能源實驗室(NREL)設計的水平軸風力機專用翼型之一.根據Bose[3]對風力機葉片覆冰形態的研究,筆者建立了S809翼型以及S809不同覆冰形態翼型的數值模型.本文的數值模型在Gambit2.0中建立,整體計算域如圖1所示,總長度為45c,高度為30c,c為翼型的弦長.圖2為S809翼型及其周圍網格的局部放大圖,翼型周圍采用結構化C型網格.在正式計算之前,對計算網格進行了無關性分析,翼型表面分布了300點,邊界層第一層高度為翼型弦長的 1×10-4,總網格數為18萬.根據文獻[3]對不同溫度及環境條件下覆冰形態數值仿真的結果建立了翼型可能的覆冰模型,并給出了周圍的網格,見圖3.S809翼型未覆洋情況定義為R0.

圖1 整體計算域Fig.1 The region of integ ral computation

圖2 S809翼型及其周圍的網格局部放大圖Fig.2 S809 airfoil and local enlargement of the surrounding g rid

圖3 翼型不同的覆冰形態及其周圍的網格局部放大圖Fig.3 Different icing patterns on airfoil and local enlargement of the surrounding grid

數值計算中對連續方程、二維動量方程、能量方程、S-A單方程湍流模型四個方程進行耦合求解,方程采用隱式格式,當流場中所有計算點上的最大殘差小于1×10-6時,認為計算收斂.

2 靜態翼型數據的求解和討論

靜態翼型數值計算的邊界條件:來流馬赫數為0.3,雷諾數為1×106.計算不同攻角下S809翼型的升力系數,并分析不同覆冰情況對升力系數的影響.計算得到S809翼型的升力系數曲線如圖4所示.

圖4 S809翼型不同覆冰情況下升力系數Fig.4 T he lift coefficient under different icing conditions of S809 airfoil

在翼型覆冰朝吸力面方向生長的情況下,如圖3(a)和圖3(b)的R1和R2覆冰狀態下,覆冰降低了翼型的彎曲,造成了翼型升力系數的下降.當攻角較小時(＜6°),翼型的升力系數改變并不明顯,但在較大的攻角范圍內,無論是在線性區域還是失速區域,R1和R2的覆冰情況都帶來了升力系數的較大下降.當攻角為0°時,在覆冰形態 R1的情況下,覆冰層后就已經出現了分離渦;且隨著攻角的增加,分離渦不斷地向尾緣方向生長(見圖5),造成了升力系數的較大下降.對于覆冰形態R2,覆冰層后分離渦的生長慢于R1,但同樣也造成了升力系數的下降,見圖6.

圖7為攻角同為8°時,不同的覆冰形態對葉片周圍分離渦的影響.從圖7可以看出:和R0覆冰形態相比,R1和R2的覆冰形態不但造成了前緣的分離渦,并且也加速了尾部分離渦的生成,造成了翼型失速的提前,引起升力系數的下降;對于R3覆冰形態,覆冰沿弦向生長,當攻角較小時,翼型的升力系數改變并不明顯,但隨著攻角的不斷增加,R3覆冰形態使失速提前了,造成了升力系數的下降;對于R4覆冰形態,覆冰主要朝壓力面方向生長,覆冰增加了翼型彎度,在一定的程度下反而提高了升力系數,并且R4覆冰形態推遲了翼型尾緣分離渦的形成,一定程度上推遲了失速,增加了升力系數,可能造成風力機運行的過載.

圖5 R1覆冰形態下覆冰層后的局部流場Fig.5 Local flow field following the icing layer of R1 pattern

圖6 R2覆冰形態下覆冰層后的局部流場Fig.6 Local flow field following the icing layer of R2 pattern

圖7 攻角為8°時,S809翼型不同覆冰情況下翼型周圍的流場Fig.7 The surrounding flow field under different icing conditions of S809 airfoil at an attack angle of 8

3 動態翼型數據的求解和討論

動態翼型的數值計算邊界條件同靜態翼型數據:來流馬赫數為0.3,雷諾數為1×106.根據文獻[12],翼型繞 x/c=0.25(x表示該點到翼型前緣的距離)處作正弦周期振蕩,正弦變化規律為:

式中 :α0=20°;Δα=11°;取衰減頻率為 k=0.4,則振蕩頻率ω=83.132 5;時間步長取0.001 26 s.

經過4個周期的計算后得到了基本一致的滯回曲線,則認為流場趨于穩定.分別計算S809翼型及其覆冰情況下的動態升力系數曲線,如圖8所示.

圖8 S809翼型不同覆冰狀態下的動態升力系數曲線Fig.8 The dy namic lift coefficient curve under different icing conditions of S809 airfoil

從翼型的升力系數曲線來看,其動態失速特性和靜態失速特性有著明顯的區別.在靜態流場中,翼型的升力系數和攻角是一一對應的.然而在動態流場中,當翼型突然從低于失速攻角增加到高于失速攻角時,則可以看到翼型升力系數相對于靜態升力系數的過調,在翼型的振蕩過程中,翼型的升力系數曲線形成了一條滯回曲線,這種動態失速特性也是失速顫振氣動負阻尼的來源.

結合圖8和圖9分析覆冰對動態失速特性的影響.對于S809翼型,和未覆冰狀態(見圖9(a))相比,R1和R2覆冰狀態提前了翼型失速過程的發生,而R4覆冰狀態(圖9(d))則推遲了翼型失速過程的發生.這和靜態翼型數據的分析結論是一致的,也可以看出覆冰對翼型失速的影響.對于R4覆冰狀態,其升力系數曲線的斜率要大于 R0無覆冰狀態(見圖8),這從氣動彈性穩定的角度考慮是有益的.

圖9 S809翼型不同覆冰狀態下翼型周圍流場圖Fig.9 The surrounding flow field under different icing conditions of S809 airfoil

4 結 論

(1)向吸力面生長的覆冰形態會造成覆冰層后的分離渦,并且隨著攻角的增加分離渦會向后緣生長,造成升力系數的較大下降,同時也會造成阻力系數的增加;

(2)翼型的前緣覆冰會對翼型尾緣分離渦的生成產生影響,對于向吸力面生長的覆冰,會促進尾緣分離渦的生成,而向壓力面生長的覆冰則能抑制尾緣分離渦的生長,從而造成對翼型升力系數的影響;沿弦向生長的覆冰情況,則對尾緣分離渦的生成影響較小;

(3)在動態失速情況下,翼型周圍流場比較復雜,覆冰的形態對翼型升力系統造成影響規律也比較復雜,但其分離渦生長的基本規律同靜態流場時;

(4)從翼型覆冰情況下的靜態流場特性分析,覆冰破壞了翼型的流線,從而直接影響了翼型的氣動性能.當翼型升力系數降低時,風力機的出力將受到影響,經濟性降低;當翼型失速被推遲,又容易造成風力機的過載,影響安全性;

(5)從翼型覆冰情況下的動態流場特性分析,覆冰改變了翼型動態升力系數曲線的斜率;而升力系數曲線的斜率將影響風力機的氣動彈性穩定性;

(6)由于風力機覆冰情況復雜,不同的覆冰厚度、不同的覆冰形態等都會對覆冰翼型的氣動性能造成影響.本文的工作僅是對該問題從數值模擬的角度開展前期的研究.為了了解風力機在覆冰情況下的運行狀況,提高風力機在覆冰情況下運行的安全性和穩定性,還需要長期的探索.

[1]李俊峰,蔡豐波,唐文倩,等.風光無限:中國風電發展報告 2011[M].北京:中國環境科學出版社,2011.

[2]DALILI N,EDRISY A,CARRIVEAU R.A review of surface engineering issues critical to wind turbine performance[J].Renewable and Sustainable Energy Reviews,2009,13(2):428-438.

[3]FORTIN G,PERRON J.Wind turbine icing and deicing[C]//47th AIAA Aerospace Sciences Meeting Including The NewHorizons Forumand Aerospace Exposition.Orlando,USA:AIAA,2009.

[4]BOSE N.Icing on a small horizontal-axis wind turbine-part 1:glaze ice profiles[J].Journal of Wind Engineering and Industrial Aerodynamics,1992,45(1):75-85.

[5]BOSE N.Icing on a small horizontal-axis wind turbine-part 2:three dimensional ice and wet snow formations[J].Journal of Wind Engineering and Industrial Aerodynamics,1992,45(1):87-96.

[6]IMAM URA Y,TODA K,YAMAMOTO M.Numerical simulation of ice accretion on wind turbine blade[C]//Abstracts of the Papers Presented at the Minisymposia Sessions of WCCM VI in Conjunction with APCOM.Beijing:[s.n.],2004.

[7]李聲茂,李巖.垂直軸風力機葉片翼型結冰的數值計算[J].動力工程學報,2011,31(3):214-219.LI Shengmao,LI Yan.Numerical simulation on icing of a blade aerofoil for vertical axis wind turbines[J].Journal of Chinese Society ofPower Engineering,2011,31(3):214-219.

[8]HANSEN M.Aerodynamics of wind turbines[M].[S.1.]:Earthscan,2008.

[9]SPALA RT P R,ALLM ARAS S R.A one-equation turbulence model for aerodynamic flows[J].La Recherche Aerospatiale,1994,1(1):5-21.

[10]李杰.風力機葉片翼型氣動特性數值研究[D].北京:華北電力大學能源動力與機械工程學院,2009.

[11]趙旭,肖俊,席德科.水平軸風力機翼型設計與動態失速數值模擬[J].太陽能學報,2009,30(3):348-354.ZHAO Xu,XIAO Jun,XI Deke.The design of airfoils and the simulation of dynamic stall of horizontal axis wind turbines[J].Acta Energiae Solaris Sinica,2009,30(3):348-354.

[12]陳旭,郝輝,田杰,等.水平軸風力機翼型動態失速特性的數值研究[J].太陽能學報,2003,24(6):735-740.CHEN Xu,HAO Hui,TIAN Jie,et al.Investigation on airfoil dynamic stall of horizontal axis wind turbine[J].Acta Energiae Solaris Sinica,2003,24(6):735-740.