趙一霖，嚴 立，來霄毅，馬祿創(chuàng)，侯凌霄，葉哲霄

（上海宇航系統(tǒng)工程研究所，上海 201109）

0 引 言

隨著航天產(chǎn)業(yè)的快速發(fā)展，對運載火箭運載能力的要求越來越高，而采用多噴管火箭發(fā)動機是提高運載火箭運載能力的主要途徑之一。多噴管火箭發(fā)動機工作時會產(chǎn)生多股噴流相互干擾的現(xiàn)象，同時噴流的引射和相互干擾的共同作用形成回流，底部復(fù)雜的熱環(huán)境及流動特性與火箭的飛行參數(shù)、推進劑種類、發(fā)動機出口的燃氣參數(shù)、發(fā)動機噴管結(jié)構(gòu)參數(shù)及其排列方式等因素相關(guān)，多因素的影響對運載火箭底部結(jié)構(gòu)防熱設(shè)計帶來風(fēng)險。

目前國內(nèi)的現(xiàn)役常溫推進劑火箭發(fā)動機推力較小，其噴流溫度和速度也較小，同時發(fā)動機噴管尺寸較小，噴管之間的距離在500 mm 左右，這使得現(xiàn)役常溫推進劑的火箭底部最大熱流密度在200 kW/m2左右。現(xiàn)役新一代捆綁運載火箭助推器與芯級的間距較大，噴流的交互作用較弱，噴流輻射作用對底部加熱影響更大，火箭底部最大熱流密度在500 kW/m2左右。對比現(xiàn)役火箭，新一代四機并聯(lián)大推力液氧煤油運載火箭的特點為：發(fā)動機噴管出口的噴流溫度和速度高，噴流中CO2和H2O 的含量高，發(fā)動機噴管尺寸大且噴管的間距小于400 mm。這些特點導(dǎo)致噴流之間形成的相互干擾作用更強，噴流的輻射加熱和對流加熱的共同作用使得火箭底部熱環(huán)境更加惡劣。因此，開展對新一代四機并聯(lián)大推力液氧煤油發(fā)動機噴流特性及噴流干擾對底部熱環(huán)境影響的研究是必要的。

多噴管火箭發(fā)動機是由單噴管火箭發(fā)動機并聯(lián)組成的，研究單噴管噴流流場可為火箭底部熱環(huán)境的研究打下堅實的基礎(chǔ)。文獻［1?4］針對單噴管火箭發(fā)動機，主要研究了不考慮噴流的二次燃燒時，不同飛行馬赫數(shù)、飛行高度對尾流場的膨脹壓縮波系個數(shù)、膨脹半徑、核心區(qū)長度的影響。文獻［5?8］研究了尾焰流場的復(fù)燃下不同飛行條件對單噴管尾焰軸向和徑向溫度及組分的影響。文獻［9］采用LU 隱式算法數(shù)值模擬了有無噴流干擾對軸對稱超聲速底部流場結(jié)構(gòu)的影響。文獻［10］采用RANS/LES 混合方法研究了超聲速底部流場的流動機理，得出IDDES方法預(yù)測的下游渦結(jié)構(gòu)更加清晰。

單噴管火箭發(fā)動機的研究結(jié)果不能簡單地疊加之后運用于多噴管火箭發(fā)動機，因此文獻［11?12］通過試驗測量研究了單機和雙機氫氧發(fā)動機噴流的激波總壓力場，結(jié)果表明在相同的軸向位置，雙機干擾噴流的最高壓力是單機噴流的5 倍，噴流的干擾范圍隨軸向距離的增大而增大。文獻［13］數(shù)值仿真研究了雙機并聯(lián)氫氧火箭發(fā)動機尾焰流場，發(fā)現(xiàn)兩臺發(fā)動機尾焰流場的相互擠壓，形成撞擊激波，使內(nèi)區(qū)燃氣溫度、壓力升高，速度降低。文獻［14?15］對四噴管火箭發(fā)動機進行三維數(shù)值模擬，研究了不同飛行高度對尾流流場和火箭底部輻射熱流密度的影響，結(jié)果表明：隨飛行高度的增加，噴流擴張角變大，四噴管噴流之間形成干擾激波，并產(chǎn)生回流，在飛行高度為20 km左右，底板的中心位置出現(xiàn)最大熱流密度。

