固液捆綁火箭上升過程噴流干擾特性及底部熱環境研究

嚴立，孫培杰，趙一霖，包軼穎，秦川，江澤敏，劉寬

（1.上海宇航系統工程研究所，上海 201109；2.上海航天技術研究院，上海 201109）

0 引言

固體推進系統和液體推進系統的動力組合是國外主要采用的動力組合方式，該方式擁有固體推力大、加速快，液體工作時間長、推力調節容易等技術優點，可實現運載火箭的系列化、模塊化、通用化和組合化生產制造。目前，國外捆綁運載火箭中固體助推器約占81%，大型、重型運載火箭幾乎都將固體火箭發動機作為助推級的首選動力［1］，如美國的德爾它系列（除重型德爾它4）、大力神系列、宇宙神系列以及航天飛機［2-4］，西歐的阿里安系列3 和5火箭，日本的M 系列、N 系列及H 系列火箭等［5］都是如此。

在固液捆綁火箭的發展歷程中，底部熱環境預示及防護一直作為關乎成敗的問題受到關注。固體火箭發動機排出熾熱高速氣流含有的Al2O3固體顆粒，除了大大增加火箭底部熱輻射外，還影響火箭底部氣流的流場分布，在助推與主發動機火焰間造成排氣交叉效應和附加干擾，加劇固體發動機羽流和液體芯級羽流間的相互作用。

美國在研制德爾它Ⅰ［6］和大力神Ⅲ［7-8］的過程中，開展了底部熱環境地面風洞試驗研究，探究底部不同區域加熱熱流呈現的變化趨勢。荷蘭的REIJASSE等［9］參照大力神III 的地面試驗方法，對阿里安5 芯級及固體助推3 個噴管狀態底部噴流流態進行了縮比風洞噴流的可視化試驗，詳細分析了不同噴流壓力與外界壓力比情況下，噴流與噴管噴流之間、來流與噴流之間的干擾情況，直觀了解噴流狀態變化的機理。文獻［10］采用Fluent 軟件，對H-IIA 固液捆綁火箭發動機底部加熱問題進行了數值模擬研究，分析了芯級液體發動機與4 枚固體助推器噴流相互干擾、回流現象。文獻［11］對雙發動機噴流干擾進行了實驗研究，拍攝了噴流干擾紋影圖片，測量了2 股噴流之間中心軸線的皮托壓力，研究了高空噴流干擾所引起的回流現象。

我國在固液捆綁火箭研究方面起步較晚，目前主要針對單一固體發動機噴流及熱環境開展了研究。文獻［12-16］采用Fluent 軟件對單一固體發動機的兩相噴流流場進行了數值仿真，其中：烏岳［12］分析了不同飛行馬赫數、飛行高度和不同顆粒相粒徑對流場的影響；黃顆［13］分析了后燃反應模型和顆粒相對流場的影響；武利敏［14］對顆粒軌道模型和擬流體模型進行了分析。文獻［17-18］采用數值仿真和試驗研究的方法研究了固體火箭發動機噴流反流對底部的影響，楊學軍［17］利用線性化熱學參數單一介質數值模擬方法分析了真實飛行軌道條件下的熱環境預示結果，并與地面試驗及飛行試驗測量結果進行了比對；李國良［18］對冷噴與熱噴方法進行對比分析，結果表明2 組分熱噴方法獲得的數據與試驗值一致性較好，較冷噴方法精度提升。文獻［19-20］分析了捆綁火箭底部噴流及熱環境，史亞男［19］數值模擬的方法分析了固體芯級捆綁固體助推火箭，不同助推器個數、不同助推器與芯級間距離對底部熱環境的影響；戴欣怡［20］采用數值模擬方法對助推器內部點火、燃燒和燃面推移等瞬態過程以及助推器尾部超聲速氣-固兩相流場分布規律進行了研究。綜上所述，至今未見含芯級液體火箭捆綁多枚固體助推器發動機底部噴流交互作用及底部熱環境研究的分析報道。

