液體火箭發(fā)動機(jī)多級跨聲速渦輪流面迭代計算

雷宗琪 梁國柱

(北京航空航天大學(xué) 宇航學(xué)院，北京100191)

渦輪是液體火箭發(fā)動機(jī)渦輪泵的核心部件，對渦輪內(nèi)部三維流動過程的深刻認(rèn)識對于提高液體火箭發(fā)動機(jī)渦輪泵的性能具有重要的意義.目前國內(nèi)外對于葉輪機(jī)械內(nèi)部三維流動的計算主要分為直接求解三維控制方程組的CFD(Computational Fluid Dynamics)方法和相對流面法兩大類基本方法.其中CFD方法可以直接求解三維Euler方程[1－2]和Navier-Stokes方程[3－4]，具有較高的計算精度，能夠比較真實地反映流動情況，但對計算資源和時間的消耗也過于巨大.而相對流面法的思路是基于吳仲華教授提出的相對流面理論[5]，這一理論把三維流場分解為S1流面和S2流面兩組二維流面，通過流面迭代得到準(zhǔn)三維或全三維的解.流面內(nèi)二維流場的求解常采用流線曲率法.與CFD方法相比，相對流面法和流線曲率法對計算資源和時間的消耗小，可以同時應(yīng)用于流場計算和葉型設(shè)計，而且配合適當(dāng)?shù)膿p失模型也能夠得到高精度結(jié)果[6]，從而在工程計算中具有很高的實用價值.國內(nèi)外有大量采用多個S1流面和一個S2流面迭代求解葉輪機(jī)械內(nèi)部準(zhǔn)三維流場的研究，但多個S1流面和多個S2流面迭代的全三維方法卻由于復(fù)雜的計算流程而沒有得到推廣[7－8]，僅有的一些研究報道也僅局限于亞聲速葉輪機(jī)械流場的求解[9].目前關(guān)于渦輪的公開研究主要集中于民用和航空領(lǐng)域，針對液體火箭發(fā)動機(jī)渦輪泵中渦輪部分的研究較為少見.相對于民用渦輪機(jī)和航空發(fā)動機(jī)內(nèi)的燃?xì)鉁u輪而言，液體火箭發(fā)動機(jī)渦輪泵使用的渦輪具有尺寸小、負(fù)荷高、轉(zhuǎn)速快等特點.馬冬英等[10]曾采用基于正交線的流線曲率法對某型液體火箭發(fā)動機(jī)跨聲速多級渦輪流場進(jìn)行了計算并得到了不同流面內(nèi)的二維無激波結(jié)果，但未能將研究拓展至三維流場，也沒有考慮激波對流場的影響.

本文以某型液體火箭發(fā)動機(jī)二級跨聲速渦輪為例，使用Fortran語言編制計算程序，采用基于任定準(zhǔn)正交線的流線曲率法求解二維流面內(nèi)的流場，并通過流面迭代的方法進(jìn)一步得到了三維流場.在計算過程中采用混合平面法處理動、靜葉間的非定常干涉，通過Miller正激波模型計算激波損失，用總壓恢復(fù)系數(shù)計算流動過程中的黏性損失，并對結(jié)果以及計算方法的精確性進(jìn)行了分析.

1 理論模型

相對流面的基本思想是把葉輪機(jī)械內(nèi)部的三維流場分解在兩類相對流面上，分別稱為S1流面和S2流面[5].在軸流葉輪機(jī)械中S1流面和S2流面的示意圖和相對關(guān)系如圖1所示.

圖1 軸流式葉輪中的流面

1.1 基本假設(shè)

本文采用的理論模型基于以下假設(shè):

1)渦輪內(nèi)的燃?xì)馐抢硐霘怏w;

2)渦輪內(nèi)的流動是相對定常的;

3)渦輪內(nèi)的流動是絕熱的;

4)忽略渦輪內(nèi)燃?xì)獾捏w積力.

