劉慶海，張紹華，鄭慧奇，王俊峰，任瓊英，趙 華，易 忠，黃建國

（北京衛(wèi)星環(huán)境工程研究所，北京 100094）

0 引言

高比沖、高效率、長壽命的電推進器在衛(wèi)星姿態(tài)控制、軌道保持及深空探測等領(lǐng)域具有性能優(yōu)勢和競爭力[1-2]。其中螺旋波等離子體推進器尤為引人關(guān)注。美國、法國、澳大利亞、荷蘭、意大利、日本等國家的相關(guān)研究機構(gòu)正在對螺旋波等離子體的產(chǎn)生機制、電雙層加速機制等開展理論研究，并通過實驗手段診斷螺旋波等離子體的特性參數(shù)[2-6]。北京衛(wèi)星環(huán)境工程研究所正致力于螺旋波等離子體推進器的研制，已取得了初步實驗成果[7-10]。

本文主要針對螺旋波等離子體推進器系統(tǒng)對磁場的要求，采取數(shù)值仿真的手段開展推進器磁場位形的研究，建立三維磁場仿真模型。

1 螺旋波等離子體推進器基本原理

圖1所示為螺旋波等離子體推進器基本原理。圖中白色虛線是磁力線，彩色點線是由推進工質(zhì)電離的離子之運動軌跡，右端噴口附近白色霧狀區(qū)域顯示的是電雙層加速區(qū)。推進工質(zhì)從左端送入帶有螺旋波天線的放電管被電離，電離后的離子順著磁力線作螺旋運動，而后進入出口處的磁場中，由電雙層加速為高速離子束流，從右端噴出而產(chǎn)生反推力。電雙層由正電荷層和負電荷層組成，盡管整體顯電中性，但在電雙層中存在電勢梯度，因而可以用于離子加速[11-13]。

圖1 螺旋波等離子體推進器的基本原理示意圖 Fig.1 The working principle of the HPT

在螺旋波等離子體推進器中，推進工質(zhì)的電離、加速、噴出過程以及電羽流的污染控制等都與磁場密切相關(guān)：1）螺旋波是在徑向被約束的磁化等離子體中傳播的右旋極化電磁波，在螺旋波等離子體產(chǎn)生機制中，磁場束縛自由電子的徑向運動；2）通道內(nèi)強磁場的約束作用可減少高溫等離子體與推進器壁面的接觸，削弱離子的壁面濺射侵蝕，延長推進器壽命；3）在噴口附近利用喇叭口形狀的磁場結(jié)構(gòu)產(chǎn)生穩(wěn)定的電雙層結(jié)構(gòu)，是加速離子和產(chǎn)生超高速噴流的關(guān)鍵部位，喇叭口狀的磁場結(jié)構(gòu)也稱為磁噴嘴。De Laval 化學(xué)噴嘴把燃燒室內(nèi)高溫、高壓、低速的化學(xué)氣體變成低溫、低壓、高速的氣流噴出，產(chǎn)生反推力[14-15]。因此，De Laval 噴嘴的尺寸設(shè)計必須非常精確，才能有效完成等離子體能量和動量的轉(zhuǎn)換，提高電推進器效率?；诨瘜W(xué)噴嘴尺寸設(shè)計要求，應(yīng)進行磁噴嘴結(jié)構(gòu)精確設(shè)計。磁場位形的設(shè)計對于控制螺旋波等離子體的產(chǎn)生、通道內(nèi)螺旋波的傳播與吸收、磁噴嘴出口處離子的加速和逃逸，以及提高電推進器的性能至關(guān)重要。

電推進器所需的磁場可用銅線圈、永磁體或者超導(dǎo)磁體技術(shù)來實現(xiàn)[6]。銅線圈方案可通過電流調(diào)節(jié)來控制磁場，但發(fā)熱嚴(yán)重，需要配置冷卻系統(tǒng)；永磁體方案不需要電力和冷卻系統(tǒng)，但磁場位形設(shè)計難度大，磁場強度不易調(diào)節(jié)；超導(dǎo)線圈方案不發(fā)熱、體積小、重量輕，產(chǎn)生的磁場強度高，但需要低溫設(shè)施，技術(shù)實現(xiàn)難度較大。綜合考慮后采用永磁體產(chǎn)生螺旋波等離子體推進器所需的磁場。

