楊鵬飛, 方洋旺, 柴棟, 高翔
(空軍工程大學(xué) 航空航天工程學(xué)院, 陜西 西安 710038)
高超聲速鈍錐體等離子體流場數(shù)值模擬
楊鵬飛, 方洋旺, 柴棟, 高翔
(空軍工程大學(xué) 航空航天工程學(xué)院, 陜西 西安 710038)
以小鈍錐體為對象建立二維幾何模型,利用COMSOL Multiphysics軟件對等離子體流場控制方程組進(jìn)行求解,對高超聲速氣流流過模型所產(chǎn)生的等離子體流場進(jìn)行數(shù)值模擬。仿真得到了在模型頭部激波、尾部膨脹波和尾后激波的共同作用下,考慮燒蝕反應(yīng)時(shí)流場的速度、溫度、壓強(qiáng)和等離子體電子數(shù)密度分布,重點(diǎn)分析了鈍錐體模型等離子體包覆層和尾后電子分布的特征和成因。可以為高超聲速飛行器的探測和控制提供重要的借鑒和參考。
高超聲速飛行器; 等離子體流場; 流場模擬; 數(shù)值模擬
由于強(qiáng)烈的氣動(dòng)加熱,高超聲速飛行器表面具有很高的溫度,因此飛行器表面的燒蝕材料通過高溫氣體與物面材料發(fā)生的氣固熱化學(xué)反應(yīng)和質(zhì)量交換,吸收大部分氣動(dòng)加熱能量,從而使飛行器表面保持在允許的溫度,有效地保護(hù)基體。在高溫、高熱流的作用下,燒蝕材料將燒蝕產(chǎn)物釋放到邊界層,又與流場中的高溫空氣組分進(jìn)一步進(jìn)行復(fù)雜化學(xué)反應(yīng)。燒蝕產(chǎn)物中含有少量的堿金屬雜質(zhì),其電離勢比空氣中任何組分的電離勢低得多,極易電離產(chǎn)生大量電子,進(jìn)一步影響了流場的電離度。
當(dāng)繞流流場的電離度達(dá)到一定程度時(shí),在飛行器周圍會(huì)形成含有大量自由電子的等離子體包覆流場,即等離子體鞘。通常把這種現(xiàn)象稱為真實(shí)氣體效應(yīng)[1]。上述這些現(xiàn)象使氣體的物理特性和完全氣體有明顯的差別,針對該特殊物理化學(xué)現(xiàn)象,通過對高速飛行器等離子體流場進(jìn)行建模仿真來研究其特性,對于高超聲速飛行器的探測和控制具有重要的意義。
1.1 流場控制方程
根據(jù)文獻(xiàn)[1],在直角坐標(biāo)系下,二維鈍體繞流的流動(dòng)方程如下:
連續(xù)方程:
(1)
x方向動(dòng)量方程:
(2)
y方向動(dòng)量方程:
(3)
能量方程:
(4)
式中,u,v分別為氣流速度在x,y方向上的分量;γ為氣體比熱比,空氣的比熱比取為1.4。
上述方程組共有p,ρ,u,v四個(gè)未知量,求解該方程組還必須利用狀態(tài)方程,對于完全氣體,其狀態(tài)方程為:
p=ρRT
(5)
式中,R為氣體常數(shù)。引入新的變量T,因此壓強(qiáng)比與溫度比的關(guān)系為:
dT/T=(γ-1)dp/p
(6)
1.2 邊界層流動(dòng)控制方程
要考慮燒蝕作用,必須單獨(dú)分析邊界層的流動(dòng)方程,二維定常粘性流動(dòng)的邊界層流動(dòng)方程有量綱形式如下:
連續(xù)方程:
?(ρu)/?x+?(ρv)/?y=0
(7)
x方向動(dòng)量方程:
(8)
y方向動(dòng)量方程:
?p/?y=0
(9)
能量方程:
(10)
式中,h為焓;pe為邊界層外緣壓力分布,由于在邊界層內(nèi),沿物體表面的法線方向,包括物面在內(nèi),壓強(qiáng)都是常數(shù)[2],所以可以采用邊界層外緣壓力分布pe(x)來作為當(dāng)?shù)氐膲毫Ψ植肌?/p>
邊界層流動(dòng)控制方程組中,由于p=pe(x)是已知的,未知數(shù)有u,v,h,ρ,對于完全氣體,有:
p=ρRT
(11)
h=cpT
(12)
則式(7)~式(12)可以求解,此時(shí)的邊界條件為:
在壁面:
u=0,v=vw, (?T/?n)w=0
(13)
在邊界層外緣:
u→ue,T→Te
(14)
式中,由于燒蝕作用,氣流沿壁面法線方向速度v為給定速度vw;壁面采用絕熱壁邊界條件;下標(biāo)w,e分別表示壁面和邊界層外緣。
1.3 真實(shí)氣體效應(yīng)和尾跡計(jì)算
真實(shí)氣體效應(yīng)對于正確地預(yù)測高超聲速飛行器的氣動(dòng)特性具有重要作用,文獻(xiàn)[3]中給出了常用的燒蝕材料中堿金屬的含量,文獻(xiàn)[4]給出了不同含量的Na對電子數(shù)密度分布的影響。從中都可以看出堿金屬燒蝕材料對高超聲速等離子體流場的影響??諝饣瘜W(xué)模型采用11組分26個(gè)化學(xué)反應(yīng)的Dunn-Kang模型,模型見文獻(xiàn)[5]。燒蝕反應(yīng)仿真采用文獻(xiàn)[6]中碳基復(fù)合材料燒蝕化學(xué)反應(yīng)。
高超聲速飛行器在大氣中飛行時(shí)由于粘性摩擦,再入體表面形成粘性附面層,由附面層在再入體后部形成一個(gè)粘性尾跡的內(nèi)核[5],尾跡結(jié)構(gòu)及電子密度的計(jì)算公式見文獻(xiàn)[7]。
2.1 計(jì)算模型與網(wǎng)格劃分
利用COMSOL Multiphysics建立仿真模型。高超聲速飛行器為了減少頭部熱流,一般采用鈍的形狀[1],因此采用小鈍錐體為對象建立二維幾何模型。
以現(xiàn)有資料的高超聲速飛行器尺寸為依據(jù)[8],設(shè)定錐體長3.6 m,底部直徑為1.6 m。計(jì)算網(wǎng)格采用非結(jié)構(gòu)網(wǎng)格,在邊界層及邊界點(diǎn)處網(wǎng)格劃分更加精細(xì)。仿真模型網(wǎng)格劃分如圖1所示,網(wǎng)格數(shù)為346 582個(gè)。
圖1 流場計(jì)算非結(jié)構(gòu)網(wǎng)格Fig.1 Non-structured grids of flow calculation
2.2 物理設(shè)定
等離子體流場仿真的物理設(shè)定包括三個(gè)部分:高馬赫數(shù)流體模塊、流固耦合模塊和等離子體模塊。
高馬赫數(shù)流場計(jì)算首先在模型計(jì)算中設(shè)定流場計(jì)算的控制方程、氣體狀態(tài)方程和相關(guān)關(guān)系式。然后要設(shè)定邊界層的控制方程、狀態(tài)方程和邊界條件。在氣流經(jīng)過模型后的坐標(biāo)范圍還需加入尾跡以及尾跡中電子密度的計(jì)算公式。
計(jì)算初始條件和邊界條件為:高度H=40 km,Ma=7.5,Re=5×105;流場進(jìn)口靜壓Pin=105 Pa,進(jìn)口靜溫Tin=216 K,流場出口靜壓為Pout=95 Pa。
流固耦合模塊設(shè)定邊界條件和進(jìn)出口初始條件與高馬赫數(shù)模塊一致,還需選定流體與固體范圍。
在等離子體模塊中,用1.3節(jié)中所述的空氣反應(yīng)模型以及燒蝕反應(yīng)來代替其中默認(rèn)的氣體反應(yīng)。其他設(shè)定與高馬赫數(shù)模塊設(shè)定保持一致。
流場馬赫數(shù)分布如圖2所示。圖中x為仿真區(qū)域長度,H為仿真區(qū)域高度。仿真結(jié)果顯示,高超聲速氣流流過模型時(shí),產(chǎn)生了一道附體激波,氣流經(jīng)過激波后速度產(chǎn)生明顯的下降。在模型尾部,由于氣流通道突然增大,在模型兩側(cè)產(chǎn)生對稱的兩道膨脹波。氣流經(jīng)過膨脹波后速度增大,流場氣壓減小。流場氣壓的減小會(huì)造成模型尾后軸線氣流氣壓的壓縮作用,在模型尾后會(huì)產(chǎn)生兩道對稱的尾后激波。尾后氣流在經(jīng)過此激波后,氣流速度再次產(chǎn)生下降。所以在圖2中尾后氣流速度會(huì)呈現(xiàn)兩側(cè)低、中間高的三個(gè)區(qū)域。緊挨鈍錐體的后部產(chǎn)生回流區(qū),回流區(qū)中氣流速度很低,與外側(cè)氣流共同作用使得軸線上氣流速度低于來流速度。仿真結(jié)果與文獻(xiàn)[6]的試驗(yàn)結(jié)果一致,并且通過仿真得到了尾后遠(yuǎn)流場的馬赫數(shù)分布情況。
圖2 流場馬赫數(shù)分布圖Fig.2 The distribution of flow mach number
圖3為流場溫度分布仿真結(jié)果??梢钥闯?鈍錐體表面溫度很高,這是由于在仿真中加入的燒蝕化學(xué)反應(yīng)燒蝕造成,通過高溫氣體與物面材料發(fā)生的氣固熱化學(xué)反應(yīng)和質(zhì)量交換,吸收大部分氣動(dòng)加熱能量,從而使飛行器內(nèi)部保持在允許的溫度,可以有效的保護(hù)基體。氣流經(jīng)過附體激波后,溫度急劇上升。流過附體激波的氣流在流出模型時(shí)經(jīng)過膨脹波的作用,溫度轉(zhuǎn)而下降。隨后氣流流經(jīng)尾后激波,氣流溫度突躍上升。仿真結(jié)果的變化規(guī)律與文獻(xiàn)[9]的計(jì)算結(jié)果基本吻合,溫度值略高,這主要是因?yàn)榧尤肓藟A金屬的燒蝕作用,使得飛行器本身的溫度降低,而周圍大氣溫度略有升高。
圖3 流場溫度分布圖Fig.3 The distribution of flow temperature
流場的氣壓分布結(jié)果如圖4所示,從圖中可見,流場氣壓分布與流場速度、溫度的分布相似,在附體激波、尾部膨脹波和尾后激波的共同作用下,使得流場的氣壓分布對稱于軸線。氣流經(jīng)過附體激波氣壓明顯增大,而在經(jīng)過尾部膨脹波后氣壓略有減小,經(jīng)過尾后激波后,氣壓得到小范圍升高。
圖4 流場氣壓分布圖Fig.4 The distribution of flow pressure
流場的等離子體電子數(shù)密度分布如圖5所示。從中可以看出,高超聲速飛行器周圍形成了等離子體包覆層,且在弓形激波后尤為明顯,這是因?yàn)榧げê蟮母邷睾惋w行器表面邊界層中的高溫可以激發(fā)氣體分子的振動(dòng)引起離解和電離,同時(shí)燒蝕產(chǎn)物會(huì)與高溫氣體組分進(jìn)行反應(yīng),進(jìn)一步增加了這一區(qū)域的等離子體電子數(shù)密度。
繞流流場產(chǎn)生的等離子體是非均勻的,物面附近、邊界層內(nèi)電子數(shù)密度高,激波層外電子數(shù)密度很小,在物面法向方向等離子體參數(shù)變化梯度明顯比沿流向的等離子體參數(shù)梯度大。在經(jīng)過尾部激波后,形成的等離子體電子數(shù)密度比較高,這是因?yàn)闊g產(chǎn)物在尾后氣流的高溫、高熱流作用下,與流場中的高溫空氣組分進(jìn)一步進(jìn)行復(fù)雜化學(xué)反應(yīng)。其中由于燒蝕產(chǎn)物中含有少量的堿金屬雜質(zhì),其電離勢比空氣中任何組分的電離勢低的多,極易電離產(chǎn)生大量電子。尾后氣流軸線方向電子數(shù)密度較低,則是因?yàn)檫@片區(qū)域燒蝕產(chǎn)物含量較少,僅僅依靠空氣的離解很難產(chǎn)生密度較高的電子。
圖5 流場等離子體電子數(shù)密度分布圖Fig.5 The distribution of plasma electron number density
尾跡軸心電子數(shù)密度如圖6所示,可見電子數(shù)密度的最高點(diǎn)出現(xiàn)在模型頭部,范圍非常小,尾跡軸心電子數(shù)密度隨著尾跡長度基本呈指數(shù)規(guī)律變化,整個(gè)過程中的電子數(shù)密度的范圍在106~1012cm-3數(shù)量級之內(nèi)。仿真結(jié)果與文獻(xiàn)[10]中試驗(yàn)結(jié)果一致,且數(shù)值更為詳細(xì),并且是加入燒蝕作用以后的結(jié)果,更符合飛行器實(shí)際飛行過程。
圖6 尾跡軸心電子數(shù)密度分布Fig.6 The distribution of electron number density of wake axis
(1)高超聲速氣流在鈍錐體模型頭部形成弓形脫體激波,波后氣體溫度、壓強(qiáng)都急劇升高。
(2)高超聲速氣流流經(jīng)鈍錐體尾部時(shí),由于氣流通道突然加寬形成膨脹波,尾后氣流由于氣壓原因形成尾后激波,在膨脹波和尾后激波的共同作用下,流場的溫度、壓強(qiáng)高低呈現(xiàn)沿尾流軸線對稱的相間分布。
(3)高超聲速飛行器周圍形成了等離子體包覆層,并且在弓形激波后更加明顯。這是因?yàn)榧げê蟮母邷睾惋w行器表面邊界層中的高溫激發(fā)了氣體分子的振動(dòng),引起離解和電離,同時(shí)燒蝕產(chǎn)物與流場中的高溫空氣組分進(jìn)行反應(yīng)使該處的等離子體電子數(shù)密度高于其他區(qū)域。
(4)繞流流場產(chǎn)生的等離子體是非均勻的,物面附近、邊界層內(nèi)電子數(shù)密度高,激波層外電子數(shù)密度很小,在物面法向方向等離子體參數(shù)變化梯度明顯比沿流向的等離子體參數(shù)梯度大。
(5)在經(jīng)過尾部激波后,形成的等離子體電子數(shù)密度比較高,這是因?yàn)闊g產(chǎn)物中所含有少量的堿金屬雜質(zhì)的電離勢比空氣中任何組分的電離勢低得多,極易電離產(chǎn)生大量電子。尾后氣流軸線方向電子數(shù)密度較低,這是因?yàn)檫@片區(qū)域燒蝕產(chǎn)物含量較少,僅僅依靠空氣的離解很難產(chǎn)生密度較高的電子。
[1] 黃志澄.高超聲速飛行器空氣動(dòng)力學(xué)[M].北京:國防工業(yè)出版社,1995:3-4,293-294.
[2] 錢翼稷.空氣動(dòng)力學(xué)[M].北京:北京航空航天大學(xué)出版社,2004:81.
[3] 李森茂.燒蝕材料中堿金屬和堿土金屬含量的原子吸收測試方法的研究[R].北京:航天材料及工藝研究所,1991.
[4] Bhutta B A,Lewis C H.A new technique for low-to-high altitude predictions of ablative hypersonic flowfields[R].AIAA-91-1392,1991.
[5] 樂嘉陵,高鐵鎖,曾學(xué)軍,等.再入物理[M].北京:國防工業(yè)出版社,2005:31-33,259-267.
[6] 張志成,潘梅林,劉初平,等.高超聲速氣動(dòng)熱和熱防護(hù)[M].北京:國防工業(yè)出版社,2003:234-237.
[7] 吳建明,高本慶.再入飛行物等離子體尾跡RCS的計(jì)算[J].電波科學(xué)學(xué)報(bào),1997,12(1):26-32.
[8] 郭朝邦,李文杰.高超聲速飛行器結(jié)構(gòu)材料與熱防護(hù)系統(tǒng)[J].飛航導(dǎo)彈,2010,35(4):88-94.
[9] 劉青云,安冬梅,張平.超音速再入體遠(yuǎn)場湍流尾跡電子密度計(jì)算[J].北京理工大學(xué)學(xué)報(bào),1999,19(S1):115-119.
[10] 馬平,曾學(xué)軍,柳森,等.高超聲速球模型尾跡電子密度試驗(yàn)研究[J].實(shí)驗(yàn)流體力學(xué),2010,24(1):20-25.
Numericalsimulationofhypersonicplasmaflowfieldsaroundbluntconemodel
YANG Peng-fei, FANG Yang-wang, CHAI Dong, GAO Xiang
(Aeronautics and Astronautics Engineering College, Air Force Engineering University, Xi’an 710038, China)
The control equations of plasma flow fields are solved and the flow fields around the blunt cone model are simulated numerically with COMSOL Multiphysics. Considering the joint influences of head shock wave, the tail expansion wave and the tail shock wave, the distribution of flow fields, pressure, temperature and the electron number density of plasma are simulated under ablation effect. It mainly analyzes the characteristics and the causes of plasma sheath and the electron distribution of blunt cone model. By doing this research, an important reference can be provided for the detection and control of hypersonic vehicles.
hypersonic vehicle; plasma flow fields; flow fields simulation; numerical simulation
V211.7
A
1002-0853(2013)04-0313-04
2012-11-29;
2013-03-21; < class="emphasis_bold">網(wǎng)絡(luò)出版時(shí)間
時(shí)間:2013-06-06 13:21
國家自然科學(xué)基金資助(60874040)
楊鵬飛(1988-),男,山西侯馬人,碩士研究生,研究方向?yàn)楦叱曀倌繕?biāo)的分析與探測。
(編輯:方春玲)