本文用連續相模型模擬氣相、離散相模型（Discrete Phase Model，DPM）模擬固體粒子相、離散坐標模型（Discrete-Ordinates Model，DOM）模擬含有固體粒子的介質輻射，對固液捆綁火箭在上升飛行過程中的氣固兩相噴流的演變過程進行了仿真研究，分析了不同飛行高度下固體發動機與液體發動機交會噴流干擾特性、高溫固體顆粒影響作用以及箭體底部表面熱環境，并將仿真結果與飛行結果進行了比較。研究成果對固液捆綁火箭底部熱防護設計具有重要指導意義。

1 仿真計算模型

1.1 數學計算模型

分析過程考慮固體粒子對噴流流場的影響，模擬過程分為3 個步驟：①氣相噴流場模擬；②加入固相粒子DPM 模型，計算氣固兩相噴流；③加入熱輻射模型，計算底部熱流。

1.1.1 連續相模型

連續性方程為

式中：ρ為密度；xi為i方向矢量；ui為i方向的速度。

動量方程為

式中：P為靜壓；μ為動力黏度；xj為j方向矢量。

能量方程為

式中：xk為k方向矢量。

采用標準k-ε方程模型，控制方程如下：

湍動能k方程為

湍動耗散率ε方程為

式中：K為湍流動能；ε為湍流能量耗散率，常數系數σμ=0.09；c1=1.44，c2=1.92。

1.1.2 離散相模型

采用拉格朗日離散相模型考慮顆粒相與連續相之間的相互作用力。首先是作用在顆粒物上的兩相間的拖拽力，表達式為

其中

式中：FD為單位質量阻尼系數；CD為阻力系數；u為流體相對速度；up為顆粒速度；μ為流體動力黏度；ρp為顆粒密度；dp為顆粒直徑；Re為相對雷諾數；a1、a2、a3為常數。

1.1.3 離散坐標方法

DOM 計算含有固體顆粒的噴流介質輻射作用。離散坐標法基于對輻射強度的方向變化進行離散，通過求解覆蓋整個4π 空間立體角上一系列離散方向上的輻射傳遞方程而得到問題的解。它可以計算所有光學厚度的輻射問題，并且計算范圍涵蓋了從表面輻射、半透明介質輻射到燃燒問題中出現的介入輻射在內的各種輻射問題。基于譜帶模型的吸收、發射、散射性介質內輻射傳遞方程的表達式為

在三維直角坐標系(x，y，z)下，采用離散坐標法，上式右端積分項近似由一數值積分代替，并在離散的方向上對輻射傳遞方程求解如下：

式中：ξm、ηm、μm為輻射傳輸方向的方向余弦；ωl為積分系數；l、m為空間方向離散的第l個和第m個立體角（l，m=1，2，…，NΩ），NΩ為4π 空間方向離散的立體角總數；=Φk(Ωm，Ωl)為離散后的輻射相函數。

介質的吸收系數采用基于網格的加權求和模型，散射相函數選擇各項同性，散射吸收選擇常數。

1.2 物理計算模型

某固液捆綁火箭芯級采用2 臺液氧/煤油發動機，捆綁4 個固體助推器，如圖1 所示。由于箭體外部結構輪廓和發動機噴管布置位置的對稱性，仿真模型采用箭體周向1/4 進行建模，計算區域包含1/4芯級以及1 個助推器。

圖1 固液捆綁火箭幾何模型Fig.1 Geometrical model of solid-liquid bound rocket

模型網格采用結構化六面體網格，網格數量1.2×107。綜合考慮計算資源和計算精度，針對不同的計算區域采用不同大小的網格，對發動機推力室內部及噴管出口附近進行網格加密處理。計算網格如圖2 所示。

圖2 火箭模型計算區域及網格Fig.2 Rocket model calculation area and grid

