基于CFD的對轉(zhuǎn)槳扇發(fā)動機性能確定技術(shù)研究

王定奇， 李秋鋒， 于洋， 黎森

(中國飛行試驗研究院發(fā)動機所， 西安 710089)

開式轉(zhuǎn)子發(fā)動機，也稱為對轉(zhuǎn)槳扇發(fā)動機，這一概念美國國家航空航天局(National Aeronautics and Space Administration，NASA)于20世紀(jì)70年代首先提出來[1]，是一種介于渦槳與渦扇之間的一種發(fā)動機。對轉(zhuǎn)槳扇發(fā)動機結(jié)合了渦槳發(fā)動機和渦扇發(fā)動機的優(yōu)勢，相對于單排槳葉的渦槳發(fā)動機，對轉(zhuǎn)槳扇發(fā)動機可以在高空大馬赫數(shù)(Ma=0.8、10 km)保持相對較低的耗油率[2]；與同等技術(shù)的渦扇發(fā)動機相比，由于其超高的涵道比(通常為25～60)，槳扇發(fā)動機具有推進(jìn)效率高，燃油消耗率低，使用費用低，循環(huán)效率高，起飛性能好的技術(shù)優(yōu)勢；與渦槳發(fā)動機相比，其巡航速度更高，后排反轉(zhuǎn)槳葉減小了前排槳葉的渦流，流場畸變較小，后排轉(zhuǎn)子抵消了前排槳葉的力矩。但由于高速旋轉(zhuǎn)的槳葉造成的噪音和震動問題難以解決；且由于槳葉高轉(zhuǎn)速，無外涵道殼體保護(hù)，一旦葉片斷裂會機體產(chǎn)生嚴(yán)重的破壞，限制了槳扇發(fā)動機的應(yīng)用。近年來各航空公司為了壓縮運營成本及降低污染物排放，美國通用電氣公司、Rolls-Royce公司等主要航空發(fā)動機制造商重新開始重視槳扇發(fā)動機[3-4]。美國通用電氣公司曾在MD-80飛機上，使用兩種發(fā)動機進(jìn)行對比試飛，結(jié)果使用槳扇發(fā)動機時油耗下降了49%。裝配D-27槳扇發(fā)動機[5]的安-70，其以最大巡航速度飛行時，油耗比使用同等推力渦扇發(fā)動機要少20%～30%。

對于槳扇發(fā)動機的部件性能確定通常有風(fēng)洞試驗和數(shù)值模擬兩種方法。美國曾開展587DX槳扇發(fā)動機風(fēng)洞試驗，并形成了相關(guān)推阻力劃分體系[6]。在數(shù)值計算方面，F(xiàn)erraro等[7]在Euler/Navier-Stokes(N-S)方程求解中嵌入動量理論作用盤模型，應(yīng)用于螺旋槳的數(shù)值模擬，但計算中螺旋槳近似為作動盤，僅能得到平均計算結(jié)果，無法考慮葉片的外形，及其周圍流場細(xì)節(jié)。Stuermer[8]利用計算流體力學(xué)(computational fluid dynamics,CFD)開展了對轉(zhuǎn)槳扇非定常三維流場數(shù)值仿真，分析了轉(zhuǎn)子葉尖渦系結(jié)構(gòu)變化、雙排葉片干擾流動及不同攻角下槳扇的氣動力。Node-Langlois等[9]建立了噪音預(yù)測模型，分析了單獨槳扇發(fā)動機和裝機后槳扇發(fā)動機噪音。近年來，學(xué)者們針對螺旋槳流場的數(shù)值分析開展了大量研究。Xu等[10]采用基于非結(jié)構(gòu)重疊網(wǎng)格方法研究了螺旋槳與機身的氣動干擾；段中喆等[11]通過對螺旋槳滑流三維流場進(jìn)行數(shù)值模擬，對比了不同工況下滑流區(qū)渦系結(jié)構(gòu)；楊小川等[12]研發(fā)了TRIP3.0軟件平臺，運用動態(tài)拼接結(jié)構(gòu)網(wǎng)格技術(shù)，模擬了螺旋槳非定?；鳌ａ槍劝l(fā)動機的分析，王逸維等[13]建立了三軸拉力式對轉(zhuǎn)槳扇發(fā)動機仿真模型，基于相似理論，評估了槳扇發(fā)動機的性能；張帥等[14]建立了一種用于總體方案論證階段的開始轉(zhuǎn)子發(fā)動機分析模型，能基于少量參數(shù)估算發(fā)動機推力與耗油率特性。

