磁控熱防護系統在天地往返運載器上的應用仿真

丁明松,劉慶宗,江濤,董維中,高鐵鎖,傅楊奧驍

中國空氣動力研究與發展中心 計算空氣動力研究所, 綿陽 621000

從1981年4月12日美國“哥倫比亞”號(OV-102)航天飛機進行首次載人軌道飛行到2011年7月8號“阿特蘭蒂斯”號航天飛機第135次升空，航天飛機完成了一次又一次航天壯舉[1]。盡管持續30年的航天飛機時代已經宣告結束，但人們從未放棄滑翔返回式天地往返運載器這一尖端領域，相關研究正不斷推動航天科學關鍵技術發展[2]。

技術瓶頸是制約航天項目跨越式發展的關鍵。如氣動熱防護這類技術瓶頸問題得不到有效解決，航天型號的使用依然代價巨大。天地往返系統高速飛行時，面臨的熱障問題尤為突出，從未來發展趨勢來看，運載器性能要求越來越高，氣動熱環境越來越嚴酷，給飛行器熱防護系統設計帶來極大的挑戰[2]。

通過配置合適的機載電磁場系統，對高超聲速流動進行控制，改善飛行器氣動熱環境，這就是高超聲速飛行器磁控熱防護技術[3]。對于高超聲速流動，速度較高時，波后氣體熱電離產生等離子體具有足夠的電導率，可直接應用于磁流體力學控制[4]；同時機載磁場通過電流響應控制，并安裝于飛行器內部，無需改變飛行外部結構，因而具備響應快速、不改變氣動外形、不影響其他被動熱防護措施等優點[5]，滿足了未來高超聲速飛行器高效穩定、可重復使用的要求。

自20世紀90年代以來，得益于稀土磁工業的發展，易于制造便攜的機載磁系統[4,6]，磁控熱防護技術掀起了新的研究熱潮。美國空軍研究實驗室Bisek團隊基于低磁雷諾假設開展了帶尾舵鈍頭橢圓錐體磁流體數值研究[7]，分析了外磁場對尾舵前緣熱流的影響，但由于前緣局部電導率較低(<0.1 S/m)，磁控熱防護效果不明顯；俄羅斯高溫研究所Bityurin和Bocharov[8]基于完全氣體模型和偶極子磁場配置方案，開展返回艙簡化外形磁控熱防護數值模擬，發現通過配置磁場可有效降低駐點區域熱流；日本Otsu[9]、Fujino[10]等開展了多種條件下簡單鈍體外形磁流體力學控制數值模擬研究，發現絕緣壁面條件霍爾效應影響幾乎可以忽略，磁控熱防護效果較為顯著；德國宇航中心和歐洲航天局共同開展了簡單鈍體外形磁控熱防護試驗[11]，發現采用螺線管磁場，可使平頭圓柱熱流下降85%；國防科技大學李開和劉偉強[12]針對簡單鈍體外形開展了磁控熱防護系統建模研究，發現增大螺線管半徑有利于提高磁控熱防護效果，而螺線管長度對磁控熱防護效果影響相對較小。

盡管國內外對磁控熱防護系統進行了不少研究，熱防護效果得到了有效驗證，但很少見到磁控熱防護系統在實際復雜外形飛行器的工程應用。這主要是由于高超聲速流動磁流體力學控制涉及高溫氣體化學反應效應、熱力學效應、等離子體碰撞與遷移機制、帶電粒子與電磁場相互作用、電磁場交叉感應等多種物理效應或現象[13-16]，每一種效應或現象包含復雜的變化規律和多種影響因素，對磁控熱防護效果造成不容忽視的影響；各種物理機制的特征尺度或特征時間存在的差異較大，多種物理現象耦合時會出現強烈的剛性問題和收斂性問題，極大地增加了仿真模擬難度。基于此，國內外高超聲速流動磁流體力學控制研究大多針對簡單外形，如球頭或鈍體外形的一維/二維或軸對稱磁流體數值模擬；物理模型方面也多以簡化模型為主，如完全氣體或平衡氣體模型等；磁場方面多采用簡化磁場，如偶極子或均勻磁場等，極少見到復雜外形飛行器磁控熱防護系統仿真研究。

本文綜合考慮高超聲速流動磁流體力學控制涉及的多種物理效應或現象，通過耦合求解電磁場泊松方程和帶電磁源項的高溫熱化學非平衡流動控制方程組，搭建磁控熱防護仿真模擬平臺。結合美國航天飛機 “哥倫比亞”號近似外形和多種磁場配置方案(均勻磁場、偶極子磁場、螺線管磁場及多個磁場組合)，較為系統地開展磁控熱防護系統在滑翔返回式天地往返運載器上的應用仿真研究。

