基于CFD的旋轉(zhuǎn)風(fēng)輪氣動性能分析

李少華, 岳巍澎, 王東華, 岳 征

(1.東北電力大學(xué)能源與動力工程學(xué)院,吉林 132012;2.張家口出入境檢驗檢疫局,張家口 075000)

風(fēng)力發(fā)電機是將風(fēng)能轉(zhuǎn)換為電能的機械裝置,葉片以及風(fēng)輪的氣動性能直接決定了風(fēng)力發(fā)電機組的效率.由葛勞渥旋渦理論可知:空氣流過旋轉(zhuǎn)葉片時會形成中心渦、附著渦和螺旋后緣尾渦[1].由于渦系對風(fēng)輪氣動性能有較大干擾,因此風(fēng)輪中葉片所處的流場特性分布與單葉片有明顯區(qū)別,旋轉(zhuǎn)葉片渦的脫落對旋轉(zhuǎn)下游葉片造成很大影響,導(dǎo)致風(fēng)輪下風(fēng)向流場周圍環(huán)境發(fā)生變化.

目前,國內(nèi)外有關(guān)大型風(fēng)力發(fā)電機風(fēng)輪流場的研究有試驗和數(shù)值模擬兩種方法:試驗方面,有關(guān)研究者[2-4]多采用PIV粒子圖像等技術(shù)測量風(fēng)力發(fā)電機葉尖及尾流區(qū)的流場工況,但試驗成本較高,對環(huán)境要求嚴(yán)格,且試驗研究本身針對性強、局限性大、不適合觀察和分析旋轉(zhuǎn)風(fēng)輪的整體流場特性;數(shù)值模擬方面,郭婷婷等[5]和吳殿文[6]分析了千瓦級風(fēng)力發(fā)電機二維翼型失速特性,并獲得了旋渦產(chǎn)生、發(fā)展和脫落規(guī)律;包飛[7]采用SST k-ω模擬單葉片周圍流場,對葉片截面所在平面的流場域進行了分析,并系統(tǒng)研究了壓力面與吸力面的壓差以及升力系數(shù);S?RENSEN 等[8]和SEZER-UZOL等[9]對千瓦級風(fēng)力發(fā)電機風(fēng)輪進行了模擬,并重點描述了風(fēng)輪旋轉(zhuǎn)平面的流動特性.

到目前為止,國內(nèi)外系統(tǒng)分析兆瓦級風(fēng)力發(fā)電機旋轉(zhuǎn)風(fēng)輪流場特性以及介紹旋轉(zhuǎn)葉片間相互影響方面的文獻較少.筆者對1.2 MW風(fēng)力發(fā)電機旋轉(zhuǎn)風(fēng)輪和單葉片流場進行了數(shù)值模擬,經(jīng)迭代計算收斂后,截取葉片周向曲面,等長度拉伸成平面,并按周期性邊界進行了排列,重點研究了旋轉(zhuǎn)葉片運動軌跡所在平面的流場分布,并通過觀察湍動能場渦脫落云圖分析了葉片間的相互影響.

1 數(shù)學(xué)模型

假設(shè)葉片為剛體且模擬過程不考慮葉片表面的變形.基于穩(wěn)態(tài)不可壓縮流動三維定常雷諾時均NS方程(RANS)進行了數(shù)值模擬,采用Segregated隱式求解器三維穩(wěn)態(tài)算法,紊流模型采用SST k-ω模型,壓力-速度耦合使用Simple算法,對流項差分格式則采用二階迎風(fēng)格式[10].

式中:φ為通用變量;u為速度;ρ為流質(zhì)密度;Γ為廣義擴散系數(shù);S為廣義源項.

2 數(shù)值計算方法

2.1 幾何模型建立與網(wǎng)格劃分

采用Fluent軟件的前處理軟件Gambit進行了幾何建模,并對基于Matlab自主編程設(shè)計的1.2 MW風(fēng)力發(fā)電機葉片進行了模擬[11],其設(shè)計風(fēng)速為11.26 m/s,風(fēng)輪額定轉(zhuǎn)速為18.44 r/min.4組葉片均采用相同的流場區(qū)域.圖1為風(fēng)輪流場分區(qū)的主要尺寸及其網(wǎng)格劃分.在圖1中,葉片長為35 m,取輪轂直徑為3.25 m.

