石安華， 李海燕， 石衛(wèi)波， 梁世昌

(中國(guó)空氣動(dòng)力研究與發(fā)展中心 超高速空氣動(dòng)力研究所， 四川 綿陽(yáng) 621000)

0 引言

臨近空間高超聲速飛行器近年來受到各軍事強(qiáng)國(guó)的熱捧，在臨近空間高超聲速武器的攻防對(duì)抗研究中，目標(biāo)輻射特性是重要的研究?jī)?nèi)容，除了試驗(yàn)研究外，建模研究作為主要方法之一也得到大力發(fā)展。臨近空間高超聲速飛行器可分為滑翔類和帶動(dòng)力類，滑翔類高超聲速飛行器輻射特性主要受氣動(dòng)熱和高溫繞流氣體影響，而帶動(dòng)力類高超聲速飛行器輻射特性主要受氣動(dòng)熱和發(fā)動(dòng)機(jī)噴焰輻射影響。針對(duì)帶動(dòng)力類高超聲速飛行器輻射特性，盡管已有一些學(xué)者對(duì)此開展了相關(guān)研究工作，但鮮有同時(shí)考慮時(shí)變氣動(dòng)加熱、殼體結(jié)構(gòu)與熱傳導(dǎo)、噴焰中氣體分子能級(jí)躍遷輻射對(duì)帶動(dòng)力類高超聲速飛行器輻射特征影響的研究，使得分析結(jié)果較為粗糙。

本文針對(duì)上述問題，從考慮流動(dòng)過程化學(xué)反應(yīng)、燃燒化學(xué)反應(yīng)、氣動(dòng)熱、熱傳導(dǎo)、本體輻射、分子躍遷輻射等輻射產(chǎn)生過程及機(jī)制出發(fā)，開展以超燃沖壓發(fā)動(dòng)機(jī)為動(dòng)力的類X-51A(見圖1)高超聲速巡航飛行器在典型飛行狀態(tài)下的紅外輻射特征分析研究，分析結(jié)果可為臨近空間高超聲速導(dǎo)彈的相關(guān)探測(cè)技術(shù)研究提供數(shù)據(jù)參考。

圖1 類X-51A高超聲速巡航飛行器數(shù)值模型Fig.1 Digital model for hypersonic cruise vehicle like X-51A

1 計(jì)算方法

1.1 流場(chǎng)及氣動(dòng)熱特性計(jì)算方法

流場(chǎng)包括飛行器繞流場(chǎng)和繞流場(chǎng)與吸氣式超燃沖壓發(fā)動(dòng)機(jī)燃?xì)怦詈系膰娧媪鲌?chǎng)。

臨近空間高超聲速巡航飛行器所在飛行區(qū)域的流動(dòng)是連續(xù)流動(dòng)，飛行器燃料通常采用航空煤油，由于處于高超聲速飛行條件，氣體溫度高、化學(xué)反應(yīng)強(qiáng)烈，繞流和噴焰流場(chǎng)中考慮O、O、N、NO、OH、NH、H、CO、NO、NO、HO、HO、CO、HO、CH、H和N 17種組分，化學(xué)動(dòng)力學(xué)模型使用包括46種反應(yīng)的燃燒化學(xué)動(dòng)力學(xué)模型。假設(shè)飛行器外部繞流和超燃沖壓發(fā)動(dòng)機(jī)內(nèi)部的流動(dòng)是可壓縮湍流流動(dòng)。對(duì)于飛行器外部繞流、超燃發(fā)動(dòng)機(jī)尾噴流與飛行器外部繞流相互作用下的流場(chǎng)特性模擬，采用三維可壓縮化學(xué)反應(yīng)湍流的Favre平均方程計(jì)算流體動(dòng)力學(xué)方法模擬。在笛卡爾坐標(biāo)系下，三維守恒形式的Favre平均方程為

(1)

式中：為守恒變量向量；為時(shí)間(s)；、和分別為笛卡爾坐標(biāo)系中軸、軸和軸方向的無黏通量向量；、和分別為笛卡爾坐標(biāo)系中軸、軸和軸方向的黏性通量向量；為化學(xué)反應(yīng)源項(xiàng)向量。

為封閉平均方程，引入Boussinesq渦黏性假設(shè)，并采用-兩方程模型來模擬湍流效應(yīng)，為湍動(dòng)能(J/kg)，為湍動(dòng)能耗散率(1s)。首先對(duì)三維化學(xué)反應(yīng)湍流控制方程進(jìn)行無量綱化，然后采用有限體積法對(duì)無量綱化后的控制方程進(jìn)行離散。對(duì)離散后的無量綱控制方程組，采用上下對(duì)稱Gauss-Seidel (LU-SGS)隱式方法進(jìn)行迭代求解。氣動(dòng)熱采用表面熱流公式計(jì)算分析：