固體火箭發動機采用鋁丁羥復合HTPB 推進劑（Hydroxyl-Terminated Polybutadiene，HTPB），鋁粉質量分數17%，燃燒產物中Al2O3顆粒質量分數約32%，固體顆粒尺寸分布對固體發動機燃燒室效率和噴管效率影響較大，是發動機性能分析及仿真分析的重要輸入參數。但推進劑燃燒過程復雜，影響因素眾多，顆粒項直徑分布在幾微米到幾百微米之間，目前還沒有成熟的理論來對顆粒項直徑分布定量描述。為了綜合考慮不同直徑粒子對噴流的影響，采用Rosin-Rammler 分布來描述高溫顆粒分布，其函數關系為

式中：Yd為比指定粒徑d大的顆粒的質量分數；dm為平均粒徑；n為傳播系數。

顆粒平均直徑dm采用質量平均直徑［21-22］，公式為

中國農業大學水利土木工程學院黨委書記楊培嶺以《節水灌溉技術的未來發展方向和趨勢》為題進行了精彩演講，他呼吁要深入基礎理論研究，加快節水灌溉科研成果的轉化，實現節水灌溉技術的創新。要推廣自動化控制系統，加強節水灌溉設備質量的監管控制，加強水資源管理，合理確定水價，建立健全節水灌溉體系服務。

式中：D為噴管喉部直徑；Cm為100 g 推進劑中鋁粉的摩爾質量；Pm為燃燒室壓強；τ為顆粒在燃燒室中的駐留時間。

邊界條件設置如下：

1）入口邊界，設定壓力入口：固體發動機燃燒室總壓6.3 MPa，喉部壓力3.6 MPa，總溫3 236 K；液體發動機燃燒室總壓17.7 MPa，喉部壓力10.2 MPa，總溫3 810 K。

2）出口邊界，在噴管的遠場處設定為壓力出口：地面常壓環境設置為1 atm，溫度300 K；高空設置為當地的環境壓力和溫度。

3）壁面邊界，即設定壁溫邊界條件。

4）對稱邊界，即設定對稱邊界條件。

本文選取火箭飛行過程中6 個不同飛行高度下的工況進行分析，各工況來流參數見表1。

表1 火箭飛行計算工況來流參數Tab.1 Parameters of the coming flow under different calculation conditions

2 計算結果與分析

2.1 流場分析

提取2 個芯級軸線和芯級固發軸線截面上的參數，在不同飛行高度上外流場速度分布云如圖3所示。由圖可知，飛行高度10 km 以下時，火箭飛行馬赫數較小，環境壓力較高，發動機噴流擴張受到抑制，形成桶鼓狀的激波，即膨脹波和壓縮波交替的激波形式，發動機噴流集中在中心區域，芯級和助推噴流之間相互干擾作用較小，如圖3（a）和圖3（b）所示。隨著飛行高度增高和速度增大，箭體逐漸進入超音速飛行段，外界氣壓降低，從工況3 開始，環境壓力小于發動機出口壓力，箭體周圍逐漸產生不同類型的激波。箭體頂端形成弓形激波（Bow Shock），噴管噴流外緣與自由來流相互干擾，因噴流羽流對自由來的流阻礙形成羽流干擾激波（Plume Induced Shock），如圖3（c）所示。且隨飛行高度增加，壓力比逐漸增大，噴管噴流擴張角逐漸擴大，自由來流與發動機羽流干擾效應逐漸加強，產生噴流交互作用，噴流向底部回流，對底部形成熱、力沖刷作用，如圖3（d）～圖3（f）所示。

圖3 不同工況下火箭飛行外流場速度分布云Fig.3 Contours of the velocity distribution under different working conditions

以上噴流流場分析可知，在低空狀態下，發動機之間無直接噴流交互作用，而在高空狀態下，噴流擴張，發生空氣來流、芯級噴流和固發噴流之間復雜的交互現象，在交互過程中，產生噴流回流，造成對底部的力、熱影響。噴流干擾回流高溫區域如圖4 所示。

圖4 底部噴流交互作用高溫區分布示意圖Fig.4 Distribution diagram of high-temperature areas caused by the bottom jet interaction

