張立軍, 米玉霞, 趙昕輝, 馬東辰, 張松, 王旱祥

(1.中國石油大學(xué)(華東)機(jī)電工程學(xué)院 山東青島 266580; 2.中國石油大學(xué)(華東)化學(xué)工程學(xué)院 山東青島 266580)

現(xiàn)今清潔能源是世界能源的發(fā)展趨勢,而風(fēng)能作為一種清潔可再生能源,其在地球上的可利用的資源比水資源更豐富,共約200億kW,發(fā)電量可達(dá)到13 PW·h,且具有地域限制小等優(yōu)點(diǎn)[1]。隨著世界能源需求不斷增長,大量的陸地和海上風(fēng)電場不斷開發(fā),風(fēng)力發(fā)電機(jī)的功率也不斷增大,但風(fēng)力機(jī)受外在因素的影響也在增大,其最直觀的結(jié)果是使風(fēng)力機(jī)運(yùn)行不平穩(wěn)[2]。風(fēng)力發(fā)電機(jī)主要是靠葉片進(jìn)行動(dòng)力轉(zhuǎn)化,可從風(fēng)力機(jī)葉片研究風(fēng)力機(jī)不平穩(wěn)運(yùn)行問題[3]。Thomsen等[4]針對某一特定風(fēng)力機(jī)通過對葉片進(jìn)行加載直至其破壞,并記錄整個(gè)過程中位移的變化情況研究葉片在揮舞方向的抗載能力;Sutherland[5]對玻璃纖維的風(fēng)力機(jī)葉片的疲勞性能進(jìn)行了分析;Veldkamp等[6]對風(fēng)力發(fā)電機(jī)的各個(gè)部件進(jìn)行了使用壽命估計(jì)和疲勞分析。趙杰等[7]驗(yàn)證了用多個(gè)小型風(fēng)力機(jī)取代單個(gè)大型風(fēng)力機(jī)可有效減小風(fēng)剪切的影響;孔屹剛等[8]研究發(fā)現(xiàn)隨著風(fēng)力機(jī)容量增大,風(fēng)力機(jī)的整體尺寸相應(yīng)變大,風(fēng)剪切問題對風(fēng)力機(jī)載荷和功率的影響也不斷增大。雖然現(xiàn)今大型風(fēng)力機(jī)多為水平軸風(fēng)力機(jī),但垂直軸風(fēng)力機(jī)在建造大機(jī)型上并不遜于水平軸風(fēng)力機(jī),甚至更有發(fā)展空間[9]。畢長飛[10]提出通過翼型后部加厚的方法來增加風(fēng)輪葉片的強(qiáng)度,提高其工作穩(wěn)定性;江成生[11]基于有限元軟件ANSYS和Miner準(zhǔn)則提出了兩種預(yù)測葉片疲勞壽命的方法;張婷婷[12]針對某特定風(fēng)力機(jī),對其葉片采用風(fēng)力機(jī)載荷譜和線性累積損傷法則的葉片疲勞壽命估計(jì)算法進(jìn)行了疲勞分析。而對于風(fēng)剪切效應(yīng)的研究,尹偉[13]分析了風(fēng)剪效應(yīng)對風(fēng)力機(jī)功率系數(shù)的影響;張衡[14]等驗(yàn)證了風(fēng)剪效應(yīng)對大型風(fēng)力機(jī)葉片載荷的影響,但在計(jì)算載荷時(shí),未考慮來流風(fēng)速不等于作用在葉片上的誘導(dǎo)風(fēng)速,上下風(fēng)區(qū)的誘導(dǎo)速度不一致的問題;且計(jì)算傾角時(shí),僅考慮了葉片兩端所受的切向力,忽略了法向力載荷的問題。筆者研究風(fēng)剪切問題對葉片載荷的影響,將對風(fēng)剪切效應(yīng)下垂直軸風(fēng)力機(jī)運(yùn)行一周時(shí)葉片載荷進(jìn)行分析,研究傾角優(yōu)化方法以改善大型垂直軸風(fēng)力機(jī)的葉片受力情況。

1 風(fēng)剪效應(yīng)對葉片載荷影響

1.1 風(fēng)剪效應(yīng)

風(fēng)剪效應(yīng)是指隨著高度增加,風(fēng)速也不斷發(fā)生變化的現(xiàn)象。風(fēng)剪效應(yīng)的主要起因有兩個(gè):動(dòng)力因素和熱力因素,第一個(gè)是由于地面的摩擦效應(yīng),第二個(gè)主要與近地層大氣垂直穩(wěn)定度有關(guān)[15]。

