沖壓空氣渦輪葉片設(shè)計和氣動性能數(shù)值模擬

姬芬竹， 張夢杰， 王瑞， 王巖， 杜發(fā)榮

(1. 北京航空航天大學(xué) 交通科學(xué)與工程學(xué)院 北京市清潔能源與高效動力工程中心， 北京 100083；2. 北京航空航天大學(xué) 能源與動力工程學(xué)院， 北京 100083)

沖壓空氣渦輪(Ram Air Turbine, RAT)系統(tǒng)是現(xiàn)代飛機(jī)的應(yīng)急動力裝置，其能夠在飛機(jī)失去主動力和輔助動力的緊急情況下提供應(yīng)急能源，保證飛行安全[1]。通常，RAT系統(tǒng)由渦輪、傳動系統(tǒng)、電氣系統(tǒng)和控制系統(tǒng)等組成。其中，渦輪把氣流的沖壓能轉(zhuǎn)變?yōu)閯幽埽缓笸ㄟ^電氣系統(tǒng)(發(fā)電機(jī))將動能轉(zhuǎn)換為電能。可見，渦輪是提取氣流能量的核心部件。渦輪輸出功率和葉片風(fēng)能利用系數(shù)是渦輪氣動性能研究的關(guān)鍵，也是RAT研究的核心[2-3]。渦輪氣動性能研究主要包括理論研究、數(shù)值模擬和實驗驗證等。理論研究發(fā)展較早，技術(shù)相對成熟，著名的BETZ理論、動量理論、葉素理論以及葉素-動量理論等都是葉片設(shè)計的理論基礎(chǔ)。然而，理論研究的計算模型、求解條件都比較復(fù)雜，且實際工作環(huán)境中存在許多不可預(yù)測因素，通常需要對模型簡化以得到相關(guān)解。如BETZ理論假定來流方向與葉輪旋轉(zhuǎn)軸線方向一致并在整個葉輪掃掠面上均勻分布，風(fēng)能利用系數(shù)最大值約為0.593[4]。葉素理論把三維流動簡化為二維流動，忽略葉素間作用，但由于葉片旋轉(zhuǎn)使各葉素間實際存在氣體力，因此必須考慮干擾系數(shù)而使求解變得十分復(fù)雜[5]。數(shù)值模擬是隨著計算機(jī)技術(shù)而快速發(fā)展的一種研究方法，借助強大的計算能力能夠?qū)崿F(xiàn)RAT全三維混合流場計算，分析葉片氣動特性和載荷分布，研究三維旋轉(zhuǎn)效應(yīng)對渦輪氣動性能的影響，較理論研究更接近實際情況[6]。實驗驗證主要借助地面風(fēng)洞對渦輪葉片氣動性能進(jìn)行測試，以驗證理論研究和數(shù)值模擬的準(zhǔn)確性，但實驗費用高，研究周期長，且具有一定局限性[7]。

目前，渦輪氣動性能研究基本都采用螺旋槳/旋翼理論，由于沒有考慮三維旋轉(zhuǎn)效應(yīng)對RAT氣動特性的影響，因此不能滿足高性能RAT系統(tǒng)渦輪氣動性能計算需要[8]。劉勇[9]模擬了靜止和旋轉(zhuǎn)2個狀態(tài)下的風(fēng)力機(jī)渦輪流場，提出三維旋轉(zhuǎn)效應(yīng)對計算結(jié)果影響較大的結(jié)論。方祥軍等[10]對RAT進(jìn)行了流場數(shù)值模擬，研究了其全三維混合流場，計算結(jié)果具有一定的參考價值。目前，RAT氣動性能數(shù)值仿真大都針對地面工況進(jìn)行計算。然而，RAT作為飛機(jī)的應(yīng)急動力裝置，氣動性能隨載機(jī)的飛行高度而變化，這方面的研究國內(nèi)文獻(xiàn)還少見報道。雖然國外在RAT方面的研究已比較成熟，但國內(nèi)基本還停留在逆向設(shè)計與計算分析階段。本文采用正向研制方法設(shè)計RAT葉片，具有明顯的先進(jìn)性。