筆者對本文中風(fēng)力發(fā)電機小域中葉片與輪轂進行了網(wǎng)格局部加密處理,并采用非結(jié)構(gòu)網(wǎng)格劃分,而對流場的其他部分則采用了相對較稀疏的結(jié)構(gòu)網(wǎng)格.為了保證網(wǎng)格質(zhì)量,筆者還進行了分區(qū)劃分.圖2為不同葉片數(shù)量的風(fēng)輪網(wǎng)格示意圖.由于在對比不同數(shù)量葉片的葉輪時,若保持每個葉片近壁面貼體網(wǎng)格數(shù)量一致,則會導(dǎo)致整體網(wǎng)格的數(shù)量相差很多.在模擬過程中,將速度和湍流度變化梯度較大的部分采用自適應(yīng)網(wǎng)格處理,每次提高6～8萬網(wǎng)格后繼續(xù)迭代計算,直至輸出功率不再明顯提高為止.

圖1 風(fēng)輪流場分區(qū)的主要尺寸及其網(wǎng)格劃分(單位:m)Fig.1 Partition,main dimension and g rid of flow field(unit:m)

圖2 風(fēng)輪的網(wǎng)格劃分Fig.2 Grid division of wind wheel

2.2 邊界條件的設(shè)置

進口:選擇速度進口邊界條件,假定進口處具有相同的風(fēng)速,不考慮風(fēng)速變化的影響,進氣方向垂直于進口面,給出進口風(fēng)速.

出口:假定出口流動充分發(fā)展,此處選出口邊界為自由出流.

葉片及輪轂壁面無滑移且絕熱,設(shè)定為旋轉(zhuǎn)壁面,相對附近區(qū)域的旋轉(zhuǎn)速度為零.

3 性能計算與誤差分析

根據(jù)下式計算風(fēng)力發(fā)電機轉(zhuǎn)子的輸出功率.

式中:P為輸出功率,W;T為轉(zhuǎn)矩,N?m;n為葉輪的轉(zhuǎn)速,r/min;b為葉片數(shù).

在得到各葉片的輸出功率后,根據(jù)式(3)可獲得相對誤差,再根據(jù)相對誤差的大小分析風(fēng)輪各葉片及單葉片對風(fēng)力發(fā)電機出力的影響.

式中:γ為相對誤差,%;P0=1.2 MW,為設(shè)計功率.

根據(jù)式(2)和式(3),結(jié)合數(shù)值模擬結(jié)果可獲得轉(zhuǎn)矩.表1為不同葉片數(shù)的模擬結(jié)果.從表1可看出:隨著葉片數(shù)的增多,風(fēng)力發(fā)電機轉(zhuǎn)子的輸出功率逐漸減小.這是由于旋轉(zhuǎn)上游的葉片發(fā)生渦脫落,使得旋轉(zhuǎn)下游葉片周圍流場的環(huán)境發(fā)生改變,導(dǎo)致其輸出功率減小.可見,葉片數(shù)越多,葉片間距離越近,葉片間的相互影響就越大.其中,三葉片風(fēng)輪的輸出功率最接近設(shè)計輸出功率,這是由于葉片翼型的所有參數(shù)均是按目前普遍使用的兆瓦級三葉片風(fēng)輪機設(shè)計的.四葉片風(fēng)輪的總輸出功率最大,但是由于其成本較昂貴,因此在實際中不考慮使用.

表1 不同葉片數(shù)的模擬結(jié)果Tab.1 Simulated results for wind wheels with different numbers of blade

4 后處理與流場分析

4.1 后處理技術(shù)

