摘要:电子回旋共振离子推力器(electron cyclotron resonance ion thruster, ECRIT)离子源内等离子体分布会影响束流引出, 而磁场结构决定的ECR区与天线的相对位置共同影响了等离子体分布. 在鞘层作用下, 等离子体中的离子或电子被加速对壁面产生溅射, 形成壁面离子或电子电流, 造成壁面磨损和等离子体损失, 因此研究壁面电流与等离子体特征十分重要. 为此本文建立2 cm ECRIT的粒子PIC/MCC (particle-in-cell with Monte Carlo collision)仿真模型, 数值模拟研究磁场结构对离子源内等离子体与壁面电流特性的影响. 计算表明, 当ECR区位于天线上游时, 等离子体集中在天线上游和内外磁环间, 栅极前离子密度最低, 故离子源引出束流、磁环端面电流和天线壁面电流较低. ECR区位于天线下游时, 天线和栅极上游附近的等离子体密度较高, 故离子源引出束流、天线壁面电流和磁环端面电流较高. 腔体壁面等离子体分布与电流受磁场影响最小.
关键词: 电子回旋共振离子推力器/
磁场结构/
粒子模拟/
壁面电流

English Abstract

全文HTML

--> --> -->

ECRIT具有无阴极烧蚀特点, 除了用于空间飞行器控制^[1,2], 还能用于材料处理等方面^[3,4]. 目前国外研究出10 cm和2 cm ECRIT, 并分别应用于HAYABUSA^[2,5]和50与70 kg微小卫星^[6,7]. 目前国内外热衷于微小卫星推进装置的研究, 因此开展2 cm ECRIT研究具有重要意义. 离子源是2 cm ECRIT重要的部件, 其内部物理过程包括电离、碰撞以及鞘层加速带电粒子对壁面轰击, 这些过程决定了等离子体和壁面电流特征. 等离子体特征影响推力器引出束流大小, 同时推力器内壁面溅射电流是约束工作寿命和引起等离子体损失的重要因素, 而且等离子体和壁面电流存在密切关系并同时受磁场结构和工作参数等的影响, 因此开展不同磁场结构对等离子体和壁面电流特征的影响研究十分重要. 由于实验研究^[8-10]很难观测到离子源内瞬态物理过程变化, 需要借助数值模拟手段研究离子源内部等离子体过程. 虽然PIC模拟已经揭示了ECRIT的等离子体特征^[11-14], 但2 cm ECRIT离子源壁面电流数值模拟少有报道. 为此本文采用PIC/MCC方法研究磁场结构对离子源内等离子体与壁面电流特性的影响, 为推力器的性能预估提供参考.

2

-->

2.1.物理模型

如图1所示, 2 cm ECRIT离子源由圆柱腔体、环形天线、磁轭、内外环形永磁体和双栅极组成, 其中两磁环分布在离子源底部, 其间环向均布8个进气孔; 4.2 GHz的微波能量通过环形天线馈入放电室, 在0.15 T的磁场区域形成ECR区^[15]. 在ECR区, 电子绕磁力线的回旋运动和微波电场变化同步产生谐振, 由此微波电场加热电子形成高能电子, 高能电子和原子发生碰撞激发和电离, 从而形成ECR等离子体, 等离子体中的离子经由双栅极被高速引出. 图1中的H₁, W₁, H₂和W₂分别表示外磁环高度、外磁环宽度、内磁环高度和内磁环宽度. L₁为天线环形段下表面与屏栅上表面之间的距离, L₂为内磁环下表面与天线环形段上表面之间的距离, 本文中不同磁场结构中L₁和L₂的值相同. 定义A, B, C和D分别表示天线环形段的内表面, 下表面, 外表面和上表面; E, F分别表示内磁环端面、外磁环端面; G表示腔体内表面.
图 1 2 cm ECRIT离子源结构示意图
Figure1. Schematic diagram of the 2 cm ECRIT ion source internal structure.

2 -->

2.2.计算模型

图2为2 cm ECRIT离子源的气体放电图像. 从图2中可以看出, z轴周围有一个圆形的等离子体区域, 其亮度较均匀, 环形天线的阴影区域是可见的. 由图可见, 离子源内的等离子体分布近似轴对称分布. 文献[9]也表明2 cm ECRIT离子源中静磁场轴对称, 高频场接近于轴对称性. 基于以上分析, 计算域等离子体参数分布具有轴对称性, 因此将放电室内的三维的环形电离区简化为二维轴对称区域.
图 2 2 cm ECRIT离子源的放电形貌
Figure2. Discharge image of 2 cm ECRIT ion source.

2 -->

2.3.全粒子PIC/MCC模型

3 -->

2.3.1.电磁场