本文從目標(biāo)功率出發(fā)，首先通過對不同翼型升阻比進(jìn)行比較，選用升阻比較大的翼型，采用葉素-動量理論編程計算渦輪葉片沿展向的弦長分布、扭轉(zhuǎn)角以及誘導(dǎo)因子等，生成葉片各截面參數(shù)；然后以計算流體力學(xué)(Computational Fluid Dynamics, CFD)方法分析葉片氣動性能以滿足設(shè)計要求，研究渦輪輸出功率和風(fēng)能利用系數(shù)隨飛行高度和來流速度變化特性。本文為飛機(jī)應(yīng)急能源系統(tǒng)研究奠定了基礎(chǔ)。

1 RAT系統(tǒng)渦輪葉片設(shè)計

翼型是渦輪葉片的基本要素，影響葉片主體結(jié)構(gòu)形式，翼型氣動性能直接影響葉片和渦輪的動力特性[11]。自20世紀(jì)20年代開始，世界各主要航空發(fā)達(dá)國家建立了不同的翼型系列，如美國的NACA系列、德國的DVL系列以及英國的RAF系列等[12]。其中，NACA系列的航空翼型升阻比大，非常適用于高雷諾數(shù)流動。

通常，RAT的來流速度在100～600 km/h之間變化，以10 km為參考計算高度，依據(jù)當(dāng)?shù)丨h(huán)境參數(shù)可知雷諾數(shù)在2.7×105～1.7×106之間變化，屬于高雷諾數(shù)流動，可選擇NACA翼型設(shè)計RAT葉片。利用Profili軟件計算不同翼型在迎角為-10°～15°下的升阻比，圖1為NACA系列不同翼型的升阻比Cl/Cd隨葉片迎角α變化曲線。可以看出，0°～5°迎角時NACA4412翼型升阻比最大，作為本文RAT葉片設(shè)計的基本翼型。

圖1 NACA系列不同翼型升阻比隨葉片迎角變化曲線Fig.1 Variation curves of lift drag ratio of different NACA airfoils with attack angle of blade

渦輪直徑是RAT葉片設(shè)計時最基本的設(shè)計參數(shù)，由式(1)進(jìn)行計算[13]：

(1)

式中:D為渦輪直徑，m；P為RAT功率，W；ρ為空氣密度，kg/m3；V為來流速度，m/s；CP為風(fēng)能利用系數(shù)；η1和η2分別為發(fā)電機(jī)和傳動系統(tǒng)效率。

按照某型飛機(jī)對應(yīng)急能源系統(tǒng)的功率需求，本文設(shè)計的RAT額定功率為15 kW。以10 km為參考計算高度，當(dāng)?shù)乜諝饷芏葹?.412 7 kg/m3，取發(fā)電機(jī)和傳動系統(tǒng)效率為0.95，風(fēng)能利用系數(shù)為0.40，來流速度為450 km/h，計算可得渦輪直徑為353 mm。

依據(jù)葉素-動量理論，在已知葉片迎角和該迎角下升力系數(shù)Cl、阻力系數(shù)Cd的前提下，可由式(2)和式(3)計算徑向半徑r處的葉素入流角和扭角：

(2)

θ=φ-α

(3)

式中：φ為入流角，(°)；a和a′分別為軸向誘導(dǎo)因子和切向誘導(dǎo)因子；ω為渦輪旋轉(zhuǎn)角速度，rad/s；r為不同葉素的徑向半徑，mm；θ為扭角，(°)。

半徑r處的法向力系數(shù)Cn和切向力系數(shù)Ct分別由式(4)和式(5)計算：

Cn=Clcosφ+Cdsinφ

(4)

Ct=Clsinφ-Cdcosφ

(5)

假設(shè)軸向誘導(dǎo)因子和切向誘導(dǎo)因子的初值為零，可由式(2)得到入流角φ，由式(4)和式(5)得到法向力系數(shù)和切向力系數(shù)，然后由式(6)和式(7)重新計算誘導(dǎo)因子a和a′：

(6)

(7)

式中：σ為葉片在半徑r處的實度；F為普朗特修正因子。

(8)

(9)

