飛行器熱解炭化材料燒蝕產物對等離子體流場的影響規律

聶春生,楊光,聶亮,周禹,趙良

(中國運載火箭技術研究院空間物理重點實驗室,北京 100076)

0 引言

臨近空間飛行器高馬赫數飛行時,高溫空氣發生電離,形成一個非均勻等離子體包裹流場[1-2],同時飛行器外表面的防熱材料會發生燒蝕,燒蝕產物通過引射和擴散效應進入流場,影響等離子體流場特性,因此準確獲取等離子流場參數,必須考慮表面防熱材料燒蝕產物對等離子體流場的影響[3-4]。

國外在高溫等離子體流場方面的研究較早,針對球錐再入體開展了一系列研究,最有代表的是電磁衰減測量項目(RAM)[5],整個項目進行了8次飛行試驗,從理論和試驗兩方面研究了等離子鞘套對再入體通信系統的影響,其中第1次和第2次飛行試驗端頭采用碳酚醛燒蝕材料,開展了表面材料燒蝕產物對等離子體流場的影響研究,但相關結果公開發表的很少。數值研究方面,Bhutta等[6]最先將數值方法引入燒蝕產物對高溫非平衡流場的影響研究;Lewis等[7]利用拋物化Narier-Stokes(N-S)方程數值研究了泰氟隆材料燒蝕產物對高溫流場的影響;Bianchi等[8-9]對電弧風洞中的球錐模型石墨燒蝕流場進行了計算分析;Keenan等[10]針對半球和鈍錐體外形,數值研究了碳碳材料燒蝕與流場的耦合的影響;Thomas等[11]建立了高超聲速平衡燒蝕再入三維模擬程序。

國內還沒有開展專門針對高溫等離子體流場研究的飛行試驗,但在地面風洞環境中開展了一些研究工作。袁野等[12]在高頻等離子體風洞中研究了碳/碳和碳/炭化鋯兩種復合防熱材料燒蝕情況下對下游流場的影響;聶春生等[13]在高頻等離子體風洞中開展了純碳-碳材料燒蝕對球錐體等離子體流場的影響研究。數值研究方面,高鐵鎖等[14]進行了包含和不包含碳-酚醛燒蝕產物的再入小鈍錐流場計算,分析了燒蝕產物對流場電子數密度、溫度等流動參數的影響;魏叔如等[15]求解邊界層方程,開展了碳-碳材料燒蝕對電離邊界層影響的研究;董維中等[16]建立了碳-碳材料燒蝕流場的數值計算方法,研究了碳-碳材料燒蝕對再入體頭部區域溫度和熱流的影響,沒有對等離子體參數進行分析。

總之,國內外針對高溫等離子體流場開展的試驗研究相對較少,主要通過數值模擬方法研究表面材料燒蝕產物對等離子體流場的影響,但數值方法存在多個難點:1)不同材料的燒蝕反應過程不同,燒蝕產物對等離子體流場的影響也有較大差別,需要針對不同材料建立相應的燒蝕反應模型;2)由于表面催化、氧化、燒蝕與熱解等復雜的氣體表面相互作用,氣體表面邊界條件會非常復雜,表面邊界方程在流場中的求解算法十分復雜,難以建立高效、魯棒的耦合計算方法,使得針對復雜外形應用難度較大;3)由于進入流場的燒蝕組分與高溫空氣組分差別較大,流場中的化學反應過程會變得更加復雜,剛性問題和求解難度大大增加[17]。因此,目前文獻中的研究對象主要針對簡單的球錐再入體,材料主要針對傳統的泰氟隆材料和碳-碳材料。與球錐再入體不同,臨近空間飛行器一般會采用復雜的升力體外形,并長時間高馬赫數飛行,其熱環境特點為中高焓、低熱流,防隔熱設計往往會采用低密度的熱解炭化燒蝕材料[4],其燒蝕機理主要靠材料熱解吸熱和熱解氣體注入邊界層的熱阻塞效應。調研發現,目前針對升力體外形的等離子體流場,以及熱解炭化材料的燒蝕產物注入對等離子體流場的影響研究鮮有報道。

