EAST中性束注入器束中性化区气流模拟分析

　　中性束注入(Neutral Beam Injection)加热是磁约束核聚变装置上主要的辅助加热手段。强流离子束中性化是中性束注入的关键步骤，它直接决定了中性束的注入功率，进而影响对等离子体的加热效果。本文分析了EAST 中性束注入器束中性化区的气流特性，发展了基于直接模拟蒙特卡洛(DSMC)方法的气流模拟程序，并用于EAST 中性束注入器束中性化区的模拟分析和优化设计。

　　中性束注入器(Neutral Beam Injector ，NBI)是将高能中性粒子束注入到聚变装置内用于加热等离子体、驱动等离子体电流的装置，也是EAST 托卡马克实现其科学目标所必须的辅助加热系统。EAST-NBI 工作原理如图1 所示。EAST-NBI 的离子源利用灯丝发射的电子和离子源供气系统提供的工作气体，在弧室内通过弧电场的作用形成均匀高密度源等离子体，源等离子体中的离子在离子源束引出系统静电场的作用下被引出并加速形成高能离子束。与此同时，离子源内未被电离的工作气体扩散进入中性化室，与中性化室的补充充气一起形成高能粒子束中性化所需要的电荷转换靶。在中性化室内，高能离子通过与电荷转换靶的气体靶分子碰撞，从而实现束流的中性化。由于离子束的中性化效率不能达到100%，EAST-NBI采用磁偏转模式利用高能离子通过偏转磁场时所受到的洛仑兹力的作用，将其偏转180° 而脱离原来的传输方向，进入离子吞食器。经过偏转系统分离出离子后的高能中性粒子束流经漂移管道注入到EAST 等离子体中。可见，真空技术网(http://www.chvacuum.com/)认为离子束中性化过程对EAST-NBI 的加热效果十分关键。

图1 EAST-NBI 的示意图

1、EAST-NBI 束中性化区

　　EAST-NBI 的离子束从电极系统引出加速后即开始中性化过程，因此束中性化区由加速器空腔、插板阀空腔、中性化室以及连接管道构成。根据各个部件的结构特点，将其简化为如图2 所示的二维几何模型。特别地，由于装配需要，真空室前端盖板上的束通道开口要大于中性化室的外径，这使得束中性化区在该连接部分形成与真空室相通的狭缝。束中性化区的机械结构确定以后，气体靶厚就由进气量来控制。如图2 所示，EAST-NBI 束中性化区的气源来自两部分，一部分是来自等离子体发生器的前端进气，另一部分是中性化室的补充气体。

　　离子束的中性化效率与气体含量(定义为气体靶厚)密切相关，中性效率随着靶厚的增加而增加，但增长率逐渐变小，且最终趋于饱和。在实验运行时气体靶厚需要较为准确地控制在最佳靶厚附近。对于EAST-NBI，如果进气量不足，中性化效率不足比较明显;相反地，如果进气量过大的话，不但中性化效率的变化不大，多余的靶气体还会增加真空系统负荷和真空室气体压强，进而导致再电离损失的增大。此外，束与气体靶碰撞会发生一定的散射，靶厚的增加使得散射次数变多，束因此会更为发散。因此，为了实现束中性化区的合理设计和优化，需要一套研究气体流动状态的较为准确的计算方法。

图2 EAST-NBI 束中性化区内气体流动状态

2、束中性化区气流模拟程序

　　由EAST-NBI 束中性化区内气流的稀薄程度和系统的特征长度可得，该气流属于过渡流区域和分子流区域。为了满足跨流域模拟的需要，采用直接模拟蒙特卡洛(DSMC)方法发展NBI束中性化区气流模拟程序。本文参考标准DSMC程序，其流程如图3 所示，并根据束中性化区的气流特点进行修改。

2.1、程序初始化

　　准备模拟所需的常数(例如，玻尔兹曼常数)，设定气体和流场基本属性，赋予各变量初始值等。其中气体和流场基本属性主要包括：坐标系(例如，直角坐标系或柱坐标系)、计算区域和边界条件(包括固壁边界、入口边界、出口边界)、划分的网格和子网格、模拟分子基数、模拟时间步长、入口的进气量、进气温度、进气速度、氘分子半径、氘分子质量、气体粘滞- 温度幂指数、气体分子碰撞模型、气体的自由度。

图3 标准DSMC 程序流程图

2.2、气源模拟分子注入

　　根据各入口的进气量Qin 和进气温度Tin，计算单位步长Δt 内进入中性化室的模拟分子个数Nin。然后，根据进气温度的Maxwell 速度分布函数，运用MC 方法确定入射模拟分子的各个速度分量的大小。再根据入口边界条件中对进入分子的空间分布和时间分布假设，利用蒙特卡洛方法确定每个入射模拟分子在入口平面上的入射点坐标xin 和入射时序tin(<Δt)。最终，一般假设分子时间分布均匀，则由匀速直线运动方程求得每个入射模拟分子的位置坐标x=xin+u·tin=xin+u·Δt。

2.3、模拟分子位置推进

　　先由匀速直线运动方程求出每个模拟分子以各自速度在Δt 内运动的距离，从而确定模拟分子新的位置坐标x=xo+u·Δt。模拟分子在经历迁移运动后有可能与边界发生碰撞。因此需要将每个模拟分子新的坐标与所有的边界坐标进行比较，判断是否与边界发生碰撞。如果发生碰撞的边界是真空，则将模拟分子作为逸出处理，不再跟踪计算。如果发生相互作用的边界是壁面，则先用模拟分子的直线运动方程和边界的面方程联立求出碰撞点xc，再根据碰撞前模拟分子运动距离和速度求出碰撞前的运动耗时tc。然后由对应的分子壁面碰撞模型来确定碰撞后的分子的各个速度分量ur。最终，经过碰撞的模拟分子新的位置坐标为x=xc+ur·(Δt-tc)，实现了模拟分子位置的推进。

2.4、分子排序

　　根据模拟分子新的位置坐标调整模拟分子所属的网格编号，统计各个网格内的分子总数，并按照网格的顺序(在每个网格中又按照子网格的顺序)对模拟分子进行排序，最终分子编号将被置于一组交叉引用数组中。这是为了之后选取“碰撞对”做准备，据此计算出每个网格所需要模拟的碰撞次数，并保证只有同一网格内的分子才会发生碰撞。

2.5、分子间碰撞

　　本文采用NTC 方法进行碰撞对的取样，先计算出Δt 内每个网格中需要模拟的碰撞次数Ncol。然后利用蒙特卡洛方法在分子编号的交叉引数组中，随机选取同一个子网格(或网格)中的两个模拟分子作为“碰撞对”。再计算两个分子的相对速度cr 和碰撞截面σ，本文采用可变径硬球(VHS)模型，其碰撞截面随着相对速度增加而减小。接着判断“碰撞对”是否被选中发碰撞。如果判断“碰撞对”发生碰撞，则依据VHS 模型的碰撞法则计算分子碰撞后的速度和内能。如果判断没有碰撞则继续选取“碰撞对”，依次完成每个网格碰撞次数。

2.6、流动性质采样

　　包括采样总数以及每个网格累加的模拟分子数的总和、各速度分量的总和、速度平方和的总和。每个Δt 的模拟的时间、模拟分子数、碰撞总次数、推进模拟分子次数等参数实时输出显示，可以初步判断模型是否存在错误。每经过jnis 个Δt 进行一次流动性质采样，从而得到速度、温度、密度等宏观量随时间的变化规律。重复nloop个Δt 直到统计误差满足要求，特别是非定常流动的模拟，需要获得足够多的总计样本。最后求出截止时刻的各物理量的系综平均。

2.7、模拟结果的输出

　　将每个网格累加的模拟分子数的总和、各速度分量的总和、速度平方和的总和分别除以采样总数，得到每个网格相关变量的统计平均值。最后，由这些变量统计平均值计算得到流场的性质并输出到文本。流场性质包括每个网格的坐标、气体密度、平动温度、旋转温度、速度分量。

3、模拟结果分析

3.1、EAST-NBI 束中性化区气流特性

　　在模拟计算中，工作气体D2 被分别从弧室和中性化室进入EAST-NBI 的束中性化区，为方便比较，两个进气量均设为15 Torr·l/s。计算区域中网格的大小设为Δx=0.02 m、Δy=0.01 m，每个网格又划分为2×2 个子网格。每个模拟分子代表的真实分子数为1014。设时间步长Δt=5×10-6 s，每5 个Δt 进行一次采样，总的采样次数为40000，总的模拟时间为1 s，并且气流在4000 次采样后达到稳定。

　　本文仅进行未引出束时的气流模拟，因此气体温度和壁面温度均设定为室温。进入分子在时间和空间上分布均匀， 初始速度遵循Maxwell 分布。为了模拟分子间碰撞，故采用VHS模型，而分子与壁面碰撞采用漫反射模型。特别地，设模拟分子与弧室侧的边界碰撞时，77%的模拟分子在碰撞点反射，而23%从边界其他点反射。图4 给出了两种进气模式计算得到的速度矢量分布。如图4(a)所示，氘气穿过整个电极系统的大截面，以层流的状态缓慢地进入束中性化区。但是气体进入沟槽区域后随即迅速向两侧扩散，同时其轴向速度停止增加反而急剧减小。由于沟槽区域的横截面较大，壁面对气流的粘滞作用将减弱，引起质量流速率的增加。与此同时，沟槽的广大空间使得其中的气体接近热平衡状态，这在气体密度分布中也能体现出来(如图5(a)所示)。此外在沟槽区域内的小狭缝与计算区域外的真空环境相连接，从而在较高的压力梯度下将气体推出。气体继续向前运动，进入相对十分收缩的中性化室桶体，导致在桶体入口处现了堵塞。从图5(a)可以看出，两个对称的滞止点出现在中性化室入口附近的壁面上，这是源于亚声速气流掠过粘滞平面的结果。

图4 两种进气模式的速度矢量分布

图5 两种进气模式的气体密度分布

　　而在中性化室进气模式中，氘气以高压缩比从侧面的进气口喷出并在对面的壁面上形成滞止点。不过由于实际的进气口在上面，进气口距离对面的壁面较远，很难形成滞止点。两个宏观量的分布图中最明显的是，在束中性化区前端的大部分区域内，气体速度矢量几乎为零，而且气体密度大小十分接近。其原因是电极系统侧的边界条件采用了完全漫反射模型，反射气流和来流相互碰撞而使得气流趋于各向同性。流场的另一个特征是在进气口和狭缝之间具有两个剪切区，把束中性化区内的气体向外挤压，特别是在靠近进气口一侧还出现了一个环形涡流。两种进气模式有一个相同的流场特征，径向压强梯度在束中性化区的中部形成，并推动气体向出口流动。

　　EAST-NBI 束中性化区气流特性的定量分析结果如图6 所示，两种进气模式的区别很明显。Kn 数(0.2 ~ 2.4)表明了束中性化区的气流属于过渡流和自由分子流领域，与定性分析相符。在一定体积内，气体温度会随着气体分子的碰撞频率增加而上升，所以两条温度曲线都是在进气口附近的区域较高。而在出口附近的区域各个分子趋于同向，碰撞频率相应地减小，温度也相应地骤减。考虑到模型在出口边界采用真空边界条件，而在实际情况中，出口并不是完全真空甚至还有返流，因此出口附近的温度会高于模拟结果。气体压强曲线可以跟日后的测量值进行比较，在真空边界条件下模拟得到的出口压强(< 0.03 Pa)满足EAST-NBI 的工程要求。

图6 EAST-NBI 束中性化区内各个中心线宏观量的分布：(a)Knudsen 数;(b)气体温度;(c)气体压强;(d)轴向速度

3.2、EAST-NBI 束中性化区的优化设计

　　在实验运行中，总是希望用最小的进气量来获得要求的气体靶厚，换句话说，就是尽可能地提高气体利用率，从而降低前文多次提及的剩余气体带来的一系列副作用。在实际装置中，如几何模型所示的狭缝在短边的宽度约为5 cm，在长边的宽度约为2 cm。由于这些狭缝直接与真空室的真空环境相连接，如同束中性化区内气体的又一个出口，很大程度上减少了其中的气体密度。为此，在未来EAST-NBI 的设计装配中，将借鉴国外其他NBI 装置的设计把狭缝尽可能的缩小，甚至是完全封闭起来。

　　和上述的模拟条件相同，工作气体D2 被分别从弧室和中性化室进入EAST-NBI 的束中性化区，为方便比较，两个进气量均设为15 Torr·l/s。图7 给出了从狭缝流出气体所占比例与气体靶厚随狭缝宽度的变化情况。在两种进气模式中，随着狭缝一点点增大，中性化区域的靶厚下降得十分迅速。当狭缝增大到4 cm 时，靶厚已经降到了没有狭缝时的一半，而且由于狭缝距离气源更近，此时从狭缝流出的气体已占主导，超过所有逸出气体分子的50%。根据模拟结果，在以后的装配中至少要控制狭缝的宽度小于1 cm。

图7 EAST-NBI 中性化室的安装缝隙对气体靶厚的影响

　　在EAST-NBI 综合测试平台上进行中性束引出实验时，对中性化室补气会引起离子源电极打火次数增多，而且弧室内气体放电的稳定性变差。中性化室补充进气的目的是增加靶厚，提高中性化效率。但由中性化室进入的气体会向离子源方向扩散，并通过加速器并进入弧室中。这些气体会增加沉积在电极表面的高能粒子数量，甚至影响弧室中的气体放电过程。为了解决这一问题，需要在气体靶厚不变(或减少很小)的前提下，尽可能的减少流向电极系统的气体。由于EAST-NBI 真空室内各部件的尺寸和位置很难改变，改变中性化室进气的方向是一个较为简单可行的解决方案。为此，在EAST-NBI 束中性化区模型中模拟了不同进气角度对气体靶厚和离子源运行的影响。中性化室进气角度的定义如图8 所示。在计算中进气角度范围从45°到135°，进气量均设为15 Torr·l/s。靠近弧室侧的气体数密度用于反映流向弧室的气流大小，在图8 中用圆点表示。目前装置的进气角度为90°，为了方便比较，以进气角度为90°时的气体靶厚和弧室侧气体数密度作归一化处理。可以预见，随着进气角度增加，靠近弧室侧的气体数密度也将变大， 并且在60° 和120°之间的增大趋势最明显。图8 还显示了不同进气角度对应的气体靶厚(用方点表示)。可以看出，不同进气角度的靶厚变化范围要小于弧室侧气体数密度的变化范围。这意味着，流向弧室的气体量可以通过减小进气角度实现，而同时气体靶厚的减小量更小。而且当进气角度小于60°时，两个量的减小趋势都变缓。因此提议将中性化室的进气角度从90°改为60°，弧室侧的气体数密度和束中性化区气体靶厚分别减小90%和97%。

图8 EAST-NBI 中性化室进气角度对离子源运行和靶厚的影响

4、结论

　　离子束中性化过程直接决定了中性束注入器的性能参数，而气体靶作为中性化媒介最为重要，考虑到束中性化区的气流属于过渡流领域，本文将该流域的标准方法———直接模拟蒙特卡洛(DSMC)方法引入NBI 的相关研究，建立了NBI 束中性化区气流模拟程序。利用该程序EAST-NBI束中性化区内的气流特性，计算结果显示：由于考虑了气体分子间碰撞，DSMC 方法的结果较符合气流的实际情况，并且根据气流的速度、温度、压强等宏观性质的细节分布，分析了气体靶的形成机制。此外，还将该程序对EAST-NBI 束中性化区的优化设计进行模拟计算，结果显示：中性化室和真空室前端盖板的束通道之间的狭缝宽度必须控制小于1 cm，才能有效的提高气体利用率;将中性化室的进气角度从90°改为60°，能够在气体靶厚不变(或减少很小)的前提下，减少流向电极系统的气体，增加离子源运行的稳定性。