1 等離子體生成機制及其電導率計算

(1)

基于電子/離子/中性粒子擴散、碰撞與遷移的原理，多電離組分等離子體電導率可寫為[13]

(2)

(3)

當s為中性組分O2、N2、NO、O和N時，其表達式為[10,14]

(4)

式(3)和式(4)中電子溫度Te的模擬，需要考慮等離子體流動過程中的熱力學效應，氣體分子平轉動溫度和振動-電子溫度之間的能量輸運可寫為[19]

(5)

2 電磁場對流場動量和能量輸運過程的模擬

具有弱導電性的等離子體在磁場中流動時，會產生感應電流，使流體微元受到洛倫茲力作用，這會改變氣體的動量和能量。與此同時，感應電流還會產生焦耳熱作用，進而改變氣體內能等狀態參數。由于在絕緣壁面條件下，霍爾效應對磁控效果的削弱影響相對較小[9-10,20]，因此為了增強磁控效果，飛行器表面可采用絕緣材料(或涂層)，此時忽略霍爾效應，電磁場對流場動量和能量輸運可寫為

(6)

式中：F為洛倫茲力；J為感應電流；B為磁感應強度；U為速度矢量；F·U為洛倫茲力做功；σ-1J2為焦耳熱。由于電磁焦耳熱的部分能量作用于振動-電子溫度模態，因此電磁場對振動-電子溫度模態的能量輸運可寫為

EMV=γσ-1J2

(7)

式中：γ為焦耳熱振動-電子能量配比。

感應電流J可通過離散求解各向同性電場泊松方程計算[16]：

(8)

磁感應強度按照磁場矢量疊加原理計算：

B=Bext+Bin

(9)

式中：Bext為外加磁場，由機載磁場發生裝置產生；Bin為感應磁場(或稱誘導磁場)，由感應電流誘導產生。對于磁雷諾數較低的狀態有|Bin|?|Bext|，因此可近似認為B=Bext。而對于磁雷諾數較大的狀態，Bin通過離散求解磁矢量勢泊松方程計算[10]：

(10)

3 機載磁場Bext的模擬

機載磁場的模擬，最常見的有均勻磁場、偶極子磁場和螺線管磁場。其中均勻磁場模擬較為簡單，一般是在飛行器流場某一區域，直接附加強度和方向相同的磁場。偶極子磁場是螺線管磁場的簡化模型，其表達式一般可寫為[15]

(11)

式中：(r,α)為極坐標單位矢量；B0為極軸上距離偶極子中心r0處磁感應強度。

對于實際的螺線管線圈磁場，其計算較為復雜，計算量較大。但如果忽略電流分布的非均勻性和線圈軸向的螺旋性，不僅可以極大地減少計算量，而且計算與實際測量結果之間的差異極小[21]。此時，通電電流I、線圈密度n、直徑2a、長度2b的螺線管磁場可寫為[21-22]

(12)

式中：坐標原點為螺線管中心；Br、Bφ和Bz為磁場徑向、周向和軸向分量，部分變量表達式為

其中：A為特征強度；k、k1、k2和h均為擬合函數中間變量；f(k)、fz(r,z)分別為r方向和z方向磁場擬合函數；K(k)、F(h,k)為第一類和第三類橢圓積分。

4 高超聲速電磁流動控制方程及離散

高超聲速電磁流動控制方程與一般的高超聲速熱化學非平衡Navier-Stokes(N-S)方程的主要區別在于包含電磁作用源項，其無量綱形式為[14]

(13)

式中：守恒變量Q=[ρi,ρ,ρu,ρv,ρw,ρEt,ρEv]T,ρi為組分i氣體密度,u、v、w分別為x、y、z3個方向速度分量，Et為氣體內能，Ev為氣體振動-電子能；Re為雷諾數；F、G、H與FV、GV、HV分別為x、y、z方向無黏向量與黏性向量；W和WMHD分別為非平衡源項和電磁作用源項，其形式為

(14)

其中：fx、fy、fz為洛倫茲力F3個方向分量。

對于式(13)，采用AUSMPW+(Advection Upstream Splitting Method by Pressure-based Weight functions)格式離散無黏項，采用中心差分格式離散黏性項，采用全隱式全耦合LU-SGS(Lower-Upper Symmetric Gauss Seidel)進行時間離散。采用MPI(Message Passing Interface)并行技術，應用負載平衡、非阻塞發送/阻塞接收、整合打包/解包以及對等/主從模式相組合等方法，建立流場分區和功能分解結合的并行算法，提升計算效率。

