基于渦旋升力線模型的風(fēng)機后掠葉片性能優(yōu)化*

李德源， 劉佳宇， 夏鴻建， 張海波

(廣東工業(yè)大學(xué) 機電工程學(xué)院， 廣州 510006)

目前，大功率兆瓦級風(fēng)力機的空氣動力學(xué)和結(jié)構(gòu)特性優(yōu)化設(shè)計已經(jīng)成為了研究重點.而風(fēng)輪作為風(fēng)力機的主要組成部分，其成本占到整個風(fēng)力機的20%左右，同時風(fēng)輪也是風(fēng)力機捕獲風(fēng)能的關(guān)鍵部分，因此，國內(nèi)外學(xué)者從不同的方面，對于風(fēng)力機葉片進(jìn)行優(yōu)化設(shè)計研究，以獲取良好的氣動特性和結(jié)構(gòu)特性.研究者們提出在保持風(fēng)輪直徑不變或增加較少的前提下，將葉片適當(dāng)后掠的概念，也稱為幾何自適應(yīng)[1].研究表明：作為一種自適應(yīng)的氣動降載布局設(shè)計，后掠型葉片在相同結(jié)構(gòu)載荷下，可增加風(fēng)輪直徑，從而增大捕風(fēng)面積，增加發(fā)電量，降低度電成本；在保證相同發(fā)電量下，可有效減少結(jié)構(gòu)載荷，降低風(fēng)輪與相關(guān)零部件的制造成本，提高結(jié)構(gòu)可靠性，降低維護(hù)成本，因此，這種葉片具備較好的使用價值和應(yīng)用前景[2].

后掠葉片在氣動載荷作用下發(fā)生沿順槳方向的扭轉(zhuǎn)，雖然降低了葉片上的氣動力，但是也可能降低葉片的功率輸出，因此，在設(shè)計后掠自適應(yīng)葉片時，進(jìn)行一體化設(shè)計優(yōu)化是必需的[3]，這其中就包括各運行風(fēng)況下，后掠葉片槳距角的優(yōu)化設(shè)計.丁勇鋼等[4]通過對風(fēng)機控制系統(tǒng)的優(yōu)化，實現(xiàn)發(fā)電量的最大化和葉片軸向推力的最小化.

后掠葉片在結(jié)構(gòu)和空氣動力學(xué)方面表現(xiàn)出較強非線性，需要采用恰當(dāng)?shù)姆蔷€性氣彈耦合模型表達(dá)空氣動力與葉片彈性變形以及振動時的相互作用.為此，本文通過將柔性葉片離散為由運動副和力元聯(lián)接的多體系統(tǒng)[5-6]，應(yīng)用計算多體系統(tǒng)動力學(xué)理論建立葉片非線性氣彈耦合方程[7].在葉片氣動載荷分析方面，后掠葉片表現(xiàn)出顯著的三維流動，螺旋尾渦升力線模型將葉片用環(huán)量連續(xù)變化的升力線代替，具有較好的三維流場分析能力[8].

采用“超級單元”方法[9]將后掠葉片離散成若干個由運動副和力元聯(lián)接的剛體，應(yīng)用Roberson-Wittenburg遞推建模方法，結(jié)合螺旋尾渦升力線模型建立葉片非線性氣彈耦合方程.以變槳角為優(yōu)化設(shè)計變量，最大年發(fā)電量作為目標(biāo)函數(shù)，將風(fēng)力機額定功率、最大葉尖揮舞位移和最大軸向推力設(shè)為約束條件，對變槳角進(jìn)行尋優(yōu)并對結(jié)果進(jìn)行對比分析.

1 后掠葉片非線性氣彈耦合模型

后掠葉片具有與直葉片不同的氣彈特性，有必要建立恰當(dāng)?shù)暮舐尤~片氣彈耦合模型.本文在后掠葉片結(jié)構(gòu)非線性模型的基礎(chǔ)上，結(jié)合升力線理論建立葉片氣彈耦合模型.