(2)

式中：為壁面熱流(W/m)；為熱傳導(dǎo)系數(shù)(W/(m·K))；為溫度(K)；為飛行器表面局部法向距離(m)；為化學(xué)組分種類數(shù)量；為當(dāng)?shù)貧怏w密度(kg/m)；和分別為第種組分的擴(kuò)散系數(shù)(m/s)和焓值(J/kg)；為第種組分的質(zhì)量分?jǐn)?shù)。

1.2 燃燒室流動(dòng)特性計(jì)算方法

超燃沖壓發(fā)動(dòng)機(jī)尾噴管前的燃燒室流場(chǎng)采用準(zhǔn)一維流動(dòng)模型，將準(zhǔn)一維Euler方程與有限速率化學(xué)動(dòng)力學(xué)方程耦合，在方程中考慮壁面摩擦、面積變化和燃料質(zhì)量添加的影響。然后對(duì)無量綱化的準(zhǔn)一維Euler方程應(yīng)用有限體積法進(jìn)行離散，采用隱式方法求解離散后的方程。

在笛卡爾坐標(biāo)系下，準(zhǔn)一維Euler方程和組分連續(xù)方程的強(qiáng)耦合形式為

(3)

式中：為燃燒室截面積(m)。

在準(zhǔn)一維Euler方程和組分連續(xù)方程中，質(zhì)量添加項(xiàng)出現(xiàn)在方程的右端源項(xiàng)中，此外，右端源項(xiàng)也考慮了壁面摩擦、面積變化以及化學(xué)反應(yīng)的貢獻(xiàn)。燃料采用以CH為簡(jiǎn)化分子式的JP7航空煤油，其燃燒化學(xué)動(dòng)力學(xué)模型使用CH燃燒的9種組分(H、H、O、O、OH、HO、CO、CO和CH)、13個(gè)化學(xué)反應(yīng)方程模型，高溫條件下空氣與煤油反應(yīng)混合氣體中各組分熱力學(xué)參數(shù)由Janaff表獲得。CH的當(dāng)量油氣比為，化學(xué)恰當(dāng)比為，根據(jù)油氣比關(guān)系求得噴射的煤油流量：

(4)

考慮到燃料的混合過程，質(zhì)量添加項(xiàng)采用指數(shù)分布進(jìn)行處理。采用參考溫度法來計(jì)算摩擦系數(shù)，通過雷諾比擬可以計(jì)算傳熱系數(shù)。對(duì)于燃燒室前方的隔離段而言，準(zhǔn)一維方法無法準(zhǔn)確捕捉到其中的預(yù)燃激波串這種復(fù)雜非一維結(jié)構(gòu)。對(duì)于激波串的處理，采用工程預(yù)燃激波串模型進(jìn)行計(jì)算。將預(yù)燃激波串模型耦合到準(zhǔn)一維非穩(wěn)態(tài)方法中，得到一種改進(jìn)的準(zhǔn)一維非穩(wěn)態(tài)計(jì)算方法，結(jié)合相應(yīng)的碳?xì)淙剂匣瘜W(xué)反應(yīng)動(dòng)力學(xué)模型和燃料混合模型，來獲得燃燒室出口的參數(shù)。

1.3 流場(chǎng)輻射特性計(jì)算方法

流場(chǎng)組分光譜輻射使用分子振轉(zhuǎn)光譜模型，以振轉(zhuǎn)譜線標(biāo)準(zhǔn)吸收強(qiáng)度為基礎(chǔ)，根據(jù)愛因斯坦輻射與吸收關(guān)系計(jì)算分子光譜輻射強(qiáng)度及吸收強(qiáng)度。

輻射傳輸使用窄帶吸收模型，對(duì)傳輸路徑上的輻射按組分種類分別計(jì)算再求和，計(jì)算方向從觀測(cè)點(diǎn)起至視線到達(dá)飛行器表面或穿出視場(chǎng)。

流場(chǎng)輻射使用光學(xué)薄傳輸模型計(jì)算會(huì)帶來較大誤差，在流場(chǎng)光輻射強(qiáng)度計(jì)算中使用吸收模型，在無粉塵顆粒的氣體輻射與傳輸環(huán)境中忽略光散射影響。假設(shè)流場(chǎng)中第個(gè)體積元d內(nèi)存在的化學(xué)組分共有種，d中第種氣體組分的數(shù)密度為,(個(gè)m)、輻射強(qiáng)度為,、吸收系數(shù)為,，則