5 數值驗證

采用鈍錐體外形[23]對數值計算進行校驗。頭部半徑為0.09 m，椎體的錐角為15°，全長為0.122 4 m；計算飛行高度40.0 km、馬赫數15.0的狀態，等溫壁面條件。采用與文獻[23]相同氣體模型和磁偶極子配置，偶極子中心為頭部圓心，B0=2.0 T,r0=0.09 m。M1偶極子方向為X軸負方向；M2偶極子方向為Y軸正方向，見圖1(a)。

圖1(b)給出了磁場配置對鈍錐體表面熱流Q影響，S為表面弧長?？梢?，3種磁場配置條件下，計算結果均與文獻[23]符合，這說明本文數值模擬方法可以較為準確地預測和評估高超聲速磁控氣動熱防護效果。

圖1 磁場矢量和鈍錐表面熱流

6 磁控熱防護數值仿真

6.1 計算狀態、網格及氣動熱環境模擬

以航天飛機(OV-102)近似外形為研究對象，以飛行器頭部尖點為坐標原點，飛行器長為X軸方向，飛行器高為Y軸方向，飛行器橫向為Z軸方向，建立直角坐標系。采用多塊結構網格，在壁面附近以及形變較大區域，如豎直舵、水平翼等區域，進行局部加密，如圖2(a)所示。計算狀態采用典型軌道參數條件，見表1，這些狀態均有飛行試驗熱流結果。

表1 航天飛機飛行狀態

由于航天飛機使用表面隔熱瓦，相對于原子復合和有關的離解能量調節是非催化的[24]，因此本文數值模擬時飛行器表面采用輻射平衡溫度條件和壁面非催化條件。

圖2(b)給出了兩套不同疏密度網格計算得到表面熱流分布與試驗結果的比較，飛行器全長L=30.24 m。Grid1為稀網格，壁面法向第1層網格間距為10-5m，網格總量約為1 200萬；Grid2為密網格，第1層間距為10-6mm，網格總量約為5 000萬?？梢姡嬎憬Y果與飛行試驗[24]符合較好，這說明本文方法能較好模擬航天飛機的氣動熱環境；兩套網格計算差異較小，滿足網格無關系要求。為了保證流場的分辨率，采用Grid2進行計算。

圖2 網格示意圖和航天飛機迎風中心線熱流

6.2 數值模擬結果分析

由于航天飛機頭部迎風面和翼前緣等表面區域熱流較高，為了降低局部區域熱流，可在這些區域加載磁場發生裝置，包括偶極子磁場、螺線管磁場和均勻磁場以及多個磁場組合。這里將位于X-Y平面內與X軸負方向的逆時針夾角定義為θ；位于Y=Constant平面與X軸負方向的逆時針夾角定義為φ。外加磁場配置如下：

Case1采用單個螺線管，主要用于頭部駐點區表面熱防護。螺線管中心點位于坐標(0.4,-0.273 8,0)m，該點與頭部頂點的橫坐標相差0.4 m，與對應橫坐標的下表面點距離也為0.4 m，利于螺線管磁場物理上調整方向。螺線管軸線正方向位于X-Y平面內，與X軸負方向的逆時針夾角θ=0°,20°,40°,60°。由于增大螺線管直徑有利于提高磁控熱防護效果，而螺線管長度對磁控熱防護效果影響相對較小[12]，因此這里盡量增大螺線管直徑，結合防熱層設計和螺線管磁場靈活轉向要求，這里取螺線管直徑為0.4 m、長度為0.4 m。通過調節螺線管通電電流和匝數，可調節磁感應強度，使螺線管軸線正方向距離中心點0.4 m處，磁感應強度為0.3 T，這在工程上可以實現[6,11]。

Case2單個偶極子磁場，主要用于與Case 1螺線管磁場對比。其中心點和方向設置均與Case 1相同。特征長度r0為0.4 m，特征磁場感應強度B0為0.3 T。

Case3多個螺線管磁場組合，主要用于大面積防熱。螺線管1～3中心點分別位于坐標(0.4,-0.273 8,0)m、(0.8,-0.445 6,0)m和(1.6,-0.703 8,0)m，與對應橫坐標的下表面點的距離為0.4 m；其軸線正方向均位于X-Y平面內，與X軸負方向的逆時針夾角分別為θ1～θ3。

Case4采用均勻磁場，在迎風面X=0.4～1.0 m、Y<0 m、Z=0～0.6 m的區間內，設置均勻磁場，磁感應強度方向位于X-Y平面內，與X軸負方向的逆時針夾角θ=0°～180°。需指出的是，該方案磁場邊界上不滿足磁場無源、無散度條件，這里僅作洛倫茲力的對比分析之用。

