夏旭 楊涓? 付瑜亮 吳先明 耿海 胡展

1) (西北工業(yè)大學(xué) 航天學(xué)院 西安 710072)

2) (蘭州空間技術(shù)物理研究所 蘭州 730000)

電子回旋共振離子推力器(electron cyclotron resonance ion thruster, ECRIT)離子源內(nèi)等離子體分布會(huì)影響束流引出, 而磁場(chǎng)結(jié)構(gòu)決定的ECR 區(qū)與天線的相對(duì)位置共同影響了等離子體分布.在鞘層作用下, 等離子體中的離子或電子被加速對(duì)壁面產(chǎn)生濺射, 形成壁面離子或電子電流, 造成壁面磨損和等離子體損失, 因此研究壁面電流與等離子體特征十分重要.為此本文建立2 cm ECRIT 的粒子PIC/MCC (particle-in-cell with Monte Carlo collision)仿真模型, 數(shù)值模擬研究磁場(chǎng)結(jié)構(gòu)對(duì)離子源內(nèi)等離子體與壁面電流特性的影響.計(jì)算表明, 當(dāng)ECR 區(qū)位于天線上游時(shí), 等離子體集中在天線上游和內(nèi)外磁環(huán)間, 柵極前離子密度最低, 故離子源引出束流、磁環(huán)端面電流和天線壁面電流較低.ECR 區(qū)位于天線下游時(shí), 天線和柵極上游附近的等離子體密度較高, 故離子源引出束流、天線壁面電流和磁環(huán)端面電流較高.腔體壁面等離子體分布與電流受磁場(chǎng)影響最小.

1 引 言

ECRIT 具有無(wú)陰極燒蝕特點(diǎn), 除了用于空間飛行器控制[1,2], 還能用于材料處理等方面[3,4].目前國(guó)外研究出10 cm 和2 cm ECRIT, 并分別應(yīng)用于HAYABUSA[2,5]和50 與70 kg 微小衛(wèi)星[6,7].目前國(guó)內(nèi)外熱衷于微小衛(wèi)星推進(jìn)裝置的研究, 因此開展2 cm ECRIT 研究具有重要意義.離子源是2 cm ECRIT 重要的部件, 其內(nèi)部物理過程包括電離、碰撞以及鞘層加速帶電粒子對(duì)壁面轟擊, 這些過程決定了等離子體和壁面電流特征.等離子體特征影響推力器引出束流大小, 同時(shí)推力器內(nèi)壁面濺射電流是約束工作壽命和引起等離子體損失的重要因素,而且等離子體和壁面電流存在密切關(guān)系并同時(shí)受磁場(chǎng)結(jié)構(gòu)和工作參數(shù)等的影響, 因此開展不同磁場(chǎng)結(jié)構(gòu)對(duì)等離子體和壁面電流特征的影響研究十分重要.由于實(shí)驗(yàn)研究[8?10]很難觀測(cè)到離子源內(nèi)瞬態(tài)物理過程變化, 需要借助數(shù)值模擬手段研究離子源內(nèi)部等離子體過程.雖然PIC 模擬已經(jīng)揭示了ECRIT 的等離子體特征[11?14], 但2 cm ECRIT 離子源壁面電流數(shù)值模擬少有報(bào)道.為此本文采用PIC/MCC 方法研究磁場(chǎng)結(jié)構(gòu)對(duì)離子源內(nèi)等離子體與壁面電流特性的影響, 為推力器的性能預(yù)估提供參考.

2 物理和計(jì)算模型

2.1 物理模型