風(fēng)速計(jì)算方法主要有對數(shù)律公式和指數(shù)律公式兩種。實(shí)測表明,這兩種計(jì)算方法均能較好地反映風(fēng)剪切效應(yīng),但根據(jù)對數(shù)律計(jì)算的風(fēng)速值與實(shí)測值的偏差比用指數(shù)律公式計(jì)算的風(fēng)速值與實(shí)測值的偏差比大[9]。對風(fēng)速沿高度的分布用指數(shù)律公式表示為

(1)

式中,vs為離地高度Z處的平均風(fēng)速,m/s;v0為離地參考高度Z0(一般取10 m)處的平均風(fēng)速, m/s;?為風(fēng)速廓線指數(shù)。結(jié)合某地的氣候條件,選取Z0=10 m[13]。

陸地普遍適用的風(fēng)速廓線指數(shù)[16]為?=0.20,得到風(fēng)速隨高度的變化曲線,如圖1所示。由圖1可知,在高度大于50 m時(shí),風(fēng)速增長趨勢減緩,在高度10～50 m范圍內(nèi)風(fēng)速增長趨勢很明顯。

圖1 不同高度的風(fēng)速分布Fig.1 Wind speed distribution at different heights

1.2 模型建立

對葉片進(jìn)行載荷分析,采用雙致動(dòng)盤多流管模型[17],如圖2所示?；趧?dòng)量方程,假設(shè)多個(gè)相同的流管通過風(fēng)機(jī)轉(zhuǎn)子,且風(fēng)通過每個(gè)流管作用到葉片翼型上的流向力是相等的,根據(jù)上下風(fēng)區(qū)將風(fēng)力機(jī)的轉(zhuǎn)子分為兩個(gè)致動(dòng)盤,即第一個(gè)致動(dòng)盤表示上風(fēng)區(qū)中的半個(gè)轉(zhuǎn)子掃掠面,第二個(gè)致動(dòng)盤表示下風(fēng)區(qū)中的半個(gè)轉(zhuǎn)子掃掠面。采用的流管數(shù)為8,流管在垂直方向上的高度為Δh。

對于垂直軸風(fēng)輪,通過矢量圖解法得到一個(gè)葉片分別在上風(fēng)區(qū)和下風(fēng)區(qū)的受力情況,如圖3所示。圖3中,vu和va分別代表垂直軸風(fēng)輪上、下風(fēng)區(qū)的誘導(dǎo)速度,兩者的值均小于來流風(fēng)速v∞;ω為風(fēng)輪旋轉(zhuǎn)角速度,ωR為切向風(fēng)速;W為誘導(dǎo)速度和切向速度的合成風(fēng)速;α為葉片攻角,是合成風(fēng)速方向與葉片弦長方向的夾角;θ為葉片方位角,當(dāng)方位角θ位于0°～180°時(shí),該區(qū)域稱為風(fēng)輪的上風(fēng)區(qū);當(dāng)θ位于180°～360°時(shí),稱為風(fēng)輪的下風(fēng)區(qū)。

以葉片在上風(fēng)區(qū)為例,分析圖3中參數(shù)間的幾何關(guān)系,可以得出垂直軸風(fēng)輪的葉片攻角α和合成風(fēng)速W的表達(dá)式為

(2)

(3)

葉片所受的升力L和阻力D分別為

(4)

式中,ρ為空氣密度,在標(biāo)準(zhǔn)大氣壓下取1.247 kg/m3;c為葉片弦長,m;Δh為所研究翼型的某段展向長度,m;CL和CD分別為翼型的升力系數(shù)和阻力系數(shù),其值可由NACA系列翼型空氣動(dòng)力學(xué)特性數(shù)據(jù)庫查得。

圖2 雙致動(dòng)盤多流管模型Fig.2 Double-disk multiple stream tube model

利用圖3所示的幾何關(guān)系,得到葉片切向運(yùn)動(dòng)的合力T和垂直葉片切向運(yùn)動(dòng)的合力N分別為

T=Lsinα-Dcosα,

(5)

N=Lcosα+Dsinα.

(6)

圖3 葉片受力分析Fig.3 Blade stress analysis

1.3 載荷分析

研究對象為1.5 MWH型垂直軸風(fēng)力發(fā)電機(jī)[13],采用葉片翼型為NACA0015對稱翼型,葉片弦長c=3 m,葉片數(shù)N=3,風(fēng)輪半徑R=34 m,風(fēng)輪高度H=30 m。