Case5采用單個螺線管，主要用于水平翼前緣熱防護。中心點位于坐標(21.228,-0.685 9, 5.84)m。螺線管軸線正方向單位矢量為(-0.707,0,0.707)；考慮翼前緣尺寸，這里螺線管直徑和長度均設為0.1 m；通過調節螺線管通電電流和匝數，使螺線管軸線正方距離中心點0.1 m處，磁感應強度分別為0.3 T和1.0 T。

1)不同飛行狀態和磁場方向差異

圖3(a)給出了有/無螺線管磁場條件下不同飛行高度迎風面中心線熱流比較。這里飛行器全長取為L=32.35 m；采用Case 1的螺線管磁場配置，螺線管軸線正方向與X軸負方向的逆時針夾角θ=0°?？梢钥闯?，在頭部設置螺線管磁場可以有效減低駐點附近的熱流；磁場還會使飛行器迎風面靠近底部區域的熱流略有上升；飛行高度越高，磁場使頭部熱流下降幅度越大，高度H=70.2 km 時駐點熱流下將約26.7%。

圖3(b)～圖3(d)給出了不同飛行高度電離氣體電導率分布?？梢钥闯觯妼瘦^高區域均在頭部迎風區域附近，這說明將Case 1磁場設置于頭部迎風區域是合理的；高度H=57.2 km流場中氣體電導率明顯低于其他狀態，這與圖3(a)該狀態的磁控熱防護效果較弱相符，這可能是57.2 km時飛行速度相對較低造成的(表1)。

圖3 不同飛行高度表面熱流和流場電導率分布

圖4(a)給出了Case 1螺線管磁場不同磁場方向時迎風面中心線熱流?？梢钥闯觯瑔我宦菥€管磁場方向大幅變化，對磁控熱防護效果存在一定影響；隨θ增大，熱流凹陷低點后移，但整體的定性規律變化不大，均表現為駐點區熱流顯著下降，下游迎風面局部區域出現不同程度的熱流反沖現象[16]。

圖4(b)給出了螺線管磁場Case 1不同磁場方向時Z方向電流(云圖顏色)和洛倫茲力矢量(黑色箭頭，長短表征洛倫茲力大小)分布，紅色箭頭為來流方向?？梢钥闯?，螺線管磁場方向大幅變化，環形電流和洛倫茲力分布存在一定差異；隨θ增大，洛倫茲力交匯點略微后移，但其整體特征變化不大，這與圖4(a)的熱流變化規律相符。

圖4 不同螺線管方向表面熱流和流場洛倫茲力

2)偶極子和螺線管磁場差異

圖5(a)給出了螺線管磁場(Case 1)和偶極子磁場(Case 2)磁場矢量圖和等值線圖。可以看出，當特征長度和特征磁場強度相同時，兩者存在一定差異，但差異幅度較小。圖5(b)進一步比較了螺線管磁場(Case 1)和偶極子磁場(Case 2)的磁控熱防護效果差異。可以看出，兩種磁場使表面熱流下降的整體規律基本一致，在局部細節上存在一定差異。這說明采用偶極子磁場，可在一定程度上近似模擬螺線管磁場的分布及其影響；但如果要精細地分析磁流體力學控制的效率，則需采用接近物理實際的磁場模型。

圖5 磁場差異及其對熱防護影響

3)多個螺線管磁場組合

由圖3和圖4可以看出單個螺線管磁場可以使頭部熱流有效下降，但下游區域熱流下降幅度不大，甚至出現熱流反沖現象。為了進一步減小迎風面較大區域熱流。可采用多個螺線管組合磁場的方式。

圖6(a)給出了采用螺線管組合磁場Case 3得到的迎風面中心線熱流。這里θ1=0°，由圖4可知此時局部熱流反沖較??；θ2和θ3的值相等(即圖中θ2，3)，分別取為45°、90°和135°。可以看出，當θ2，3=135°時，可降低頭部迎風面較大區域的熱流，且熱流反沖現象較小。

圖6 組合磁場條件下的表面熱流和洛倫茲力

圖6(b)給出了采用螺線管組合磁場Case 3時壓力和洛倫茲力(黑色箭頭，長短表征洛倫茲力大小)分布?？梢钥闯?，當θ2，3=45°時，在X=0.5 m 附近，洛倫茲力使高溫氣流向壁面附著，形成高壓、高熱流區；而當θ2，3=135°時，由于磁場分布變化，洛倫茲力使高溫氣流向壁面附著的區域前移至X=0.4 m附近，且該洛倫茲力矢量相對較小，因此熱流反沖區域略微前移，且反沖幅度較小。