需要注意的是，本文的研究方法為全三維方法，認(rèn)為S1流面是可以翹曲的.

1.2 流線曲率法

流線曲率法在求解時主要涉及以下方程:

1)速度梯度方程

經(jīng)過一系列變換，連續(xù)方程、動量方程和能量方程被整合為可以直接進(jìn)行數(shù)值計算的速度梯度方程:

式中，q表示連接橫向計算站的曲線的弧長;m表示子午流線，即流線在子午面上投影的弧長;β表示相對速度與子午面的夾角;φ表示子午流線與z軸的夾角;Wm表示相對速度在子午流線方向的分量;表示相對滯止轉(zhuǎn)焓;s表示熵.相關(guān)參數(shù)的關(guān)系如圖2所示.

圖2 流線相關(guān)參數(shù)示意圖

2)流量方程

S1流面和S2流面內(nèi)的流量方程具有相似的形式，均為沿計算站進(jìn)行的流量積分.

S1流面流量方程:

S2流面流量方程:

其中，q表示計算站方向的微元長度;κ，λ表示流線與S1流面和S2流面內(nèi)計算站的夾角;τ，σ表示S1流面和S2流面的厚度.

2 求解過程

2.1 邊界條件

渦輪葉片流道的邊界由若干曲面組成，如圖3所示.其中下表面A1為輪轂面，上表面A2為輪蓋面，左面A3為入口截面，右面A4為出口截面.流道側(cè)表面的中間部分A5，A6為葉片表面，前后兩段A7，A8為周期性邊界.

圖3 單級渦輪的三維流道

1)固體壁面A1，A2，A5，A6

A1，A2，A5，A6處的邊界均為固體壁面，應(yīng)當(dāng)滿足不可滲透條件，即相對速度W沿壁面法向n的分量為0:

2)進(jìn)出口截面A3，A4

由于動靜葉間存在非定常干涉，上一級葉片的總溫和總壓分布不能直接傳遞給下一級葉片，本文采用了混合平面法進(jìn)行處理.

根據(jù)混合平面法的思路，上一級葉片的出口截面A4和下一級葉片的入口截面A3相互重疊，該界面即為混合平面，如圖4所示.對于S1流面的流動計算，由圖4a可見各葉柵排的流道寬度不同，各條流線不能直接連接，所以計算過程中對周向的總溫和總壓進(jìn)行平均處理.對于S2流面的流動計算，由圖4b可見在混合平面前后沿徑向各條流線可以連接，所以混合平面后一區(qū)域S2流面各條流線的入口總溫和總壓分別取前一區(qū)域?qū)?yīng)各條流線的出口總溫和總壓.

3)周期性邊界A7，A8

圖4 混合平面示意圖

A7和A8是周期性邊界，應(yīng)當(dāng)滿足周期性條件.由于A7和A8的初始形狀和位置都是未知的，所以采用迭代修正的方法求解.

2.2 損失模型

渦輪內(nèi)部的真實流動過程非常復(fù)雜，存在各種損失，例如葉型損失、二次流損失、激波損失等，一般采用經(jīng)驗公式和模型對葉輪機(jī)械內(nèi)部的損失進(jìn)行計算.

就本文而言，由于缺乏相關(guān)數(shù)據(jù)，在計算黏性損失時無法直接套用已公開的損失模型.目前采用的方法是將這些損失等效為總壓損失，以在每一級計算站上施加總壓恢復(fù)系數(shù)的方式進(jìn)行計算.

在入口條件為超聲速的渦輪級中，激波的影響是不可忽略的.本文在計算激波損失時采用了Miller激波模型，該模型在各類跨聲速葉輪機(jī)械的流場計算中被廣泛使用，其主要特點是將葉輪流道內(nèi)的激波系簡化為入口處的一道正激波，并以波前平均馬赫數(shù)M'為依據(jù)計算通過激波的總壓損失[11]，如圖5所示.

