王金金, 查柏林,*, 張煒, 惠哲, 蘇慶東, 何齊
(1. 火箭軍工程大學導彈工程學院, 西安 710025; 2. 南昌航空大學測試與光電工程學院, 南昌 330063)
固體沖壓發動機具有結構簡單、比沖高、質量輕等優點[1],是現代火箭技術中應用非常廣泛的動力裝置。在固體沖壓發動機的研究過程中,補燃室綜合性能的研究對提高固體沖壓發動機性能具有重要意義。
為研究補燃室綜合性能,1989年,Cherng等[2]利用SIMPLE算法研究了不同進氣道角度和進氣道位置對補燃室內燃燒的影響,結果表明,通過優化圓頂高度和入口流動角,可以提高混合效率和燃燒效率,且仿真結果和試驗結果吻合。1995年,董巖等[3]采用二維k-ε模型對固體沖壓發動機補燃室二次燃燒進行了數值模擬,結果表明,經數值模擬設計的二次燃燒補燃室構型比普通補燃室構型燃燒效率有明顯提高。2001年,Stowe[4]對補燃室內的兩相流化學反應進行了數值模擬,結果表明,采用兩相化學燃燒時,燃燒效率誤差為16%,而純氣相燃燒模擬時的誤差會超過27%。2011年,劉杰等[5]將旋轉射流技術引入固體沖壓發動機設計中,開展了旋轉進氣和補燃室含硼一次燃氣三維反應流場仿真研究,研究發現存在一個最佳進氣角,可以使燃燒效率最大,且一次旋流數增加可以促進二次燃燒。2014年,胡旭等[6]對不同進氣道結構下的固體沖壓發動機補燃室二次燃燒性能進行了仿真研究,得到進氣結構對補燃室性能的影響規律,其中,對稱進氣結構有助于燃氣與空氣的摻混燃燒。2015年,王洪遠等[7]研究了空氣旋轉進氣對含硼固體沖壓發動機二次燃燒性能的影響,研究表明旋流進氣方式減小了硼顆粒點火時間,在旋流數為0.385時,點火時間最小。2016年,鞏倫昆等[8]研究了結構尺寸對固體燃料沖壓發動機燃速影響,結果表明,結構尺寸中突擴比對燃速起著關鍵的作用,燃氣燃燒速率與突擴比近似呈線性關系。2017年,王洪遠和唐田田[9]對含鋁顆粒固體燃料沖壓發動機燃燒速率特性進行了分析,仿真結果表明,在補燃室燃燒過程中,影響平均燃速的因素由強到弱依次為來流空氣質量流量、鋁顆粒含量、來流空氣總溫以及鋁顆粒直徑。2018年,李唯暄等[10]研究了旋流燃燒室構型對固體燃料沖壓發動機自持燃燒性能的影響,仿真以及實驗結果表明,在旋流工況下,相對臺階高度對火焰穩定以及燃燒特性有顯著影響,相對臺階高度的改變對特征速度與推力的作用不大,而在無旋工況下,特征速度和推力則隨相對臺階高度的增加而增加。
綜上所述,目前針對固體沖壓發動機補燃室相關的研究,主要集中于提高補燃室燃燒效率上,包括在補燃室結構設計[11-13]、硼顆粒點火燃燒[14-16]、推進劑配方[17]等,其中補燃室結構設計包括改變一次燃氣出口的角度和數目,空氣進氣道的壓力、射流的速度、方向、進氣角度[10,12,18],空氣與燃料的比值[19-20]等。但是,關于固體沖壓發動機振動環境[21]、熱防護[22]和燒蝕條件[23-24]等研究相對較少,而補燃室絕熱層耐燒蝕性能是影響補燃室綜合性能的重要指標之一。
因此,為研究不同進氣道結構對固沖發動機二次燃燒和內壁燒蝕環境的影響,本文采用標準k-ε湍流模型、單步渦耗散燃燒模型以及KING硼粒子點火燃燒模型,對目前導彈應用比較廣泛的雙下側90°進氣結構和雙側180°進氣結構沖壓發動機補燃室摻混燃燒進行數值模擬,對比分析兩種進氣道結構對補燃室燃氣摻混燃燒和內壁燒蝕的影響,為固體沖壓發動機結構設計提供參考。
兩種不同進氣道的補燃室結構簡圖及主要尺寸如圖1所示。Case 1為雙下側90°進氣結構,進氣道間的夾角為90°;Case 2為雙側180°進氣結構,進氣道間夾角為180°。補燃室結構尺寸參考文獻[6],一次燃氣由直徑為40 mm的圓形噴口沿軸線方向進入補燃室。
為簡化模型,對補燃室中的內流場作如下假設:

圖1 兩種不同進氣道結構的補燃室結構簡圖Fig.1 Structure diagram of secondary combustion chamber under two air-inlet structures
1) 補燃室內氣相射流為準定常流動,忽略燃氣與補燃室壁面的熱交換作用。
2) 補燃室燃氣為理想氣體,服從理想氣體狀態方程pV=nRT。
本文采用CFD數值模擬方法分析進氣道結構對補燃室摻混燃燒的影響。為降低計算難度,補燃室內的湍流流場和組分燃燒采用三維雷諾平均Navier-Stokes方程求解[24],即
(1)
式中:W為守恒向量;V為流體微元;n為法向方向;F(W)和G(W)為無黏和黏性通量;dS為矢量面積元;H為化學反應模型或體積力的源向量。
目前雷諾平均Navier-Stokes方程求解中,常用的湍流模型有k-ε湍流模型、k-ω湍流模型、雷諾應力模型等。由于k-ε湍流模型具有收斂穩定、精度適當、計算量小等優點,成為CFD分析中應用最廣的湍流模型。在計算精度要求相對較小時,k-ε湍流模型滿足大部分計算要求[23],故本文通過自編程CFD數值模擬方法,采用基于密度的k-ε湍流模型進行計算,并考慮輻射對流場的影響。在該模型中,湍流輸運方程可表示成湍流能量輸運方程和能量耗散方程:
Gb+Gk-ρε-Yk+Sk
(2)
Gε+Dε-Yε+Sε
(3)
式中:μ為流體動力黏度;μt為湍流黏性系數;uj為速度在j上的分量;k為湍動能;ε為湍動能耗散率;σk和σε為k和ε引起的正應力;ρ為流體密度;Gk為湍動能的產生項;Gε為湍動能耗散率ε的產生項;Dε為交叉擴散項;YM為可壓湍流中脈動擴張的貢獻;Yk和Yε分別為k和ε由于湍流產生的耗散項;Gb為浮力引起的湍流動能k的產生項;Sk和Sε為自定義的源項。
由于補燃室中流場流動比較復雜,為方便計算其中的化學反應過程,實現各組分凈反應速率的有效控制,采用單步渦耗散燃燒模型[18]。主要的氣相化學燃燒反應為
(4)
凝聚相顆粒點火燃燒模型主要有L-W模型和KING硼顆粒點火模型。目前由于KING模型的數學表達更容易表征一些關鍵的物理化學參數,因而在硼顆粒點火計算中被廣泛應用。而發動機實際工作過程中,高速射流對硼顆粒表面液態層產生氣動剝蝕效應,促進了硼顆粒的點火燃燒。結合文獻[6]中推導的氣動剝蝕計算,在原有KING模型的基礎上,加入高速射流對硼顆粒表面氧化層厚度的影響,硼顆粒氧化層厚度隨時間的變化關系為
(5)
式中:rB為硼顆粒的半徑;x為硼顆粒氧化層厚度;τ為黏性應力分量;μB2O3為液態B2O3黏性;ρB2O3為B2O3密度;MB2O3為B2O3摩爾質量;RB、RE和RH分別為硼消耗速率、B2O3蒸發速率和B2O3與水的反應速率;φ1和φ2為剝離角度。
硼的燃燒速率為
vB=4πrPDln(1+0.667ωO2)/t
(6)
式中:ωO2為顆粒周圍環境氣體中O2的質量分數;擴散系數D=1.5×10-5m2/s。
一次燃氣入口和空氣入口均采用質量流量入口,為簡化計算,一次燃氣入口由H2、H2O、CO、CO2和N2組成,各成分的質量百分比分別為10%、1%、47%、1%和41%,總質量流量為2.4 kg/s,硼顆粒含量為燃氣的35%,燃氣溫度為1 800 K;空氣入口的質量流量為4 kg/s、總溫為573 K,氧含量為23%。根據文獻[6],凝聚相粒子的粒徑設為5 μm;補燃室內壁采用絕熱無滑移壁面條件(見圖2);各組分質量梯度和壓力梯度為0;補燃室出口采用壓力出口邊界,壓強和溫度分別為1 atm(1 atm=1.01 kPa)和300 K。
噴管出口截面燃燒效率反映了補燃室結構、一次燃氣進氣、沖壓空氣進氣對燃燒的綜合影響。截面燃燒效率通過截面的組分燃燒完成率來表示。任意截面處的氣相組分二次燃燒效率為
(7)