4)磁場方向對洛倫茲力分布和流場結構的影響

為了進一步分析磁場方向對流動的影響，圖7(a)給出了均勻磁場條件下(Case 4)的磁場矢量(紅色箭頭)、壓力云圖和洛倫茲力矢量(黑色箭頭)分布。可以看出，這種情況下，無論磁場方向怎么變化，洛倫茲力矢量永遠是垂直于磁場矢量方向的，其作用效果永遠是阻礙流動方向的。當θ∈[0°, 60°)時，洛倫茲力使高溫氣流附著，隨θ增大，附著的角度減??；當θ∈(60°, 160°)時，洛倫茲力使高溫氣體遠離壁面外推，隨θ增大，外推角度增大。

進一步結合洛倫茲力強度分析，在忽略感應電場影響時，JZ∝(U×B)Z，因此洛倫茲力|F|=|J×B|∝JZ∝(U×B)Z。當θ∈[0°，60°)時，隨θ增大，速度矢量U與磁場矢量B夾角增大，因此JZ增大，洛倫茲力強度增大。但由于附著角度減小，洛倫茲力使其流向壁面的壓縮程度反而減小。當θ∈(60°，160°)時，隨θ增大，外推的角度增大，洛倫茲力越來越傾向于使高溫氣流遠離壁面方向。同時，速度矢量與磁場矢量之間的夾角逐漸減小，會使相同條件下的洛倫茲力強度減小，因此在θ≈120°時，磁場的外推作用最強。

圖7 磁場方向對流動結構的影響

5)水平翼前緣熱環境的磁控

圖8(a)給出了高度H=68.2 km時的表面熱流分布。可以看出頭部駐點區域并不是表面熱環境最嚴酷的區域，水平翼部分區域的表面熱流更高，因此針對飛行水平翼的高熱流區設計了磁場配置方案Case 5。水平翼的厚度較小，螺線管磁場發生裝置尺寸不能太大，其長度和直徑均取為0.1 m；螺線管中心位于高熱流峰值所在的區域附近Z= 5.84 m平面內，見圖8(b)。

圖8 全表面熱流分布和Case 5磁場設置

圖9給出了磁場對水平翼前緣熱流的影響。可以看出，當磁場特征強度為 1.0 T時，表面熱流下降很明顯，這說明采用合適的磁場可以有效地降低水平翼前緣的熱流；而當磁場特征強度為0.3 T時，熱流受磁場影響不明顯。

圖9 不同磁場條件下水平翼前緣熱流

為了進一步分析圖9的結論，圖10分別給出了磁場配置Case 5時3種情況下水平翼附近的電導率分布。這里顯示的平面經過方案Case 5螺線管中心點，平面方向為(0.707，0，0.707)m。

由圖10可以看出，當B0= 0.3 T時，水平翼前緣附近氣體電導率較低；當B0=1.0 T時，電導率顯著升高。這可能是由于隨磁感應強度增大，磁場使氣體減速的效果增強，氣體速度下降，更多的氣體動能轉化為內能，從而使氣體溫度上升，高溫區增大，氣體離解和電離程度增加，電導率增大；而電導率增大，又反過來增強了磁場對氣體滯留和減速等磁控效果，由此形成“正反饋”。當磁場強度較小時，這種正反饋機制強度較小，表面熱流變化不明顯；當B0=1.0 T，正反饋效應較強，電導率明顯升高，因而使表面熱流顯著下降。由此可見，磁場的作用效率與流動中等離子體生成過程緊密相關，準確模擬時必須予以綜合考慮。

圖10 水平翼前緣附近電導率分布

7 結 論

1)本文搭建的磁控熱防護數值仿真模擬平臺，能用于計算評估航天飛機等復雜外形飛行器氣動熱環境的磁流體力學控制效果，典型算例數值計算結果與文獻或飛行試驗數據符合較好。

2)采用合適的磁場能有效降低航天飛機的表面熱流，改善航天飛機的氣動熱環境；在飛行器頭部或翼前緣加載一定強度的螺線管磁場，能使典型狀態的表面熱流下降25%以上；通過多個螺線管磁場組合，能有效降低航天飛行表面大面積的熱流。

3)采用偶極子磁場，可在一定程度上近似模擬同等尺度的螺線管磁場及其磁控熱防護效果；局部磁場方向與流動方向的夾角，在一定程度上決定了洛倫茲力的強度和方向，對磁控效果的影響顯著；磁場的作用效率還與流動中等離子體生成過程緊密相關，準確模擬時必須予以考慮。