圖5 Miller激波模型

圖中為葉輪流道入口處的來流馬赫數(shù)，為激波前吸力面處的馬赫數(shù)，波前平均馬赫數(shù)的計算方法為

對于激波與附面層的相互作用，本文的處理方法是在激波后的流場計算中加入流量系數(shù).

2.3 跨聲速問題的解決方法

流線曲率法采用流量作為計算站內(nèi)的流速分布是否正確的判據(jù)，在求解跨聲速流場時會遇到雙值問題，即對于同樣的流量m1，同時存在亞聲速w1和超聲速w'1兩種速度與之對應(yīng)，見圖6.

圖6 雙值問題

針對雙值問題的解決方法是在計算過程中逐個計算站檢查計算流量以確定臨界截面的位置，之后將流場分成亞聲速區(qū)和超聲速區(qū)分別求解.

2.4 流面迭代過程

流面間的迭代是逐級進(jìn)行的，求出一級葉片流道的全三維流場之后才開始進(jìn)行下一級的計算.流面迭代是否收斂的判斷方法是相同空間位置上不同流面的求出的流動參數(shù)是否一致.本文中采用的判據(jù)是絕對馬赫數(shù)的相對誤差ε，其定義為

式中，MS1和MS2分別表示同一計算節(jié)點在S1流面和S2流面內(nèi)的馬赫數(shù).

3 計算示例和結(jié)果分析

3.1 搭建模型

以某型液氫液氧火箭發(fā)動機(jī)速度分級的兩級沖擊式氧渦輪為算例，如圖7所示.其中第2排和第4排是動葉，其余是靜葉.

圖7 渦輪葉型示意圖

在計算過程中認(rèn)為氣體常數(shù)和絕熱指數(shù)為常量，計算工況如表1所示.

表1 渦輪工作參數(shù)

根據(jù)渦輪的流道形狀建立了計算模型，每級渦輪葉片采用了7個S1流面和7個S2流面，如圖8所示.

圖8 算例渦輪流面

3.2 主要參數(shù)計算結(jié)果

表2列出了各級葉片流道在流面迭代中達(dá)到的最小誤差和相應(yīng)的迭代次數(shù).計算結(jié)果和實驗數(shù)據(jù)見表3，其中實驗數(shù)據(jù)由工程單位提供，出口壓力和溫度的實驗值是在表1工況下的測量結(jié)果，渦輪效率的實驗值則是多種工況下測量結(jié)果的范圍.

表2 流面迭代的收斂程度

表3 計算結(jié)果和實驗數(shù)據(jù)

3.3 馬赫數(shù)分布

圖9展示了各級流道中的馬赫數(shù)分布.從中可以看出氣流在第1級靜葉中從亞聲速變?yōu)槌曀伲诘?級動葉入口經(jīng)過一道激波降為亞聲速，之后在流道內(nèi)重新加速為至馬赫數(shù)較低的超聲速.在第2級靜葉入口存在一道弱激波，氣流經(jīng)過弱激波后降為亞聲速，并在之后的流動過程中始終保持亞聲速的狀態(tài).

圖9 各級流道中的馬赫數(shù)分布

3.4 誤差分析

從表3中可以看出計算得到的出口平均壓強(qiáng)比實驗數(shù)據(jù)高9.6%，出口總溫比實驗數(shù)據(jù)低2.2%，效率處于實驗數(shù)據(jù)的范圍之內(nèi).產(chǎn)生誤差的原因主要有以下幾點:

首先，流面迭代的最終計算結(jié)果是兩組流面上結(jié)果的平均值，所以流面迭代的收斂情況會對計算精度產(chǎn)生較大的影響.

其次，在計算過程中采用的控制方程并未包含黏性項，黏性損失是通過總壓恢復(fù)系數(shù)加入計算過程的.而且在計算激波后的流場時為了考慮附面層的影響引入了流量系數(shù).這在一定程度上影響到了流動參數(shù)和渦輪效率的計算精度.