(5)

(6)

式中：、分別為流場(chǎng)體積元d的輻射強(qiáng)度(W/Sr)、吸收系數(shù)(1m)。

在光譜分辨率要求不高的情況下，分子的光譜吸收系數(shù)通常采用近似算法，在本文中采用窄帶平均吸收系數(shù)模型。窄帶平均吸收系數(shù)模型是從吸收系數(shù)的嚴(yán)格公式出發(fā)，計(jì)算出每條單獨(dú)譜線的吸收系數(shù)后，將研究的波段劃分為多個(gè)窄帶Δ，計(jì)算落于此Δ內(nèi)所有譜線的平均吸收系數(shù)，采用(7)式計(jì)算：

(7)

體積元d照射到觀測(cè)點(diǎn)處的輻照度為

(8)

(8)式中的求和項(xiàng)是針對(duì)輻射傳輸路徑上所有體積元d進(jìn)行的。

1.4 飛行器本體輻射特性計(jì)算方法

141 表面溫度求解方法

飛行器表面溫度是關(guān)系到飛行器輻射計(jì)算的關(guān)鍵物理量，其受到熱環(huán)境、飛行器結(jié)構(gòu)、材料物性參數(shù)等多種因素的影響。為精細(xì)描述包含表面溫度在內(nèi)的飛行器結(jié)構(gòu)溫度場(chǎng)變化，一般需要求解三維瞬態(tài)熱傳導(dǎo)方程。本文采用有限元方法對(duì)熱傳導(dǎo)過程進(jìn)行模擬，求得飛行器表面溫度。

瞬態(tài)溫度場(chǎng)的溫度在直角坐標(biāo)中應(yīng)滿足的微分方程為

(9)

式中：為飛行器殼體材料密度(kg/m)；為材料比熱容(J/(kg·K))；、、分別為材料沿軸、軸、軸方向的熱導(dǎo)率(W/(m·K))；為物體內(nèi)部熱源的熱流密度(W/m)。

求解域內(nèi)溫度場(chǎng)采用的邊界條件為

(10)

式中：、、為邊界外法線方向余弦；為邊界上給定熱流量或熱流密度(W/m)。

求解瞬態(tài)溫度場(chǎng)問題是求解在一定初始條件下滿足瞬態(tài)熱傳導(dǎo)方程及邊界條件的溫度場(chǎng)，為空間坐標(biāo)與時(shí)間的函數(shù)。采用四節(jié)點(diǎn)四面體單元對(duì)空間進(jìn)行離散，依據(jù)泛函理論，將傳熱控制方程離散為代數(shù)方程組，寫成矩陣形式為

(11)

142 本體輻射特性計(jì)算方法

假定飛行器是朗伯輻射體，半球發(fā)射率、方向發(fā)射率和光譜發(fā)射率相等，且與方向無關(guān)，則本體紅外輻射可采用灰體輻射模型計(jì)算，表面網(wǎng)格使用溫度場(chǎng)網(wǎng)格的外表面部分。

利用普朗克公式，對(duì)波長(zhǎng)從到積分，得到黑體面源在～波段范圍內(nèi)的輻射出射度為

(12)

式中：為黑體輻射出射度(W/m)；()為黑體的光譜輻射出射度(W/m)；為第一輻射常數(shù)，=3741 8×10W·m；為第二輻射常數(shù)，=1438 8×10m·K。

飛行器表面材料發(fā)射率為、溫度為的小面源的輻射出射度=·，光譜輻射亮度()=()π。

由于投影到觀察方向的面源輻射強(qiáng)度具有獨(dú)立性，相互之間不干擾，并具有可疊加性；分析中考慮飛行器離觀察點(diǎn)很遠(yuǎn)，面元中心至觀察點(diǎn)連線的矢量幾乎相互平行，對(duì)飛行器表面各點(diǎn)輻射積分即可得到飛行器本體的總輻射。

1.5 總輻射計(jì)算方法

考慮本體輻射穿過流場(chǎng)的吸收效應(yīng)、本體對(duì)部分流場(chǎng)輻射的遮擋效應(yīng)以及本體之間的遮擋效應(yīng)，即可計(jì)算得到目標(biāo)總輻射。

2 計(jì)算結(jié)果

2.1 計(jì)算條件