氣相燃氣的總燃燒效率為
(8)
式中:QCO和QH2分別為CO和H2的燃燒熱值;ηH2和ηCO分別為H2和CO的燃燒效率;λ為CO的在氣相組分中占的質量百分比。
任意截面處硼顆粒的二次燃燒效率ηB為
(9)
(10)
式中:QB為硼的燃燒熱值;β為氣相組分占的質量比。
圖3為固體沖壓發動機補燃室的網格分布。根據固體沖壓發動機補燃室的結構特點,補燃室縱向采用雙重O型網格劃分。由于燃氣發生器出口射流速度組分相對復雜,對燃燒室出口、補燃室內壁等區域進行加密處理。
為驗證計算的準確性和可靠性,需對網格進行驗證。圖4給出了網格無關性驗證結果,當網格節點數大于120萬時,網格節點數對壁面平均氣膜有效度幾乎沒有影響,因此選用150萬網格計算。

圖3 網格劃分Fig.3 Mesh partition

圖4 網格無關性驗證Fig.4 Grid independency verification
圖5為補燃室內不同截面處的O2濃度分布,截面間相距100 mm。從圖中可以看出,兩種結構中,截面內O2濃度分布不均,呈軸對稱分布,富氧區緊靠進氣道一側。其中,Case 1富氧區呈“C”字形分布,至補燃室出口位置,富氧區O2濃度依然比較高,與貧氧區形成較大濃度差。Case 2中,進氣道出口位置富氧區O2濃度較高,隨著向補燃室下游發展,富氧區O2濃度逐漸降低,至補燃出口位置,富氧區O2濃度降低,貧氧區O2濃度增加。結果說明Case 2有利于空氣中的O2與一次燃氣的有效混合。

圖5 補燃室內O2濃度分布Fig.5 O2 concentration distribution in secondary combustion chamber

圖6 補燃室內溫度分布Fig.6 Temperature distribution in secondary combustion chamber
圖6為補燃室內溫度分布。從圖中可以看出,兩種結構補燃室頭部溫度約為2 200 K。整個補燃室中,O2濃度分布與溫度分布形成互補,呈現出高溫區O2濃度低,低溫區O2濃度高的特點。其中,Case 1中高溫區保持單側分布,最高溫度超過2 500 K;Case 2從空氣入口至補燃室出口,溫度逐步混合均勻,僅小部分區域溫度稍高,與Case 1相比較,溫度相對均勻。Case 2有利于補燃室溫度的擴散。
綜上,Case 1中,一次燃氣與空氣相互分離,O2濃度分布和溫度分布互補;Case 2中,燃氣的溫度分布和O2濃度分布逐漸趨于均勻。結果表明,與Case 1相比,Case 2更有利于空氣與一次燃氣的摻混,且有利于溫度的擴散。
圖7為兩種補燃室不同截面位置的燃燒效率對比。從圖中可以看出,至距補燃室頭部300 mm位置,兩種結構中的氣相組分燃燒效率超過90%,說明氣相化學反應主要集中在補燃室頭部位置。另外,從圖7還可以看出,至補燃室出口截面,Case 1的總燃燒效率為74%,Case 2的總燃燒效率超過95%,Case 2的燃燒效率高于Case 1。結果說明雙側180°進氣結構沖壓發動機補燃室中,一次燃氣與空氣得到了更加充分的摻混和燃燒,總燃燒效率更高。