由于對轉(zhuǎn)槳扇發(fā)動機前/后排槳葉角、轉(zhuǎn)速均可調(diào)，較渦槳發(fā)動機控制參數(shù)更多，目前尚無成熟的型號提供參考，因此槳扇相關(guān)流場氣動分析研究較少。參考587DX槳扇發(fā)動機的相關(guān)氣動分析方法，建立三維槳扇模型，并利用CFD數(shù)值仿真，定量計算出裝機條件下槳扇發(fā)動機部件特性，為后續(xù)通過CFD方法建立槳扇發(fā)動機的螺旋槳部件的安裝性能及前后排槳葉的氣動耦合性能性提供參考。

1 槳扇發(fā)動機受力分析

所研究的某型槳扇發(fā)動機采用“采用后置推進(jìn)式螺旋槳”形式，兩排列直徑為3.8 m的后掠式葉片。前排槳扇為6葉槳(順航向逆時針旋轉(zhuǎn))，后排槳扇也為6個槳葉(順航向順時針旋轉(zhuǎn))。

參考Propeller/PropfanIn-FlightThrustDetermination(SAE-AIR 4065A—2012)可以建立動力裝置總凈推力F′net的計算表達(dá)式，分別取上下游無窮遠(yuǎn)截面，取向右為受力正方向，則控制體(圖1)受力表達(dá)式為

F′net=W10(V0+w)+W10V∞-W10V0-W0V0

(1)

其中：Wi為控制體流量，i=0,1,2,…,∞；Vi為速度,i=0,1,2,…,∞；w為槳扇后方氣流速度增量。

將整個控制體劃分為內(nèi)流控制體和外流控制體。外流控制體為0-∞截面，通過短艙、吊掛和槳葉表面氣流；內(nèi)流控制體分為0-1截面的進(jìn)氣控制體，0-∞的排氣系統(tǒng)控制體，可得

F′net=W10w+φpost+φplug-W0V0-FG9

(2)

或者

F′net=TS-Dnac-Dstr+FG9-φpre

(3)

式(3)中：TS為槳葉拉力；Dnac為短艙阻力；Dstr為吊掛阻力；FG9為排氣推力；φpre附加前體力；φplug后體中心錐阻力。

由于內(nèi)部流管附加阻力和內(nèi)部流管凈推力，可以在發(fā)動機關(guān)鍵截面處加裝測量耙，通過燃?xì)獍l(fā)生器法[15]計算確定，故此兩項無需研究。主要研究目標(biāo)為確定槳扇的部件特性，通過數(shù)值積分獲取槳葉表面的壓差力和摩擦力，從而得到槳葉的拉力及功率特性。

P0為大氣靜壓；PS為槳扇后方靜壓圖1 推進(jìn)式槳扇發(fā)動機受力分析示意圖Fig.1 Schematic diagram of engine stress analysis

2 計算模型及網(wǎng)格生成

2.1 模型分析

通過CATIA(computer-graphics aided three-dimensional interactive application)建立三維實體模型，在UG(Unigraphics)中對槳扇進(jìn)行光順和模型的優(yōu)化，導(dǎo)出可用于網(wǎng)格分析的模型。參考SAE AIR 4065A中的587DX槳扇發(fā)動機(圖2)相關(guān)參數(shù)，設(shè)計前后排對轉(zhuǎn)，帶后掠形式的推進(jìn)式槳扇發(fā)動機。

圖2 587DX三維模型Fig.2 Three dimensional model of 587DX

2.2 模型驗證及網(wǎng)格劃分

槳扇的轉(zhuǎn)速為1 200～1 600 r/min，飛行條件中，螺旋槳前進(jìn)比范圍在0.6～2.8，對應(yīng)的雷諾數(shù)達(dá)到106量級。首先進(jìn)行網(wǎng)格無關(guān)性驗證，單排槳扇內(nèi)網(wǎng)格分別設(shè)置為2×106、3×106、4×106，計算中設(shè)置馬赫數(shù)Ma為0.6、高度Hp為8 km、轉(zhuǎn)速n=1 400 r/min工況點下，3×106網(wǎng)格與4×106網(wǎng)格計算槳扇拉力相差為0.5%。由于中國目前還沒有成熟的槳扇發(fā)動機，因此槳扇的建模方法參考同等尺寸及轉(zhuǎn)速范圍內(nèi)的螺旋槳建模方法[16-19]，開展螺旋槳特性的計算。