圖3為風(fēng)輪截面的后處理示意圖.在圖3中,筆者在單葉片以及風(fēng)輪流場中分別截取r/R為46%截面的柱面,并利用 Tecplot后處理軟件將圓周曲面拉伸成長度與曲面周長相等的平面,然后將單葉片圓周平面按周期性邊界條件復(fù)制,得到不考慮葉片間損失的三葉片翼型后流場,再將相同半徑的單葉片和風(fēng)輪流場導(dǎo)入Tecplot后處理軟件進行比較和分析.這種后處理方法的優(yōu)點是:沿葉片的運動軌跡取截面研究,可清晰地看到葉片之間的相互影響.

圖3 風(fēng)輪截面的后處理示意圖Fig.3 Wind turbine section for post processing

4.2 結(jié)果與分析

4.2.1 不同數(shù)目葉片的風(fēng)輪速度比較

圖4為不同數(shù)目葉片風(fēng)輪的旋轉(zhuǎn)平面速度比較云圖.從圖4可看出:單葉片和雙葉片風(fēng)輪的葉片周圍速度分布基本相似,其葉片吸力面均存在同樣的低速區(qū)和增長區(qū),且壓力面附近隨轉(zhuǎn)動軌跡的速度變化不明顯.與其他3種風(fēng)輪葉片相比,三葉片風(fēng)輪葉片間的速度分布最理想,其兩葉片間的速度大小分布均勻,隨著速度由快變慢,速度帶逐漸變寬.四葉片風(fēng)輪的速度分布不均勻,靠近葉片部分的速度變化梯度較大,而兩葉片的中間部分速度變化梯度較小.

圖4 不同葉片數(shù)風(fēng)輪的旋轉(zhuǎn)平面速度比較云圖Fig.4 Velocity comparison of rotating flat among wind wheels with different numbers of blade

4.2.2 不同葉片風(fēng)輪的湍動能比較

圖5為r/R等于46%曲面處不同葉片數(shù)目風(fēng)輪旋轉(zhuǎn)截面的湍動能對比,其截面生成方式及截面范圍見圖3,即取風(fēng)輪上風(fēng)向5 m到下風(fēng)向15 m之間部分進行分析.將圖5中三維環(huán)面區(qū)域流場簡化為二維平面問題進行分析比較,可清楚觀察到不同情況下湍流強度分布規(guī)律及差別.

圖5 不同葉片數(shù)風(fēng)輪的旋轉(zhuǎn)截面湍動能對比Fig.5 Comparison of the turbulent kinetic energy on circular cut of wind wheels with different numbers of blade

單葉片流場域獨立,無其他葉片尾流影響.尾流湍動能在2 m2/s2以上,其長度為10 m,與旋轉(zhuǎn)平面約成25°,之后出現(xiàn)湍動能較穩(wěn)定的帶狀區(qū)域.在(15,5)處,開始出現(xiàn)變化梯度較大的低湍動能尾跡渦,并向下風(fēng)向擴散延伸.3種多葉片風(fēng)輪的翼型流場相似,但它們的湍流強度與影響區(qū)域的面積不同.雙葉片風(fēng)力發(fā)電機葉片周圍高湍區(qū)域的湍動能大于單葉片,其影響面積略小于單葉片,葉片避開了另一葉片的尾流區(qū),即各葉片間的影響較小,尾流遠場出現(xiàn)低湍動能的尾跡渦.三葉片風(fēng)輪葉片間的相互影響很大,葉片尾緣的湍動能帶出現(xiàn)分離,即附著渦平穩(wěn)脫落,壓力面脫落湍動能沿x軸切向延伸至旋轉(zhuǎn)下游的葉片吸力面,湍動能在2～3 m2/s2,吸力面脫落湍動能與旋轉(zhuǎn)上游葉片切向延伸的湍動能匯合后偏y軸延伸,此段湍動能變化明顯,由5 m2/s2降為3 m2/s2,影響范圍在y=0到y(tǒng)=10 m間.四葉片風(fēng)輪葉片均處于旋轉(zhuǎn)上游的葉片尾流陰影中,葉輪間影響作用最大,但由于間隔較小,風(fēng)輪上游及下風(fēng)向y＞14 m流場區(qū)域流動平穩(wěn),沒有出現(xiàn)明顯的旋渦.

