固體火箭發(fā)動機(jī)噴管瞬態(tài)流場特性分析

徐學(xué)文，牟俊林，任建存，單 鑫

（海軍航空工程學(xué)院,山東 煙臺264001）

0 引言

固體火箭發(fā)動機(jī)噴管是將燃燒室內(nèi)高溫高壓燃?xì)馑械囊徊糠值哪芰哭D(zhuǎn)變?yōu)闅怏w動能，從而使飛行器獲得推力的機(jī)械裝置，在導(dǎo)彈飛行條件和發(fā)動機(jī)工作狀態(tài)一定的條件下，噴管的結(jié)構(gòu)參數(shù)直接決定噴管的工作狀態(tài)（欠膨脹狀態(tài)、完全膨脹狀態(tài)或過膨脹狀態(tài)），從而決定了發(fā)動機(jī)推力大小和推進(jìn)效率的高低。因此，開展固體火箭發(fā)動機(jī)噴管流場特性分析，使噴管更好、更有效地發(fā)揮其能量轉(zhuǎn)換作用，這對噴管工程應(yīng)用實(shí)踐和進(jìn)一步的優(yōu)化結(jié)構(gòu)設(shè)計(jì)具有指導(dǎo)意義。

當(dāng)前，噴管流場特性分析一般采用實(shí)驗(yàn)研究和計(jì)算機(jī)仿真手段，而實(shí)驗(yàn)研究大多采用冷場研究[1-2]，這與噴管實(shí)際工作中的流場特征還有很大差別，因此計(jì)算機(jī)仿真就成為主要的研究手段，特別是許多大型商用CFD軟件如FLUENT和PHOENICS等的出現(xiàn)。在計(jì)算機(jī)仿真中，國內(nèi)大多計(jì)算模型僅以噴管流場為研究對象，噴管進(jìn)口邊界條件也只由熱力學(xué)滯止?fàn)顟B(tài)參數(shù)（總溫T0和總壓p0）確定[3-5]，不考慮噴管實(shí)際工作過程中上游流場－燃燒室對噴管流場的影響。本文以某型固體火箭發(fā)動機(jī)噴管為研究對象，為準(zhǔn)確預(yù)測發(fā)動機(jī)噴管的流場特性，考慮發(fā)動機(jī)工作過程中燃燒室流場對噴管流場的影響，因此，建立了燃燒室-噴管一體化三維流場模型[6-9]，仿真計(jì)算發(fā)動機(jī)點(diǎn)火啟動過程中噴管內(nèi)瞬態(tài)流場特性。

1 數(shù)值計(jì)算模型

1.1 流場模型

所研究的固體火箭發(fā)動機(jī)噴管為拉瓦爾噴管，燃燒室形狀為五角星型，其頭部安裝有簍式點(diǎn)火器，燃燒室僅作為噴管流場進(jìn)口邊界條件的產(chǎn)生部件，在發(fā)動機(jī)點(diǎn)火啟動過程仿真計(jì)算中，將燃燒室與噴管作為一個整體，所建立的噴管-燃燒室一體化三維流場物理模型如圖1所示。模型中發(fā)動機(jī)燃燒室軸線為直角坐標(biāo)系的x軸，方向指向發(fā)動機(jī)噴管，中心原點(diǎn)位于發(fā)動機(jī)頭部軸線的起點(diǎn)上，由此表示流場橫截面的位置。

圖1 某固體火箭發(fā)動機(jī)燃燒室-噴管一體化流場模型Fig.1 The integrated flow field model of chamber and nozzle of the solid rocket motor

為便于計(jì)算，仿真中不考慮流場和固相之間的耦合，燃燒產(chǎn)物遵循完全氣體狀態(tài)方程，燃燒室藥柱壁面按平行層規(guī)律燃燒，不考慮藥柱的侵蝕燃燒效應(yīng)和燃?xì)廨椛鋫鳠岬挠绊憽?/p>

1.2 計(jì)算方程

流場控制方程采用三維可壓粘性方程，湍流模型采用計(jì)算精度比較高的應(yīng)用比較廣泛的k-ε二方程模型，另外，為使方程組封閉，還要補(bǔ)充氣體狀態(tài)方程。

1.3 初始和邊界條件

1.3.1 初始條件