從表1可以看出，在發(fā)動機小狀態(tài)下仿真計算誤差為4.16%，而大狀態(tài)下計算誤差為2.71%，結(jié)算結(jié)果滿足工程精度要求，驗證了仿真模型計算的準(zhǔn)確度。

葉片表面邊界層網(wǎng)格高度初始值為0.3 mm，近壁面y+控制約為1，增長率設(shè)置1.3，邊界層網(wǎng)格共15層。計算域流場交界面Interface網(wǎng)格如圖3所示，槳扇表面網(wǎng)格如圖4所示。網(wǎng)格量為3×106，總模型網(wǎng)格為12×106。

采用非結(jié)構(gòu)四面體網(wǎng)格對前/后排槳扇網(wǎng)格劃分，槳扇流場整體網(wǎng)格如圖5所示。由于槳扇做旋轉(zhuǎn)運動，故需將整個計算域分為3個區(qū)域，采用動態(tài)拼接網(wǎng)格技術(shù)：包含槳扇旋轉(zhuǎn)域，外流場靜止域。采用塊分界面上的網(wǎng)格面直接相連的方式，通過交界面進(jìn)行流場信息的傳遞。采用以上方法網(wǎng)格生成難度降低，且節(jié)點分布靈活，用UG建立槳扇模型整體計算域，區(qū)域a、b為前后排槳葉計算域，其幾何尺寸為：半徑R=2 m，長度L=0.5 m；計算域c為半模計算域：長300 m，寬50 m，高 60 m。模型計算域如圖6所示。

表1 功率系數(shù)對比Table 1 Comparison of Power factor

圖3 交界面處網(wǎng)格分布Fig.3 Grid distribution of interface

圖4 槳扇表面網(wǎng)格分布Fig.4 Grid distribution of the blade

圖5 計算域網(wǎng)格Fig.5 Grid of computational zone

在流場求解中，采用有限體積法對控制方程進(jìn)行空間離散黏性項采用中心差分格式離散，無粘項采用二階Rose迎風(fēng)偏置通量差分方法離散。時間推進(jìn)采用隱式近似因子分解法，且采用多重網(wǎng)格加速收斂算法。湍流模型選擇壁面加強型RNGk-ε(k為湍動能，ε為耗散率)湍流模型。

圖6 模型計算域Fig.6 Model of computational zone

3 計算結(jié)果及分析

3.1 槳扇特性參數(shù)確定

在槳扇轉(zhuǎn)速1 400 r/min，高度4 km條件下，槳扇的部件特性隨來流前進(jìn)比、槳葉角度的變化如圖7、圖8所示。可以看出，拉力系數(shù)CT和功率系數(shù)CP都是隨著前進(jìn)比J的增大而降低，其變化規(guī)律與渦槳發(fā)動機的槳葉變化規(guī)律一致，但由于槳扇發(fā)動機的前后排槳葉旋轉(zhuǎn)方向相反，使氣流在通過后排槳葉時壓力增大，且氣流方向與后排槳葉相反，葉盆對氣流的作用力更強。相同槳葉角下，前排槳葉(PropFan1)拉力系數(shù)隨著前進(jìn)比的降低幅度大于后排槳葉(PropFan2)；相同前進(jìn)比條件下隨著槳葉角的增大拉力系數(shù)增大，在槳葉角30°變化到35°，槳扇的拉力變化范圍可達(dá)37.8%，功率變化可達(dá)28.4%。

圖7 拉力系數(shù)隨前進(jìn)比變化曲線Fig.7 The curve for thrust coefficient followed the velocity

圖8 功率系數(shù)隨前進(jìn)比變化曲線Fig.8 The curve for power coefficient followed the velocity

3.2 滑流給機翼表面流場影響

圖9為槳葉角(θ=35°)狀態(tài)下，轉(zhuǎn)速1 400 r/min，來流Ma由0.3增大至0.6后，槳扇表面壓力云圖?？梢钥闯?槳扇葉背處壓力小于葉盆處壓力，因而使槳扇產(chǎn)生沿飛行方向的拉力。隨著飛行馬赫數(shù)的增大，前排葉片葉背處整體壓力增大，且壓力低壓區(qū)域(槳葉表面藍(lán)色部分)在葉片0.75R葉高以上增大的更加明顯(R為槳扇半徑)；槳扇葉盆處的高壓區(qū)域(槳葉表面紅色部分)隨著Ma增大而減小，且向葉尖移動，此時葉片前后壓差減小，使得槳扇的拉力降低。

圖9 不同Ma下葉片表面壓力分布(θ=35°)Fig.9 The pressure distribution of blade under the Ma(θ=35°)