本文針對采用低密度熱解炭化防熱材料的臨近空間升力體飛行器,建立表面材料燒蝕與高溫等離子體流場數值計算方法,研究不同馬赫數下燒蝕產物注入流場對等離子體的影響規律。

1 計算方法

1.1 控制方程與求解方法

求解三維熱化學非平衡N-S方程[18],守恒形式的方程如下:

式中:Q為守恒變量;F、G、H與Fv、Gv、Hv分別對應不同方向的對流項與黏性項;W為非平衡源項。

采用有限體積法求解控制方程,無黏通量離散求解采用Roe平均方法與Yee的對稱TVD限制器,并對特征值采用了Harten的熵修正[19];黏性通量離散求解采用中心差分格式;時間推進采用全隱式LU-SGS方法,非平衡源項采用隱式處理[20]。

1.2 化學模型

當飛行器表面溫度高于表面材料的燒蝕溫度門檻時就會發生燒蝕,此時流場中要考慮從壁面注入的熱解燒蝕產物的氣體組分;本文等離子體流場計算考慮了32種組分,包括N、O、N2、O2、NO、NO+、C2H、C2H2、CN、H、H2、HCN、CH、Si、SiO、SiO2、SiC、Si2、SiN,化學反應模型見表1,其中組分的熱力學參數見文獻[21],輸運參數見文獻[22-23],化學反應速率常數見文獻[24-25]。

表1 高溫空氣/材料燒蝕組分的化學反應模型Tab.1 Chemical reaction model of material ablative components

1.3 燒蝕壁面邊界條件

本文研究定常的燒蝕流場,對應于壁面溫度最高,固體材料內部溫度分布達到穩態的情況,避免了耦合熱傳導方程的求解,并且不考慮燒蝕引起的壁面外形變化對流場的影響。

假設固體碳處于與表面氣體達到平衡的狀態,存在如下平衡約束條件[26]:

式中:ρw為壁面處氣體密度;cw,C為壁面碳原子的質量濃度;MC為碳的摩爾質量;KC,C為異相反應C(固相)→C(氣相)的化學平衡常數;Tw為壁面溫度。(2)式的求解需要知道壁面碳原子的質量濃度cw,C,此碳原子質量分數假設在壁面溫度、壓力和給定的元素比例條件下處于平衡。

化學平衡常數KC,C采用擬合關系式計算:

式中:T為氣體溫度。

對于某一元素k壁面處的質量濃度可以從如下表面元素質量平衡方程中獲得元素質量分數:

壁面壓力pw采用如下壁面法向動量獲得:

式中:vw為壁面法向速度,vw來自如下質量連續方程,

式中:qC為壁面熱流;α為輻射吸收系數;qr為壁面向氣體的輻射熱流;ε為材料輻射系數;σ為斯蒂芬-玻爾茲曼常數;hw為壁面焓值。

式中:ρv、ρC分別為原始材料密度和炭化層密度。

2 方法驗證

采用公開發表的電磁衰減測量項目C方案(RAM-C)飛行試驗等離子體測試結果,對本文計算方法進行驗證。RAM-C模型為鈍球錐體,球頭半徑Rn=0.154 2 m,半錐角θ=9°,模型長度L=1.295 m。

飛行條件:速度U∞=7 650.0 m/s,高度H=71 km,壁面溫度Tw=1 500 K。本文主要考核流場中電子數密度分布,空氣化學反應模型采用Park模型,考慮了N2、O2、NO、NO+、N、O、e-共計7種組分,熱力學模型采用兩溫度模型,壁面考慮完全非催化和完全催化兩種條件。