最后，二次流等損失在計算時也沒有被考慮在內(nèi)，這也影響到了渦輪效率的計算.

4 結(jié)論

本文采用流面迭代的方法實現(xiàn)了液體火箭發(fā)動機(jī)渦輪泵多級渦輪內(nèi)跨聲速流場的全三維計算，得到如下結(jié)論:

1)計算結(jié)果較好地描述了多級跨聲速渦輪內(nèi)部的流動情況，主要性能參數(shù)的計算值與實驗值基本吻合，可為進(jìn)一步改進(jìn)渦輪設(shè)計、提高渦輪性能提供理論依據(jù);

2)采用混合平面法處理動靜葉間的非定常干涉，能夠解決多級葉輪機(jī)械葉柵進(jìn)出口邊界的參數(shù)傳遞問題;通過引入簡化的激波模型和由黏性引起的總壓恢復(fù)系數(shù)，使相對流面法可以用來近似計算渦輪內(nèi)部的實際流動.

采用無奇點的速度梯度方程可以正確求解跨聲速流場;通過流量以及出口參數(shù)的對比可以較為準(zhǔn)確地判斷多級跨聲速渦輪內(nèi)不同區(qū)域的流動狀態(tài).

References)

[1]Denton J D.An improved time-marching method for turbomachinery flow calculation[J].Journal of Engineering for Power，1983，105:514－521

[2]Denton J D.Calculation of three dimensional flow through any type of turbomachine[R].AGARD，Lecture Series 140，1985

[3]Hah C A.A navier-stokes analysis of 3d turbulent flows inside turbine blade rows[J].Journal of Engeering for Power，1984(106):421－429

[4]Dawes W N.Development of a three-dimensional navier-stokes solver for application to all types of turbomachinery[R].ASME.Paper No.88-GT-70，1988

[5]Wu Chunghua.A general theory of three-dimensional flow in subsonic and supersonic turbomachines of axial-，radial-，and mixedflow types[R].NACA，TN-2604，1952

[6]Boyer K M，O'Brien W F.Application of an improved streamline curvature approach to a modern，two-stage transonic fan:comparison with data and cfd[R].ASME Turbo Expo 2002.GT2002-30383，2002

[7]華鑫，喬渭陽，盧蕊.基于流線曲率法的航空軸流渦輪損失模型研究[J].機(jī)械設(shè)計與制造，2005(12):12－14

Hua Xin，Qiao Weiyang，Lu Rui，et al.A study on different loss models in axial-flow turbines based on the methed of streamlinecurvature[J].Machinery Design＆Manufacture，2005(12):12－14(in Chinese)

[8]Denton J D，Dawes W N.Computational fluid dynamics for turbomachinery design[J].Journal of Mechanical Engineering Science，1998，213(107):107－124

[9]汪慶桓，吳文權(quán)，朱根興，等.軸流式壓氣機(jī)兩類流面迭代三元完全解[J].工程熱物理學(xué)報，1983，4(4):336－346

Wang Qinghuan，Wu Wenquan，Zhu Genxing，et al.Solution of three-dimensional compressible flow in an axial-flow compressor stator by the use of s1 and s2 relative surfaces[J].Journal of Engineering Thermophysics，1983，4(4):336－346(in Chinese)

[10]馬冬英，梁國柱.基于S1，S2流面理論的液體火箭發(fā)動機(jī)渦輪內(nèi)部流場計算[J].推進(jìn)技術(shù)，2008，29(4):401－406

Ma Dongying，Liang Guozhu.Calculation of the internal now for turbine of liquid rocket engine based on the theory of S1，S2 stream surface[J].JournalofPropulsion Technology，2008，29(4):401－406(in Chinese)

[11]Bloch G S.Flow losses in supersonic compressor cascades[D].Blacksburg，Virgina:Virginia Polytechnic Institute and State U-niversity，1996