3.3 不同槳葉角對螺旋槳部件特性的影響

圖10為槳葉角35°，來流Ma=0.4狀態(tài)下，不同發(fā)動機轉(zhuǎn)速情況下對應(yīng)的截面流場分布，可以看出隨著轉(zhuǎn)速的增大在葉尖區(qū)域的Ma增大；且槳葉對于氣流的做功能力明顯增強，槳扇后方的高壓氣流區(qū)域范圍更大，表明槳扇轉(zhuǎn)速增大，飛行拉力隨之增大。由于葉根處半徑小，線速度低，葉根后擾流強度大，Ma較小。

圖10 發(fā)動機截面Ma分布Fig.10 Distribution of Ma for section of engine

圖11為槳葉角35°，Ma=0.4不同轉(zhuǎn)速(NPR)下槳葉的拉力系數(shù)和功率系數(shù)變化曲線，可以看出，當(dāng)轉(zhuǎn)速小于1 600 r/min時，隨著轉(zhuǎn)速的增大拉力系數(shù)、功率系數(shù)基本保持線性增大；在轉(zhuǎn)速1 600 r/min下，前后排槳葉的拉力系數(shù)最大相差33%，功率系數(shù)相差17%。當(dāng)轉(zhuǎn)速達(dá)到1 800 r/min后，拉力系數(shù)增加幅度減小，而功率系數(shù)開始降低。通過分析主要是由于當(dāng)來流Ma和槳葉角不變時，轉(zhuǎn)速增大，使得氣流的和速度增大，氣流通過前后排葉片氣流偏轉(zhuǎn)角度大，在高轉(zhuǎn)速下氣流分離較嚴(yán)重。

4 結(jié)論

結(jié)合SAE AIR 4065A中槳扇發(fā)動機的氣動特征，自行設(shè)計了一對轉(zhuǎn)槳扇發(fā)動機的三維模型，并基于該模型進(jìn)行了槳扇發(fā)動機的槳葉部件特性的CFD計算，研究了不同工況下槳葉特性的變化規(guī)律，得到如下結(jié)論。

(1)對轉(zhuǎn)槳扇的前排葉片的拉力系數(shù)、功率系數(shù)均小于后排槳葉，主要是由于氣流通過前排槳葉后受到葉片做功，使氣流在通過后排槳葉時壓力增大，且氣流方向與后排槳葉相反，葉盆對氣流的作用力更強。

圖11 不同NPR下的槳扇特性參數(shù)變化曲線Fig.11 Variation curve of propfan characteristic parameters under different NPR

(2)隨著轉(zhuǎn)速的增大，拉力系數(shù)功率系數(shù)隨之增大，在轉(zhuǎn)速1 600 r/min下，前后排槳葉的拉力系數(shù)最大相差33%，功率系數(shù)相差17.6%。

(3)對轉(zhuǎn)槳扇的前后排槳葉氣動耦合特性明顯，且由于前后排槳葉旋轉(zhuǎn)方向不同，在高轉(zhuǎn)速情況下，會發(fā)生氣流較大的分離，槳葉的最大效率為82.6%。

(4)槳扇的槳葉角對于發(fā)動機的總推力影響非常大，是關(guān)鍵的性能調(diào)節(jié)參數(shù)。在相同的飛行Ma、轉(zhuǎn)速和高度條件下，槳葉角30°～35°的變化范圍，槳扇的拉力變化范圍可達(dá)37.8%，功率變化可達(dá)28.4%。

(5)參考SAE AIR 4065A中槳扇發(fā)動機的推/阻力劃分體系，結(jié)合CFD獲得槳葉的部件特性，為后續(xù)槳扇發(fā)動機飛行中部件特性的獲取提供參考。

由于對轉(zhuǎn)槳扇發(fā)動機，采用前后排槳葉對轉(zhuǎn)設(shè)計，轉(zhuǎn)速較高，槳葉角度，轉(zhuǎn)速都可調(diào)，因此槳扇發(fā)動機的推力相對于渦槳發(fā)動機受到更多變量影響。本文限于篇幅，計算中限定前后排槳葉角相同、轉(zhuǎn)速均相同，后續(xù)可以在本文的研究基礎(chǔ)上進(jìn)行更多維度的轉(zhuǎn)速及槳葉角組合，進(jìn)行前后排槳葉氣動優(yōu)化的計算；以及結(jié)合槳扇發(fā)動機內(nèi)流的狀態(tài)，進(jìn)行槳葉與進(jìn)排氣系統(tǒng)的氣動耦合分析。