依照相同功率的水平軸風(fēng)力機(jī)的相關(guān)研究,假定葉片底部距離地面為15 m,即對于垂直軸風(fēng)力機(jī)而言,葉片風(fēng)速、載荷等參數(shù)的計(jì)算從h=15 m開始。為了計(jì)算葉片在風(fēng)剪切效應(yīng)作用下受到的載荷,首先計(jì)算Δh=0.1 m葉片截段在高度為h=45 m(葉片頂端離地高度)和高度h=15 m(葉片底端離地高度)所受載荷,以分析葉片由于風(fēng)剪效應(yīng)引起的在高度方向上載荷不均勻現(xiàn)象。葉片截段在一周范圍內(nèi)所受切向力T和法向力N對比如圖4所示,豎線上端和下端分別對應(yīng)葉片截段的較大載荷和較小載荷,豎線長度代表載荷差值。

圖4 不同高度葉片所受切向力和法向力Fig.4 Blade tangential force and normal force at different heights

從圖4可知,每單位長度的葉片在高度為45和15 m處所受切向力相差較小,最大差值僅為20 N;而每單位長度葉片受到的法向力相差值較大,最大差值高達(dá)60 N。由圖4還發(fā)現(xiàn),考慮風(fēng)剪效應(yīng)下,葉片載荷不均勻現(xiàn)象主要體現(xiàn)在所受法向力載荷。

2 方 案

為改善葉片受力不均勻情況,提出了3種垂直軸風(fēng)輪類型(圖5),葉片繞運(yùn)動(dòng)切向方向旋轉(zhuǎn),定義葉片與豎直方向的夾角為傾角β,以逆時(shí)針為正,分別對這3種風(fēng)輪所受的法向力進(jìn)行對比。以一周范圍內(nèi)葉片上下兩端受到的法向力的差值之和為比較依據(jù),研究葉片受力情況隨傾角的變化。

圖5 不同風(fēng)輪類型Fig.5 Different types of wind wheels

對1.5 MW的風(fēng)力機(jī)進(jìn)行分析,不同傾角下葉片的受力情況如圖6所示?？梢园l(fā)現(xiàn),相對于H型風(fēng)力機(jī)β=0°,傾角β>0°時(shí),傾角越大,葉片上下、兩端受力差值之和越大,即葉片受力越不均勻;傾角β<0°時(shí),傾角越小,葉片上下、兩端受力差值之和先減小,后增大,即葉片受力不均勻情況先得到改善,后又加劇。

圖6 隨傾角變化葉片載荷不均勻情況Fig.6 Non-uniform load of with dip angle change

由圖6可以判定,針對這一1.5 MW垂直軸風(fēng)力機(jī),對應(yīng)的最佳葉片傾角β應(yīng)小于零,即對應(yīng)的風(fēng)輪形狀如圖5(β<0°)所示。由文獻(xiàn)[18]可知,葉片傾角小于零的對應(yīng)形狀的風(fēng)輪還可有效地改善啟動(dòng)特性。對該風(fēng)力機(jī)整個(gè)葉片進(jìn)行傾角優(yōu)化,尋找最優(yōu)葉片傾角,以進(jìn)一步改善其葉片受力不均勻情況。

3 傾角優(yōu)化

3.1 優(yōu)化思路

計(jì)算葉片傾角的主要思路如圖7所示。以高度h=45 m開始,假定半徑R,葉片長度Δh=0.1 m,求得葉片截段受到的載荷;然后在高度為(h-Δh)位置,半徑r在(R-Range,R+Range)變化,計(jì)算過程中Range=1,在不同半徑、高度組合下求出一系列葉片截段的受力N,與上一高度葉片截段所受載荷N0求差值之和的平方,在這一系列中最小差值之和所對應(yīng)的半徑r,即為高度上的最佳半徑。再以此為基準(zhǔn),進(jìn)行下一高度的最佳半徑計(jì)算,直至高度達(dá)到葉片最低端(h=15 m)。對這一系列最佳半徑與高度的組合進(jìn)行數(shù)據(jù)擬合,求出斜率,即葉片傾角。

圖7 葉片傾角優(yōu)化流程Fig.7 Flow chart of blade angle optimization

保證傾斜后的風(fēng)力機(jī)與H型風(fēng)力機(jī)的風(fēng)功率相等,通過不斷改變初始半徑R(高度h=45 m),計(jì)算變傾角后風(fēng)力機(jī)的功率P,與H型風(fēng)力機(jī)功率P0進(jìn)行比較,直至相等為止,得到葉片的最優(yōu)傾角。具體計(jì)算過程如下。

目標(biāo)函數(shù):

ΔP=|P(β,r)=-P0(β0,r0)|=0.

(7)

約束條件:15≤h≤45，min(N-N0)2.