以固體火箭發(fā)動機(jī)點(diǎn)火器點(diǎn)火向燃燒室噴射燃?xì)獾臅r刻為流場仿真的計(jì)算起點(diǎn)，此時，燃燒室內(nèi)和噴管內(nèi)氣體的壓力、溫度為周圍大氣的壓力和溫度。

1.3.2 流場壁面邊界條件

噴管壁面和燃燒室內(nèi)在被點(diǎn)燃之前的藥柱面為流場的壁面邊界，由于在近壁區(qū)域，燃?xì)獾牧鲃忧闆r受壁面影響比較大，燃?xì)庹承粤屯牧髑袘?yīng)力在每個子層所主導(dǎo)的地位不一樣，求解規(guī)律也不一樣，本文采用壁面函數(shù)法[10]來求解近壁區(qū)域物理量。

1.3.3 質(zhì)量進(jìn)口邊界條件

把燃燒室頭部點(diǎn)火器的殼體表面作為燃燒室質(zhì)量進(jìn)口邊界條件來確定，并假設(shè)點(diǎn)火器的單位面積質(zhì)量流率沿整個殼體表面均勻分布，燃?xì)饬飨虼怪庇跉んw表面。當(dāng)燃燒室裝藥表面的某點(diǎn)達(dá)到著火溫度被點(diǎn)燃后，燃?xì)忾_始從該點(diǎn)處向燃燒室內(nèi)噴射，此處的壁面邊界條件變?yōu)橘|(zhì)量進(jìn)口邊界條件，給定入口單位面積的質(zhì)量流率，即：

式中：ρprop為推進(jìn)劑的密度；r為推進(jìn)劑的燃速。

1.3.4 固相傳熱邊界條件

在噴管殼體和燃燒室藥柱壁面上均按第三類傳熱邊界條件給出，即：固相壁面導(dǎo)熱量等于流場表面?zhèn)鳠崃俊?/p>

1.3.5 噴管出口邊界條件

噴管出口邊界條件按壓力出口邊界條件給出，當(dāng)噴管出口處燃?xì)饬魉贋閬喴羲贂r，給定噴管出口上的背壓力；當(dāng)噴管出口處燃?xì)饬魉龠_(dá)到超音速時，無須給定噴管出口上的任何參數(shù)，出口處的參數(shù)由內(nèi)部流場參數(shù)推得。

2 計(jì)算結(jié)果及分析

采用有限體積法來離散噴管-燃燒室一體化三維流場模型，應(yīng)用隱式的一階迎風(fēng)格式來離散控制方程，采用simple算法求解差分方程，仿真計(jì)算出發(fā)動機(jī)點(diǎn)火啟動過程中噴管內(nèi)燃?xì)饬鲌龇植肌?/p>

2.1 燃?xì)鈦喴羲倭鲃与A段

發(fā)動機(jī)點(diǎn)火啟動后，從燃燒室頭部的點(diǎn)火器噴出高溫高壓燃?xì)?，燃?xì)饣鹧驿h和壓力波開始向燃燒室四周及下游傳播。在傳播過程中，燃?xì)鈮毫Σㄤh明顯超前火焰鋒，燃燒室內(nèi)的空氣受到壓力波的擾動而壓力升高，而溫度還未升高。0.3 ms時燃?xì)饣鹧驿h到達(dá)燃燒室星根，開始加熱星根處推進(jìn)劑，藥柱壁面受熱面積逐漸擴(kuò)大。0.9 ms時在燃燒室軸向0.08～0.14 m處的星根壁面上某些點(diǎn)處推進(jìn)劑經(jīng)過熱量的積累達(dá)到著火溫度開始燃燒，燃燒室頭部壓力升高，燃?xì)鈮毫Σㄏ蛳掠蝹鞑ニ俣燃涌臁?/p>

1.5 ms時，燃?xì)鈮毫Σǖ竭_(dá)噴管入口處（0.774 m），如圖2所示，燃燒室內(nèi)氣體壓力從頭部到噴管呈單調(diào)下降趨勢，而逐漸接近外界氣壓；燃?xì)饣鹧驿h在燃燒室軸向傳播的距離為0.22 m，在火焰鋒處突然降溫；燃?xì)饣鹧驿h向前傳播過程中先增速，在0.12 m處達(dá)到最高速度，然后速度迅速下降；在火焰鋒前方，燃燒室內(nèi)氣體受火焰鋒加熱溫度的影響越來越低，根據(jù)速度公式，壓力波沿軸向傳播速度也逐漸降低，但下降幅度比較平緩；燃?xì)鈮毫?、溫度及速度的這種變化趨勢直至2.2 ms時燃?xì)鈮毫Σㄤh到達(dá)噴管出口。

