燃?xì)饬骷訜崞髯杂缮淞髁鲌鼋Y(jié)構(gòu)數(shù)值分析

田 寧，曹知紅，姜一通，李文浩

（北京航天長征飛行器研究所，北京，100076）

0 引 言

燃?xì)饬鞯孛嬖囼?yàn)設(shè)備是高超聲速飛行器防熱方案考核不可替代的手段之一，有功率大、熱流高等特點(diǎn)。試驗(yàn)系統(tǒng)一般有自由射流試驗(yàn)系統(tǒng)和風(fēng)洞試驗(yàn)系統(tǒng) 2種，其核心設(shè)備燃?xì)饬骷訜崞鞯墓ぷ髟砼c液體火箭發(fā)動機(jī)相同，利用氧氣和煤油燃燒并通過拉瓦爾噴管加速產(chǎn)生高溫超聲速燃?xì)饬鲌觯瑢︼w行器模型進(jìn)行考核。自由射流式燃?xì)饬骷訜崞鞯奈擦髁鲌鍪且环N自由射流的射流流場，具有明顯的自由邊界，與大外場來流相比，射流流場有波系多、結(jié)構(gòu)復(fù)雜等特點(diǎn)，其噴管類型一般有錐形噴管和型面噴管兩類。

國外相關(guān)機(jī)構(gòu)對液體火箭發(fā)動機(jī)內(nèi)的燃燒進(jìn)行了建模和研究[1,2]，中國進(jìn)行的燃?xì)饬鲌鲅芯恐饕ㄡ槍鸺l(fā)動機(jī)壓力沖擊效應(yīng)的影響[3]、推力室參數(shù)[4]和燃燒室噴管內(nèi)流場特性研究[5～10]，超燃沖壓發(fā)動機(jī)流場研究[11～13]，以及航空發(fā)動機(jī)尾流流場研究[14,15]等。此類研究關(guān)注點(diǎn)大多在發(fā)動機(jī)內(nèi)部流場以及發(fā)動機(jī)性能直接相關(guān)的方向。國內(nèi)外對燃?xì)饬髟囼?yàn)系統(tǒng)關(guān)心的噴管出口流場結(jié)構(gòu)及燃?xì)饬鲌鰺岘h(huán)境的研究尚未深入展開。

本文利用燃燒反應(yīng)方程和CFD-FASTRAN相結(jié)合的方法進(jìn)行了燃?xì)饬髁鲌龅臄?shù)值模擬。由于燃?xì)鈿怏w成份復(fù)雜，熱物理性能參數(shù)基于化學(xué)反應(yīng)的理論模型進(jìn)行計(jì)算，沒有成熟的數(shù)值模擬氣體模型，因此首先利用燃燒平衡方程求得加熱器噴管內(nèi)部關(guān)鍵截面的氣體成分、物性參數(shù)等，然后利用數(shù)值模擬軟件計(jì)算外流場。分析錐形噴管和型面噴管兩種噴管結(jié)構(gòu)的外流場特點(diǎn)，模擬了燃?xì)饬鳝h(huán)境下端頭試驗(yàn)?zāi)Ｐ偷臒岘h(huán)境分布，并與試驗(yàn)進(jìn)行了對比。

1 計(jì)算模型

1.1 加熱器介紹

燃?xì)饬骷訜崞饔扇紵液统曀賴姽芙M成，采用氧氣和煤油作為燃料，在燃燒室內(nèi)燃燒，利用拉瓦爾噴管產(chǎn)生高溫超聲速燃?xì)猓瑢υ囼?yàn)?zāi)Ｐ瓦M(jìn)行防隔熱材料的篩選和部段考核試驗(yàn)。

錐形噴管和型面噴管加熱器分別是北 YA6804和YA6814燃?xì)饬骷訜崞鳌ｅF形噴管加熱器的燃燒室總壓為1.2 MPa和1.5 MPa兩種狀態(tài)，噴管出口馬赫數(shù)為2.4，燃?xì)饪倻貫?300 K左右。型面噴管加熱器的燃燒室總壓為5 MPa，噴管出口馬赫數(shù)為3.0，燃?xì)饪倻貫?650 K左右。

