聶春生,楊光,聶亮,周禹,趙良
(中國運載火箭技術(shù)研究院空間物理重點實驗室,北京 100076)
臨近空間飛行器高馬赫數(shù)飛行時,高溫空氣發(fā)生電離,形成一個非均勻等離子體包裹流場[1-2],同時飛行器外表面的防熱材料會發(fā)生燒蝕,燒蝕產(chǎn)物通過引射和擴散效應(yīng)進入流場,影響等離子體流場特性,因此準確獲取等離子流場參數(shù),必須考慮表面防熱材料燒蝕產(chǎn)物對等離子體流場的影響[3-4]。
國外在高溫等離子體流場方面的研究較早,針對球錐再入體開展了一系列研究,最有代表的是電磁衰減測量項目(RAM)[5],整個項目進行了8次飛行試驗,從理論和試驗兩方面研究了等離子鞘套對再入體通信系統(tǒng)的影響,其中第1次和第2次飛行試驗端頭采用碳酚醛燒蝕材料,開展了表面材料燒蝕產(chǎn)物對等離子體流場的影響研究,但相關(guān)結(jié)果公開發(fā)表的很少。數(shù)值研究方面,Bhutta等[6]最先將數(shù)值方法引入燒蝕產(chǎn)物對高溫非平衡流場的影響研究;Lewis等[7]利用拋物化Narier-Stokes(N-S)方程數(shù)值研究了泰氟隆材料燒蝕產(chǎn)物對高溫流場的影響;Bianchi等[8-9]對電弧風(fēng)洞中的球錐模型石墨燒蝕流場進行了計算分析;Keenan等[10]針對半球和鈍錐體外形,數(shù)值研究了碳碳材料燒蝕與流場的耦合的影響;Thomas等[11]建立了高超聲速平衡燒蝕再入三維模擬程序。
國內(nèi)還沒有開展專門針對高溫等離子體流場研究的飛行試驗,但在地面風(fēng)洞環(huán)境中開展了一些研究工作。袁野等[12]在高頻等離子體風(fēng)洞中研究了碳/碳和碳/炭化鋯兩種復(fù)合防熱材料燒蝕情況下對下游流場的影響;聶春生等[13]在高頻等離子體風(fēng)洞中開展了純碳-碳材料燒蝕對球錐體等離子體流場的影響研究。數(shù)值研究方面,高鐵鎖等[14]進行了包含和不包含碳-酚醛燒蝕產(chǎn)物的再入小鈍錐流場計算,分析了燒蝕產(chǎn)物對流場電子數(shù)密度、溫度等流動參數(shù)的影響;魏叔如等[15]求解邊界層方程,開展了碳-碳材料燒蝕對電離邊界層影響的研究;董維中等[16]建立了碳-碳材料燒蝕流場的數(shù)值計算方法,研究了碳-碳材料燒蝕對再入體頭部區(qū)域溫度和熱流的影響,沒有對等離子體參數(shù)進行分析。
總之,國內(nèi)外針對高溫等離子體流場開展的試驗研究相對較少,主要通過數(shù)值模擬方法研究表面材料燒蝕產(chǎn)物對等離子體流場的影響,但數(shù)值方法存在多個難點:1)不同材料的燒蝕反應(yīng)過程不同,燒蝕產(chǎn)物對等離子體流場的影響也有較大差別,需要針對不同材料建立相應(yīng)的燒蝕反應(yīng)模型;2)由于表面催化、氧化、燒蝕與熱解等復(fù)雜的氣體表面相互作用,氣體表面邊界條件會非常復(fù)雜,表面邊界方程在流場中的求解算法十分復(fù)雜,難以建立高效、魯棒的耦合計算方法,使得針對復(fù)雜外形應(yīng)用難度較大;3)由于進入流場的燒蝕組分與高溫空氣組分差別較大,流場中的化學(xué)反應(yīng)過程會變得更加復(fù)雜,剛性問題和求解難度大大增加[17]。因此,目前文獻中的研究對象主要針對簡單的球錐再入體,材料主要針對傳統(tǒng)的泰氟隆材料和碳-碳材料。與球錐再入體不同,臨近空間飛行器一般會采用復(fù)雜的升力體外形,并長時間高馬赫數(shù)飛行,其熱環(huán)境特點為中高焓、低熱流,防隔熱設(shè)計往往會采用低密度的熱解炭化燒蝕材料[4],其燒蝕機理主要靠材料熱解吸熱和熱解氣體注入邊界層的熱阻塞效應(yīng)。調(diào)研發(fā)現(xiàn),目前針對升力體外形的等離子體流場,以及熱解炭化材料的燒蝕產(chǎn)物注入對等離子體流場的影響研究鮮有報道。
本文針對采用低密度熱解炭化防熱材料的臨近空間升力體飛行器,建立表面材料燒蝕與高溫等離子體流場數(shù)值計算方法,研究不同馬赫數(shù)下燒蝕產(chǎn)物注入流場對等離子體的影響規(guī)律。
求解三維熱化學(xué)非平衡N-S方程[18],守恒形式的方程如下:
式中:Q為守恒變量;F、G、H與Fv、Gv、Hv分別對應(yīng)不同方向的對流項與黏性項;W為非平衡源項。
采用有限體積法求解控制方程,無黏通量離散求解采用Roe平均方法與Yee的對稱TVD限制器,并對特征值采用了Harten的熵修正[19];黏性通量離散求解采用中心差分格式;時間推進采用全隱式LU-SGS方法,非平衡源項采用隱式處理[20]。
當飛行器表面溫度高于表面材料的燒蝕溫度門檻時就會發(fā)生燒蝕,此時流場中要考慮從壁面注入的熱解燒蝕產(chǎn)物的氣體組分;本文等離子體流場計算考慮了32種組分,包括N、O、N2、O2、NO、NO+、C2H、C2H2、CN、H、H2、HCN、CH、Si、SiO、SiO2、SiC、Si2、SiN,化學(xué)反應(yīng)模型見表1,其中組分的熱力學(xué)參數(shù)見文獻[21],輸運參數(shù)見文獻[22-23],化學(xué)反應(yīng)速率常數(shù)見文獻[24-25]。
表1 高溫空氣/材料燒蝕組分的化學(xué)反應(yīng)模型Tab.1 Chemical reaction model of material ablative components
本文研究定常的燒蝕流場,對應(yīng)于壁面溫度最高,固體材料內(nèi)部溫度分布達到穩(wěn)態(tài)的情況,避免了耦合熱傳導(dǎo)方程的求解,并且不考慮燒蝕引起的壁面外形變化對流場的影響。
假設(shè)固體碳處于與表面氣體達到平衡的狀態(tài),存在如下平衡約束條件[26]:
式中:ρw為壁面處氣體密度;cw,C為壁面碳原子的質(zhì)量濃度;MC為碳的摩爾質(zhì)量;KC,C為異相反應(yīng)C(固相)→C(氣相)的化學(xué)平衡常數(shù);Tw為壁面溫度。(2)式的求解需要知道壁面碳原子的質(zhì)量濃度cw,C,此碳原子質(zhì)量分數(shù)假設(shè)在壁面溫度、壓力和給定的元素比例條件下處于平衡。
化學(xué)平衡常數(shù)KC,C采用擬合關(guān)系式計算:
式中:T為氣體溫度。
對于某一元素k壁面處的質(zhì)量濃度可以從如下表面元素質(zhì)量平衡方程中獲得元素質(zhì)量分數(shù):
壁面壓力pw采用如下壁面法向動量獲得:
式中:vw為壁面法向速度,vw來自如下質(zhì)量連續(xù)方程,
式中:qC為壁面熱流;α為輻射吸收系數(shù);qr為壁面向氣體的輻射熱流;ε為材料輻射系數(shù);σ為斯蒂芬-玻爾茲曼常數(shù);hw為壁面焓值。
式中:ρv、ρC分別為原始材料密度和炭化層密度。
采用公開發(fā)表的電磁衰減測量項目C方案(RAM-C)飛行試驗等離子體測試結(jié)果,對本文計算方法進行驗證。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。本文主要考核流場中電子數(shù)密度分布,空氣化學(xué)反應(yīng)模型采用Park模型,考慮了N2、O2、NO、NO+、N、O、e-共計7種組分,熱力學(xué)模型采用兩溫度模型,壁面考慮完全非催化和完全催化兩種條件。
圖1給出了沿軸向位置電子數(shù)密度N峰值分布與飛行試驗結(jié)果對比,橫軸為軸向位置與頭部半徑的無量綱距離,縱軸為電子數(shù)密度。由圖1可以看出,本文電子數(shù)密度峰值分布計算結(jié)果與飛行數(shù)據(jù)符合較好,誤差在半個數(shù)量級以內(nèi)。圖2給出了典型截面處法向電子數(shù)密度分布與飛行試驗結(jié)果對比。由圖2可以看出:完全非催化和完全催化壁面條件會影響壁面附近的電子數(shù)密度,完全非催化壁結(jié)果高于完全催化壁結(jié)果,遠離壁面的空間等離子體流場受影響較小??傮w上,完全非催化壁結(jié)果更接近飛行試驗測量值,表明本文等離子體流場計算方法在電子數(shù)密度分布計算中有較好的準確性。
圖1 電子數(shù)密度峰值分布比較Fig.1 Comparison of peak distributions of electron number density
圖2 法向電子數(shù)密度分布比較Fig.2 Comparison of normal electron number density distributions
針對某升力體氣動布局的臨近空間飛行器前體部分,開展燒蝕產(chǎn)物對等離子體流場的影響分析。圖3給出了計算外形和網(wǎng)格示意圖,采用六面體結(jié)構(gòu)化網(wǎng)格,在近壁面附近進行了加密處理,網(wǎng)格總量150萬。計算來流工況為高度60 km、馬赫數(shù)分別為20、25、28,飛行攻角0°。
圖3 計算模型和網(wǎng)格Fig.3 Computational model and grid
為研究不同熱解炭化材料燒蝕產(chǎn)物對等離子體流場的影響,選取表2所示成分存在差異的同一類輕質(zhì)燒蝕防熱材料[27-28],其中材料1為文獻[27-28]中給出的阿波羅飛船采用的輕質(zhì)燒蝕材料;材料2和材料3分別在材料1基礎(chǔ)上進行了工藝改進,分別去除材料中的氫元素和硅元素,分別記為輕質(zhì)燒蝕材料2和輕質(zhì)燒蝕材料3。3種材料密度均取544.63 kg/m3,熱解燒蝕后的炭化層密度取256.29 kg/m3。
表2 材料參數(shù)Tab.2 Material parameters
3.2.1 純空氣等離子體流場規(guī)律分析
圖4給出了高度60 km、馬赫數(shù)27狀態(tài)的流場溫度和電子數(shù)密度分布。由圖4可以看出,高馬赫數(shù)來流在飛行器周圍形成了脫體激波,激波對等離子體流場的影響主要取決于激波后的流場溫度;對于頭部區(qū)域,激波強度大、壓縮性強,波后溫度超過了10 000 K,空氣發(fā)生大規(guī)模離解和電離反應(yīng),流場中的O2分子幾乎完全離解,O、N原子的電離反應(yīng)是整個流場中等離子體產(chǎn)生的主要源頭,該區(qū)域的電子數(shù)密度最高,峰值超過了1×1014cm-3;飛行器身部區(qū)域的激波強度較弱,激波后的流場溫度約為5 000 K,該區(qū)域的空氣離解和電離程度較弱,飛行器身部流場的等離子體主要由頭部流場輸運擴散而來;受激波層影響,飛行器身部區(qū)域的等離子體層主要貼近飛行器壁面,其中前緣部位的等離子體層厚度相對迎風(fēng)面較薄,隨著流動向下游發(fā)展,流場溫度逐漸降低,飛行器身部電子數(shù)密度比頭部低2~3個量級。
圖4 流場溫度和電子數(shù)密度分布Fig.4 Temperature and electron number density distribution of flow field
圖5給出了頭部駐點線流場溫度分布。由圖5可以看出,駐點線上的振動溫度與平動溫度差別較大,流場的熱力學(xué)非平衡效應(yīng)顯著,由于激波層內(nèi)高溫區(qū)的氣體電離反應(yīng)是等離子體層形成的主要原因,需要采用兩溫度模型進行等離子體流場計算。
圖5 沿駐點線溫度分布Fig.5 Temperature distribution along stagnation line
圖6給出了駐點線流場氣體電離組元的數(shù)密度分布。由圖6可以看出:馬赫數(shù)為20時,NO+離子是流場中電子數(shù)密度的主要貢獻源;馬赫數(shù)增大為25、27時,除了NO+離子,N+對電子數(shù)密度也有重要貢獻,其中NO+主要由O與N原子的締合電離反應(yīng)機制產(chǎn)生,N+主要由N與其他組分的碰撞電離反應(yīng)機制產(chǎn)生。綜合上述分析可知,對流場中電子數(shù)密度貢獻最重要的離子是NO+和N+,其次依次為這些離子對電子數(shù)密度的貢獻很小。
圖6 沿駐點線電離組元的數(shù)密度分布Fig.6 Number density distribution of ionized components along stagnation line
3.2.2 燒蝕產(chǎn)物對等離子流場影響分析
圖7給出了材料3在高度60 km、馬赫數(shù)為20狀態(tài)的表面溫度和燒蝕質(zhì)量流率分布。由圖7可以看出:飛行器頭部和前緣溫度超過了2 000 K,該區(qū)域發(fā)生明顯燒蝕,其中頭部燒蝕最嚴重,燒蝕質(zhì)量流率超過了0.2 kg/(m2·s);飛行器身部迎風(fēng)面的溫度約為1 000 K,燒蝕質(zhì)量流率相對較小,沿流向燒蝕效應(yīng)逐漸減弱。
圖7 表面溫度和燒蝕質(zhì)量流率(高度60 km、馬赫數(shù)20,材料3)Fig.7 Surface temperature and ablation mass flow rate(height of 60 km,Mach number of 20,and Material 3)
圖8給出了流場中主要燒蝕組分的數(shù)密度。由圖8可見:注入流場的主要熱解燒蝕產(chǎn)物是CO和H2,其他燒蝕氣體組分的質(zhì)量分數(shù)相對較小,并主要集中在頭部區(qū)域;下游飛行器身部的燒蝕效應(yīng)逐漸減弱,但隨著流動向下游發(fā)展,頭部燒蝕產(chǎn)物會隨流動的擴散效應(yīng)輸運到下游,燒蝕效應(yīng)對流場等離子分布的影響范圍擴大。
圖8 對稱面燒蝕產(chǎn)物分布(高度60 km、馬赫數(shù)20,材料3)Fig.8 Distribution of ablated products on symmetry plane(height of 60 km,Mach number of 20,and Material 3)
圖9給出了頭部駐點線的流場參數(shù)對比。由圖9可以看出:燒蝕產(chǎn)物注入流場后,流場基本結(jié)構(gòu)沒有發(fā)生明顯變化,溫度和電子數(shù)密度峰值未受影響;由于壁面燒蝕的引射效應(yīng),使得頭部激波脫體距離增大,駐點線上的溫度、電子數(shù)密度和組分質(zhì)量分數(shù)的分布也發(fā)生了一定變化,流場中高溫區(qū)域變大,O2和N2分子更多地分解成O和N原子,O原子參加氧化反應(yīng)被消耗而質(zhì)量分數(shù)明顯減小,消耗的O原子主要生成了CO。
圖9 頭部軸線流場參數(shù)分布比較(高度60 km、馬赫數(shù)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在不同馬赫數(shù)下駐點線的溫度和電子數(shù)密度曲線。由圖10和圖11可以看出,對于材料3,隨著飛行馬赫數(shù)變大,表面燒蝕速率增強,燒蝕產(chǎn)物進入流場使得激波脫體距離增大,但駐點線上的電子密度峰值變化很小,這主要是因為流場中電子數(shù)密度貢獻最重要的離子是NO+。圖12給出了材料3在不同馬赫數(shù)駐點線NO+離子數(shù)密度分布。從圖12中可以看出,不同馬赫數(shù)下燒蝕效應(yīng)對頭部駐點線上的NO+離子數(shù)密度影響相對較小導(dǎo)致的。
圖10 不同馬赫數(shù)下駐點線溫度分布(材料3)Fig.10 Distribution of temperature along stagnation line at different Mach numbers(Material 3)
圖11 不同馬赫數(shù)下駐點線電子數(shù)密度分布(材料3)Fig.11 Distribution of electron number density along stagnation line at different Mach numbers(Material 3)
圖12 不同馬赫數(shù)駐點線NO+離子數(shù)密度分布(材料3)Fig.12 Distribution of NO+ion number density along stagnation line at different Mach numbers(Material 3)
圖13給出了材料3的對稱面電子數(shù)密度峰值對比曲線。由圖13可以看出:隨著馬赫數(shù)增大,表面燒蝕速率增大,由于進入流場的燒蝕產(chǎn)物增多,使得流場中的電子數(shù)密度有一定的增大;隨著流動從頭部向下游發(fā)展,飛行器身部的熱流會逐漸降低,表面材料的燒蝕速率逐漸減小,但燒蝕產(chǎn)物對飛行器身部流場電子數(shù)密度的影響逐漸增強;下游0.5 m的位置,材料3的燒蝕產(chǎn)物對電子數(shù)密度的影響增量約為20%~30%。
圖13 對稱面電子數(shù)密度峰值對比(材料3)Fig.13 Comparison of electron number density peaks on symmetry plane(Material 3)
圖14給出了表2中不同成分的低密度防熱材料在高度60 km、馬赫數(shù)25狀態(tài)下飛行器表面中心線的燒蝕質(zhì)量流率對比曲線。由圖14可以看出:由于飛行器頭部的氣動加熱最為嚴酷,使得燒蝕質(zhì)量流率最大,材料1和材料3接近,約為0.55 kg/(m2·s),材料2為0.33 kg/(m2·s);與頭部相比,飛行器身部的燒蝕質(zhì)量流率會迅速降低約2~3個量級,表2中3種材料的燒蝕質(zhì)量流率均低于0.013 kg/(m2·s),其中材料2要略高于材料1和材料3。
圖14 不同材料飛行器表面對稱線的燒蝕質(zhì)量流率Fig.14 Ablation mass flow rate of surface center lines of aircraft with different materials
圖15、圖16分別為高度60 km、馬赫數(shù)25狀態(tài)下駐點線溫度和電子數(shù)密度對比曲線。由圖15和圖16可以看出:該狀態(tài)不同材料的燒蝕效應(yīng)存在差異,從而會影響頭激波的脫體距離,但駐點線的峰值溫度和峰值電子數(shù)密度幾乎不受影響;駐點線脫體距離從大到小依次為材料1、材料2和材料3;材料1和材料2的燒蝕質(zhì)量流率接近,但由于材料的成分不同,注入流場的燒蝕產(chǎn)物不同,使得兩種材料燒蝕后對流場激波脫體距離的影響程度存在差異。
圖15 不同材料條件下駐點線溫度分布Fig.15 Distribution of temperatures of different materials along stagnation line
圖16 不同材料條件下駐點線電子數(shù)密度分布Fig.16 Distribution of electron number density along stagnation line under different material conditions
圖17給出了高度60 km、馬赫數(shù)25狀態(tài)下對稱面電子數(shù)密度峰值對比曲線。由圖17可以看出:對稱面材料的燒蝕產(chǎn)物對流場峰值電子數(shù)密度影響沿流動向下游發(fā)展逐漸增強,并且不同材料燒蝕產(chǎn)物對流場中的電子數(shù)密度峰值影響量不同,其中材料2的影響相對最大;在0.5 m位置,電子數(shù)密度的影響增量約為300%~320%。結(jié)合圖15表面燒蝕質(zhì)量流率分析可知,飛行器身部流場的電子數(shù)密度一方面受頭部區(qū)域燒蝕產(chǎn)物流動擴散的影響,另一方面還要受當?shù)氐奈⒘繜g產(chǎn)物的影響,二者均會使得流場中的電子數(shù)密度增大,其影響機理為燒蝕產(chǎn)物中CH、CO會發(fā)生反應(yīng)產(chǎn)生O原子,O原子與N原子發(fā)生締合電離反應(yīng),產(chǎn)生更多NO+。圖18所示為對稱面x=0.4 m位置流場中的NO+離子對比。由圖18可以看出,材料燒蝕產(chǎn)物注入后流場中的NO+離子顯著增大,從而使得流場中的電子數(shù)密度均比純空氣狀態(tài)增大。
圖17 不同材料條件下對稱面電子數(shù)密度峰值對比Fig.17 Comparison of peak values of electron number density on symmetry planes of different materials
圖18 x=0.4 m位置流場NO+離子數(shù)密度對比曲線Fig.18 Comparison curves of NO+ion number density in flow field for x=0.4 m
本文建立了考慮燒蝕產(chǎn)物影響的等離子體流場數(shù)值計算方法,針對阿波羅飛船采用的輕質(zhì)燒蝕材料,開展了升力體外形的燒蝕產(chǎn)物對等離子體流場的影響規(guī)律研究。得到如下主要結(jié)論:
1)對于升力體飛行器的等離子體流場,影響最重要的離子是NO+和N+。
2)熱解炭化材料最主要的熱解燒蝕產(chǎn)物是CO和H2;燒蝕產(chǎn)物進入流場會使得激波脫體距離增加,等離子體厚度增加。
3)燒蝕產(chǎn)物進入流場,對頭部的峰值電子數(shù)密度影響很小,隨著流動向下游發(fā)展影響逐漸增大;影響機理為,燒蝕產(chǎn)物中CH、CO組分會發(fā)生反應(yīng),產(chǎn)生O原子,O原子與N原子發(fā)生締合電離反應(yīng),產(chǎn)生更多NO+,使得流場中電子數(shù)密度增大。
4)成分存在差異的同一類防熱材料,燒蝕產(chǎn)物對流場電子數(shù)密度的影響增量存在差異;本文研究的一類防熱材料,高度60 km、馬赫數(shù)為25狀態(tài)下飛行器下游0.5 m位置對稱面峰值電子數(shù)密度,材料2的增量約為300%~320%。