如圖1 所示, 2 cm ECRIT 離子源由圓柱腔體、環(huán)形天線、磁軛、內(nèi)外環(huán)形永磁體和雙柵極組成, 其中兩磁環(huán)分布在離子源底部, 其間環(huán)向均布8 個(gè)進(jìn)氣孔; 4.2 GHz 的微波能量通過環(huán)形天線饋入放電室, 在0.15 T 的磁場(chǎng)區(qū)域形成ECR 區(qū)[15].在ECR 區(qū), 電子繞磁力線的回旋運(yùn)動(dòng)和微波電場(chǎng)變化同步產(chǎn)生諧振, 由此微波電場(chǎng)加熱電子形成高能電子, 高能電子和原子發(fā)生碰撞激發(fā)和電離, 從而形成ECR 等離子體, 等離子體中的離子經(jīng)由雙柵極被高速引出.圖1 中的H1,W1,H2和W2分別表示外磁環(huán)高度、外磁環(huán)寬度、內(nèi)磁環(huán)高度和內(nèi)磁環(huán)寬度.L1為天線環(huán)形段下表面與屏柵上表面之間的距離,L2為內(nèi)磁環(huán)下表面與天線環(huán)形段上表面之間的距離, 本文中不同磁場(chǎng)結(jié)構(gòu)中L1和L2的值相同.定義A, B, C 和D 分別表示天線環(huán)形段的內(nèi)表面, 下表面, 外表面和上表面; E, F 分別表示內(nèi)磁環(huán)端面、外磁環(huán)端面; G 表示腔體內(nèi)表面.

圖1 2 cm ECRIT 離子源結(jié)構(gòu)示意圖Fig.1.Schematic diagram of the 2 cm ECRIT ion source internal structure.

2.2 計(jì)算模型

圖2 2 cm ECRIT 離子源的放電形貌Fig.2.Discharge image of 2 cm ECRIT ion source.

圖2 為2 cm ECRIT 離子源的氣體放電圖像.從圖2 中可以看出,z軸周圍有一個(gè)圓形的等離子體區(qū)域, 其亮度較均勻, 環(huán)形天線的陰影區(qū)域是可見的.由圖可見, 離子源內(nèi)的等離子體分布近似軸對(duì)稱分布.文獻(xiàn)[9]也表明2 cm ECRIT 離子源中靜磁場(chǎng)軸對(duì)稱, 高頻場(chǎng)接近于軸對(duì)稱性.基于以上分析, 計(jì)算域等離子體參數(shù)分布具有軸對(duì)稱性, 因此將放電室內(nèi)的三維的環(huán)形電離區(qū)簡(jiǎn)化為二維軸對(duì)稱區(qū)域.

2.3 全粒子PIC/MCC 模型

2.3.1 電磁場(chǎng)

本文將ECRIT 離子源內(nèi)部的電磁場(chǎng)分解為3 個(gè)模塊: 靜電場(chǎng)、靜磁場(chǎng)和微波高頻電場(chǎng), 這里忽略高頻磁場(chǎng), 因?yàn)槠浔褥o磁場(chǎng)小得多.靜電場(chǎng)通過離散泊松方程?·?φ=?ρ/ε0直接求解得到, 其中φ,ρ和ε0分別為電勢(shì), 電荷密度和真空介電常數(shù).除了對(duì)稱邊界條件為其余邊界電勢(shì)均為0 V; 靜磁場(chǎng)作為背景場(chǎng), 由COMSOL 多物理場(chǎng)軟件計(jì)算, 靜磁場(chǎng)的計(jì)算參數(shù)設(shè)置: 永磁體(釤鈷合金)的相對(duì)磁導(dǎo)率為1.36; 磁軛(軟鐵)的相對(duì)磁導(dǎo)率為4000; 空氣的相對(duì)磁導(dǎo)率為1.靜磁場(chǎng)的求解域?yàn)槲锢砟Ｐ偷乃袇^(qū)域, 邊界條件僅需設(shè)置外邊界為零磁標(biāo)量勢(shì)邊界.微波高頻電場(chǎng)采用COMSOL 計(jì)算, 計(jì)算時(shí)不考慮等離子體電流, 這是由于微波功率較低時(shí)等離子體電流對(duì)微波電場(chǎng)影響甚微[16].微波電場(chǎng)的計(jì)算參數(shù)設(shè)置: 壁面(不銹鋼)的相對(duì)介電常數(shù)為1, 電導(dǎo)率為 4×106S/m ;天線(銅)的相對(duì)介電常數(shù)為1, 電導(dǎo)率為 6×107S/m ;空氣的相對(duì)介電常數(shù)為1, 電導(dǎo)率為0 S/m.微波電場(chǎng)的求解域?yàn)榉烹娛覅^(qū)域, 邊界條件僅需設(shè)置微波輸入端口邊界, 輸入端口類型設(shè)置為同軸, 輸入功率設(shè)為2 W.PIC 程序計(jì)算時(shí), 每個(gè)時(shí)間步上高頻場(chǎng)的瞬時(shí)值與靜電場(chǎng)疊加共同推進(jìn)電子運(yùn)動(dòng); 而離子質(zhì)量較大, 故高頻場(chǎng)與離子作用弱, 認(rèn)為離子不受高頻場(chǎng)的影響, 只有靜電場(chǎng)推動(dòng)離子運(yùn)動(dòng).