1.1 后掠葉片結(jié)構(gòu)非線性模型

為了準(zhǔn)確描述葉片的結(jié)構(gòu)變形，采用“超級單元”建立葉片模型.其中，每個超級單元由4個剛體(B1～B4)組成，各剛體通過運動副、彈簧和阻尼器連接[5-6]，用來描述葉片的揮舞、擺振和扭轉(zhuǎn)變形，其模型如圖1所示.

圖1 超級單元模型Fig.1 Super-element model

應(yīng)用多體系統(tǒng)動力學(xué)中的R-W方法，取各鉸的轉(zhuǎn)動方向為廣義坐標(biāo)，則從柔性葉片的拓?fù)錁?gòu)形可得其廣義坐標(biāo)陣為

q(t)=(q1(t)，q2(t)，…，qs(t))

(1)

(2)

1.2 基于升力線理論的后掠葉片氣動載荷計算

由于后掠葉片具有顯著的展向流動和積疊線彎曲特點，采用升力線模型模擬后掠葉片流場的三維流動.

升力線模型的核心思想是將風(fēng)力機葉片三維流場中的渦量分布簡化為一根渦量沿葉片展向變化的線渦，表示葉片與流場之間的相互作用.該線渦代表了葉片三維流場中的環(huán)量分布，因而該線渦也稱作附著渦.由Helmholtz第二定理可知，在無粘流體環(huán)境中，渦量并不能在流體內(nèi)部終止，而只能延伸至流體邊界或構(gòu)成封閉區(qū)間，附著渦沿葉展方向變化的渦量會從附著渦上脫落形成尾渦，渦方法模型如圖2所示.

圖2 渦方法模型示意圖Fig.2 Schematic model diagram of vortex method

從附著渦上脫落的螺旋尾渦流向葉片下游無限遠(yuǎn)處，形成誘導(dǎo)速度場，本文采用剛性尾渦模型，尾渦坐標(biāo)系如圖3所示.

圖3 剛性尾渦及其坐標(biāo)系示意圖Fig.3 Schematic diagram of rigid wake vortex and its coordinate system

文獻(xiàn)[10]指出，誘導(dǎo)速度Wn垂直于未受擾動的相對風(fēng)速W′，并給出了附著渦任意位置e′(已知慣性坐標(biāo)(x′，y′，z′))處由所有尾渦引起的誘導(dǎo)速度，即

(3)

式中：Γ(ξ)為螺旋尾渦的環(huán)量；R為風(fēng)輪半徑；ξ、ξ′為葉展方向無量綱位置參數(shù)；ξhub為葉展方向輪轂位置處的無量綱參數(shù)；λ0=V0/(RΩcosΛ)；V0為來流風(fēng)速；Λ為當(dāng)?shù)睾舐咏?；Ω為風(fēng)輪轉(zhuǎn)速；ψ為風(fēng)輪錐角；θ為螺旋尾渦的方位角；θk為控制點所在的葉片與發(fā)出尾渦的葉片之間的夾角；N為葉片數(shù)；K為控制點所在葉片序號.

附著渦線為空間曲線，渦量從輪轂沿著葉片徑向流向葉尖，形成誘導(dǎo)速度場，附著渦坐標(biāo)系如圖4所示.

葉片展向位置處受到附著渦上某一直線渦量微段dn′的誘導(dǎo)速度為

(4)

圖4 附著渦線及其坐標(biāo)系示意圖Fig.4 Schematic diagram of attached vortex line and its coordinate system

式中：dWi2為附著渦誘導(dǎo)速度；dwx2、dwy2、dwz2分別為dWi2在x、y、z軸的分量；dn′為附著渦上的微段向量；s′為輪轂中心指向微段dn′的向量；dn′和s′與積疊線形狀有關(guān).附著渦誘導(dǎo)速度垂直于相對風(fēng)速，經(jīng)過葉片積分得到所有葉片附著渦對控制點的誘導(dǎo)速度為