其中：B為渦輪葉片數(shù)目，本文取2片葉槳；c為葉片弦長，mm；R為渦輪半徑。

把計算得到的誘導(dǎo)因子a和a′與前次計算結(jié)果進(jìn)行比較，直至誤差滿足要求時迭代終止，本文的計算誤差取0.001。通過編程計算可得不同半徑r處的葉片弦長和扭角，如圖2所示。根據(jù)NACA4412翼型數(shù)據(jù)，結(jié)合計算所得的弦長和扭角，基于點的坐標(biāo)變換求出各葉素截面在三維坐標(biāo)系中的坐標(biāo)，完成RAT的三維建模。

圖2 葉片弦長和扭角隨徑向半徑變化曲線Fig.2 Variation curves of blade chord length and twist angle with radial radius

2 RAT系統(tǒng)渦輪氣動性能理論分析

RAT系統(tǒng)中，渦輪以迎面氣流作為動力源。由空氣動力學(xué)理論可知，迎面氣流的沖壓能量與速度的三次方成線性關(guān)系，即

E=0.5ρAV3

(10)

式中：E為迎面氣流的沖壓能量，W；A為空氣流過的橫截面積, m2。

RAT在迎面氣流作用下旋轉(zhuǎn)并把氣流的沖壓能量轉(zhuǎn)換為渦輪功率而輸出，渦輪輸出功率可通過理論分析或數(shù)值模擬得到，由式(11)計算：

(11)

式中:M為渦輪葉片氣動升力沿展向的積分，N·m；n為渦輪轉(zhuǎn)速，r/min。

RAT系統(tǒng)渦輪的功率提取系數(shù)，即風(fēng)能利用系數(shù)為

(12)

葉尖速比是表征RAT系統(tǒng)渦輪運轉(zhuǎn)性能的一個重要參數(shù)，數(shù)值上等于渦輪葉片圓周速度與來流速度之比，即

(13)

式中：λ為葉尖速比。

通常，飛機(jī)處于失去動力的應(yīng)急狀態(tài)時飛行速度很低，而RAT正常工作時渦輪轉(zhuǎn)速由調(diào)節(jié)機(jī)構(gòu)控制并可保持恒定，因此葉尖速比較大，有助于提升風(fēng)能利用系數(shù)；而小的葉尖速比表示來流速度較高，此時雖然風(fēng)能利用系數(shù)不大，但由于來流速度較大，渦輪提取的功率也能滿足要求。

3 RAT系統(tǒng)渦輪氣動性能數(shù)值模擬

3.1 計算模型

本文研究的RAT裝有槳距角自動調(diào)節(jié)裝置，可根據(jù)來流狀態(tài)和工作轉(zhuǎn)速調(diào)整槳距角到最佳迎風(fēng)角度。為此，把渦輪流場劃分為3個區(qū)域，分別為最外層的靜止區(qū)域1、隨渦輪一起旋轉(zhuǎn)的旋轉(zhuǎn)域2以及包圍葉片的2個圓柱形區(qū)域31和32，各區(qū)域間以交界面連接，如圖3所示。其中，交界面1為旋轉(zhuǎn)域和靜止區(qū)域的交界，交界面2和3分別為2個葉片與旋轉(zhuǎn)域的交界。這里劃分的2個圓柱形區(qū)域31和32可以方便地改變?nèi)~片槳距角，且由于2個葉片各自獨立，因此分別建立交界面2和3以便于數(shù)值模擬時旋轉(zhuǎn)坐標(biāo)系方法的運用。

RAT旋轉(zhuǎn)對周圍流場產(chǎn)生影響，為了把RAT對流場的擾動降低到最小，理論上外邊界應(yīng)為無限大。實際上，在數(shù)值仿真時設(shè)置合適的外邊界以提高計算效率。若外邊界設(shè)置太小，氣體流動受到干擾，進(jìn)而影響計算精度。通常，取外邊界尺寸為模型尺寸的10倍或以上。考慮到本文所研究RAT直徑為353 mm，計算時取外邊界(區(qū)域1)直徑為4400mm，入口端長度為1000mm，出口端長度為3 600 mm；旋轉(zhuǎn)域2半徑為210 mm，長為160 mm；區(qū)域3為2個半徑為40 mm、長為155 mm的圓柱，如圖4所示。此外，區(qū)域劃分時必須滿足2個要求：①2個圓柱區(qū)域在靠近輪轂的端面與葉片根部相切；②2個圓柱區(qū)域和葉片共軸線，其目的是保證每次轉(zhuǎn)動區(qū)域31和32時RAT具有完整性。