2.3.2 粒子碰撞

離子源以氙氣為工質(zhì), 假設(shè)中性原子為背景粒子, 均勻分布在計(jì)算域中.考慮的碰撞類型有: 電子與原子的彈性碰撞、激發(fā)碰撞和一級(jí)電離碰撞,以及帶電粒子間的庫(kù)倫碰撞.電子與中性原子之間的碰撞是通過蒙特卡洛(MCC)方法處理, 碰撞截面采用Szabo[17]的截面數(shù)據(jù); 采用直接蒙特卡洛(DSMC)方法來(lái)描述帶電粒子間的庫(kù)倫碰撞[18].帶電粒子的運(yùn)動(dòng)采用Boris“蛙跳”算法求解[17].粒子運(yùn)動(dòng)的邊界條件: 對(duì)稱邊界時(shí), 粒子在該邊界上反射, 其余邊界上, 粒子完全吸收.

圖3 不同磁場(chǎng)結(jié)構(gòu)下磁場(chǎng)分布: (a) 1 號(hào)源; (b) 2 號(hào)源; (c) 3 號(hào)源Fig.3.Distribution of magnetic flux density inside of the discharge chamber: (a) 1; (b) 2; (c) 3.

2.3.3 磁場(chǎng)結(jié)構(gòu)

如表1 所示, 計(jì)算中模擬了3 種不同磁場(chǎng)結(jié)構(gòu).通過COMSOL 軟件計(jì)算離子源中磁場(chǎng)強(qiáng)度,其ECR 區(qū)如圖3 所示.可以看出, 1, 2 和3 號(hào)源的ECR 區(qū)分別位于天線上游、與天線相交、天線下游.隨著外磁環(huán)的高度和寬度增大, ECR 區(qū)逐漸向屏柵移動(dòng)并遠(yuǎn)離天線.

表1 磁場(chǎng)結(jié)構(gòu)參數(shù)Table 1.Magnetic circuit structure parameters.

2.4 計(jì)算條件

計(jì)算條件為氙氣流量0.3 sccm (1 sccm =98.3 μg/s), 微波功率2 W, 微波頻率4.2 GHz.中性粒子密度當(dāng)作背景氣體代入PIC 程序, 0.3 sccm流量時(shí), 估算得到的中性原子密度范圍為1.5×1020m?3[19].

根據(jù)文獻(xiàn)[12]給出的離子源內(nèi)等離子體密度為 2×1017m?3, 再假設(shè)平均電子溫度為10 eV, 則德拜長(zhǎng)度約為0.053 mm.保持真實(shí)的介電常數(shù)下,為避免數(shù)值過熱, 將網(wǎng)格步長(zhǎng)設(shè)定為0.025 mm,即使等離子體密度達(dá)到5 × 1017m–3也能保證計(jì)算精度.