圖1給出了沿軸向位置電子數密度N峰值分布與飛行試驗結果對比,橫軸為軸向位置與頭部半徑的無量綱距離,縱軸為電子數密度。由圖1可以看出,本文電子數密度峰值分布計算結果與飛行數據符合較好,誤差在半個數量級以內。圖2給出了典型截面處法向電子數密度分布與飛行試驗結果對比。由圖2可以看出:完全非催化和完全催化壁面條件會影響壁面附近的電子數密度,完全非催化壁結果高于完全催化壁結果,遠離壁面的空間等離子體流場受影響較小。總體上,完全非催化壁結果更接近飛行試驗測量值,表明本文等離子體流場計算方法在電子數密度分布計算中有較好的準確性。

圖1 電子數密度峰值分布比較Fig.1 Comparison of peak distributions of electron number density

圖2 法向電子數密度分布比較Fig.2 Comparison of normal electron number density distributions

3 燒蝕產物的影響分析

3.1 計算模型及防熱材料屬性

針對某升力體氣動布局的臨近空間飛行器前體部分,開展燒蝕產物對等離子體流場的影響分析。圖3給出了計算外形和網格示意圖,采用六面體結構化網格,在近壁面附近進行了加密處理,網格總量150萬。計算來流工況為高度60 km、馬赫數分別為20、25、28,飛行攻角0°。

圖3 計算模型和網格Fig.3 Computational model and grid

為研究不同熱解炭化材料燒蝕產物對等離子體流場的影響,選取表2所示成分存在差異的同一類輕質燒蝕防熱材料[27-28],其中材料1為文獻[27-28]中給出的阿波羅飛船采用的輕質燒蝕材料;材料2和材料3分別在材料1基礎上進行了工藝改進,分別去除材料中的氫元素和硅元素,分別記為輕質燒蝕材料2和輕質燒蝕材料3。3種材料密度均取544.63 kg/m3,熱解燒蝕后的炭化層密度取256.29 kg/m3。

表2 材料參數Tab.2 Material parameters

3.2 計算結果分析

3.2.1 純空氣等離子體流場規律分析

圖4給出了高度60 km、馬赫數27狀態的流場溫度和電子數密度分布。由圖4可以看出,高馬赫數來流在飛行器周圍形成了脫體激波,激波對等離子體流場的影響主要取決于激波后的流場溫度;對于頭部區域,激波強度大、壓縮性強,波后溫度超過了10 000 K,空氣發生大規模離解和電離反應,流場中的O2分子幾乎完全離解,O、N原子的電離反應是整個流場中等離子體產生的主要源頭,該區域的電子數密度最高,峰值超過了1×1014cm-3;飛行器身部區域的激波強度較弱,激波后的流場溫度約為5 000 K,該區域的空氣離解和電離程度較弱,飛行器身部流場的等離子體主要由頭部流場輸運擴散而來;受激波層影響,飛行器身部區域的等離子體層主要貼近飛行器壁面,其中前緣部位的等離子體層厚度相對迎風面較薄,隨著流動向下游發展,流場溫度逐漸降低,飛行器身部電子數密度比頭部低2～3個量級。

圖4 流場溫度和電子數密度分布Fig.4 Temperature and electron number density distribution of flow field

圖5給出了頭部駐點線流場溫度分布。由圖5可以看出,駐點線上的振動溫度與平動溫度差別較大,流場的熱力學非平衡效應顯著,由于激波層內高溫區的氣體電離反應是等離子體層形成的主要原因,需要采用兩溫度模型進行等離子體流場計算。

圖5 沿駐點線溫度分布Fig.5 Temperature distribution along stagnation line

圖6給出了駐點線流場氣體電離組元的數密度分布。由圖6可以看出:馬赫數為20時,NO+離子是流場中電子數密度的主要貢獻源;馬赫數增大為25、27時,除了NO+離子,N+對電子數密度也有重要貢獻,其中NO+主要由O與N原子的締合電離反應機制產生,N+主要由N與其他組分的碰撞電離反應機制產生。綜合上述分析可知,對流場中電子數密度貢獻最重要的離子是NO+和N+,其次依次為這些離子對電子數密度的貢獻很小。

圖6 沿駐點線電離組元的數密度分布Fig.6 Number density distribution of ionized components along stagnation line