1.2 燃?xì)鈪?shù)

燃?xì)饬骷訜崞鞯娜剂蠟?號航空煤油和氧氣，余氧系數(shù)為0.8。煤油分子式為C7.15H14.2。

燃燒劑和氧化劑分子量為

式中 μ為分子量；c，h，o分別為分子數(shù)；下標(biāo)g為燃燒劑，o為氧化劑。

各元素質(zhì)量組成：

式中 C，H，O分別為元素質(zhì)量組成。

完全燃燒給定的單位質(zhì)量燃燒劑，所需氧化劑的理論需要量為

式中oχ為氧化劑質(zhì)量。

余氧系數(shù)為 α?xí)r燃燒單位質(zhì)量燃燒劑所需氧化劑消耗量為

給定壓力下，燃料產(chǎn)物的理論溫度和成分，通過以下方程組確定[16]。

燃料元素的物質(zhì)平衡方程：

式中 下標(biāo)t為燃料。

燃料燃燒產(chǎn)物的分壓平衡方程：

燃料及其燃燒產(chǎn)物的能量守恒方程：

式中 I為能焓；下標(biāo)k為燃燒產(chǎn)物。

式（1）～（11）適用于氣體沿噴管長度完全達(dá)到化學(xué)平衡和能量平衡的假設(shè)。噴管出口溫度及氣流成分通過等熵條件迭帶求得：

式中 S為熵；R為氣體常數(shù)。

經(jīng)過迭代求解后，各燃燒產(chǎn)物質(zhì)量成分如表 1所示。

表1 燃燒產(chǎn)物質(zhì)量成分Tab.1 Composition of Gas Flow

1.3 計(jì)算外形及條件

計(jì)算采用二維軸對稱外形，如圖1所示。錐形噴管空流場模型網(wǎng)格如圖1a所示，以加熱器喉道為流場計(jì)算入口邊界，計(jì)算外形范圍為噴管出口軸向距離約15倍噴管直徑，徑向半徑3倍噴管直徑。型面噴管由于噴管出口速度可認(rèn)為無徑向分量，因此以噴管出口作為入口邊界，其帶球錐模型的模擬外形如圖1b所示。

圖1 計(jì)算網(wǎng)格示意Fig.1 Grid of Simulation

續(xù)圖1

前處理采用CFD-GEOM劃分網(wǎng)格，計(jì)算網(wǎng)格節(jié)點(diǎn)總數(shù)為75 052，壁面網(wǎng)格加密至10-5m，模型壁面網(wǎng)格雷諾數(shù)小于10。求解器采用CFD-FASTRAN，通過基于密度的有限體積法來求解N-S方程組。氣體模型為混合氣體模型，其成分組成由燃燒方程在噴管出口條件下求得。計(jì)算模型分別采用層流模型和k-ε湍流模型，空間離散采用Roe_FDS通量差分格式，時間離散采用全隱式時間積分格式。

質(zhì)量守恒方程：

式中 ρ為混合氣體密度；ju為j方向的氣體速度。

動量守恒方程：

式中 p為氣體壓力；ijδ為克羅內(nèi)克函數(shù)；ijτ為粘性切應(yīng)力；k為比熱比。

能量守恒方程：

式中 h為焓值；jq為j方向的熱流密度。

2 空流場結(jié)構(gòu)分析

圖2為錐形噴管加熱器空流場的壓力等值線圖。由圖2可知，錐形噴管出口出現(xiàn)菱形波系。氣流在噴管出口首先產(chǎn)生錐形激波，激波相交后仍為激波，在自由邊界處經(jīng)過反射轉(zhuǎn)變?yōu)榕蛎洸ǎ蛎洸ㄏ嘟缓笤俳?jīng)過自由邊界反射變回激波，如此多次反復(fù)，形成了由激波-膨脹波組成的噴管出口菱形波系。壓力最高區(qū)域和最低區(qū)域分別出現(xiàn)在激波交點(diǎn)后和膨脹波交點(diǎn)后。