圖3 RAT渦輪流場網(wǎng)格Fig.3 Grid in turbine flow field of RAT

高質(zhì)量有限元網(wǎng)格是高精度數(shù)值模擬的前提。考慮到黏性氣流參數(shù)沿垂直于壁面方向變化劇烈，采用密集細(xì)化的網(wǎng)格。第1層流體網(wǎng)格高度為湍流邊界層，可根據(jù)布拉修斯方程估算[14]：

(14)

圖4 流體區(qū)域劃分Fig.4 Division of fluid region

取初始值y+=30，可計算得到第1層網(wǎng)格高度，隨來流速度而變化。為精確求解湍流邊界層，在獲得流場收斂解后，把求解得到的y+代入式(14)重新計算第1層網(wǎng)格高度，可以對網(wǎng)格進(jìn)一步優(yōu)化。

為提高模擬精度，對葉片和輪轂進(jìn)行表面網(wǎng)格加密，并在貼近渦輪壁面處布置6層棱柱網(wǎng)格。流場網(wǎng)格總數(shù)約200萬。

3.2 邊界條件

1) 進(jìn)出口邊界條件。RAT工作條件是大雷諾數(shù)流動，理想氣體并做湍流運動。進(jìn)出口邊界設(shè)置為速度入口、壓力出口和壓力遠(yuǎn)場條件。來流速度V垂直于入口邊界，渦輪轉(zhuǎn)速為n。

2) 壁面邊界條件。壁面滑移特性與渦輪外部流動特性有關(guān)。若流體為無黏流動，則渦輪壁面具有滑移特性，視具體情況給定粗糙度和壁面滑移等。本文研究的RAT為無滑移壁面。此外，為減小邊界干擾，設(shè)置靜止域外圍遠(yuǎn)場邊界與渦輪距離為渦輪直徑的10倍以上，取3 600 mm。

3.3 數(shù)值模擬

RAT處在一個可壓縮有黏性的非定常流場中，采用多重旋轉(zhuǎn)坐標(biāo)系(MRF)模型把問題轉(zhuǎn)化為穩(wěn)定流動。選用CFD方法對渦輪氣動性能進(jìn)行數(shù)值模擬。空氣設(shè)為理想氣體，黏性系數(shù)隨溫度的變化以Sutherland公式計算，即

(15)

式中：μ0為溫度等于288.15 K時空氣的黏性系數(shù)，μ0=1.789 4×10-5N·s/m2；C為常數(shù)，取值110.4。

湍流模型選用k-ω模型[15]，以SIMPLEC算法進(jìn)行求解，求解方式為速度與壓力的耦合。設(shè)置旋轉(zhuǎn)域2的轉(zhuǎn)速為渦輪轉(zhuǎn)速，2個圓柱形區(qū)域31和32相對于旋轉(zhuǎn)域2靜止，相對速度為零。

4 計算結(jié)果分析

飛機(jī)失去主動力后，RAT進(jìn)入工作狀態(tài)，渦輪從迎面氣流獲取能量并轉(zhuǎn)換機(jī)械能，輸出功率與渦輪葉尖速比、槳距角、氣體入流角以及環(huán)境參數(shù)等有關(guān)，在飛行包線內(nèi)環(huán)境參數(shù)隨飛行高度而變化。輸出功率和風(fēng)能利用系數(shù)是評價渦輪氣動性能的2個重要參數(shù)。

圖5 渦輪風(fēng)能利用系數(shù)隨葉尖速比變化曲線Fig.5 Variation curves of rotor power coefficient of turbine with tip speed ratios

4.1 不同葉尖速比時渦輪氣動性能