3.2.2 燒蝕產物對等離子流場影響分析

圖7給出了材料3在高度60 km、馬赫數為20狀態的表面溫度和燒蝕質量流率分布。由圖7可以看出:飛行器頭部和前緣溫度超過了2 000 K,該區域發生明顯燒蝕,其中頭部燒蝕最嚴重,燒蝕質量流率超過了0.2 kg/(m2·s);飛行器身部迎風面的溫度約為1 000 K,燒蝕質量流率相對較小,沿流向燒蝕效應逐漸減弱。

圖7 表面溫度和燒蝕質量流率(高度60 km、馬赫數20,材料3)Fig.7 Surface temperature and ablation mass flow rate(height of 60 km,Mach number of 20,and Material 3)

圖8給出了流場中主要燒蝕組分的數密度。由圖8可見:注入流場的主要熱解燒蝕產物是CO和H2,其他燒蝕氣體組分的質量分數相對較小,并主要集中在頭部區域;下游飛行器身部的燒蝕效應逐漸減弱,但隨著流動向下游發展,頭部燒蝕產物會隨流動的擴散效應輸運到下游,燒蝕效應對流場等離子分布的影響范圍擴大。

圖8 對稱面燒蝕產物分布(高度60 km、馬赫數20,材料3)Fig.8 Distribution of ablated products on symmetry plane(height of 60 km,Mach number of 20,and Material 3)

圖9給出了頭部駐點線的流場參數對比。由圖9可以看出:燒蝕產物注入流場后,流場基本結構沒有發生明顯變化,溫度和電子數密度峰值未受影響;由于壁面燒蝕的引射效應,使得頭部激波脫體距離增大,駐點線上的溫度、電子數密度和組分質量分數的分布也發生了一定變化,流場中高溫區域變大,O2和N2分子更多地分解成O和N原子,O原子參加氧化反應被消耗而質量分數明顯減小,消耗的O原子主要生成了CO。

圖9 頭部軸線流場參數分布比較(高度60 km、馬赫數20,材料3)Fig.9 Distribution of axial flow field parameters on aircraft head(height of 60 km,Mach number of 20,and Material 3)

圖10和圖11分別給出了材料3在不同馬赫數下駐點線的溫度和電子數密度曲線。由圖10和圖11可以看出,對于材料3,隨著飛行馬赫數變大,表面燒蝕速率增強,燒蝕產物進入流場使得激波脫體距離增大,但駐點線上的電子密度峰值變化很小,這主要是因為流場中電子數密度貢獻最重要的離子是NO+。圖12給出了材料3在不同馬赫數駐點線NO+離子數密度分布。從圖12中可以看出,不同馬赫數下燒蝕效應對頭部駐點線上的NO+離子數密度影響相對較小導致的。

圖10 不同馬赫數下駐點線溫度分布(材料3)Fig.10 Distribution of temperature along stagnation line at different Mach numbers(Material 3)

圖11 不同馬赫數下駐點線電子數密度分布(材料3)Fig.11 Distribution of electron number density along stagnation line at different Mach numbers(Material 3)

圖12 不同馬赫數駐點線NO+離子數密度分布(材料3)Fig.12 Distribution of NO+ion number density along stagnation line at different Mach numbers(Material 3)

圖13給出了材料3的對稱面電子數密度峰值對比曲線。由圖13可以看出:隨著馬赫數增大,表面燒蝕速率增大,由于進入流場的燒蝕產物增多,使得流場中的電子數密度有一定的增大;隨著流動從頭部向下游發展,飛行器身部的熱流會逐漸降低,表面材料的燒蝕速率逐漸減小,但燒蝕產物對飛行器身部流場電子數密度的影響逐漸增強;下游0.5 m的位置,材料3的燒蝕產物對電子數密度的影響增量約為20%～30%。

圖13 對稱面電子數密度峰值對比(材料3)Fig.13 Comparison of electron number density peaks on symmetry plane(Material 3)

