王博,祁文軍,孫文磊,姜超
(1.新疆大學(xué)機(jī)械工程學(xué)院,新疆烏魯木齊830000;2.三峽新能源有限公司,新疆烏魯木齊830000)
風(fēng)力發(fā)電機(jī)葉片是整個(gè)風(fēng)力機(jī)中最為關(guān)鍵的部分,其氣動(dòng)性能和結(jié)構(gòu)性能的優(yōu)劣直接決定風(fēng)力機(jī)的工作效率和運(yùn)行壽命。因此,風(fēng)力機(jī)葉片的氣動(dòng)設(shè)計(jì)和氣動(dòng)性能計(jì)算是非常關(guān)鍵的工作。為了不斷提高設(shè)計(jì)和計(jì)算的準(zhǔn)確性,國(guó)內(nèi)外許多專(zhuān)家和學(xué)者在不斷研究和探索的基礎(chǔ)上提出了很多有效的理論和方法,主要包括葉素動(dòng)量理論、風(fēng)洞實(shí)驗(yàn)、外場(chǎng)實(shí)驗(yàn)及CFD數(shù)值模擬。其中CFD數(shù)值模擬方法是近年逐漸發(fā)展起來(lái)的理論和方法,在計(jì)算流體力學(xué)和計(jì)算機(jī)硬件的不斷發(fā)展推動(dòng)下,逐步發(fā)展并不斷走向成熟,為認(rèn)識(shí)風(fēng)力機(jī)的流動(dòng)性能提供了新的途徑。與前3種理論方法相比,有許多方面的優(yōu)勢(shì)和更好的發(fā)展前景。CFD數(shù)值模擬技術(shù)可以方便快捷地計(jì)算出流場(chǎng)流動(dòng)的各個(gè)定量結(jié)果,便于對(duì)各種流動(dòng)參數(shù)進(jìn)行分析,對(duì)設(shè)計(jì)方案進(jìn)行評(píng)估,對(duì)后續(xù)的實(shí)驗(yàn)驗(yàn)證提供指導(dǎo),避免了設(shè)計(jì)過(guò)程的盲目重復(fù)性,提高了風(fēng)力機(jī)的整體設(shè)計(jì)水平,大幅度縮短了研制周期,降低了研制成本,代替了大量的實(shí)驗(yàn),成為風(fēng)力機(jī)設(shè)計(jì)與分析的不可替代的工具[1-6]。
CFD是在計(jì)算機(jī)上實(shí)現(xiàn)一個(gè)特定的計(jì)算,就好像是在計(jì)算機(jī)上做一次全尺寸的風(fēng)洞實(shí)驗(yàn)。目前的應(yīng)用主要集中在風(fēng)電葉片翼型氣動(dòng)特性數(shù)值模擬、風(fēng)力機(jī)葉片氣動(dòng)性能模擬及氣動(dòng)彈性計(jì)算等方面[7]。然而,要將CFD直接應(yīng)用于工程計(jì)算,并獲得可信賴(lài)的數(shù)值模擬計(jì)算結(jié)果,還存在計(jì)算域尺度與邊界條件的合理選取、能夠滿足計(jì)算精度要求的網(wǎng)格質(zhì)量與網(wǎng)格數(shù)目的選擇以及合適的湍流模型選擇等問(wèn)題。
作者應(yīng)用UG軟件完成750 kW風(fēng)力機(jī)葉輪和用于數(shù)值模擬的風(fēng)力機(jī)流場(chǎng)的三維幾何建模,利用專(zhuān)用網(wǎng)格劃分軟件GAMBIT劃分網(wǎng)格,應(yīng)用計(jì)算流體力學(xué)CFD中已被廣泛應(yīng)用的Fluent軟件,對(duì)風(fēng)力機(jī)建模葉片在額定風(fēng)速和12個(gè)非額定風(fēng)速工況下的氣動(dòng)性能進(jìn)行數(shù)值模擬,分析多個(gè)風(fēng)速工況下的風(fēng)力機(jī)建模葉片的氣動(dòng)性能特點(diǎn)和差異,模擬仿真葉片氣動(dòng)流場(chǎng),計(jì)算葉輪的受力、扭轉(zhuǎn)力矩、輸出軸功率和風(fēng)能利用效率等性能參數(shù) ,觀察葉輪表面的壓強(qiáng)分布、流速分布、湍流強(qiáng)度、流速矢量等流態(tài)圖,驗(yàn)證風(fēng)力機(jī)氣動(dòng)性能數(shù)值模擬的可靠性,并以現(xiàn)場(chǎng)750 kW葉片實(shí)際測(cè)得的功率數(shù)據(jù)曲線作為依據(jù),驗(yàn)證葉片建模的合理性及計(jì)算機(jī)仿真葉片的氣動(dòng)性能是否滿足現(xiàn)場(chǎng)的要求。
文中研究對(duì)象屬于湍流模型,通過(guò)分析模型屬于重整化k-ε模型 (RNG k-ε模型)。
RNG k-ε模型的輸運(yùn)方程為:
其中:
式中:Gk表示平均速度梯度引起的湍動(dòng)能的產(chǎn)生; Gb表示浮力影響引起的湍動(dòng)能的產(chǎn)生;YM表示可壓湍流脈動(dòng)膨脹對(duì)總的耗散率的影響。
FLUENT中,默認(rèn)常數(shù)值為,C1ε=1.44,C2ε= 1.92,Cμ=0.09,湍動(dòng)能k和耗散率ε的湍流普朗特?cái)?shù)分別為σk=1.0,σε=1.3,αk和αε分別為湍動(dòng)能k和耗散率ε的有效湍流普朗特?cái)?shù)的倒數(shù)。
湍流黏性系數(shù)計(jì)算公式為:
在FLUENT中,用RNG k-ε模型時(shí),默認(rèn)設(shè)置為高雷諾數(shù)流動(dòng)問(wèn)題。如果是低雷諾數(shù)問(wèn)題,必須進(jìn)行相應(yīng)的設(shè)置。
RNG模型在ε方程中加了一個(gè)條件,又考慮到了湍流漩渦,還提供了一個(gè)考慮低雷諾數(shù)流動(dòng)黏性的解析公式,這些特點(diǎn)使得RNG k-ε模型比標(biāo)準(zhǔn)k-ε模型有更高的可信度和精度。
風(fēng)力發(fā)電機(jī)基本參數(shù)為:額定功率為750 kW,額定風(fēng)速為15.0 m/s,啟動(dòng)風(fēng)速為4.0 m/s,切出風(fēng)速為25.0 m/s,極大風(fēng)速為70.0 m/s,功率調(diào)節(jié)方式為失速,風(fēng)輪直徑為48.4 m,輪轂高度為50 m,葉片材料為玻璃纖維增強(qiáng)樹(shù)脂,葉片端線速度為57.4 m/s,風(fēng)力機(jī)尖速比為3.83,葉輪額定轉(zhuǎn)速為22.5 r/min,葉輪直徑為48.4 m,掃風(fēng)面積為1 840.0 m2,三葉片,選用NACA 4412翼型進(jìn)行葉片建模。
風(fēng)力機(jī)幾何建模包括:葉片建模,3葉片復(fù)制和旋轉(zhuǎn)體建模,輪轂建模,葉片周?chē)鲌?chǎng)建模,風(fēng)力機(jī)流場(chǎng)建模和整機(jī)流場(chǎng)建模,見(jiàn)圖1、圖2。
圖1 葉輪幾何建模
圖2 整機(jī)流場(chǎng)幾何建模
葉輪旋轉(zhuǎn)流場(chǎng)是包羅葉輪的扁圓盤(pán),旋轉(zhuǎn)中心指向進(jìn)風(fēng)口方向4.0 m,指向出風(fēng)口方向6.0 m,半徑30.0 m。整機(jī)流場(chǎng)是包羅整個(gè)葉輪旋轉(zhuǎn)流場(chǎng)的流場(chǎng),整個(gè)風(fēng)場(chǎng)由一個(gè)長(zhǎng)方體和一個(gè)半圓體構(gòu)成,長(zhǎng)方體的長(zhǎng)、寬、高分別為350、250、50 m,半圓體的半徑為125 m。旋轉(zhuǎn)中心到進(jìn)風(fēng)口 150 m,到出風(fēng)口200 m。
網(wǎng)格劃分采用4面體網(wǎng)格。其中:旋轉(zhuǎn)流場(chǎng)節(jié)點(diǎn)209 568個(gè),網(wǎng)格1 150 464個(gè);靜止流場(chǎng)節(jié)點(diǎn)54 464個(gè),網(wǎng)格282 124個(gè);總共節(jié)點(diǎn)264 032個(gè),網(wǎng)格1 432 588個(gè),見(jiàn)圖3、圖4。
圖3 葉輪旋轉(zhuǎn)流場(chǎng)網(wǎng)格劃分
圖4 整機(jī)流場(chǎng)網(wǎng)格劃分
此次三維模型的氣動(dòng)數(shù)值模擬計(jì)算采用分離式求解器 (Segregated),三維穩(wěn)態(tài)流動(dòng) (3D,Steady),速度是絕對(duì)速度 (Absolute),隱式求解 (Implicit),RNG k-ε模型,流體材料設(shè)為空氣。動(dòng)量采用二階迎風(fēng)離散格式 (Second Order Up-wind)。壓力和速度耦合采用SIMPLE算法。
速度進(jìn)口為15.0 m/s,同時(shí)計(jì)算變工況下風(fēng)力機(jī)氣動(dòng)數(shù)值模擬。風(fēng)速變化為:4.0,6.0,8.0,10.0,12.0,14.0,15.0,16.0,18.0,20.0,22,24.0,25 m/s。自由出口。進(jìn)出口之間壁面模擬為風(fēng)洞壁面,靜止無(wú)滑移。葉輪旋轉(zhuǎn)流體采用MRF模型,額定轉(zhuǎn)速22.55 r/min,方向右手定則確定。
對(duì)于葉片壁面和近壁面區(qū)域采用動(dòng)參考系(MRF)模型 (Multiple Reference Frame Model)。其中葉片壁面為旋轉(zhuǎn)壁面,相對(duì)旋轉(zhuǎn)流場(chǎng)旋轉(zhuǎn)速度為0,方向右手定則確定。對(duì)于近壁面區(qū)域,采用壁面函數(shù)法。在黏性流體流場(chǎng)中壁面處默認(rèn)為非滑移邊界條件。
在額定工況和變工況下對(duì)風(fēng)力機(jī)進(jìn)行氣動(dòng)數(shù)值模擬,計(jì)算收斂結(jié)束后,輸出各個(gè)工況下風(fēng)力機(jī)的繞軸旋轉(zhuǎn)扭矩及軸力,并對(duì)記錄的旋轉(zhuǎn)數(shù)據(jù)和葉片已知數(shù)據(jù)進(jìn)行整理計(jì)算,得到對(duì)應(yīng)的轉(zhuǎn)矩功率,將其與現(xiàn)場(chǎng)測(cè)試的750 kW風(fēng)力機(jī)實(shí)驗(yàn)功率曲線對(duì)比,二者的功率曲線對(duì)比如圖5所示。
圖5 750 kW風(fēng)力機(jī)計(jì)算機(jī)模擬與實(shí)測(cè)輸出功率曲線對(duì)比
從圖5中可以看出:
(1)利用FLUENT軟件對(duì)750 kW建模葉片在額定風(fēng)速和啟動(dòng)風(fēng)速到切出風(fēng)速的各個(gè)工況下的氣動(dòng)性能進(jìn)行模擬計(jì)算,因?yàn)槟M時(shí)不考慮風(fēng)力機(jī)變槳距以及風(fēng)力機(jī)內(nèi)部的功率損耗等影響,可以說(shuō)是在理想狀態(tài)下對(duì)風(fēng)力機(jī)進(jìn)行的數(shù)值模擬,風(fēng)力機(jī)輸出功率與風(fēng)速基本上成3次方曲線增長(zhǎng),這和理論上的風(fēng)力機(jī)輸出功率計(jì)算公式相符。通過(guò)對(duì)風(fēng)力機(jī)建模葉片進(jìn)行氣動(dòng)數(shù)值模擬計(jì)算,得出風(fēng)力機(jī)氣動(dòng)旋轉(zhuǎn)扭矩功率,在15 m/s額定風(fēng)速情況下,模擬得到的轉(zhuǎn)矩功率為768.4 kW,與額定功率750 kW誤差只有2.5%,基本達(dá)到預(yù)期的數(shù)值,說(shuō)明基于FLUENT的風(fēng)力發(fā)電機(jī)氣動(dòng)性能數(shù)值模擬計(jì)算和仿真是可靠的,也是可行的。
(2)在低于額定風(fēng)速15 m/s時(shí),風(fēng)力發(fā)電機(jī)的實(shí)際輸出功率和建模葉片氣動(dòng)性能數(shù)值模擬理想轉(zhuǎn)矩功率相差比較小,主要原因是在低于額定風(fēng)速時(shí)風(fēng)力機(jī)變槳距不對(duì)風(fēng)機(jī)進(jìn)行調(diào)節(jié),而風(fēng)速較小時(shí)對(duì)應(yīng)的功率損耗也較小,這就使風(fēng)力機(jī)實(shí)際工作狀態(tài)接近于數(shù)值模擬時(shí)的邊界條件,即理想狀態(tài)。
(3)在風(fēng)速大于15 m/s時(shí),實(shí)際輸出功率和模擬計(jì)算的理想輸出功率相差較遠(yuǎn),原因是風(fēng)速大于額定風(fēng)速時(shí),風(fēng)力機(jī)變槳距開(kāi)始工作,會(huì)通過(guò)槳距調(diào)節(jié),限制風(fēng)輪的轉(zhuǎn)速,以保護(hù)風(fēng)機(jī),延長(zhǎng)使用壽命,提高風(fēng)力機(jī)工作穩(wěn)定性,其結(jié)果必然降低風(fēng)力機(jī)的功率和風(fēng)能利用率;同時(shí)也因風(fēng)速加大,變速箱及傳動(dòng)軸間的摩擦等功率損耗會(huì)加大,這些都導(dǎo)致風(fēng)速大于額定風(fēng)速時(shí)實(shí)際輸出功率幾乎不增大,甚至減小。而數(shù)值模擬中并不會(huì)考慮這些因素,還按照理論上的功率輸出關(guān)系模擬計(jì)算其模擬輸出功率,按照3次方曲線方式增長(zhǎng)。
(4)結(jié)合風(fēng)力機(jī)的實(shí)際工作狀態(tài),不可能讓風(fēng)力機(jī)功率隨著風(fēng)速不斷增大,因此數(shù)值模擬計(jì)算的理想輸出功率在小于額定輸出功率時(shí)會(huì)更接近于實(shí)際,而風(fēng)力機(jī)絕大部分時(shí)間工作在小于或接近額定風(fēng)速的工作狀態(tài)下,這樣,對(duì)于風(fēng)力機(jī)氣動(dòng)性能的研究在小于額定風(fēng)速時(shí)會(huì)更有實(shí)際意義。而從圖3中可以看到,小于額定風(fēng)速時(shí)模擬結(jié)果與實(shí)際結(jié)果非常接近,從而證明了建模葉片氣動(dòng)性能較好地滿足了實(shí)際的需要,也證明了利用FLUENT軟件研究風(fēng)力機(jī)的氣動(dòng)性能振動(dòng)特性是完全可行的。
圖6—17分別為葉輪、流場(chǎng)出入口、葉輪周?chē)D(zhuǎn)流場(chǎng)這3個(gè)區(qū)域的壓強(qiáng)與速度的分布圖與跡線圖,能直觀地、感性地認(rèn)識(shí)和觀察葉輪及周?chē)鲌?chǎng)的氣動(dòng)特性和流動(dòng)現(xiàn)象,以及壓力、流速等重要?dú)鈩?dòng)參數(shù)的分布規(guī)律,并能很好地觀察葉輪以及風(fēng)力機(jī)周?chē)牧鲌?chǎng)分布。
從圖6、圖7可以看到:葉輪正面為迎風(fēng)面,在遇到葉輪阻礙之前,沒(méi)有發(fā)生流場(chǎng)的較大變化,因而氣動(dòng)壓強(qiáng)分布較為均勻;而葉輪背面卻因?yàn)闅饬魍ㄟ^(guò)旋轉(zhuǎn)葉輪時(shí),紊流現(xiàn)象急劇加強(qiáng),導(dǎo)致葉片表面流場(chǎng)發(fā)生了較大變化,故氣動(dòng)壓強(qiáng)分布表現(xiàn)出明顯的不均勻性。
圖6 葉輪正面氣動(dòng)壓強(qiáng)分布圖
圖7 葉輪背面氣動(dòng)壓強(qiáng)分布圖
從圖8、圖9可以看到:葉輪正面和背面的速度分布幾乎沒(méi)有變化,但是由于風(fēng)輪的轉(zhuǎn)動(dòng)帶動(dòng)了周?chē)鷼饬鞯男D(zhuǎn),從而在接近葉輪周?chē)膱A柱形區(qū)域,空氣流場(chǎng)幾乎和葉輪同步旋轉(zhuǎn)。因而在選擇風(fēng)力機(jī)氣動(dòng)數(shù)值模擬邊界條件時(shí),設(shè)定這一區(qū)域是和風(fēng)輪同步旋轉(zhuǎn)的動(dòng)流場(chǎng),計(jì)算結(jié)果表明是正確合理的。
圖8 葉輪正面氣動(dòng)速度分布圖
圖9 葉輪背面氣動(dòng)速度分布圖
從圖10、圖11可以看到:在風(fēng)力機(jī)葉輪周?chē)諝馐茏?,氣?dòng)壓強(qiáng)大,遠(yuǎn)離風(fēng)力機(jī)區(qū)域空氣受阻力小,氣動(dòng)壓強(qiáng)小,符合實(shí)際流場(chǎng)的空氣流動(dòng)特性。
從圖12、圖13可以看到:進(jìn)出口的空氣速度分布有較大的不同,這主要是因?yàn)榭諝鈿饬髟谕ㄟ^(guò)葉片風(fēng)輪前紊流較弱,而通過(guò)風(fēng)力機(jī)葉輪后,空氣氣流發(fā)生了明顯的紊流現(xiàn)象,流場(chǎng)中速度分布也隨之發(fā)生明顯變化。
圖10 流場(chǎng)入口氣動(dòng)壓強(qiáng)分布圖
圖11 流場(chǎng)出口氣動(dòng)壓強(qiáng)分布圖
圖12 流場(chǎng)入口氣動(dòng)速度分布圖
圖13 流場(chǎng)出口的氣動(dòng)速度分布圖
從圖14、圖15可以看到:風(fēng)力機(jī)葉輪及周?chē)D(zhuǎn)流場(chǎng)中壓強(qiáng)和速度的分布規(guī)律,葉輪處于旋轉(zhuǎn)中,故靠近葉輪旋轉(zhuǎn)區(qū)域速度大,遠(yuǎn)離葉輪旋轉(zhuǎn)區(qū)域的速度減小;同樣,由于空氣氣流速度的影響,葉輪及周?chē)鷧^(qū)域的氣流壓強(qiáng)也隨之發(fā)生變化,呈現(xiàn)和速度相呼應(yīng)的分布特點(diǎn)。
圖14 葉輪周?chē)D(zhuǎn)流場(chǎng)氣動(dòng)壓強(qiáng)圖
圖15 葉輪周?chē)D(zhuǎn)流場(chǎng)速度分布圖
從圖16、圖17中可以看到:通過(guò)FLUENT軟件 建立風(fēng)力機(jī)葉輪旋轉(zhuǎn)流場(chǎng),并進(jìn)行數(shù)值模擬計(jì)算,得到的葉輪表面的壓強(qiáng)和速度跡線圖非常直觀、形象地描述了葉輪表面的壓強(qiáng)和速度流場(chǎng)的動(dòng)態(tài)變化和分布,總的來(lái)看,按照旋轉(zhuǎn)方向呈流線型分布。
圖16 葉輪旋轉(zhuǎn)流場(chǎng)氣動(dòng)壓強(qiáng)跡線圖
圖17 葉輪旋轉(zhuǎn)流場(chǎng)速度跡線圖
(1)以NACA63XXX翼型葉片為基礎(chǔ),通過(guò)UG軟件完成750 kW風(fēng)力機(jī)葉輪和用于數(shù)值模擬的風(fēng)力機(jī)流場(chǎng)的三維幾何建模,利用專(zhuān)用網(wǎng)格劃分軟件GAMBIT劃分網(wǎng)格,最后在FLUENT軟件中針對(duì)額定風(fēng)速和12個(gè)非額定風(fēng)速的工況進(jìn)行氣動(dòng)數(shù)值模擬計(jì)算,分析不同風(fēng)速工況下的風(fēng)力機(jī)氣動(dòng)性能的特點(diǎn)和差異。通過(guò)風(fēng)力機(jī)氣動(dòng)數(shù)值模擬計(jì)算發(fā)現(xiàn):風(fēng)力機(jī)轉(zhuǎn)矩輸出功率與風(fēng)速成3次方曲線增長(zhǎng),符合風(fēng)力機(jī)輸出功率公式;得出的風(fēng)力機(jī)氣動(dòng)旋轉(zhuǎn)扭矩功率和750 kW風(fēng)力機(jī)額定功率基本一致,得到的功率曲線在低于額定功率時(shí)也與現(xiàn)場(chǎng)實(shí)測(cè)750 kW葉輪得到的功率曲線相符,說(shuō)明建模葉片的氣動(dòng)性能基本符合750 kW機(jī)組的技術(shù)要求,基于FLUENT的風(fēng)力發(fā)電機(jī)氣動(dòng)性能數(shù)值模擬計(jì)算和仿真是可靠的、可行的。
(2)觀察葉輪表面以及整機(jī)周?chē)鲌?chǎng)的壓強(qiáng)分布、流速分布、湍流強(qiáng)度、速度矢量等動(dòng)態(tài)云圖,可以直觀和感性地認(rèn)識(shí)和觀察葉輪和整機(jī)周?chē)鲌?chǎng)的氣動(dòng)特性和流動(dòng)現(xiàn)象,以及其壓力、流速等氣動(dòng)參數(shù)的分布特點(diǎn),并且能很好地觀察葉輪以及整機(jī)周?chē)牧鲌?chǎng)分布。
(3)將模擬數(shù)據(jù)與實(shí)驗(yàn)數(shù)據(jù)結(jié)合起來(lái)進(jìn)行分析,有利于研究風(fēng)力機(jī)葉輪及整機(jī)氣動(dòng)性能和運(yùn)行情況,數(shù)值模擬仿真研究在一定程度上可以替代實(shí)際風(fēng)力機(jī)的運(yùn)行實(shí)驗(yàn),在成本、效率方面有不可替代的優(yōu)勢(shì),也豐富了對(duì)風(fēng)力機(jī)進(jìn)行氣動(dòng)分析的手段和方法,為風(fēng)力機(jī)葉片及整機(jī)的設(shè)計(jì)、改型和研發(fā)工作提供技術(shù)參數(shù)和指導(dǎo)意見(jiàn)。
【1】SIMMS D,SCHREEK S,HAND M,et al.NREL Unsteady Aerodynamics Experiment in the NASA-Ames Wind Tunnel:A Comparison of Predictions to Measurements[R].NREL/TP-500-29494,NREL,2001.
【2】BENJANIRAT S,SANKAR L N,XU G.Evaluation of turbulence models for the prediction of wind turbine aerodynamics[C]//41st AIAA Aerospace Sciences Meeting and Exhibit,Reno,American,2003.
【3】DUQUE E P N,BURKLUND M D,JOHNSON W.Navier-Stokes and comprehensive analysis performance predictions of the NREL Phase VI experiment[C]//41st ALAA Aerospace Sciences Meeting and Exhibit,Reno,America,2003.
【4】SEZER-UZOL Nilay,LONG Lyle N.3-D Time-Accurate CFD Simulations of Wind Turbine Rotor Flow Fields[C]//44th AIAA Aerospace Sciences Meeting and Exhibit,Reno,America,2006.
【5】ZAHLE Frederik,JPHANSEN Jeppe,SORENSEN Niels N.Wind TurbineAerodynamicsUsingan Incompressible Overset Grid Method[C]//Proceedings of the European Wind Energy Conference&Exhibition,Milan,Italy,2007.
【6】van ROOIJ R P J O M,ARENS E A.Analysis of the Experimental and Computational Flow characteristics with Respect to the Augmented Lift Phenomenon Caused by Blade Rotation[C]//Journal of Physics:Conference Series,Lyngby,Denmark,2007.
【7】劉丹,李軍向,薛忠民,等.風(fēng)力發(fā)電機(jī)葉片數(shù)值模擬綜述[J].電氣技術(shù),2010(7):7-13.