葉尖速比和槳距角是影響渦輪氣動性能的2個重要參數(shù)。當(dāng)渦輪轉(zhuǎn)速恒定，來流速度變化時葉尖速比隨之變化。圖5為地面性能計算時恒速工況下(n=1×104r/min)風(fēng)能利用系數(shù)隨葉尖速比變化曲線。可見，同一槳距角下，隨著葉尖速比的增加，風(fēng)能利用系數(shù)先增大后減小，每一槳距角對應(yīng)一個最大風(fēng)能利用系數(shù)和最佳葉尖速比，且槳距角越小，最大風(fēng)能利用系數(shù)越大，對應(yīng)的葉尖速比也相應(yīng)增大；槳距角為40°時，最大風(fēng)能利用系數(shù)只有0.1，對應(yīng)的葉尖速比只在小范圍變化，已不能正常運轉(zhuǎn)。因此，為提高渦輪風(fēng)能利用系數(shù)，應(yīng)根據(jù)來流狀態(tài)適當(dāng)調(diào)整槳距角。

4.2 不同槳距角時渦輪氣動性能

圖6為渦輪恒速(n=1×104r/min)旋轉(zhuǎn)工況下，地面性能計算時輸出功率和風(fēng)能利用系數(shù)隨槳距角變化曲線。來流速度相同時，隨著槳距角的增大，輸出功率和風(fēng)能利用系數(shù)都是先增大后減小，每一來流速度對應(yīng)一個最佳槳距角，且來流速度越大，最大風(fēng)能利用系數(shù)和渦輪輸出功率越大，但當(dāng)來流速度大于某一值時風(fēng)能利用系數(shù)反而下降，如圖6(b)中來流速度為450 km/h時所示。由圖6(a)看出，來流速度為340 km/h時，RAT的最大輸出功率為11.9 kW，小于設(shè)計功率15 kW；來流速度為450 km/h時，RAT的最大輸出功率為19 kW, 大于設(shè)計功率15 kW。因此，渦輪輸出功率隨來流速度增加而增大，不同來流速度時可通過調(diào)節(jié)槳距角改變輸出功率，以滿足設(shè)計要求。當(dāng)調(diào)整到最佳槳距角之前，隨著槳距角的增大，渦輪輸出功率緩慢上升，到達(dá)最佳槳距角之后，輸出功率迅速下降。

圖6 渦輪輸出功率和風(fēng)能利用系數(shù)隨槳距角變化曲線Fig.6 Variation curves of turbine output power and rotor power coefficient with pitch angles

4.3 不同入流角時渦輪氣動性能

入流角是來流方向與渦輪旋轉(zhuǎn)平面的夾角，數(shù)值模擬時可通過調(diào)節(jié)3個坐標(biāo)分量改變?nèi)肓鹘恰．?dāng)入流角變化時，氣流作用于葉片表面的升力發(fā)生變化，渦輪輸出功率和效率(風(fēng)能利用系數(shù))隨之改變。為模擬不同入流角時渦輪氣動性能，需要對外層靜止區(qū)域1重新劃分網(wǎng)格，本文把來流方向的外層靜止區(qū)域處理為橢球形，截面為長軸7 000 m、短軸5 000 m的橢圓，而內(nèi)部旋轉(zhuǎn)域2和區(qū)域3由于之前的合理劃分則保持不變。

圖7(a)為不同飛行高度下渦輪輸出功率隨入流角變化曲線。可以看出，飛行高度一定時，隨著入流角增大，輸出功率逐漸減小；入流角一定時，輸出功率隨飛行高度增加而減小。因此，當(dāng)其他條件相同時，入流角使渦輪輸出功率降低。以3 km飛行高度計算結(jié)果為例，圖7(b)給出了輸出功率隨入流角變化曲線，葉尖速比一定時，隨著入流角增大，輸出功率逐漸減小，且葉尖速比越大，輸出功率下降越快；當(dāng)入流角一定時，葉尖速比越小，輸出功率越大。因此，入流角的存在使渦輪輸出功率降低，最大功率為14.2 kW, 小于額定功率15 kW。渦輪轉(zhuǎn)速一定時，為獲得較大輸出功率，應(yīng)盡可能減小葉尖速比，即提高飛行速度，同時保持較小的入流角。