圖14給出了表2中不同成分的低密度防熱材料在高度60 km、馬赫數25狀態下飛行器表面中心線的燒蝕質量流率對比曲線。由圖14可以看出:由于飛行器頭部的氣動加熱最為嚴酷,使得燒蝕質量流率最大,材料1和材料3接近,約為0.55 kg/(m2·s),材料2為0.33 kg/(m2·s);與頭部相比,飛行器身部的燒蝕質量流率會迅速降低約2～3個量級,表2中3種材料的燒蝕質量流率均低于0.013 kg/(m2·s),其中材料2要略高于材料1和材料3。

圖14 不同材料飛行器表面對稱線的燒蝕質量流率Fig.14 Ablation mass flow rate of surface center lines of aircraft with different materials

圖15、圖16分別為高度60 km、馬赫數25狀態下駐點線溫度和電子數密度對比曲線。由圖15和圖16可以看出:該狀態不同材料的燒蝕效應存在差異,從而會影響頭激波的脫體距離,但駐點線的峰值溫度和峰值電子數密度幾乎不受影響;駐點線脫體距離從大到小依次為材料1、材料2和材料3;材料1和材料2的燒蝕質量流率接近,但由于材料的成分不同,注入流場的燒蝕產物不同,使得兩種材料燒蝕后對流場激波脫體距離的影響程度存在差異。

圖15 不同材料條件下駐點線溫度分布Fig.15 Distribution of temperatures of different materials along stagnation line

圖16 不同材料條件下駐點線電子數密度分布Fig.16 Distribution of electron number density along stagnation line under different material conditions

圖17給出了高度60 km、馬赫數25狀態下對稱面電子數密度峰值對比曲線。由圖17可以看出:對稱面材料的燒蝕產物對流場峰值電子數密度影響沿流動向下游發展逐漸增強,并且不同材料燒蝕產物對流場中的電子數密度峰值影響量不同,其中材料2的影響相對最大;在0.5 m位置,電子數密度的影響增量約為300%～320%。結合圖15表面燒蝕質量流率分析可知,飛行器身部流場的電子數密度一方面受頭部區域燒蝕產物流動擴散的影響,另一方面還要受當地的微量燒蝕產物的影響,二者均會使得流場中的電子數密度增大,其影響機理為燒蝕產物中CH、CO會發生反應產生O原子,O原子與N原子發生締合電離反應,產生更多NO+。圖18所示為對稱面x=0.4 m位置流場中的NO+離子對比。由圖18可以看出,材料燒蝕產物注入后流場中的NO+離子顯著增大,從而使得流場中的電子數密度均比純空氣狀態增大。

圖17 不同材料條件下對稱面電子數密度峰值對比Fig.17 Comparison of peak values of electron number density on symmetry planes of different materials

圖18 x=0.4 m位置流場NO+離子數密度對比曲線Fig.18 Comparison curves of NO+ion number density in flow field for x=0.4 m

3 結論

本文建立了考慮燒蝕產物影響的等離子體流場數值計算方法,針對阿波羅飛船采用的輕質燒蝕材料,開展了升力體外形的燒蝕產物對等離子體流場的影響規律研究。得到如下主要結論:

1)對于升力體飛行器的等離子體流場,影響最重要的離子是NO+和N+。

2)熱解炭化材料最主要的熱解燒蝕產物是CO和H2;燒蝕產物進入流場會使得激波脫體距離增加,等離子體厚度增加。

3)燒蝕產物進入流場,對頭部的峰值電子數密度影響很小,隨著流動向下游發展影響逐漸增大;影響機理為,燒蝕產物中CH、CO組分會發生反應,產生O原子,O原子與N原子發生締合電離反應,產生更多NO+,使得流場中電子數密度增大。

4)成分存在差異的同一類防熱材料,燒蝕產物對流場電子數密度的影響增量存在差異;本文研究的一類防熱材料,高度60 km、馬赫數為25狀態下飛行器下游0.5 m位置對稱面峰值電子數密度,材料2的增量約為300%～320%。