圖2 錐形噴管出口壓力等值線圖及實(shí)際流場Fig.2 Pressure Contour Line Chart of Conical Nozzle and Real Gas Feild

加熱器燃燒室總壓為1.8 MPa、1.4 MPa和1.1 MPa 3種狀態(tài)的出口軸線壓力分布如圖3所示。軸線上的壓力分布整體呈現(xiàn)交替膨脹壓縮的過程，并隨著軸向距離的增加，膨脹和壓縮的強(qiáng)度顯著衰減。由于錐形噴管的膨脹作用，中心氣流在噴管出口首先繼續(xù)膨脹，壓力減小，在激波相交點(diǎn)前達(dá)到最低，然后開始經(jīng)過激波，壓力急劇上升，在激波相交點(diǎn)后達(dá)到最大，然后再轉(zhuǎn)為膨脹。3種總壓條件下出口波系均為激波，可見錐形噴管的波系性質(zhì)是由噴管結(jié)構(gòu)決定的，狀態(tài)變化不會改變波的性質(zhì)。加熱器總壓的降低導(dǎo)致噴管出口壓力降低，菱形波軸向尺度減小。

圖3 噴管出口軸線壓力分布曲線Fig.3 Pressure Distribution Along the Axis of Nozzle

型面噴管出口氣流方向基本平行，因此其出口波系主要由于加熱器狀態(tài)變化導(dǎo)致的出口壓力變化引起。當(dāng)加熱器余氧系數(shù)變化時，燃?xì)獬煞謱l(fā)生變化，燃?xì)獾奈镄詤?shù)也會隨之改變，而燃?xì)獾慕^熱指數(shù)直接影響噴管出口的氣流壓力。根據(jù)加熱器設(shè)計(jì)指標(biāo)，余氧系數(shù)0.7對應(yīng)出口壓力約為92 kPa, 余氧系數(shù)0.8對應(yīng)出口壓力約為108 kPa。

噴管出口壓力為92 kPa時的壓力等值線如圖4所示。

圖4 略低壓力型面噴管出口壓力等值線Fig.4 Pressure Contour Line Chart of Contoured Nozzle in Lower Pressure

錐形噴管激波的形成，主要是由于噴管的錐度導(dǎo)致出口速度方向?yàn)閲姽芙嵌壬⑸洌毕虻某曀贇饬饔龅酵饨绲退龠吔绠a(chǎn)生反射，形成較強(qiáng)的斜激波，類似于超聲速氣流作用在斜板上的激波效果。型面噴管出口處噴口角度為水平，且內(nèi)部圓滑過渡避免產(chǎn)生過強(qiáng)波系，所以出口處激波或膨脹波的形成主要是由于壓差產(chǎn)生。由于噴口超聲速氣流靜壓與環(huán)境壓力差別較小，故激波或膨脹波強(qiáng)度較小。

噴管出口壓力為92 kPa時噴管出口壓力沿軸線分布如圖5所示。氣流在噴管出口形成菱形激波，由圖5可知，激波的強(qiáng)度比錐形噴管形成的要小很多，壓力經(jīng)過激波上升了20%，而錐形噴管上升了約3～4倍。可見出口壓力先略抬升，氣流邊界反射后下降，后面按照此規(guī)律周期性衰減。

圖5 型面噴管出口軸線壓力分布曲線Fig.5 Pressure Distribution Along the Axis of Contoured Nozzle

噴管出口壓力為108 kPa時的壓力分布如圖6所示。噴管出口形成弱膨脹波。膨脹波前氣流壓力略有升高，因此在膨脹波狀態(tài)下要提高熱環(huán)境應(yīng)該將模型盡量放在膨脹波交點(diǎn)之前。

圖6 略高壓力型面噴管出口壓力等值線Fig.6 Pressure Contour Line Chart of Contoured Nozzle in Higher Pressure