圖7 補燃室不同截面燃氣成分燃燒效率Fig.7 Combustion efficiency of gas components with different cross sections
固體沖壓發動機補燃室內的流動存在頭部回流和軸向渦流。頭部回流是由一次燃燒軸向流動過程中形成的引射抽吸作用形成的,頭部回流在補燃室摻混燃燒中發揮重要作用。圖8為補燃室頭部對稱截面上的射流速度矢量圖。由圖可知,空氣射流和燃氣對擊后,一部分空氣在進氣道出口位置形成小渦旋,大部分空氣與中心一次燃氣沖擊后分流。其中一部分空氣形成逆向回流,回到補燃室頭部,形成較強的回流區漩渦,并與軸向運動的一次燃氣發生激烈的摩擦卷吸作用。對比圖8中Case 1和Case 2兩種補燃室頭部渦旋可以看出,Case 2補燃室頭部形成2個大的渦旋,空氣與一次燃氣反生劇烈摻混;Case 1空氣入口在補燃室單側形成較大沖擊,頭部回流相對較小,導致補燃室前段的燃燒效率偏低。
圖9為x=400 mm位置截面的射流速度矢量分布。空氣射流和燃氣對擊后,一部分空氣沿著斜向進氣道向補燃室下游發展。如圖9(a)所示,Case 1中,空氣射流單向壓制一次燃氣,空氣與燃氣間形成小旋渦,結合圖5和圖6可知,燃氣的燃燒面積相對較小,燃燒效率低。而Case 2結構中,空氣與一次燃氣撞擊后,沿著補燃室內壁向兩側運動,形成4個大的對稱漩渦(見圖9(b)),促進了燃氣與空氣的有效摻混。
圖10為補燃室中粒子的運動軌跡和速度分布。由圖可知,Case 1中,補燃室凝聚相粒子與一次燃氣運動軌跡一致,燃氣中的氧含量低,不利于凝聚相粒子的點火燃燒。Case2中,燃氣與空氣摻混的過程中,一次燃氣中顆粒相發生擴散,增強了燃氣與空氣的摻混,形成強烈的摻混擴散區域和擴散火焰峰面。在火焰面上,燃氣與空氣中的O2發生劇烈的化學反應,并釋放熱能,提高燃氣溫度,達到硼顆粒的點火條件,促進硼顆粒燃燒。隨著燃燒向下游發展,反應基本完成,O2質量濃度和燃氣溫度逐漸均勻并降低。對比Case 1和Case 2可以看出,雙側180°進氣結構可促進凝聚相粒子的擴散,實現凝聚相粒子與空氣混合燃燒。

圖8 補燃室頭部射流速度矢量分布Fig.8 Jet velocity vector distribution on head of secondary combustion chamber

圖9 截面400 mm位置的射流速度矢量分布Fig.9 Jet velocity vector distribution on cross section of 400 mm

圖10 凝聚相粒子運動軌跡和速度分布Fig.10 Motion path line of condensed phase particles and their velocity distribution
因此,Case 1中燃氣與空氣間的摻混漩渦小,凝聚相粒子沿高溫低氧區運動,燃氣無法充分地混合燃燒,燃燒效率低;Case 2有利于補燃室中燃氣與空氣間形成較大的摻混漩渦,且有利與凝聚相粒子與空氣混合,促進了燃氣的二次點火,總燃燒效率相對更高。
固體沖壓發動機補燃室燒蝕工況是影響補燃室綜合性能的重要指標之一,在分析補燃室燃燒效率的同時,需考慮補燃室內壁附近的射流環境對補燃室結構完整性的影響。目前固沖補燃室內壁主要采用硅橡膠復合材料等作為絕熱材料。影響絕熱材料耐燒蝕性能的主要因素有溫度、O2含量、氣流沖刷和凝聚相粒子的侵蝕。
圖11為兩種進氣道結構補燃室內壁溫度分布,從圖中可以看出,Case 1補燃室頭部溫度約為2 200 K,補燃室中段射流混合區溫度超過2 500 K。在該溫度環境下,絕熱材料表面會發生劇烈的熱分解作用,導致絕熱層燒蝕率增加,并形成局部熱應力集中,損傷材料結構完整性,進而影響補燃室的工作壽命。因此,補燃室中段是熱防護的重點位置。結合圖5還可以看出,Case 1補燃室壁面高溫區具有單側分布、面積大、O2濃度低的特點;Case 2高溫區對稱分布,整體受熱相對均勻,且壁面O2濃度梯度相對較小。為進一步分析補燃室內壁燒蝕工況特征,提取Case 1和Case 2壁面縱向溫度最高和最低方向的O2濃度、溫度和射流速度分布特征進行分析(沿圖11黑色標記線方向提取數據)。
圖12~圖14為兩種補燃室結構內壁面高溫區和低溫區附近溫度、O2濃度及速度分布曲線。從圖中可以看出,Case 1和Case 2高溫區具有溫度高、O2濃度低、射流速度快的特點;Case 1低溫區的溫度約為600 K,高溫區溫度超過2 600 K,溫差較大,容易造成材料熱應變不均,導致材料龜裂。另外,Case1中補燃室單側高溫區表面氣流速度超過200 m/s,凝聚相粒子在氣流的帶動下沖刷在補燃室內壁(見圖10(a)),形成高溫熱分解、高速射流沖刷、高濃度粒子侵蝕和熱應力集中的綜合作用,對補燃室內壁產生嚴重破壞。目前,大部分絕熱材料的抗粒子侵蝕性能較差,在高溫熱作用和粒子沖擊作用下,表面熱防護層脫落,引起絕熱材料局部損傷變形,破壞熱防護結構的完整性。