圖2 1.5 ms時壓力沿軸向變化曲線Fig.2 Change curves of pressure along axial direction at 1.5 ms

2.7 ms時燃?xì)饣鹧驿h到達(dá)燃燒室中部，燃燒室內(nèi)燃?xì)鈮毫ι?，并沿軸向分布趨向平緩，在噴管內(nèi)，燃?xì)饬鲃铀俣仍龃螅詾閬喴羲?，在喉部之前燃?xì)饧铀伲诤聿狂R赫數(shù)達(dá)到最高值，但小于1，在喉部以后燃?xì)鉁p速，如圖3所示；燃?xì)鈮毫υ趪姽軆?nèi)沿軸向仍呈單調(diào)下降趨勢，噴管內(nèi)還未產(chǎn)生激波，這種變化過程一直持續(xù)到燃?xì)庠诤聿窟_(dá)到音速之前。

圖3 2.7 ms時軸向速度沿軸向變化曲線Fig.3 Change curves of axial velocity along axial direction at 2.7ms

2.2 燃?xì)廑杖鲃与A段

隨著燃燒室內(nèi)燃?xì)鈮毫υ絹碓礁?，噴管的進(jìn)出口壓力比也越來越高。3.4 ms時在噴管喉部燃?xì)獾牧鲃铀俣冗_(dá)到臨界速度，此時馬赫數(shù)Ma=1，但速度在喉部之后又迅速下降；此時噴管進(jìn)出口壓力比為2.165，燃?xì)鈮毫υ趪姽芎聿恐俺氏陆第厔荩诤聿恐髩毫ι陨陨咴仝呄蛴趪姽艹隹诒硥骸饨绱髿鈮骸?/p>

3.8 ms時，燃?xì)庠趪姽芎聿恐皝喴羲偌铀倭鲃?，在喉部達(dá)到音速，在喉部之后燃?xì)獬羲偌铀僖欢尉嚯x之后突然減速流動，如圖4(a)所示，這是說明在噴管擴(kuò)張段出現(xiàn)了一道激波，正是由于激波的存在，使得燃?xì)庠趪姽芰鲃又袎毫υ诩げㄖ耙恢碧幱谙陆第厔?，遇到激波后壓力上升，而后又趨向出口處外界氣壓，如圖4(b)所示，此時，噴管進(jìn)出口壓力比為9.19，噴管內(nèi)燃?xì)饬鲃映霈F(xiàn)壅塞。

圖4 3.8 ms時馬赫數(shù)和壓力沿軸向變化曲線Fig.4 Change curves of Mach number and pressure along axial direction at 3.8ms

同樣，根據(jù)噴管一維等熵流函數(shù)也可證明噴管內(nèi)此時出現(xiàn)了一道激波。已知噴管出口截面積A2和喉部截面積A*之比為：A2/A*=7.3 382；出口背壓為：p2=1.01e5 Pa；噴管喉部處燃?xì)獾臏箟毫椋簆o=1.06e6 Pa。根據(jù)連續(xù)方程，對于噴管內(nèi)等熵壅塞流動有：

根據(jù)喉部Ma=1，查等熵函數(shù)表[11]并計(jì)算可得：當(dāng)噴管出口截面上為亞音速時（即文特利管流）p2=1.057e6 Pa；當(dāng)噴管出口截面上為超音速時p2=1.245e4 Pa，此時出口截面上的氣體壓力均不等于出口背壓1.01e5 Pa，并且相差很大，因此在噴管喉部和出口截面之間必定有一道正激波存在。

2.3 燃?xì)獬羲倭鲃与A段

隨著燃燒室內(nèi)燃?xì)鈮毫Φ纳撸細(xì)庠诤聿恐蟪羲偌铀俚木嚯x不斷擴(kuò)大，激波在噴管內(nèi)不斷向外移動，燃?xì)鈮毫υ趪姽軆?nèi)仍呈現(xiàn)出先下降，然后突升，再下降至外界氣壓的趨勢。