火箭底部防熱問題是火箭研制階段應(yīng)重點關(guān)注的技術(shù)問題。研制階段對火箭底部熱環(huán)境評估過低會導(dǎo)致防熱措施不足，降低火箭的可靠性；評估過高則使防熱結(jié)構(gòu)質(zhì)量增加，降低運載能力。因此研究火箭發(fā)動機噴流對箭體底部的加熱和箭體底部流場結(jié)構(gòu)有利于了解火箭底部熱環(huán)境。本文基于標準κ-ε湍流模型與離散坐標輻射（discrete ordinates，DO）模型，采用FLUENT 軟件對四機并聯(lián)火箭發(fā)動機整體建模計算，綜合比較分析了不同飛行高度下噴流流場及噴流干擾形成的回流對火箭底部熱流密度的影響，可為合理評估火箭底部熱環(huán)境，進而設(shè)計出合理的熱防護結(jié)構(gòu)提供理論依據(jù)。

1 物理模型及數(shù)值計算方法

1.1 物理模型及網(wǎng)格化分

火箭一級采用四個液氧煤油的液體發(fā)動機，火箭幾何模型如圖1 所示。模型網(wǎng)格采用ICEM 軟件完成，采用結(jié)構(gòu)化網(wǎng)格，針對不同的計算區(qū)域采用不同大小的網(wǎng)格，對發(fā)動機推力室及噴管出口附近進行網(wǎng)格加密處理，如圖2所示，共生成1 500萬六面體網(wǎng)格。發(fā)動機的參數(shù)如表1所示。

表1 發(fā)動機主要參數(shù)Tab.1 Parameters of rocket Engine

圖1 運載火箭幾何模型Fig.1 Geometry model of launch vehicle

圖2 模型網(wǎng)格和尾部網(wǎng)格劃分Fig.2 Model meshes and afterbody meshes

1.2 控制方程

1.2.1多組分N-S方程

基于液體運載火箭多組分燃氣滿足連續(xù)介質(zhì)假設(shè)和理想氣體狀態(tài)方程的特點，建立三維多組分可壓縮Navier-Stokes 方程，包括質(zhì)量、動量與能量輸運的控制方程，形式如下。

多組分輸運方程為

式中：t為時間；ρ為密度；ν為速度向量；l為多組分中的一組分；Yl為組分l的質(zhì)量分數(shù)；Rl為組分l在化學(xué)反應(yīng)后的凈生成率；Sl為自定義源項的離散相引起的生成率。

組分擴散通量Jl可表示為

式中：T為環(huán)境溫度；Dl，m為組份l的質(zhì)量耗散系數(shù)；DT，l為組份l的熱擴散系數(shù)。

在直角坐標系中，單一組分l的可壓縮N-S 方程的守恒形式為

式（3）～（6）中：U為流動變量；F、G、H為氣流通量矢量，F(xiàn)v、Gv、Hv為黏性通量矢量；K為熱傳導(dǎo)系數(shù)；T為環(huán)境溫度；p、ρ、e、τ、μ分別為壓力、密度、比動能、應(yīng)力、黏性系數(shù)；u、v、w分別為速度在x、y、z方向上的分量。

注：式（6）中，ux為u在x方向上的偏導(dǎo)。

1.2.2湍流模型

標準κ-ε湍流模型，近壁面為標準壁面函數(shù)，此湍流模型本身具有較好的穩(wěn)定性、經(jīng)濟性和比較高的計算精度。

湍流輸運方程為

式（7）～（9）中：k為湍流動能；ε為湍流耗散率；Gk為由平均速度梯度引起的湍動能k的產(chǎn)生項；Gb為由浮力引起的湍動能k的產(chǎn)生項；YM為可壓湍流中的脈動擴張項；σk、σε為與湍動能k和耗散率ε二者相對應(yīng)的Prandtl數(shù)；Sk、Sε為用戶定義的源項；C1ε、C2ε、C3ε、Cμ是經(jīng)驗系數(shù)，通常為常數(shù)。

1.2.3輻射傳輸模型

離散坐標輻射（DO）模型，DO 模型使用范圍最大的輻射模型，可以計算所有光學(xué)厚度的輻射問題。輻射模型中計算氣體的吸收系數(shù)時采用了FLUENT 軟件提供的灰氣體加權(quán)平均模型（WSGGM）。

輻射方程為

式中：r為位置向量；s為方向向量；a為氣體吸收系數(shù)；σs為氣體散射系數(shù)；n為氣體折射率；σ為黑體輻射常數(shù)；I（r，s）為位置r處的微元體在s方向上的輻射出射度；s′為散射方向向量；Φ（s，s′）為散射相函數(shù)；Ω′為空間立體角。

1.2 邊界條件