3 結(jié)果和討論

3.1 計(jì)算結(jié)果的收斂性

計(jì)算開始時(shí), 設(shè)初始電子和離子溫度分別為2,0.5 eV, 初始等離子體密度為 1×1016m?3.圖4 記錄了不同磁場(chǎng)結(jié)構(gòu)離子源的帶電粒子數(shù)隨計(jì)算時(shí)間的變化.如圖4 所示, 隨著計(jì)算時(shí)間的進(jìn)行,在50 萬(wàn)步以后, 電子和離子數(shù)都逐漸達(dá)到了線性穩(wěn)態(tài).

圖4 帶電粒子數(shù)隨計(jì)算步的變化Fig.4.Quantities of charged particles versus computation step.

圖5 不同磁場(chǎng)結(jié)構(gòu)下等離子體分布: (a)電子密度; (b)離子密度Fig.5.Plasma distribution of ion source with different magnetic circuit: (a) Electron density; (b) ion density.

3.2 等離子體特性分布

不同磁場(chǎng)結(jié)構(gòu)下離子源內(nèi)部的電子和離子密度分布計(jì)算結(jié)果如圖5 所示.從圖中可以看出不同離子源內(nèi)的電子和離子密度峰值都出現(xiàn)在ECR 區(qū)附近, 且電子密度和離子密度峰值分別是3.5×1017m?3/3.9×1017m?3,2.6×1017m?3/3.0×1017m?3, 4.1×1017m?3/4.2×1017m?3.由此可見, 3 種離子源內(nèi)部的等離子體密度峰值均大于4.2 GHz微波頻率的截止密度ncut-off=ε0meω2/e2≈2.1×1017m?3[20,21], 其中ω為微波頻率, 說(shuō)明了ECR等離子體的高密度特征.不同離子源內(nèi)的等離子體密度峰值區(qū)面積大小為: 3 號(hào)源 > 1 號(hào)源 > 2 號(hào)源.2 號(hào)源的密度峰值區(qū)域最小, 是因?yàn)镋CR 區(qū)被天線切斷, 從而影響磁鏡中的電子被ECR 區(qū)持續(xù)加熱.1 號(hào)源的等離子體密度峰值位于天線上游,并且ECR 區(qū)位于較小的磁鏡區(qū), 所以天線上游的內(nèi)外磁環(huán)間的等離子體密度較高而靠近屏柵區(qū)域的等離子體密度較低.這是因?yàn)樘炀€環(huán)形段阻礙了等離子體擴(kuò)散.3 號(hào)源的ECR 區(qū)面積大于1, 2 號(hào)源, 而且ECR 區(qū)位于較大的磁鏡區(qū), 所以其等離子體密度最高而且靠近柵極前的等離子體高于1,2 號(hào)源.柵極前的等離子體密度越高越有利于離子束流的引出, 符合文獻(xiàn)[15]的實(shí)驗(yàn)結(jié)果, 即ECR區(qū)位于天線下游的離子源引出效果優(yōu)于ECR 區(qū)位于天線左側(cè)的離子源.