計(jì)算中使用的飛行器外形見圖1，飛行器殼體各部位使用的材料和厚度取值見表1，殼體材料物性參數(shù)見表2，計(jì)算飛行狀態(tài)為飛行器以馬赫數(shù)=6、攻角=8°在海拔為30 km、40 km兩種高度下長(zhǎng)時(shí)間飛行150 s、400 s和800 s共6個(gè)狀態(tài)。計(jì)算中，燃料使用JP7航空煤油(簡(jiǎn)化分子式為CH)，油氣當(dāng)量比假設(shè)為0.8；飛行器初始溫度假定為298 K。

表1 飛行器殼體材料及厚度

表2 飛行器殼體材料物性參數(shù)

本文光譜輻射強(qiáng)度計(jì)算光譜范圍為1～14 μm，光輻射強(qiáng)度及光輻射亮度計(jì)算波段為3～5 μm、8～12 μm。計(jì)算輻射方向?yàn)閺念^部開始往背風(fēng)面方向繞到尾部共17個(gè)方向(=10°～170°，間隔10°)。

2.2 典型計(jì)算結(jié)果

圖2給出了飛行器在30 km高度飛行時(shí)的典型流場(chǎng)參數(shù)和飛行150 s時(shí)的本體溫度。

圖2 典型飛行狀態(tài)下飛行器流場(chǎng)參數(shù)、表面熱流及溫度(H=30 km、Ma=6、α=8°)Fig.2 Flow field parameters, heat flux and surface temperature (H=30 km, Ma=6, α=8°)

圖3、圖4給出了飛行器在30 km高度飛行150 s、400 s、800 s時(shí)，飛行器在90°輻射方向上波段分別為3～5 μm、8～12 μm的輻射亮度。從圖3和圖4中可以看出，飛行器噴焰輻射亮度隨飛行時(shí)間基本無變化，而本體輻射亮度隨飛行時(shí)間增加而增強(qiáng)。

圖3 飛行器在典型飛行時(shí)刻在3～5 μm波段的輻射亮度(H=30 km, Ma=6, α=8°, θ =90°)Fig.3 Luminance in 3～5 μm wave band at typical flight time (H=30 km, Ma=6, α=8°,θ =90°)

圖4 飛行器在典型飛行時(shí)刻在8～12 μm波段的輻射亮度(H=30 km, Ma=6, α=8°, θ =90°)Fig.4 Luminance in 8～12 μm wave band at typical flight time (H=30 km, Ma=6, α=8°, θ=90°)

圖5、圖6給出了飛行器在30 km高度飛行 800 s 時(shí)，飛行器在50°和130°輻射方向上波段分別為3～5 μm、8～12 μm的輻射亮度。從圖5和圖6中可以看出，不同輻射方向上的輻射圖像是不一樣的，50°方向上的輻射亮度強(qiáng)于130°方向上的輻射亮度。

圖5 飛行器飛行800 s時(shí)典型輻射方向3～5 μm波段的 輻射亮度(H=30 km、Ma=6、α=8°)Fig.5 Luminance in 3～5 μm wave band at 800 s and different θ (H=30 km, Ma=6, α=8°)

圖6 飛行器飛行800 s時(shí)典型輻射方向8～12 μm波段的輻射亮度(H=30 km、Ma=6、α=8°)Fig.6 Luminance in 8～12 μm wave band at 800 s and different θ (H=30 km, Ma=6, α=8°)

圖7給出了飛行器在30 km高度飛行800 s時(shí)飛行器典型輻射方向的光譜輻射強(qiáng)度隨波段的變化曲線。從圖7中可以看出：飛行器輻射總體體現(xiàn)為在本體熱輻射上疊加噴焰氣體的特征輻射，主要為HO、CO、CO氣體的特征輻射，特征輻射為峰值波長(zhǎng)在1.52 μm、2.68 μm、4.39 μm處的3個(gè)輻射帶。其中，1.52 μm輻射帶主要為HO輻射，2.68 μm輻射帶主要為CO、HO、CO輻射，4.39 μm輻射帶主要為CO輻射；特征譜波長(zhǎng)不隨輻射方向變化，但特征譜輻射強(qiáng)度隨輻射方向變化劇烈，在90°左右方向上特征譜輻射明顯強(qiáng)于本體輻射，這是因?yàn)檫@些方向上噴焰輻射投影面積大，強(qiáng)輻射區(qū)穿過的噴焰氣體距離短。此外，在該時(shí)刻的本體紅外輻射主要集中在短波紅外波段，主要是因?yàn)榇藭r(shí)飛行器本體大部分區(qū)域溫度已到上千度。而在飛行初始階段，本體紅外輻射應(yīng)主要集中在中長(zhǎng)波紅外波段，這是因?yàn)轱w行器在初始飛行階段表面溫度還不高造成的。