噴管出口壓力108 kPa時噴管出口壓力沿軸線分布如圖7所示。與圖5相比，2種狀態(tài)下出口壓力分布規(guī)律不同，這是由于一種是超聲速氣流出口靜壓低于環(huán)境壓力，出口產(chǎn)生壓縮；一種是出口靜壓高于環(huán)境壓力，出口產(chǎn)生膨脹。因此，針對試驗(yàn)需求，對不同狀態(tài)下尋找合適的熱環(huán)境點(diǎn)位置的規(guī)律也不同。

圖7 型面噴管出口軸線壓力分布曲線Fig.7 Pressure Distribution Along the Axis of Contoured Nozzle

3 模型熱環(huán)境分析及試驗(yàn)驗(yàn)證

球錐模型在燃?xì)饬骷訜崞髁鲌鲋旭R赫數(shù)等值線如圖8所示。選取加熱器出口壓力條件為92 kPa時，在噴管出口形成錐形壓縮波的狀態(tài)。

圖8 燃?xì)饬鲌鲋星蝈F模型馬赫等值線Fig.8 March Contour Line Chart of Ball Head Model in Gas Flow

燃?xì)馍淞髟诙祟^模型表面的波系如圖9所示，圖中黑點(diǎn)表示模型在熱流測試試驗(yàn)中的測點(diǎn)。模型頭部的弓形激波在遇到氣流自由邊界后反射，變?yōu)榕蛎洸ǎ蛎洸ń?jīng)過模型壁面的反射仍為膨脹波，在自由邊界處再度被反射為壓縮波，如此反復(fù)，在模型表面形成由激波和膨脹波交替變化的波系。

圖9 模型表面波系示意Fig.9 Wave System on the Surface of Model

球錐模型表面母線熱流密度分布的數(shù)值計(jì)算結(jié)果與試驗(yàn)結(jié)果對比如圖10所示。由圖10可知，模型駐點(diǎn)熱流與層流計(jì)算結(jié)果符合較好，錐身熱流與湍流計(jì)算結(jié)果符合較好。

圖10 模型表面熱流分布Fig.10 Heat Flux Distribution on the Surface of Mode

試驗(yàn)測得的錐面熱流分布規(guī)律與分析所得的波系分布規(guī)律符合，在壓縮波影響范圍內(nèi)，熱流密度逐漸升高，在膨脹波影響范圍內(nèi)，熱流密度開始降低。

湍流狀態(tài)下計(jì)算的錐身熱流密度分布更符合試驗(yàn)測試結(jié)果，模型表面的實(shí)際流場狀態(tài)在端頭附近為層流，之后在一定位置轉(zhuǎn)捩為湍流，數(shù)值模擬結(jié)果較好地驗(yàn)證了這種特點(diǎn)，模擬顯示的轉(zhuǎn)捩點(diǎn)在球頭與錐面的連接點(diǎn)附近。

4 結(jié) 論

通過數(shù)值模擬對氣氧煤油加熱器的燃?xì)饬鲌鲞M(jìn)行了分析，包括空流場的分析和端頭模型熱環(huán)境分析，得到以下結(jié)論：

a）錐形噴管和型面噴管的波系形成機(jī)理不同，造成波強(qiáng)度差別較大；加熱器狀態(tài)改變會影響波系角度和位置，特定的試驗(yàn)環(huán)境需根據(jù)特定狀態(tài)下的波系結(jié)構(gòu)進(jìn)行選擇；

b）燃?xì)饬魃淞髟谇蝈F模型駐點(diǎn)前形成脫體弓形激波，在模型錐身周圍形成由錐形膨脹波和錐形激波交替變化的波系，模擬顯示模型表面氣流存在層流向湍流的轉(zhuǎn)捩點(diǎn)；

c）經(jīng)過驗(yàn)證，化學(xué)平衡反應(yīng)方程與 CFD結(jié)合的數(shù)值模擬方法，可以實(shí)現(xiàn)試驗(yàn)環(huán)境的準(zhǔn)確預(yù)測，為高超聲速飛行器防隔熱試驗(yàn)的地面試驗(yàn)?zāi)M提供支持。