圖11 補燃室內壁溫度分布Fig.11 Temperature distribution on internal walls of secondary combustion chamber

圖12 溫度分布曲線Fig.12 Curves of temperature distribution

圖13 O2質量分數分布曲線Fig.13 Curves of O2 mass fraction distribution
Case 2中補燃室中段和后段低溫區溫度高達1 000 K以上,高溫區和低溫區溫差相對較小,平衡了高溫區熱膨脹引起的局部應力集中,保證補燃室工作過程中相對穩定的力學環境。但低溫區O2濃度大于15%,該處壁面以氧化燒蝕為主。另外,Case 2中,凝聚相粒子沿中心軸擴散并燃燒(如圖10(b)所示),粒子濃度逐漸降低,到達壁面的粒子相對較小,避免了高濃度粒子對絕熱層壁面的侵蝕作用,保證補燃室工作過程中相對穩定的燒蝕環境。Case 2補燃室高溫區的主要燒蝕行為是高溫熱分解作用,低溫區為熱氧化燒蝕。

圖14 速度分布曲線Fig.14 Curves of velocity distribution
綜上分析,Case 1在遠離補燃室進氣道一側形成高溫熱分解、高濃度粒子侵蝕、高速射流沖刷和熱應力集中的綜合破壞作用;Case 2有效避免了高濃度粒子侵蝕和高速射流的直接沖刷,燒蝕的重點部位為高溫區,燒蝕的主要模式為化學燒蝕和高溫熱分解作用;Case 2內壁整體射流環境優于Case 1。
本文采用k-ε湍流模型、單步渦耗散燃燒模型以及KING硼粒子點火燃燒模型,開展了雙下側90°進氣結構和雙側180°進氣結構沖壓發動機補燃室內燃燒數值模擬,對比分析了補燃室燃氣燃燒特征和內壁的燒蝕環境。獲得以下結論:
1) 兩種進氣道結構對應補燃室中溫度和O2濃度對稱分布,呈現出高氧低溫和低氧高溫分布特征。雙下側90°進氣結構中,一次燃氣與空氣相互分離,O2濃度分布和溫度分布互補;雙側180°進氣結構中,燃氣的溫度分布和O2濃度分布逐漸趨于均勻。Case 2更有利于空氣與一次燃氣的摻混,且有利于溫度的擴散。
2) 兩種進氣道結構的補燃室內均形成頭部回流和軸向渦流。其中雙側180°進氣結構中,一次燃氣與空氣之間形成大漩渦,燃氣中的氣相組分和凝聚相顆粒與空氣中的氧得到充分混合和燃燒,至補燃室出口位置,總燃燒效率超過95%;雙下側90°進氣結構補燃室總燃燒效率為74%。
3) 雙下側90°進氣結構中,補燃室遠離進氣道一側壁面受高溫熱分解、高速射流沖擊和高濃度粒子侵蝕作用,周向高溫燒蝕和粒子侵蝕不均,容易形成局部燒蝕破壞和熱應力集中。
4) 雙側180°進氣沖壓發動機補燃室結構中,補燃室中后段軸向內壁面溫度相對均勻,且凝聚相粒子集中在中心軸線方向,擴散過程中濃度降低,形成的高溫熱燒蝕作用和粒子侵蝕作用相對較低。燒蝕的重點部位為高溫區,燒蝕的主要模式為化學燒蝕和高溫熱分解作用,Case 2內壁整體射流環境優于Case 1。