本文将ECRIT离子源内部的电磁场分解为3个模块: 静电场、静磁场和微波高频电场, 这里忽略高频磁场, 因为其比静磁场小得多. 静电场通过离散泊松方程$\nabla \cdot \nabla \varphi =- {\rho }/{{\varepsilon }_{0}}$直接求解得到, 其中$ \varphi $, $ \rho $和$ {\varepsilon }_{0} $分别为电势, 电荷密度和真空介电常数. 除了对称边界条件为$ \dfrac{\partial \varphi }{\partial r}=0\;{\rm{V}}/{\rm{m}} $, 其余边界电势均为0 V; 静磁场作为背景场, 由COMSOL多物理场软件计算, 静磁场的计算参数设置: 永磁体(钐钴合金)的相对磁导率为1.36; 磁轭(软铁)的相对磁导率为4000; 空气的相对磁导率为1. 静磁场的求解域为物理模型的所有区域, 边界条件仅需设置外边界为零磁标量势边界. 微波高频电场采用COMSOL计算, 计算时不考虑等离子体电流, 这是由于微波功率较低时等离子体电流对微波电场影响甚微^[16]. 微波电场的计算参数设置: 壁面(不锈钢)的相对介电常数为1, 电导率为$ 4\times {10}^{6} \;{\rm{S}}/{\rm{m}}$; 天线(铜)的相对介电常数为1, 电导率为$ 6\times {10}^{7} \;{\rm{S}}/{\rm{m}}$; 空气的相对介电常数为1, 电导率为0 S/m. 微波电场的求解域为放电室区域, 边界条件仅需设置微波输入端口边界, 输入端口类型设置为同轴, 输入功率设为2 W. PIC程序计算时, 每个时间步上高频场的瞬时值与静电场叠加共同推进电子运动; 而离子质量较大, 故高频场与离子作用弱, 认为离子不受高频场的影响, 只有静电场推动离子运动.
3 -->

2.3.2.粒子碰撞

离子源以氙气为工质, 假设中性原子为背景粒子, 均匀分布在计算域中. 考虑的碰撞类型有: 电子与原子的弹性碰撞、激发碰撞和一级电离碰撞, 以及带电粒子间的库伦碰撞. 电子与中性原子之间的碰撞是通过蒙特卡洛(MCC)方法处理, 碰撞截面采用Szabo^[17]的截面数据; 采用直接蒙特卡洛(DSMC)方法来描述带电粒子间的库伦碰撞^[18]. 带电粒子的运动采用Boris“蛙跳”算法求解^[17]. 粒子运动的边界条件: 对称边界时, 粒子在该边界上反射, 其余边界上, 粒子完全吸收.
3 -->

2.3.3.磁场结构

如表1所示, 计算中模拟了3种不同磁场结构. 通过COMSOL软件计算离子源中磁场强度, 其ECR区如图3所示. 可以看出, 1, 2和3号源的ECR区分别位于天线上游、与天线相交、天线下游. 随着外磁环的高度和宽度增大, ECR区逐渐向屏栅移动并远离天线.

	H1/mm	W1/mm	H2/mm	W2/mm
1号源	5.4	2	5.4	1.65
2号源	5.6	2.7	5.8	1.8
3号源	5.8	3	5.6	1.8

表1磁场结构参数
Table1.Magnetic circuit structure parameters.

图 3 不同磁场结构下磁场分布:　(a) 1号源; (b) 2号源; (c) 3号源
Figure3. Distribution of magnetic flux density inside of the discharge chamber: (a) 1; (b) 2; (c) 3.

2 -->

2.4.计算条件

计算条件为氙气流量0.3 sccm (1 sccm = 98.3 μg/s), 微波功率2 W, 微波频率4.2 GHz. 中性粒子密度当作背景气体代入PIC程序, 0.3 sccm流量时, 估算得到的中性原子密度范围为$1.5\times $$ {10}^{20}\;{\rm{m}}^{-3}$^[19].
根据公式${w}_{{\rm{p}}{{e}}}=\sqrt{\dfrac{{n}_{{\rm{e}}}{\left|{{e}}\right|}^{2}}{{\varepsilon }_{0}{m}_{{\rm{e}}}}}$和${w}_{{\rm{c}}{\rm{e}}}=\dfrac{\left|{\rm{e}}\right|B}{{m}_{{\rm{e}}}}$, 其中$ {w}_{{\rm{p}}{\rm{e}}}, {w}_{{\rm{c}}{\rm{e}}} $, $ {n}_{{\rm{e}}} $, B, e和$ {m}_{{\rm{e}}} $分别表示等离子体频率、电子回旋频率、电子密度、磁场强度、元电荷量和电子质量. 假设$ {n}_{{\rm{e}}}=5\times {10}^{17}\;{\rm{m}}^{-3} $和B = 0.35 T, 估算出离子源内部分区域的等离子体周期与电子回旋周期低于10^–10 s, 为保证精度, 模拟时间步长Δt取为5 × 10^–12 s. 若以离子在推力器中迁移速度表示扰动在离子源中的传播速度, 则以扰动在离子源内的传播时间为模拟总时长. 2 cm ECRIT离子源内离子迁移速度为10⁴—10⁵ m/s, 则离子源内等离子体演化特征时间为10^–6—10^–7 s, 预计模拟总时间为7.5 × 10^–6 s以保证精确.
根据文献[12]给出的离子源内等离子体密度为$ 2\times {10}^{17}\;{\rm{m}}^{-3} $, 再假设平均电子温度为10 eV, 则德拜长度约为0.053 mm. 保持真实的介电常数下, 为避免数值过热, 将网格步长设定为0.025 mm, 即使等离子体密度达到5 × 10¹⁷ m^–3也能保证计算精度.