目標(biāo)函數(shù)中,β0=0°,r0=34 m。給定一個(gè)初始半徑r的變化范圍,如(R-Range,R+Range),變化梯度為Δr,計(jì)算每變化一個(gè)初始半徑,對應(yīng)傾斜后的風(fēng)力機(jī)與H型風(fēng)力機(jī)的風(fēng)功率之差,尋找min(P-P0)對應(yīng)初始半徑r1。啟動(dòng)程序,輸入初始半徑r,r在(r1-Δr,r1+Δr)內(nèi)變化,變化梯度為Δr1(Δr1<Δr,縮小半徑變化值,提高精度),再計(jì)算每變化一個(gè)初始半徑,對應(yīng)傾斜后的風(fēng)力機(jī)與H型風(fēng)力機(jī)的風(fēng)功率之差,尋找min(P-P0)對應(yīng)初始半徑r2,提高精度,直至ΔP=0,得到最優(yōu)葉片傾角β以及最佳初始半徑R。

3.2 計(jì)算結(jié)果

對該1.5 MW垂直軸風(fēng)力機(jī)進(jìn)行傾角優(yōu)化,不斷調(diào)整輸入的初始半徑,經(jīng)過多次迭代,得到計(jì)算結(jié)果如表1所示。計(jì)算結(jié)果中,當(dāng)初始半徑精度達(dá)到0.000 1時(shí),實(shí)際風(fēng)功率相差值ΔP為0.485 8 W,已近似為零,不必再提高精度。為安裝方便,將表1中的最優(yōu)葉片傾角計(jì)算結(jié)果圓整為整數(shù),編寫MATLAB程序,使風(fēng)力機(jī)傾斜前后的風(fēng)功率差值等于零,代入最佳傾角,計(jì)算出初始半徑,從而確定傾角優(yōu)化后的葉片安裝位置。表1中的計(jì)算結(jié)果是在離地10 m處風(fēng)速為8 m/s條件下獲得,其他風(fēng)速環(huán)境下的優(yōu)化計(jì)算方法與此類似。

表1 計(jì)算結(jié)果Table 1 Calculation results

為更清楚地說明葉片的受力情況,引入了方差σ,表示為

(8)

其中

不同方位的不同高度上Δh葉片截段受到載荷的方差可以表征數(shù)據(jù)的離散程度。不同方位對應(yīng)的方差越小,說明所受載荷越集中,即載荷相差不大,葉片受力更加均勻;對應(yīng)的方差越大,說明所受載荷越離散,即葉片在高度方向上受力更加不均勻。驗(yàn)證表1中的計(jì)算結(jié)果:當(dāng)傾角為-7°時(shí),葉片在不同方位上受力的對應(yīng)方差數(shù)值在改善前后的對比如圖8所示。

圖8 不同方位載荷的方差對比Fig.8 Variance values for different bearing loads

由圖8可知,風(fēng)力機(jī)經(jīng)傾角優(yōu)化后,其大部分方位的方差均有所下降,且全部方差均小于1;且經(jīng)過計(jì)算,一周范圍內(nèi)經(jīng)過傾角優(yōu)化后,方差平均值由0.66降到0.52。這說明經(jīng)傾角優(yōu)化后,葉片受力不均勻的情況得到了明顯改善。傾角優(yōu)化后,風(fēng)力機(jī)的形狀如圖5所示(β=-7°)。由于傾角變化值為-7°,且保證了風(fēng)功率一致,所需零部件尺寸變化較小,對風(fēng)力機(jī)的經(jīng)濟(jì)性影響不大。

4 結(jié) 論

(1)風(fēng)剪效應(yīng)對葉片載荷影響較大,且影響主要集中在法向力載荷上,所研究的1.5 MW H型垂直軸風(fēng)力機(jī)的單個(gè)葉片的法向力最大差值可達(dá)60 N。

(2)相對于H型垂直軸風(fēng)力機(jī),當(dāng)葉片傾角大于零時(shí),傾角越大,葉片受力越不均勻;當(dāng)傾角小于零時(shí),傾角越小,葉片受力不均勻情況先得到改善,后又加劇。

(3) 1.5 MW垂直軸風(fēng)力機(jī)在考慮風(fēng)剪效應(yīng)下所受載荷均勻時(shí),對應(yīng)的最優(yōu)葉片傾角為-7°,初始半徑為R=32.750 7 m。

(4)風(fēng)力機(jī)經(jīng)過傾角優(yōu)化后,葉片在一周范圍內(nèi)其方差最大為0.96,而原始風(fēng)力機(jī)最大值為1.53,且方差平均值由0.66降到0.52,葉片的受力情況明顯改善。