(5)

式中：φ為入流角；Vh和Vb分別為葉片的揮舞、擺振速度.當(dāng)N=K時，即|s′-e′|=0，式(4)中分子和分母均為0，這不符合物理實際.實際上，由于粘性的作用會使得靠近渦線的中心速度趨于有限值，形成所謂的粘性渦核[8].基于粘性渦核模型建立修正的Biot-Savart公式，在計算誘導(dǎo)速度時引入粘性修正因子，即

(6)

式中，Re為渦核的有效半徑.修正的Biot-Savart公式為

dwy2sinφ)dξ

(7)

控制點上總的誘導(dǎo)速度為

Wn=Wn1+Wn2

(8)

誘導(dǎo)速度表達(dá)式(3)和(7)中環(huán)量項Γ(ξ)未知.假設(shè)附著渦的環(huán)量分布為Fourier級數(shù)，并且滿足邊界條件Γ(ξhub)=Γ(1)=0，設(shè)環(huán)量分布函數(shù)為

(9)

一般m取9即可[10].根據(jù)Kutta-Joukovski定理可知，環(huán)量與升力系數(shù)之間的關(guān)系為

Γ=0.5cWCL

(10)

式中：c為對應(yīng)翼型的弦長；CL為升力系數(shù)；W為截面位置受擾動的相對風(fēng)速.

聯(lián)立式(3)和(7)～(10)，代入m個控制點的已知數(shù)據(jù)，即可得到m個方程，方程組中只有Am未知，求解該方程組可得到Am的初值，利用迭代法[10]可以得到數(shù)值精度更高的有效攻角αe.迭代結(jié)束得到有效攻角αe和誘導(dǎo)攻角αi，通過二維翼型數(shù)據(jù)表線性插值得到升力系數(shù)CL和阻力系數(shù)CD后，沿著葉片徑向積分求得軸向推力F和轉(zhuǎn)動功率P，即

(11)

式中：ρ為空氣密度；φ′為實際入流角.

1.3 后掠葉片氣彈耦合模型驗證算例

為驗證后掠葉片非線性氣彈耦合模型的準(zhǔn)確性和可行性，選用NREL實驗室的5 MW風(fēng)力機直葉片數(shù)據(jù)[11]和后掠4 m葉片數(shù)據(jù)[12]進(jìn)行模型驗證.

采用本文建立的氣彈耦合模型分別分析了直葉片和后掠葉片風(fēng)輪在5～25 m/s風(fēng)速下的氣動性能與結(jié)構(gòu)響應(yīng)，并與文獻(xiàn)[12]采用葉素動量理論(BEM)氣動模型的分析結(jié)果進(jìn)行對比，結(jié)果如圖5所示.

從圖5可以看出，無論直、后掠葉片，應(yīng)用本文所建立的氣彈耦合模型能夠比較準(zhǔn)確地進(jìn)行氣動性能計算，結(jié)果與文獻(xiàn)[12]的結(jié)果基本一致.兩種分析模型結(jié)果隨風(fēng)速增大相差更大一些，這主要是因為對于大型柔性葉片，存在著沿葉片展向的流動，后掠葉片這種流動更為顯著一些，使用葉素動量理論計算氣動性能時，忽略了氣流在葉片展向的流動影響，之后可通過葉根和葉尖損失等經(jīng)驗修正模型進(jìn)行彌補.而升力線模型考慮了葉片的展向流動與尾流對整個葉片的影響，無需引入經(jīng)驗修正模型，考慮的影響因素和分析效率更高.在結(jié)構(gòu)上由于本文所建立的“超級單元”模型和文獻(xiàn)[12]鐵木辛柯梁的模型差異，帶來計算模擬結(jié)果的誤差.

圖5 升力線理論和葉素動量理論的穩(wěn)態(tài)計算結(jié)果Fig.5 Steady-state calculation results for lift line and blade element momentum theories