圖6為風(fēng)輪旋轉(zhuǎn)平面湍動能的比較.由圖6可知:隨著葉片數(shù)的增加,葉片間湍動能逐漸增大.單葉片和雙葉片旋轉(zhuǎn)風(fēng)輪葉片附近的流場基本一致,即雙葉片風(fēng)輪葉片間影響可以忽略不計.三葉片風(fēng)輪葉片處于旋轉(zhuǎn)上游葉片脫落附著渦衰減后的平穩(wěn)部分,受影響較小.四葉片風(fēng)輪的流場較混亂,葉片周圍湍動能較大.

圖6 不同葉片數(shù)風(fēng)輪旋轉(zhuǎn)平面湍動能的比較Fig.6 Comparison of the turbulent kinetic energy on rotating flat of wind wheels with different numbers of blade

5 結(jié) 論

(1)由于旋轉(zhuǎn)上游的葉片發(fā)生渦脫落,導(dǎo)致旋轉(zhuǎn)下游葉片周圍流場的環(huán)境發(fā)生變化.隨著葉片數(shù)的增加,旋轉(zhuǎn)風(fēng)輪中單個葉片的輸出功率逐漸減小,風(fēng)力發(fā)電機效率與運行質(zhì)量都受到影響.

(2)多葉片旋轉(zhuǎn)風(fēng)輪背風(fēng)側(cè)流場的平均湍動能比單葉片和雙葉片風(fēng)輪大.脫落后的附著渦通過旋轉(zhuǎn)下游葉片時沒有形成類似單葉片后的較大旋渦,因而對風(fēng)輪下風(fēng)向的穩(wěn)定性有利.

[1]HANSEN M O L.Areodynamics of wind turbines[M].London,UK:Earthscan Science Publishers,2008:18-26.

[2]胡丹梅,田杰,杜朝輝.水平軸風(fēng)力機尾跡的測量與分析[J].動力工程,2006,26(5):751-760.HU Danmei,TIAN Jie,DU Zhaohui.Measurement and analysis of wake behind horizontally orientated air turbines[J].Journal of Power Engineering,2006,26(5):751-760.

[3]BART HELEMIE R J,RATHMANN O,FRANDSEN S T,et al.Modelling and measurements of wakes in large wind farms[J].Journal of Physics:Conference Series,2007,75(1):12-49.

[4]DOBREV I,MAALOUF B,TROLDBORG N,et al.Investigation of the wind turbine vortex structure[C]//14th Int Symp on Applications of Laser Techniques to Fluid Mechanics.Lisbon,Portugal:[s.n.],2008.

[5]郭婷婷,吳殿文,王成蔭,等.風(fēng)力發(fā)電機葉片預(yù)彎設(shè)計及其數(shù)值研究[J].動力工程學(xué)報,2010,30(6):450-455.GUO Tingting,WU Dianwen,WANG Chengyin,et al.Prebend design and numerical simulation of wind turbo-generator[J].Journal of Chinese Society of Power Engineering,2010,30(6):450-455.

[6]吳殿文.風(fēng)力機優(yōu)化設(shè)計及氣動性能的數(shù)值模擬[D].吉林:東北電力大學(xué)能源與動力工程學(xué)院,2010.

[7]包飛.風(fēng)力機葉片幾何設(shè)計與空氣動力學(xué)仿真[D].大連:大連理工大學(xué)機械工程學(xué)院,2008.

[8]S?RENSEN N N,MICHELSEN J A.Drag prediction for blades at high angle of attack using CFD[J].Journal of Solar Energy Engineering,2004,126(4):1011-1016.

[9]SEZER-UZOL N,LONG LN.3-D time-accurate CFD simulations of wind turbine rotor flow fields[C]//44th AIAA Aerospace Sciences Meeting and Exhibit.Reno,Nevada:American Institute of Aeronautics and Astronautics,2006:6-19.

[10]何顯富,盧霞,楊躍進,等.風(fēng)力機設(shè)計、制造與運行[M].北京:化學(xué)工業(yè)出版社,2009:189-196.

[11]GUO Tingting,WU Dianwen,XU Jihui.The method of large-scale wind turbine blades design based on Matlab programming[C]//2009 Supergen Conference.Nanjing,China:[s.n.],2009.