无限大平行板间轴对称电子束二维静电PIC模拟程序设计

黄 鳌,杨中海,黄 桃,胡 权,金晓林

摘 要:利用粒子模拟方法(PIC),用Fortran语言进行程序设计,对无限大平行板间轴对称电子束的运动进行二维模拟。在该程序设计中,宏粒子采用环宏粒子模型,用体积加权模型将每一个宏粒子电量分配到格点上,相对于面积加权模型更加精确。与通常数值模拟方法相比,PIC方法得出各个宏粒子的瞬态运动特性。最后通过以Magic仿真软件的计算结果为标准进行验证,证明了该算法的准确性。

关键词:PIC;环宏粒子模型;体积加权模型;有限差分方法

中图分类号:TP311文献标识码:B

文章编号:1004-373X(2009)12-108-03

PIC Simulated Program Design of Two-dimensional Axial Symmetry

Electron Beam in Infinite Parallel

HUANG Ao,YANG Zhonghai,HUANG Tao,HU Quan,JIN Xiaolin

(School of Physical Electronics,University of Electronic Science and Technology of China,Chengdu,610054,China)

Abstract:The movement of axial symmetry electron beam located between infinite parallel by using particle simulation method in Fortran platform.Ring structure macro particle is distributed on grid by the volume-weighted model method and can be more accurate than face-weighted method.Compared with traditional numerical simulation method,the position and velocity of macro electron is calculated at any time step in this PIC method.According to the comparison with the calculation result from MAGIC,the accuracy of such algorithm mentioned in this paper has been verified.

Keywords:PIC;ring marco particle model;volume-weighted mode;finite difference method

0 引 言

近年来计算机技术的迅速发展,在传统的解析理论分析和实验研究以外,发展了第3种强有力的重要研究方法:计算机模拟,即使用计算机模拟计算跟踪大量微观粒子的运动,再对组成物体(包括气态、液态、固态和等离子态)的大量微观粒子进行统计平均,由此得到宏观物体的物质特征和运动规律。粒子模拟的特点在于直接从最基本的电磁运动规律和粒子运动的力学规律出发,利用高速计算机直接求解完整的Maxwell方程组和Lorentz方程,追踪每一个宏粒子的运动,它能更真实的反映实际物理过程[1]。

这里将采用粒子模拟的方法对无限大平行板间轴对称电子束进行静电PIC模拟程序设计,算法的最终设计目标为了实现任意形状型号电子枪的通用软件设计。

1 物理模型

本文所模拟的模型为无限大金属平行板系统这一理想化模型。电子以某初速度从阴极板的某一圆形区域连续发射,到达阳极板后电子完全被吸收。

2 粒子模拟方法的二维静电模型

无限大平行板间轴对称电子束运动模拟采用二维静电PIC模型。

根据以上的物理模型,将数值模拟分析分为以下4个部分:Poisson方程的求解得出空间电势分布、空间电场的求解、粒子所在位置场的求解和粒子位置的更新。

2.1 求解Poisson方程

2φ=-ρ/ε0(1)

在旋转对称场的情况下,式(1)变为:

2φ祌2+1r•郸摘祌+2φ祕2=-ρε0(2)

采用有限差分法(FDM),对式(2)进行差分处理。采用如图1所示的网格,轴线(z轴)处j=1,而0点落在第j行,网格为边长为h的正方形,通过对式(2)做等间距差分化后:

图1 旋转对称场的等距网络

2.1.1 对称轴以外的点(j>1)

4φ0=φ2+φ4+[1+1/2(j-1)]φ1+

[1-1/2(j-1)]φ3+ρ/h2ε0(3)

2.1.2 对称轴上(j=1)

4φ0=φ2+φ4+2φ1+ρ/h2ε0(4)

2.2 求解网格点上的电场

E0z=(φ2-φ4)/2h;E0r=(φ3-φ1)/2h(5)

2.3 求粒子所受到的力

采用如图2所示面积加权法,粒子所受到的电场力为:

Fi(ri)=qiE=qi[E1A1+E2A2+E3A3+E4A4]/h2(6)

图2 二维面积加权法

2.4 速度和位置的更新

粒子的洛仑兹运动方程:

dvi/dt=F(rj,t)/mi(7)

dri/dt=vi(t)(8)

用蛙跳算法来推动粒子,将式(7),式(8)用时间中心差分得:

vin+1=vin+(Fin+1/2/mi)Δt(9)