算例驗證結(jié)果表明，相較于葉素動量理論，基于升力線所建立的氣彈耦合模型具有更好的三維流場分析能力，驗證了該方法的有效性與準(zhǔn)確性.

2 優(yōu)化設(shè)計模型

2.1 后掠葉片外形設(shè)計與運行條件

文獻(xiàn)[12]后掠葉片積疊線形狀計算表達(dá)式為

(12)

此處分別取a=8，b=12，即葉尖后掠4 m，R0為葉片長度，R0=61.5 m，則根據(jù)式(12)積疊線形狀為

(13)

根據(jù)文獻(xiàn)[12]取后掠葉片的切入風(fēng)速為5 m/s、切出風(fēng)速為25 m/s.在較低風(fēng)速時葉片保持尖速比8.7運行，隨著風(fēng)速的增加，葉片轉(zhuǎn)速從6 r/min開始增大，當(dāng)轉(zhuǎn)速提高到12.1 r/min時不再保持額定尖速比，在額定風(fēng)速達(dá)到11.4 m/s后，風(fēng)輪轉(zhuǎn)速保持為12.1 r/min.

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

年發(fā)電量(AEP)與風(fēng)場風(fēng)速的威布爾分布有關(guān)[13]，理論年發(fā)電量的計算公式為

(14)

式中：H為年小時，此處取8 766 h；vi和vi+1為風(fēng)速；P(vi)和P(vi+1)為風(fēng)速vi和vi+1對應(yīng)的功率；B和u分別為風(fēng)速威布爾分布的尺度因子和形狀因子.本文中計算AEP時u取2，B取7.334 5.

2.3 約束條件

為了確保傳遞到發(fā)電機的軸功率不超過其最大容量，需要約束風(fēng)力機風(fēng)輪在各風(fēng)速下的轉(zhuǎn)動功率不能超過額定功率5.3 MW[11]，即

Pmax≤5.3 MW

(15)

為保證葉片不會因為揮舞方向變形過大引起結(jié)構(gòu)破壞，變槳角優(yōu)化后的風(fēng)機葉片最大葉尖揮舞位移小于NREL 5 MW直葉片的最大葉尖揮舞位移[12]，即

δmax≤5.6 m

(16)

風(fēng)力機風(fēng)輪的工作壽命、結(jié)構(gòu)強度均與最大軸向推力有關(guān).變槳角優(yōu)化后的風(fēng)機葉片最大軸向推力小于NREL 5 MW直葉片的最大軸向推力[12]，即

Fmax≤0.779 MN

(17)

2.4 優(yōu)化變量與優(yōu)化方法

采用一維搜索優(yōu)化方法，以變槳角為優(yōu)化變量，尋求各個風(fēng)速下，變槳角與功率的擬合函數(shù)關(guān)系，并求出擬合函數(shù)最大功率所對應(yīng)的變槳角，計算氣動性能和結(jié)構(gòu)響應(yīng).

3 優(yōu)化結(jié)果分析

按照前面所建立的優(yōu)化設(shè)計方法，對后掠量為4 m的風(fēng)機葉片在各個風(fēng)速下的變槳角進(jìn)行優(yōu)化，分析中采用升力線模型.圖6為優(yōu)化的后掠葉片與直葉片、后掠葉片[12]在5～25 m/s風(fēng)速范圍內(nèi)氣動性能和結(jié)構(gòu)響應(yīng)的結(jié)果對比(仿真環(huán)境為勻流風(fēng)，不考慮塔影效應(yīng)和重力影響，設(shè)計尖速比為8.7，轉(zhuǎn)速范圍為6～12.1 r/min，優(yōu)化變量為變槳角).

表1為直葉片、后掠葉片與變槳角優(yōu)化的后掠葉片在風(fēng)速11 m/s時最大葉尖揮舞位移、最大軸向推力.變槳角優(yōu)化的后掠葉片與原后掠葉片最大葉尖揮舞位移小于5.6 m，最大軸向推力小于0.779 MN，保證了葉片不會因為揮舞方向變形過大引起結(jié)構(gòu)破壞且不會對風(fēng)輪的工作壽命和結(jié)構(gòu)強度產(chǎn)生較大影響.