2

-->

3.1.计算结果的收敛性

计算开始时, 设初始电子和离子温度分别为2, 0.5 eV, 初始等离子体密度为$ 1\times {10}^{16}\;{\rm{m}}^{-3} $. 图4记录了不同磁场结构离子源的带电粒子数随计算时间的变化. 如图4所示, 随着计算时间的进行, 在50万步以后, 电子和离子数都逐渐达到了线性稳态.
图 4 带电粒子数随计算步的变化
Figure4. Quantities of charged particles versus computation step.

2 -->

3.2.等离子体特性分布

不同磁场结构下离子源内部的电子和离子密度分布计算结果如图5所示. 从图中可以看出不同离子源内的电子和离子密度峰值都出现在ECR区附近, 且电子密度和离子密度峰值分别是$3.5\times $$ {10}^{17}\;{\rm{m}}^{-3}/3.9\times{10}^{17}\;{\rm{m}}^{-3}$, $2.6\times{10}^{17}\;{\rm{m}}^{-3}/3.0\times {10}^{17} $$ \;{\rm{m}}^{-3}$, $4.1\times {10}^{17}\;{\rm{m}}^{-3}/4.2\times {10}^{17}\;{\rm{m}}^{-3}$. 由此可见, 3种离子源内部的等离子体密度峰值均大于4.2 GHz微波频率的截止密度${n}_{{\rm{cut}}\text-{\rm{off}}}= {{\varepsilon }_{0}{m}_{{\rm{e}}}{\omega }^{2}}/{{{{e}}}^{2}}\;\approx \;2.1\;\times $$ {10}^{17}\;{\rm{m}}^{-3}$^[20,21], 其中ω为微波频率, 说明了ECR等离子体的高密度特征. 不同离子源内的等离子体密度峰值区面积大小为: 3号源 > 1号源 > 2号源. 2号源的密度峰值区域最小, 是因为ECR区被天线切断, 从而影响磁镜中的电子被ECR区持续加热. 1号源的等离子体密度峰值位于天线上游, 并且ECR区位于较小的磁镜区, 所以天线上游的内外磁环间的等离子体密度较高而靠近屏栅区域的等离子体密度较低. 这是因为天线环形段阻碍了等离子体扩散. 3号源的ECR区面积大于1, 2号源, 而且ECR区位于较大的磁镜区, 所以其等离子体密度最高而且靠近栅极前的等离子体高于1, 2号源. 栅极前的等离子体密度越高越有利于离子束流的引出, 符合文献[15]的实验结果, 即ECR区位于天线下游的离子源引出效果优于ECR区位于天线左侧的离子源.
图 5 不同磁场结构下等离子体分布:　(a)电子密度; (b)离子密度
Figure5. Plasma distribution of ion source with different magnetic circuit: (a) Electron density; (b) ion density.