4.4 ms時，燃?xì)馑俣仍趪姽軆?nèi)流動過程中持續(xù)增大，在出口處達(dá)到最大值，平均馬赫數(shù)為3.39，噴管內(nèi)無激波存在，如圖5(a)所示；燃?xì)鈮毫υ趪姽軆?nèi)呈單調(diào)下降趨勢，但愈靠近噴管出口處壓力變化愈平緩，如圖5(b)所示，噴管進(jìn)出口壓力比為15.47。此后，噴管內(nèi)燃?xì)馑俣群蛪毫Π催@種變化規(guī)律一直增加到穩(wěn)定狀態(tài)為止。

圖5 4.4 ms時馬赫數(shù)和壓力沿軸向變化曲線Fig.5 Change curves of Mach number and pressure along axial direction at 4.4 ms

根據(jù)噴管一維等熵流函數(shù)以及噴管A2和A*之比算得的噴管此時出口馬赫數(shù)為：Mat=3.57，計(jì)算結(jié)果與理論分析結(jié)果比較相符。

3 結(jié)論

1）在藥柱點(diǎn)燃初期，燃?xì)鈮毫Σㄤh先于火焰鋒到達(dá)噴管，噴管內(nèi)燃?xì)獬尸F(xiàn)亞音速流動，壓力沿軸向呈單調(diào)下降趨勢，噴管內(nèi)無激波產(chǎn)生。

2）隨著燃燒室壓力升高，噴管內(nèi)燃?xì)獬尸F(xiàn)壅塞流動，在擴(kuò)張段出現(xiàn)了一道激波，燃?xì)鈮毫υ趪姽軆?nèi)呈現(xiàn)出先下降，然后突升，再下降至外界氣壓的變化趨勢。

3）隨著時間推移，燃燒室內(nèi)燃?xì)鈮毫χ饾u達(dá)到穩(wěn)態(tài)，噴管內(nèi)燃?xì)庠诤聿恐蟪羲偌铀俚木嚯x不斷擴(kuò)大，激波逐漸移出噴管外，最終，在噴管擴(kuò)展段燃?xì)獬尸F(xiàn)超音速流動，壓力在噴管內(nèi)呈單調(diào)下降趨勢，噴管內(nèi)燃?xì)饬鲃禹槙场?/p>

[1]丁英濤,謝君堂,仲順安.微型噴管內(nèi)氣體流動的流量壅塞現(xiàn)象[J].北京理工學(xué)學(xué)報(bào),2006,26(12)：1086-1089.

[2]張修峰,楊立,宋智勇.噴管通道內(nèi)氣體流動特性的實(shí)驗(yàn)分析與數(shù)值模擬[J].機(jī)電工程技術(shù),2010,39(4)：25-28.

[3]劉君,郭健.雙斜噴管固體火箭發(fā)動機(jī)流動特性數(shù)值模擬[J].固體火箭技術(shù),2002,25(1)：8-9.

[4]黃宏艷,王強(qiáng).過膨脹狀態(tài)下軸對稱收-擴(kuò)噴管內(nèi)外流場計(jì)算及分析[J].航空動力學(xué)報(bào),2007,22(7)：1069-1073.

[5]羅靜,王強(qiáng),額日其太.軸對稱矢量噴管內(nèi)流特性的不同湍流模型計(jì)算[J].推進(jìn)技術(shù),2003,24(4)：326-329.

[6]趙汝巖,隋玉堂,周紅梅.具有潛入式噴管的翼柱形藥柱發(fā)動機(jī)火焰?zhèn)鞑ミ^程研究[J].固體火箭技術(shù),2010,33(6)：526-630.

[7]楊樂,余貞勇,何景軒.基于FLUENT的固體火箭發(fā)動機(jī)點(diǎn)火瞬態(tài)內(nèi)流場仿真影響因素分析[J].固體火箭技術(shù),2011,34(4)：474-477.

[8]WILLIAM A D,MICHAEL T H,ROBERT A F.Integrated 3-D simulation of solid propellant rockets,AIAA 2001-3949[R].USA：AIAA,2001.

[9]FIEDLER R,JIAO X,NAMAZXIFARD A.Coupled fluid-structure 3-D solid rocketmotorsimulations,AIAA2001-3549[R].USA：AIAA,2001.

[10]VERSTEEG H K,MALAASEKERAW.An introduction to computational fluid dynamics-The Finite Volume Method[M].Beijing：Beijing World Publishing Corpration.2000.

[11]蘇長蓀.高等工程熱力學(xué)[M].北京：高等教育出版社,1987.