2 磁場仿真模型

2.1 等效磁荷原理

對于簡單形狀的永磁體，可直接給出磁場的解析表達式；但對于復(fù)雜形狀的永磁體，則需要借助于數(shù)值計算來得到磁場位形分布。永磁體磁場仿真計算的主要方法有等效磁荷法、等效電流法和有限元法[16-17]。其中，等效磁荷法所用的等效磁荷模型認為磁化體內(nèi)部存在密度為ρm的體磁荷，磁化邊界處存在密度為ρms的面磁荷；等效電流法所用的等效電流模型認為存在體電流和面電流。

本研究中采用等效磁荷法。等效磁荷模型給出永磁體周圍空間磁場強度表達式為

式中：R是永磁體內(nèi)源點到磁場求解點的矢徑；R是源點到求解點的距離；V是永磁體積分區(qū)域；A是永磁體積分區(qū)域表面邊界。體磁荷密度ρm=-μ0?·M，其中μ0是真空磁導(dǎo)率，M是磁化強度。在均勻磁體內(nèi)部，由于?·M=0，所以體磁荷密度ρm=0。在磁體表面邊界，由于磁化強度M不連續(xù)，存在面磁荷密度ρms=μ0en·M，en是磁體邊界的外法向單位矢量。

2.2 模型驗證

用Fortran 程序開發(fā)出三維等效磁荷模型軟件，并進行驗證。在三維仿真計算中，只有表面磁荷對磁場有貢獻，因此，利用式(1)第二項求解，磁場在空間某一點的值是所有表面離散網(wǎng)格點的積分值[17]。對三維空間區(qū)域進行網(wǎng)格劃分，求解每個網(wǎng)格點的磁場值，就得到三維空間的磁場位形分布。

實驗中使用的永磁體是銣鐵硼永磁體，它的質(zhì)量密度是8.0 g/cm3，面磁荷密度是1.0 Wb/cm2。為了驗證模型的正確性，計算了直徑40 mm、高40 mm 的圓柱形永磁體周圍的磁力線及磁場強度分布。圓柱體充磁方向沿著z軸方向。圖2是圓柱形永磁體軸心截面的磁矢（黑色箭頭）、磁場分布及磁力線（白色帶箭頭線）?？梢钥闯?，磁場強度在永磁體附近最強（約1 T 左右），在離磁體50 mm 之外衰減到0.1 T 以下，這是永磁體磁場強度分布的典型特征。

圖2 圓柱形永磁體軸心截面（x-z）的磁場分布 Fig.2 The distribution of magnetic field in x-z section

3 磁場仿真結(jié)果

首先，針對圓柱和圓環(huán)永磁體，構(gòu)造螺旋波等離子體推進器所需的磁場位形，各磁體同軸分布，且磁體相互之間的軸向距離可調(diào)（見圖3）。磁體的具體尺寸和質(zhì)量參見表1。

圖3 磁體形狀及位置示意圖 Fig.3 The permanent magnet configuration of the HPT system

表1 磁體尺寸及質(zhì)量 Table1 The sizes and weights of the permanent magnets

圖4是圓柱磁體B1 與圓環(huán)磁體B2 之間，以及4 個圓環(huán)磁體（B2～B5）之間的軸向間距均為10 mm 時，軸截面的磁矢（白色箭頭）、磁場分布和磁力線（紅色實線）。圖中的黃色方框代表放電管所在位置。可以看出，在磁體內(nèi)側(cè)和外側(cè)與磁體相鄰的地方以及放電管的右端有磁場較弱（小于0.01 T）的區(qū)域，這是不同方向磁場交匯的結(jié)果。從追蹤出的幾條典型磁力線可以看出，磁場在放電管內(nèi)部磁場基本是軸向均勻的，在放電管右端出口處有發(fā)散的趨勢。

圖4 磁場疊加磁力線的分布 Fig.4 The contour of magnetic field and magnetic force line

為更清楚地展示電推進器工作區(qū)域的磁場強度分布，圖5給出了磁場強度沿著中心軸的分布及局部放大圖，其中黃色線段示意的是圓環(huán)磁體內(nèi)部區(qū)域，即放電管所在位置。從圖5(b)可以看出，圓柱形永磁體內(nèi)部的磁場最強，達到約1 T 量級。圓環(huán)磁體區(qū)域的磁場強度，在最左端是0.2 T，到中心區(qū)域下降至約0.03 T，在放電管右端出口處接近于0；遠離放電管出口50 mm 處上升至約0.02 T，之后往更遠處磁場強度再次緩慢衰減至0。