圖6 為不同磁場(chǎng)結(jié)構(gòu)下離子源的電子溫度分布.1, 2 和3 號(hào)源中電子在ECR 區(qū)被加熱后, 沿磁力線擴(kuò)散到ECR 區(qū)附近的磁鏡中, 電子溫度峰值分別為13.7 eV, 13.3 eV 和12.7 eV.1 號(hào)源中電子溫度最高, 因?yàn)殡娮覧CR 加熱區(qū)靠近微波電場(chǎng)幅值高的天線上游[15].與圖6 對(duì)比, 1 號(hào)源的電子溫度峰值出現(xiàn)在內(nèi)外磁環(huán)間的ECR 區(qū)和天線上表面附近的ECR 區(qū)附近, 但離子密度峰值出現(xiàn)在天線上游附近的ECR 區(qū).這是因?yàn)樘炀€上表面ECR 區(qū)的高能電子未與原子碰撞電離就損失在天線上表面, 而內(nèi)外磁環(huán)間ECR 區(qū)的高能電子會(huì)與原子發(fā)生電離碰撞.2 號(hào)源的電子溫度峰值靠近在天線下表面ECR 區(qū), 而高密度區(qū)偏移了ECR 區(qū).這也是因?yàn)槲窗l(fā)生電離碰撞的高能電子損失在天線下表面.3 號(hào)源中電子被ECR 區(qū)加熱后, 沿磁力線擴(kuò)散到ECR 區(qū)和柵極上游附近的磁鏡中.柵極上游磁鏡中的高能電子多于ECR 區(qū)附近的高能電子, 這是因?yàn)闁艠O上游的磁鏡范圍更大, 而且ECR區(qū)附近的高能電子可能損失在天線環(huán)形段表面或與原子發(fā)生電離碰撞.3 號(hào)源的高密度區(qū)位于ECR區(qū)和柵極上游之間的磁鏡中.這是ECR 區(qū)加熱產(chǎn)生的高能電子向柵極擴(kuò)散中, 部分高能電子與原子發(fā)生電離碰撞, 從而ECR 區(qū)和柵極上游之間的電子溫度較低.部分高能電子和未與原子發(fā)生電離碰撞的電子可能在磁鏡的約束下又往復(fù)穿越ECR區(qū), 被ECR 區(qū)加熱并約束在柵極上游的磁鏡區(qū).2 和3 號(hào)源在靠近柵極表面附近有少量高能電子(大于10 eV), 這是因?yàn)楸诿娓浇蕦拥臍W姆加熱效應(yīng).

圖6 不同磁場(chǎng)結(jié)構(gòu)下離子源的電子溫度Fig.6.Electron temperature of ion source with different magnetic circuit.

3.3 壁面電流

圖7 不同磁場(chǎng)結(jié)構(gòu)下離子源內(nèi)各壁面的電流密度統(tǒng)計(jì)Fig.7.Current density statistics on each surface of ion source with different magnetic circuit.

4 結(jié) 論

本文采用全粒子PIC/MCC 方法模擬了2 cm ECRIT 離子源內(nèi)部的等離子體與壁面電流.結(jié)論如下:

(1)電子被磁鏡約束在天線ECR 區(qū)附近, 使得等離子體也分布在ECR 區(qū)附近.并且ECR 區(qū)與天線的相對(duì)位置會(huì)共同影響離子源的引出束流.當(dāng)ECR 區(qū)位于天線下游時(shí), 更有利于離子引出.從等離子體分布結(jié)果中發(fā)現(xiàn)柵極前等離子體主要集中在外磁環(huán)與軸線之間, 推測(cè)當(dāng)保持孔數(shù)不變時(shí), 而減少壁面附近的柵極孔分布范圍, 將柵極孔集中分布在軸線中心附近, 能提高引出束流以及推進(jìn)劑利用效率.

(2) ECR 區(qū)位于天線下游時(shí), 等離子體主要分布在天線下游而且柵極前等離子體密度高, 但離子對(duì)天線濺射嚴(yán)重.這樣的磁場(chǎng)設(shè)計(jì)雖然能夠提高引出的離子束流, 但會(huì)使天線工作壽命縮減.不同壁面的電流分布規(guī)律為: 天線壁面 > 磁體端面 > 腔體壁面.由此可見離子源天線是壁面電流的主要來(lái)源, 有可能是影響工作壽命的主要因素.未來(lái)設(shè)計(jì)推力器時(shí), 應(yīng)權(quán)衡引出性能和壽命來(lái)選取合適的磁場(chǎng)結(jié)構(gòu).