圖6a結(jié)果表明，在6～11 m/s風(fēng)速下，變槳角優(yōu)化的后掠葉片輸出功率顯著提升；在12～25 m/s風(fēng)速下，各葉片的輸出功率已經(jīng)達(dá)到額定功率5.3 MW，功率不再增加.通過式(14)計算變槳角優(yōu)化的后掠葉片、原后掠葉片和直葉片的年發(fā)電量，其值分別為14 393 942、13 257 036、13 378 251 kW·h，變槳角優(yōu)化的后掠葉片年發(fā)電量相較于其余葉片的年發(fā)電量有明顯提高，與原后掠葉片相比提高約8.58%、與直葉片相比提高約7.59%.

圖6b、c結(jié)果表明，在6～11 m/s風(fēng)速下，后掠葉片的軸向推力、葉尖揮舞位移有所增加；在12～25 m/s風(fēng)速下，軸向推力、葉尖揮舞位移顯著減少.通過優(yōu)化后掠葉片的變槳角，使葉片在高風(fēng)速下維持5.3 MW的同時，減少了對風(fēng)輪和塔架的載荷，提高結(jié)構(gòu)可靠性.

圖6d結(jié)果表明，原后掠葉片在5～11 m/s風(fēng)速下變槳角一直減小，這是由于葉片積疊線后掠，在相同條件下后掠葉片切向力要乘上后掠角余弦值，使功率(轉(zhuǎn)動力矩)比直葉片要??；要讓后掠葉片獲得與直葉片相近的功率(轉(zhuǎn)動力矩)則需要增大氣動載荷，所以后掠葉片變槳距角比直葉片小，在11 m/s風(fēng)速下取得最小值.在12～25 m/s風(fēng)速下變槳角一直增大，使得功率維持在5.3 MW附近，變槳角優(yōu)化的后掠葉片在7～11 m/s風(fēng)速下，變槳角進(jìn)一步減小，增大了有效攻角，提高了升力系數(shù)CL，從而提高了功率，增大軸向推力和葉尖揮舞位移.

圖6 氣動性能和結(jié)構(gòu)響應(yīng)計算結(jié)果Fig.6 Calculation results for aerodynamic performance and structural response

表1 11 m/s風(fēng)速下最大葉尖揮舞位移、最大軸向推力Tab.1 Maximum flapwise tip deflection and axial thrust at wind speed of 11 m/s

4 結(jié) 論

本文建立葉片非線性氣彈耦合方程，并驗證模型的有效性和準(zhǔn)確性.以5 MW后掠4 m葉片為原型，變槳角為優(yōu)化變量，年發(fā)電量作為目標(biāo)函數(shù)，將風(fēng)力機額定功率、最大軸向推力、最大葉尖揮舞位移設(shè)為約束條件，通過后掠葉片非線性氣彈耦合程序，計算葉片的氣動性能和結(jié)構(gòu)響應(yīng)，得到結(jié)論如下：

1) 利用基于升力線的葉片非線性氣彈耦合程序，計算后掠葉片的氣動性能和結(jié)構(gòu)響應(yīng)具有較高的效率和精度，為風(fēng)力機的優(yōu)化設(shè)計提供理論依據(jù).

2) 對變槳角進(jìn)行優(yōu)化，所得到的槳距角隨風(fēng)速變化曲線可以作為風(fēng)機運行時槳距角變化的實際控制運行規(guī)律曲線.

3) 通過該優(yōu)化模型，變槳角優(yōu)化的后掠葉片相較于原后掠葉片和直葉片，在年發(fā)電量上分別提高了8.58%、7.59%；該優(yōu)化模型也為風(fēng)力機葉片的優(yōu)化設(shè)計提供了參考.