工況2、工況3 和工況5 的噴流溫度分布如圖5～圖7 所示。由圖5 可知，工況2 時，由于外界氣壓較高，各發動機噴流集中在中心區域，芯級和助推噴流之間相互干擾較小，高溫區出現在發動機的出口下方。在工況3 時，發動機之間有了明顯的噴流交互作用，噴流干擾引起了局部高溫區，如圖6 所示，首先是芯級2 個液發噴流之間發生交互作用，在5區出現高溫區。其次由于2 固發與2 液發共4 個發動機噴流之間的交互作用，在1 區和3 區也出現高溫區，最高溫達2 200 K，高于同截面噴管中心噴流溫度。在工況5 時，隨著噴流擴張角增大，噴流干擾高溫區域向底部移動，在2 區和4 區出現明顯的噴流交互高溫區，該高溫區由2 固發與1 液發共3 個發動機噴流之間的交互作用產生，如圖7 所示。

圖5 工況2 噴流溫度分布Fig.5 Temperature distribution of the jet flow under working condition 2

圖6 工況3 噴流溫度分布Fig.6 Temperature distribution of the jet flow under working condition 3

圖7 工況5 噴流溫度分布Fig.7 Temperature distribution of the jet flow under working condition 5

從以上分析可以看出，噴流交互作用受發動機擴張角、大氣環境、飛行速度和噴管相對位置等因素影響，發生交互作用的形態和時刻都有所不同。最高溫區在1、3 和5 區，次高溫度在2、4 區，且隨著噴流的擴張，高溫交互出現的軸向方向距離箭體底部越來越近。

2.2 固體顆粒影響分析

分別對小粒徑（3 μm）顆粒和大粒徑（100 μm）顆粒對固液捆綁火箭氣固兩相流流場的影響進行了數值仿真分析，數值仿真結果如圖8 所示，由圖可知，3 μm 直徑的顆粒由于粒徑較小，隨流性較好，受氣體膨脹波影響較大，分布在噴流與空氣的邊界混合層區域；而100 μm 的大粒徑，隨流性較差，受氣體膨脹波影響較小，噴出噴管后，呈自由膨脹狀態流動，對交互流場核心區的影響較大。但從2 種直徑顆粒的交互流場來看，都不存在交互流場固體顆粒反流的情況，這與交互流場的氣體反流的狀態有本質區別。

圖8 不同粒徑的兩相流流場分布Fig.8 Flow distribution of the two-phase flow with different particle sizes

2 種直徑顆粒的溫度場和速度場對比結果如圖9 所示，2 個噴流流場的膨脹波外形尺寸相當，邊界混合層的參數分布相似，固體顆粒大小對噴流流場的膨脹程度影響不大，但核心區溫度與速度參數分布存在差異，100 μm 顆粒核心區的溫度衰減滯后于3μm顆粒，且速度低于3μm顆粒，主要是由于100μm 顆粒的流場中，顆粒呈自由膨脹狀態流動，4個固體助推器噴出的固體顆粒在噴流流場的中后部產生交互作用，對核心區的流場影響較大。

圖9 不同粒徑的溫度場及速度場對比Fig.9 Comparison of the temperature field and velocity field with different particle sizes

固液捆綁火箭底部從芯級底部的防熱板中心到流場計算域邊界的軸線上的參數分布如圖10 所示。由壓強變化曲線可知，3 μm 顆粒的流場參數與100 μm 流場參數吻合較好，這與之前二維流場分析所獲得的結果一致，固體顆粒對流場的壓力影響可以忽略；由溫度變化曲線可知，3 μm 流場的軸線溫度與100 μm 流場的軸線溫度變化趨勢一致，但在軸向距離40 m 后，100 μm 顆粒流場的溫度高于3 μm顆粒流場，主要是由于4 個固體助推器的固體顆粒在流場中心處產生交互作用，使溫度進一步升高；由速度以及馬赫數軸向曲線可知，在軸向距離40 m后，100 μm 顆粒流場的速度低于3 μm 顆粒流場，與溫度場變化規律一致，這是由于固體顆粒交互作用對流場產生減速作用。

圖10 流場沿軸線的參數分布Fig.10 Parameter distribution along the axis in the flow field

2.3 底部熱環境分析