圖5 磁場沿著z 軸的分布及局部放大圖 Fig.5 The profiles of magnetic field along z-axis and the local zoom-in plot

對于螺旋波等離子體推進器，若以13.56 MHz的射頻源產(chǎn)生1018m-3的高密度等離子體，則需要磁場強度約為0.02T[3-5]。圖5(b)顯示永磁體組合產(chǎn)生的磁場在放電管內(nèi)均勻度低于20%，強度也基本滿足產(chǎn)生螺旋波等離子體所需條件。但電推進器噴口附近存在的弱磁場區(qū)導(dǎo)致磁場不連續(xù)，可能對等離子體羽流噴射造成一定的不利影響[18]，需要深入研究分析。

本研究中的三維磁場仿真模型是專門針對螺旋波等離子體推進器開發(fā)設(shè)計的。圖6給出優(yōu)化改進磁體組合之后，再仿真得到的磁場位形。由于永磁體組合是軸對稱分布的，產(chǎn)生的磁場也是三維空間軸對稱分布的，所以只給出軸截面的磁場強度和磁力線的分布。

圖6 優(yōu)化磁體組合之后的電推進器磁場位形分布 Fig.6 The magnetic field distribution for the optimized assembly of permanent magnets

基于仿真得到的三維磁場位形，進行了幾十萬個離子加速的測試模擬，給定離子初始位置和速度分布，求解離子運動控制方程，追蹤離子運動軌跡。離子加速的測試模擬中采用的是相對論離子運動方程：

式中：P是離子動量；q是離子電荷量；r是離子在三維空間的位置坐標(biāo)；γ是相對論因子；mp是靜止離子質(zhì)量；c是光速；B是離子所在位置的磁感應(yīng)強度，而B隨著離子運動在空間維度中不斷更新，由計算網(wǎng)格點上的磁場通過雙線性插值得到。

4 模擬結(jié)果與實驗結(jié)果對比

圖7是模擬結(jié)果與實驗結(jié)果的對比，兩者形狀基本吻合。觀察實驗點火結(jié)果發(fā)現(xiàn)，放電室右端出口附近的等離子體羽流周圍較亮，中心區(qū)域偏暗。而理論上講，如果右端出口處的磁場較大，離子會被磁場有效束縛[14]，等離子體羽流中心密度高，應(yīng)較周圍區(qū)域偏亮。因此，實驗點火的結(jié)果與理論預(yù)期的結(jié)果不一致。其原因可能是放電室右端出口區(qū)域的磁場較弱，使得較重的離子不能被有效約束而是向周圍甩出，造成羽流中心的離子密度低，因而亮度偏暗。

圖7 測試粒子模擬結(jié)果與實驗結(jié)果對比 Fig.7 The comparison of results between test particle simulation and discharge experiment

5 結(jié)束語

本文基于等效磁荷方法，建立了專門針對螺旋波等離子體推進器的三維磁場仿真模型。該模型的計算方法簡單，計算效率高，且使用靈活。在驗證磁場仿真模型正確性基礎(chǔ)上，模擬出圓柱形與圓環(huán)形永磁體組合產(chǎn)生的磁場位形以及磁場強度分布特征。給定不同的磁鐵構(gòu)造，通過大量仿真計算，給出放電室內(nèi)的三維磁場位形和磁場強度分布，從中優(yōu)選出滿足螺旋波等離子體推進器的磁鐵構(gòu)造。

三維磁場仿真模型和磁鐵構(gòu)造方案還有待在電推力測試和羽流特性試驗中不斷改進和優(yōu)化，其中，解決出口磁場強度偏弱導(dǎo)致羽流發(fā)散的問題，是提高電推進器工作效率的關(guān)鍵。

（References）

[1] 于達仁,武志文,鄂鵬,等.穩(wěn)態(tài)等離子體發(fā)動機磁場設(shè)計的發(fā)展及其展望[J].中國工程科學(xué),2005,7(4): 22-29 Yu Daren,Wu Zhiwen,E Peng,et al.Research status and development of the magnetic field design for stationary plasma thruster[J].Engineering Science,2005,7(4): 22-29