rin+3/2=rin+1/2+vin+1Δt(10)

初始时刻t=0开始,由给定的粒子速度和位置分布v0i和r0i,求得r1/2i=r0i+v0i(Δt/2),进而求得F1/2i。

2.5 模拟流程图

模拟流程图如图3所示。

图3 PIC模拟流程图

3 计算模型

为节省计算机内存和计算时间,在进行计算机模拟时不必跟踪每一个点粒子(微观粒子),只需计算能代表这个荷电粒子群的一个粒子就可以,从而把这个荷电粒子群处理为一个特殊粒子,即所谓“宏粒子”。

两板间的粒子束会产生空间电荷场(即自洽场),其大小由泊松方程求解,在由宏粒子的位置变化而得到电荷密度分布和求宏粒子所受到的电场力的过程中,需要借助面积加权模型[5] 来进行电荷、电场分配。在本模型中,阴极板所发射出来的电子束为一轴对称模型,可以采用二维的z-r坐标系静电模型来计算。在对电荷向各个网格点分配时,面积加权模型在z-r坐标系中的精度并不是很高,尤其在对称轴附近的区域会产生较大误差,在电子粒中心附近会产生振荡,与实际不符[6]。程序设计采用一种体积加权模型,将宏粒子视为环形体,粒子在其内均匀分布。随着粒子的向前推进,环的半径增大,但宏粒子环本身的z向宽度和r向宽度不变,它的电量也不变。如图2所示,环宏粒子的电量为qi,环横截面的面积为A,则分配到1,2,3,4个格点所在环上的电量分别为qiA1/A,qiA2/A,qiA3/A,qiA4/A。进而求得每一个格点上的电荷体密度。由泊松方程可得出每一格点的电势,从而得出格点上每一点的电场强度,采用面积加权的方法求出每个粒子上的电场,进而对粒子的速度和位置进行更新,进入下一个时间步长。

4 模拟结果

这里所选用的物理模型是:相隔距离为s=0.01 m的两无限大(只要两板间的距离远小于板的线度就可视为无限大平行板)平行板,从阴极板的一半径为r=0.001 m的圆形区域内发射一束电子,电流密度为jcathode,粒子初速度为cell(i)%vz。程序设计的目标是观察粒子在不同时刻两板间运动粒子位置状态分布。

图4和图5是在如下初始条件分别通过Magic和该程序设计输出的结果:

粒子初速度为:

cell(i)%vz= 8.381675d6 m/s

发射电流密度为:

jcathode=-2.5d3 A/m2

运行总时间为:

t=4.5 d-9 s。

图4 Magic输出图

图5 本文程序设计输出图

通过输出结果比较,其包络面形状基本一致。Magic结果输出图在阳极板的最大高度为:2.055 387 370 310 987 E-003 m,所设计程序结果在阳极板最大高度为:2.007 482 300 816 977 E-003 m。如果该程序将网格化分的更细,两者的差别将更小,证明本程序设计的正确性。

5 结 语

这里用粒子模拟的方法对无限大平行板间轴对称粒子束运动进行模拟,采用环宏粒子模型进行程序设计,在宏粒子电量向网格分配时采用比面积加权模型更加精确的体积加权模型(CIC)。相对于通常的计算机数值模拟方法,该程序设计通过追踪每一个粒子,能观察到任一时刻的粒子状态(位置、速度)分布,能输出每个宏粒子在任意时刻的位置和速度数据,直接求解完整的Maxwell方程组和Lorentz方程,它能更真实地再现整个微观过程。在相同的初始条件下,程序通过Magic通用软件进行验证是正确的。

参考文献

[1]廖平,杨中海,雷文强,等.微波管电子枪三维粒子模拟研究[J].强激光与粒子束,2004,16(3):353-357.

[2]王秉中.计算电磁学[M].北京:科学出版社,2002.

[3]Freeman J C.Preliminary Study of Electron Emission for Use in the PIC Portion of MAFIA[R].NASA/TM220012210890,Glenn Research Center,2001.

[4]Birdsall C K,Fuss D.Clouds-in-Clouds,Clouds-in-Cell′s Physics for Many Body Plasma Simulation[J].Comput.Phys.,1969(3):494-511.

[5]李永东,何锋,刘纯亮.轴对称平板二极管空间电荷限制流的二维效应[J].强激光与粒子束,2005,17(6):913-916.

[6]盛新庆.计算电磁学要论[M].北京:科学出版社,2004.