火箭底部熱流分布云如圖11～14 所示。由圖可知，底部熱流分布與發動機噴流狀態相關，底部對流熱流密度最大值最早出現在2 芯級發動機噴管之間，隨著飛行高度的增加，回流沖擊位置向徑向移動，一直到工況5 時，移動到底板邊緣，并對火箭底部側壁產生熱沖刷作用。回流沖擊底板熱對流大小受回流溫度、速度和密度等影響，在工況4 時出現峰值，最大對流熱流密度為90 kW/m2。噴流對底部的輻射作用與表面和噴流之間的相對可視距離密切相關，可以看出芯級發動機與助推發動機軸線之間的輻射熱流較大，輻射熱流變化范圍為160～400 kW/m2。

圖11 工況2 時底部熱流分布Fig.11 Heat flux distribution at the bottom under working condition 2

在火箭飛行試驗中，底部防熱板上布置了3 個熱流測點，如圖15所示。

圖15 底部測流點計算位置Fig.15 Locations for calculation

圖12 工況3 時底部熱流分布Fig.12 Heat flux distribution at the bottom under working condition 3

圖13 工況4 時底部熱流分布Fig.13 Heat flux distribution at the bottom under working condition 4

圖14 工況5 時底部熱流分布Fig.14 Heat flux distribution at the bottom under working condition 5

測流點的對流、輻射與總熱流值如圖16 所示。與測量值的比較如圖17 所示。測流點位置熱流在起飛10 s 達到最大值，最大值分別410 kW/m2和472 kW/m2，平均熱流為191 kW/m2和222.7 kW/m2。仿真的熱流最大值為434.1 kW/m2和452.4 kW/m2，平均熱流為334 kW/m2和305 kW/m2。由圖可知，仿真最大熱流與測量值吻合較好，平均值是測量值的1.5 倍左右，變化趨勢上略有不同。測量值在10 s 內快速上升，然后緩慢下降，在80～100 s 內增大，最后降低；仿真值是起飛后一直緩慢上升，到90 s 左右增大到最大值，然后開始下降。分析仿真值與測量值趨勢不同的原因在于：1）仿真中未考慮起飛過程發射臺的遮擋作用，起飛過程中的熱流偏低；2）實際飛行過程中壁面溫度會一直上升，而仿真中未考慮流場與壁面的耦合溫升過程，熱流值一直增大。

圖16 測流點位置熱流值Fig.16 Calculated heat flux results

圖17 測流點熱流仿真值與計算值比較Fig.17 Simulated and calculated heat flux results

3 結束語

采用數值分析方法，對固液捆綁火箭飛行過程中發動機噴流相互干擾特性和箭體底部熱環境進行了分析，得到以下結論：

1）固液捆綁火箭飛行過程中，低空階段噴流交互作用不明顯；隨著噴流擴張，出現噴流交互作用，首先是2 芯級發動機出現較強交互作用，隨后在高空區域出現兩固體火箭發動機和芯級發動機噴流之間的交互作用，并隨著噴流擴張，高溫區域溫度升高，且向底部移動。

2）小直徑的固體顆粒隨流性較好，受氣體膨脹波影響較大，分布在噴流與空氣的邊界混合層區域；大粒徑的固體顆粒其隨流性較差，受氣體膨脹波影響較小，噴出噴管后，呈自由膨脹狀態流動，對交互流場核心區的影響較大。

3）噴流交互反流直接沖擊箭地底部，形成較大的加熱熱流密度，不同位置受羽流影響不同所體現的熱流密度也有較大差異，熱流密度最大點出現在4 個噴管（2 芯級和2 固體）噴流干擾的底部位置，最大對流熱流密度為90 kW/m2，底部輻射熱流變化范圍為160～400 kW/m2。

4）通過與飛行測量值比較，最大熱流仿真值與測量值吻合較好，其平均值約是測量值的1.5倍，變化趨勢上略有不同。

本文對固體捆綁火箭噴流流場的分析可以有效解釋底部不同位置熱流密度變化的機理，對箭體底部熱防護設計具有較強的指導意義。