圖7 飛行器在典型輻射方向上的光譜輻射強(qiáng)度(H=30 km、Ma=6、α=8°、t=800 s)Fig.7 Spectral intensity at different θ (H=30 km,Ma=6, α=8°, t=800 s)

圖8所示為飛行器在30 km高度飛行800 s時(shí)典型特征光譜(波長(zhǎng)1.52 μm、2.68 μm、4.39 μm)的輻射強(qiáng)度與波長(zhǎng)1.6 μm處的光譜輻射強(qiáng)度(基本為本體光譜輻射強(qiáng)度)隨輻射方向的變化曲線。從圖8中可以看出，在特征譜輻射波長(zhǎng)位置，噴焰氣體特征譜輻射明顯強(qiáng)于本體輻射，如果用特征輻射波長(zhǎng)做窄帶測(cè)量，特別是在飛行器正上方附近方向，將具有最高的信噪比。

圖8 飛行器在不同方向的特征光譜輻射強(qiáng)度(H=30 km、Ma=6、α=8°、t=800 s)Fig.8 Typical spectral intensity at different θ (H=30 km,Ma=6, α=8°, t=800 s)

飛行器在40 km高度飛行時(shí)的輻射出射度與光譜輻射強(qiáng)度與在30 km高度飛行時(shí)類似，只是相同時(shí)刻輻射出射度和光譜輻射強(qiáng)度變小。表3給出了飛行器在70°輻射方向上、兩個(gè)飛行高度下在3～5 μm和8～12 μm波段本體輻射強(qiáng)度和總輻射強(qiáng)度。從表3中可以看出：飛行器噴焰輻射隨飛行時(shí)間變化不大，而本體輻射隨飛行時(shí)間增加而增加；3～5 μm波段的輻射遠(yuǎn)強(qiáng)于8～12 μm波段的輻射；在3～5 μm和8～12 μm波段，噴焰輻射占比在13%～17%左右，與在特征輻射波長(zhǎng)的窄帶內(nèi)輻射相比差異十分明顯，窄帶特征輻射波長(zhǎng)內(nèi)噴焰輻射占比可達(dá)50%～90%。

表3 典型飛行狀態(tài)下的飛行器輻射強(qiáng)度

3 結(jié)論

本文采用高超聲速空氣動(dòng)力學(xué)、燃燒學(xué)、傳熱學(xué)、氣體輻射理論等相結(jié)合的計(jì)算分析方法，開展了類X-51A高超聲速巡航飛行器在典型飛行狀態(tài)下的紅外輻射光譜及紅外輻射亮度和紅外輻射強(qiáng)度研究。得到如下主要結(jié)論：

1)類X-51A高超聲速巡航飛行器具有典型的特征譜輻射，特征輻射譜集中在HO的1.52 μm波長(zhǎng)輻射帶，CO、HO、CO的2.68 μm波長(zhǎng)輻射帶以及CO的4.39 μm波長(zhǎng)輻射帶。

2)類X-51A高超聲速巡航飛行器在強(qiáng)特征輻射譜帶內(nèi)的總輻射中，氣體輻射占主導(dǎo)；而在非特征輻射譜帶內(nèi)的總輻射中，飛行器本體輻射占主導(dǎo)。

3)在50°～130°輻射方向，給定飛行器在3～5 μm波段的輻射強(qiáng)度在數(shù)萬瓦/Sr量級(jí)，而在8～12 μm波段的輻射強(qiáng)度在數(shù)千瓦/Sr量級(jí)。

4)基于天基探測(cè)，選擇低層地球大氣(20 km以下)中CO和HO的強(qiáng)吸收帶特性，使用1.52 μm、2.68 μm、4.39 μm窄帶探測(cè)類X-51A高超聲速巡航飛行器具有強(qiáng)的探測(cè)目標(biāo)信號(hào)和高的信噪比。

開展臨近空間高超聲速巡航飛行器輻射特性研究的意義在于攻防對(duì)抗中的探測(cè)與反探測(cè)應(yīng)用。因此，下一步研究重點(diǎn)是在綜合考慮飛行器的具體結(jié)構(gòu)、防熱模式、飛行彈道基礎(chǔ)上，結(jié)合探測(cè)器特性、環(huán)境對(duì)輻射的傳輸特性及背景輻射特性進(jìn)行綜合分析，提出有利于探測(cè)臨近空間高超聲速巡航飛行器的可靠探測(cè)方案。