圖7 渦輪輸出功率隨入流角變化曲線Fig.7 Variation curves of turbine output power with inflow angle

圖8 葉片壓力面和吸力面壓力云圖Fig.8 Pressure contours of acting surface and suction surface

4.4 渦輪流場分析

氣流沿葉片分布和流動情況直接影響渦輪輸出功率和能量提取效率。圖8為地面性能計算時槳距角為8°、渦輪轉(zhuǎn)速為9 700 r/min、來流速度為210 km/h時葉片表面壓力云圖。可以看出，葉片壓力面的壓力中心靠近葉尖區(qū)域，高壓區(qū)域位于80%相對葉高處；在相對葉高2%～5%的葉尖區(qū)域等值線密布，靜壓梯度大，而葉根區(qū)域等值線比較稀疏，壓力梯度小，主要原因是氣流在葉片表面有一向心加速度，因此沿徑向存在壓力差，使得從葉根到葉尖，靜壓逐漸增大，即為負(fù)的逆壓梯度。

在葉片吸力面上從葉尖到葉根，靜壓先降低，然后維持一段低壓區(qū)域后逐漸增大，葉片根部壓力較大。可見，壓力面與吸力面的壓差也是靠近葉尖區(qū)域較大，因此葉片主要做功區(qū)域是中部靠近葉尖部分，即相對葉高40%～95%區(qū)域。此外，吸力面葉根區(qū)域靜壓分布不甚合理，原因可能是受到葉根迎角和附面層影響，因此應(yīng)設(shè)法減小葉根迎角，改善流動狀態(tài)，并采用不等迎角設(shè)計。

渦輪轉(zhuǎn)速一定時，來流速度和槳距角影響葉片周圍流線分布。圖9(a)、(b)分別給出了渦輪轉(zhuǎn)速為9 700 r/min、槳距角為8°時，不同來流速度時葉片截面流線分布。隨著來流速度的升高，在葉片吸力面靠近后緣區(qū)域，有明顯的氣體流動分離現(xiàn)象，且吸力面上葉尖與根部區(qū)域氣體流速明顯增大，與壓力面間的壓差增大，做功能力增加；圖9(c)、(d)為來流速度為240 km/h時，不同槳距角下葉片截面流線分布。隨著槳距角增大，吸力面一側(cè)靠近葉片中間區(qū)域有明顯的漩渦，導(dǎo)致葉片所受升力減小，進(jìn)而引起輸出功率下降。但從整體效應(yīng)看，小槳距角和低來流速度下，本文所設(shè)計的RAT尾緣厚度小，計算結(jié)果沒有粗大雜亂的尾跡，流動狀況比較理想，但仍有改進(jìn)的空間。

圖9 不同來流速度和槳距角下流線圖Fig.9 Streamlines at different inflow velocities and pitch angles

5 結(jié) 論

依據(jù)某型應(yīng)急能源系統(tǒng)對沖壓空氣渦輪(RAT)的功率需求，采用正向研制方法，基于葉素-動量理論設(shè)計了RAT葉片；采用CFD方法仿真計算了不同工況下RAT氣動性能，主要結(jié)論如下：

1) 渦輪轉(zhuǎn)速一定時，風(fēng)能利用系數(shù)隨葉尖速比的增加先增大后減小，每一槳距角對應(yīng)一個最大風(fēng)能利用系數(shù)和最佳葉尖速比，且槳距角越小，最大風(fēng)能利用系數(shù)越大。

2) 針對不同來流速度，適當(dāng)調(diào)整葉片槳距角不僅可以提高渦輪風(fēng)能利用系數(shù)，同時還可以實現(xiàn)渦輪恒功率輸出，滿足設(shè)計要求。

3) 從整個流場來看，葉片主要做功區(qū)域位于中部靠近葉尖部分、相對葉高40%～95%的區(qū)域，吸力面葉根區(qū)域靜壓分布不甚合理，仍有改進(jìn)空間。

4) 采用正向研究方法，基于葉素-動量理論設(shè)計的RAT葉片能夠滿足某型RAT需要，研究方法具有明顯的先進(jìn)性。