噴管喉部設(shè)置為壓力進口邊界，采用發(fā)動機穩(wěn)態(tài)工作時的總壓和總溫；箭體頭部和周圍外部進口邊界為壓力遠場邊界，計算考慮外界環(huán)境壓力和溫度的變化對噴流的影響；箭體尾部邊界為壓力出口邊界，邊界設(shè)置與飛行高度有關(guān)；箭體、噴管外壁和箭體底部為300 K定溫邊界。

根據(jù)海拔高度不同設(shè)置4 個工況，對底部熱環(huán)境進行仿真分析，工況設(shè)置如表2所示。

表2 自由來流參數(shù)Tab.2 Parameters of free flow

2 計算結(jié)果分析

2.1 不同飛行高度馬赫數(shù)分布分析

圖3 是Z=0 截面馬赫數(shù)等值線圖，飛行高度為7 km 時，在箭體頂端沒有形成激波，噴流之間幾乎無交互作用。隨著飛行時間的增加，火箭進入超音速飛行，箭體頭部形成弓形激波，噴流的擴張角變大，自由來流與噴流相互作用形成羽流激波、桶形激波和自由剪切層即膨脹邊界。在噴流內(nèi)部即噴流相互作用的區(qū)域，噴流之間形成干擾激波，從圖4 可以看出，速度會出現(xiàn)負值，這表明底部出現(xiàn)回流，回流的軸向速度隨飛行高度的增加而升高。速度為零的位置是噴流之間發(fā)生交匯并產(chǎn)生返流的位置，不同飛行高度返流點位置如表3所示。當飛行高度高于13 km，從X軸向速度線圖可以看出，開始回流點的軸向位置隨飛行高度的增加而減小。

表3 不同飛行高度下馬赫數(shù)為零的軸向位置Tab.3 The axial position where Mach number is zero at different flight altitudes

圖3 不同飛行高度馬赫數(shù)等值線圖Fig.3 Isograms of Mach number at different flight altitudes

圖4 不同飛行高度X軸向馬赫數(shù)分布Fig.4 Axial distribution curve of the Mach number at different flight altitudes

2.2 高低空底部流場對比分析

圖5（a）、（b）是飛行高度為7 km 時底部參數(shù)分布圖。從圖可以看出，在飛行高度為7 km 時，火箭飛行馬赫數(shù)為0.72，發(fā)動機噴流主要集中在噴管中心區(qū)域，噴流膨脹角較小，噴流之間因膨脹產(chǎn)生相互作用而形成一個壓力較高的區(qū)域，噴流的引射作用將空氣吸入底部區(qū)域，被吸入底部的自由流空氣向下游流動時不能通過噴流相互作用區(qū)，在發(fā)動機出口附近形成腎形渦對，增強了與噴流摻混。噴流相互作用區(qū)與底部的壓差較小，返流以亞音速沖擊底板中心，并在底板附近形成沿噴管間隙流動的壁面射流，此時壁面射流的馬赫數(shù)小于0.1。圖5（c）、（d）給出了飛行高度為21 km 時底部參數(shù)分布圖。在飛行高度為21 km 時，飛行馬赫數(shù)為1.89，環(huán)境壓力降低，噴流膨脹角增大，四噴管噴流的相互作用在四噴管出口附近形成一個高壓噴流相互作用區(qū)。噴管內(nèi)壁面的噴流邊界層隨噴流的膨脹波流動到噴流相互作用區(qū)，但低動能的噴流邊界層［16］與被吸入底部的自由流空氣不能通過高壓的噴流相互作用區(qū)，因此出現(xiàn)返流現(xiàn)象。噴流相互作用區(qū)與底部之間形成很大的壓差，返流逐漸加速到超音速，形成返流激波，最后在距離底板中心0.2 m 左右形成一個正激波，超音速返流穿過正激波后，導(dǎo)致返流降為亞音速，壓力和溫度急劇升高。相對于低空，高空的壁面射流的馬赫數(shù)增大到0.6 左右。從圖6（a）、（b）可以看出，低空時返流的溫度較低，高空時返流溫度較高，致使底板的熱環(huán)境更加惡劣。

圖5 底部截面馬赫數(shù)流線圖Fig.5 Mach number flow field at the bottom of rocket

圖6 溫度分布云圖Fig.6 Temperature contour