[2] 夏廣慶,王冬雪,薛偉華,等.螺旋波等離子體推進研究進展[J].推進技術(shù),2011,32(6): 857-863 Xia Guangqing,Wang Dongxue,Xue Weihua,et al.Progress on the research of helicon plasma thruster[J].Journal of Propulsion Technology,2011,32(6): 857-863

[3] Chen F F,Boswell R W.Helicons: the past decade[J].IEEE Transactions on Plasma Science,1997,25(6): 1245-1257

[4] 房同珍.螺旋波激發(fā)等離子體源的原理和應(yīng)用[J].物理,1999,28(3): 162-167 Fang Tongzhen.Principle and applications of helicon wave excited plasma sources[J].Physics,1999,28(3): 162-167

[5] Chen F F.Physics of helicon discharges[J].Physics Plasmas,1996,3(5): 1783-1793

[6] 李波,王一白,張普卓,等.VASIMR 中螺旋波等離子體源設(shè)計[J].北京航空航天大學(xué)學(xué)報,2012,38(6): 720-725 Li Bo,Wang Yibai,Zhang Puzhuo,et al.Design of helicon plasma source for VASIMR[J].Journal of Beijing University of Aeronautics and Astronautics,2012,38(6): 720-725

[7] 趙華,任瓊英,張紹華,等.螺旋波等離子體推力器機制的實驗研究[C]//2013年（第九屆）中國電推力器學(xué)術(shù)研討會.蘭州,2013: 260-265

[8] 任瓊英,丁永杰,彭輝,等.螺旋波電推進器氙氣工質(zhì)電離率實驗測量分析[C]//2013年（第九屆）中國電推力器學(xué)術(shù)研討會.蘭州,2013: 282-293

[9] 劉慶海,任瓊英,王俊峰,等.螺旋波電推力器氙氣、氮氣羽流比較實驗研究[C]//2013年（第九屆）中國電推力器學(xué)術(shù)研討會.蘭州,2013: 274-281

[10] 趙華,嚴(yán)世林,任瓊英,等.螺旋波電推力器氙氣工質(zhì)電離狀態(tài)轉(zhuǎn)換光譜實驗研究[C]//2013年（第九屆）中國電推力器學(xué)術(shù)研討會.蘭州,2013: 294-300

[11] Charles C,Boswell R W.Laboratory evidence of supersonic ion beam generated by a current-free “helicon” double layer[J].Physics of Plasma,2004,11(4): 1706-1714

[12] Charles C,Boswell R W.Time development of a current-free double-layer[J].Physics of Plasma,2004,11(8): 3808-3812

[13] Charles C,Boswell R W.Effect of exhaust magnetic field in a helicon double-layer thruster operating in xenon[J].IEEE Transactions on Plasma Science,2008,36(5): 2141-2146

[14] Ebersohn F H,Girimaji S S.Magnetic nozzle plasma plume: review of crucial physical phenomena,AIAA 2012-4274[R]

[15] Winglee R,Ziemba T,Giersch L.Simulation and laboratory validation of magnetic nozzle effects for the high power helicon thruster[J].Physics of Plasmas,2007,14(6): 3501

[16] 湯雙清,陳習(xí)坤,唐波.永磁體空間磁場的分析計算及其在永磁磁力軸承中的應(yīng)用[J].大學(xué)物理,2005,24(3): 32-36 Tang Shuangqing,Chen Xikun,Tang Bo.Analysis and calculation on spatial magnetic field emanated from a permanent magnet and its application to permanent magnetic bearing[J].College Physics,2005,24(3): 32-36

[17] 李金,鄭小林,侯文生,等.用于磁定位的永磁體磁場仿真與實驗研究[J].系統(tǒng)仿真學(xué)報,2009,21(18): 5919-5922 Li Jin,Zheng Xiaolin,Hou Wensheng,et al.Simulation and experimental research of magnetic field produced by permanent magnet used in magnetic localization[J].Journal of System Simulation,2009,21(18): 5919-5922

[18] Chen F F.Permanent magnet helicon source for ion propulsion[J].IEEE Transactions on Plasma Science,2008,36(5): 2095-2110