图6为不同磁场结构下离子源的电子温度分布. 1, 2和3号源中电子在ECR区被加热后, 沿磁力线扩散到ECR区附近的磁镜中, 电子温度峰值分别为13.7 eV, 13.3 eV和12.7 eV. 1号源中电子温度最高, 因为电子ECR加热区靠近微波电场幅值高的天线上游^[15]. 与图6对比, 1号源的电子温度峰值出现在内外磁环间的ECR区和天线上表面附近的ECR区附近, 但离子密度峰值出现在天线上游附近的ECR区. 这是因为天线上表面ECR区的高能电子未与原子碰撞电离就损失在天线上表面, 而内外磁环间ECR区的高能电子会与原子发生电离碰撞. 2号源的电子温度峰值靠近在天线下表面ECR区, 而高密度区偏移了ECR区. 这也是因为未发生电离碰撞的高能电子损失在天线下表面. 3号源中电子被ECR区加热后, 沿磁力线扩散到ECR区和栅极上游附近的磁镜中. 栅极上游磁镜中的高能电子多于ECR区附近的高能电子, 这是因为栅极上游的磁镜范围更大, 而且ECR区附近的高能电子可能损失在天线环形段表面或与原子发生电离碰撞. 3号源的高密度区位于ECR区和栅极上游之间的磁镜中. 这是ECR区加热产生的高能电子向栅极扩散中, 部分高能电子与原子发生电离碰撞, 从而ECR区和栅极上游之间的电子温度较低. 部分高能电子和未与原子发生电离碰撞的电子可能在磁镜的约束下又往复穿越ECR区, 被ECR区加热并约束在栅极上游的磁镜区. 2和3号源在靠近栅极表面附近有少量高能电子(大于10 eV), 这是因为壁面附近鞘层的欧姆加热效应.
图 6 不同磁场结构下离子源的电子温度
Figure6. Electron temperature of ion source with different magnetic circuit.

2 -->

3.3.壁面电流

图7为不同磁场结构下离子源各表面的电流密度统计$J=\dfrac{{Q}_{\rm{num}}\times {\rm{e}}}{S\times \Delta t}$, $ {Q}_{\rm{num}} $表示到达壁面的带电粒子数, S表示壁面面积. 横坐标A, B, C, D, E, F, G的含义如图1所示. 由图7可见天线各表面上离子电流密度大于电子电流密度. 这是因为离子没有被磁化, 所以流入天线的离子电流大于流入天线的电子电流. 另一方面, 由于磁镜的存在, 电子受到很好的束缚, 文献[13]中也证实了这一现象. 天线壁面电流密度之和的高低顺序为3号源 > 2号源 > 1号源. 因为3号源中天线附近的等离子体密度最高, 其次是2号源、1号源. 到达G表面上的离子电流密度大小为2号源 > 3号源 > 1号源, 而G表面上的电子电流密度几乎为0. 因为从图5可知2号源中等离子体分布更均匀, 其靠近腔体表面附近的等离子体密度较高, 壁面鞘层阻止电子到达腔体壁面而加速离子到达壁面. 在内磁环端面E和外磁环端面F的电子电流密度大于离子电流密度. 流向E, F表面上的电流密度之和的高低顺序为3号源 > 2号源 > 1号源. 对于1号源, 流向外磁环端面F的电子电流密度大于离子电流密度, 而流向内磁环端面E的情况正好相反. 这是因为该离子源内等离子体主要集中在内磁环端面附近, 其离子密度高于外磁环端面. 而对于2, 3号源, 流向内外磁环端面的电子电流密度大于离子电流密度. 这是因为电子质量较轻, 所以比离子运动得更快, 即使磁环表面存在鞘层, 更多的电子绕磁力线垂直运动到磁环下表面上. 由于G表面上的离子电流密度很小, 所以认为对于2 cm ECRIT离子源, 可以忽略等离子体对腔体壁面的侵蚀. 由图7可知, 各壁面附近等离子体密度高低顺序为: 天线壁面 > 磁体端面 > 腔体壁面, 其对应的壁面电流密度高低顺序也相同. 综上, 天线壁面溅射最严重, 可得出天线材料决定离子源工作寿命. 虽然ECR区位于天线下游时引出的离子束流高, 但对天线壁面磨损严重, 所以需要权衡寿命与性能二者关系.
图 7 不同磁场结构下离子源内各壁面的电流密度统计
Figure7. Current density statistics on each surface of ion source with different magnetic circuit.

本文采用全粒子PIC/MCC方法模拟了2 cm ECRIT离子源内部的等离子体与壁面电流. 结论如下:
(1)电子被磁镜约束在天线ECR区附近, 使得等离子体也分布在ECR区附近. 并且ECR区与天线的相对位置会共同影响离子源的引出束流. 当ECR区位于天线下游时, 更有利于离子引出. 从等离子体分布结果中发现栅极前等离子体主要集中在外磁环与轴线之间, 推测当保持孔数不变时, 而减少壁面附近的栅极孔分布范围, 将栅极孔集中分布在轴线中心附近, 能提高引出束流以及推进剂利用效率.
(2) ECR区位于天线下游时, 等离子体主要分布在天线下游而且栅极前等离子体密度高, 但离子对天线溅射严重. 这样的磁场设计虽然能够提高引出的离子束流, 但会使天线工作寿命缩减. 不同壁面的电流分布规律为: 天线壁面 > 磁体端面 > 腔体壁面. 由此可见离子源天线是壁面电流的主要来源, 有可能是影响工作寿命的主要因素. 未来设计推力器时, 应权衡引出性能和寿命来选取合适的磁场结构.