從圖7（a）、（b）可知低空時返流主要為空氣，高空時返流主要集中在軸線附近的區(qū)域，影響區(qū)域變大且高溫噴流氣體摩爾分數(shù)占比增大。從圖8 可以看出，在X=2.042 m 時X軸向速度為零，此位置是開始返流的位置，在X=0～2.042 m 速度會出現(xiàn)負值，這表明底部出現(xiàn)返流。從開始返流到底板，返流空氣摩爾分數(shù)從90%左右逐漸升高到99%左右，這說明返流氣體主要是被吸入底部的空氣，只有少量被卷吸的噴流氣體。返流溫度從900 K左右逐漸下降到357 K左右，如圖9 所示。飛行高度為21 km 時開始返流的位置為X=2.105 m。從開始返流到底板，返流空氣摩爾分數(shù)從19%左右逐漸升高到42%左右，如圖10所示，空氣摩爾分數(shù)逐漸升高是因為返流過程中不斷地與被吸入底部區(qū)域的空氣摻混，但返流氣體組分中噴流氣體占50%以上。返流溫度從3 200 K 左右波動下降到2 780 K左右，如圖11所示。

圖7 空氣摩爾分數(shù)流線圖Fig.7 Mole fraction of air

圖8 7 km X軸向馬赫數(shù)和空氣摩爾分數(shù)分布曲線Fig.8 Axial distribution curve of the Mach number and mole fraction of air at 7 km

圖10 21 km X軸向馬赫數(shù)和空氣摩爾分數(shù)分布曲線Fig.10 Axial distribution curve of the Mach number and mole fraction of air at 21 km

圖11 21 km X軸向壓力和溫度分布曲線Fig.11 Axial distribution curve of the pressure and temperature at 21 km

2.3 不同飛行高度底部熱流分布分析

圖12 為底部熱流密度分布云圖。從圖中可以看出，防熱底板的最大總熱流密度出現(xiàn)在中心位置，高熱流密度區(qū)主要分布在中心區(qū)域和靠近邊緣的各噴管之間狹長區(qū)域。這一是因為這兩區(qū)域受輻射的面積較大；二是因為高溫反流直接沖擊防熱底板后沿著各噴管之間狹長區(qū)域流出。噴管周圍區(qū)域的熱流密度都較小，主要是噴流作用相當于引射器，將自由流空氣吸入底部區(qū)域，起到冷卻和減弱熱輻射的作用。隨著飛行高度的升高，防熱底板熱流密度呈現(xiàn)先升高后降低的趨勢，如圖13 所示，最大總熱流密度出現(xiàn)在飛行高度21 km 處，總熱流密度達到748 kW/m2。飛行高度為30 km 到60 km 期間，噴流返流的速度和溫度增大幅度較小，但噴流返流的密度下降，導(dǎo)致總熱流密度呈下降趨勢。

圖12 熱流密度分布云圖Fig.12 Heat flux contour

圖13 防熱底板最大熱流密度隨時間變化圖Fig.13 Max.heat flux vs time on the bottom of the rocket

發(fā)動機的總熱流密度隨高度的增加呈現(xiàn)先升高后下降的趨勢，如圖14 所示，在飛行高度為7 km 時總熱流密度達到最大，為297 kW/m2。發(fā)動機大熱流密度區(qū)主要分布在發(fā)動機出口靠近箭體軸線的位置，這主要是因為四個發(fā)動機的高溫交匯區(qū)距離發(fā)動機出口較近，輻射作用對發(fā)動機的加熱占主要地位。尾段的熱流密度隨飛行高度的增大，從負熱流密度向正熱流密度變化，總熱流密度是先升高后降低，如圖15 所示。隨著飛行高度的繼續(xù)增大，尾段主要受輻射加熱的影響，尾段總熱流密度增大，總熱流密度在飛行高度42 km時達到最大，為18 kW/m2。

圖14 發(fā)動機最大熱流密度隨時間變化圖Fig.14 Max.heat flux vs time of the rocket engine

圖15 尾段最大熱流密度隨時間變化圖Fig.15 Max.heat flux vs time of the afterbody

3 結(jié)束語

本文通過對四機并聯(lián)火箭發(fā)動機噴流流場進行數(shù)值仿真，分析了不同飛行高度下噴流干擾形成的回流特性和箭體底部熱流分布。分析結(jié)果表明，四臺發(fā)動機噴流在不同飛行高度下，在發(fā)動機出口附近產(chǎn)生一個交互回流區(qū)。通過分析高低空底部流場發(fā)現(xiàn)：隨飛行高度的增加，返流向底板流動速度從亞音速加速到超音速，返流溫度也逐漸升高，防熱底板熱流密度呈現(xiàn)先升高后降低的趨勢，最大熱流密度點出現(xiàn)在飛行高度21 km 處。本文對流場的結(jié)構(gòu)、返流的溫度和組分等方面的研究，可以較為清晰地了解不同飛行高度的熱流密度變化原因，對箭體底部防熱區(qū)域細化，進而減輕防熱結(jié)構(gòu)重量具有較強的指導(dǎo)意義。