一种多粒子多碰撞微放电阈值蒙特卡罗计算方法与流程

1.本发明涉及一种多粒子多碰撞微放电阈值蒙特卡罗计算方法，属于空间微波技术领域。


背景技术：

2.微放电效应也称二次电子倍增效应，是指部件处于1
×
10
‑3pa或更低压强时，在承受大功率的情况下发生的谐振放电现象。航天器载荷中大功率微波部件如输出多工器、滤波器、开关矩阵、天线馈源等极易产生微放电效应，微放电效应一旦发生将造成严重后果：噪声电平抬高，输出功率下降；微波传输系统驻波比增大，反射功率增加，信道阻塞；微波部件表面损坏，载荷寿命缩短；航天器载荷永久性失效，因此大功率微波部件的研制过程中微放电是必须要克服的效应之一。
3.在设计阶段通过微放电敏感区域或分析方法对星载微波部件微放电风险评估是减少反复设计、避免长周期地面试验的有效手段。目前，微放电分析方法方面已经发展出了包括经典理论和统计理论的理论分析方法，以及包括蒙特卡罗(monte carlo，mc)和粒子模拟(particle
‑
in
‑
cell，pic)的数值模拟方法。由于理论方法需要电子运动轨迹的解析解，因而难以应用到复杂结构的微波部件微放电阈值计算中。pic方法模拟粒子在电磁场下的运动以及粒子与部件表面碰撞情况所作的忽略和假设少，既考虑了电子出射表面的速度、方向、相位等齐全的运动参数，又考虑了电子与电磁场间的相互作用，可以实现阻抗变换器、脊波导滤波器以及同轴腔体滤波器等复杂微波部件结构的微放电阈值仿真。然而粒子模拟方法需要在剖分的网格内更新粒子状态和网格节点上的场分布，计算硬件开销大、费时长。mc方法由于不考虑电子与电磁场间的互作用，在时间开销上相对于pic方法更有优势，然而，目前提出的单粒子
‑
多碰撞、多粒子
‑
单碰撞、伪多粒子
‑
多碰撞三种蒙特卡罗微放电阈值计算方法由于算法本身的问题，在计算精度上仍有待于提高。


技术实现要素：

4.本发明解决的技术问题是：克服现有技术的不足，提供了一种多粒子多碰撞微放电阈值蒙特卡罗计算方法，
5.本发明的技术解决方案是：一种多粒子多碰撞微放电阈值蒙特卡罗计算方法，包括如下步骤：
6.s1，按照初始电子能量、角度、相位所符合的概率分布，随机生成n个初始电子，n个初始电子状态为：位置r
i
、能量e
i
、角度(θ
i
,ф
i
)，相位t
i
；
7.s2，按照输入频率f和输入功率w，计算微波部件内部的电磁场；
8.s3，更新电子进入电磁场的状态，电子状态为：位置r
in
、能量e
in
、角度(θ
in
,ф
in
)，相位t
in
，对于初始电子两者状态相同；
9.s4，计算时间步长δt内的电子运动轨迹；
10.s5，判断时间步长δt内的电子所处位置是否到达部件边界，如果没有到达，重复
s4～s5；反之，则记录到达边界处电子碰撞的状态：位置r
p
、能量e
p
、角度(θ
p
,ф
p
)，相位t
p
；
11.s6，采用碰撞时刻的电子能量e
p
和角度(θ
p
,ф
p
)，根据微波部件表面材料特性，计算出二次电子发射系数(δ
ts
、δ
e
、δ
r
)、二次电子电子的各能谱分布(f
1,e
、f
1,r
、f
n,ts
)及产生n个电子的概率p
n
，并计算出射电子的能量e
si
、角度(θ
si
,ф
si
)；i＝0,1,
…
,n；
12.s7，判断是否达到设定的截止条件，如果达到，则判断是否发生微放电，并结束；如果没有达到，则返回s3。
13.进一步地，所述产生n个电子的概率为
[0014][0015]
其中，分别为不产生本征二次电子的概率，产生1个本征二次电子的概率，产生n个本征二次电子的概率。
[0016]
进一步地，所述计算出射电子的能量e
s
的方法包括如下步骤：
[0017]
产生0～1内的随机数rand1，比较rand1与产生n个电子的概率p
n
，判断出实际产生的电子个数n；
[0018]
若n＝0，不产生二次电子；
[0019]
若n＝1，根据x(1)＝f
1,e
f
1,r
f
1,ts
的概率分布计算出射电子的能量e
s1
，其中，f
1,e
为产生1个弹性散射电子的能谱，f
1,r
为产生1个反射电子的能谱，f
1,ts
为产生1个本征二次电子的能谱；
[0020]
若n≥2，根据x(n)＝f
n,ts
的概率分布计算出射电子的能量e
s1
、e
s2
…
e
sn
；其中，f
n,ts
为产生n个本征二次电子的能谱。
[0021]
进一步地，第i个出射电子的出射角度(θ
si
,ф
si
)按照如下步骤计算：
[0022]
sinθ
si
为[
‑
1,1]之间抽样的随机数，ф
si
为[0,2π]之间抽样的随机数。
[0023]
进一步地，所述截止条件为达到最大循环次数，或计算时间，或电子数目，或人为设定的其他条件。
[0024]
进一步地，所述判断是否发生微放电的方法包括如下步骤：
[0025]
若以循环次数collsion_num作为判断微放电阈值方法，则计算δ
equl
＝(δ1 δ2
…
δ
collsion_num
)/collsion_num；其中，δ1、δ2、
…
、δ
collsion_num
分别为第i次碰撞的二次电子发射系数平均值；如果δ
equl
>1，则该计算条件下发生微放电；如果δ
equl
<1，则该计算条件下不发生微放电；如果δ
equl
＝1，该电压为微放电阈值电压；
[0026]
若以随时间演化的电子数目变化趋势作为判断微放电阈值方法，则电子数目随时间而减少，则该计算条件下不发生微放电；如电子数目随时间而增多，则该计算条件下发生微放电；如果电子数目随时间不变，则该电压为微放电阈值。
[0027]
进一步地，所述二次电子发射系数δ＝δ
ts
δ
e
δ
r
，δ
ts
为本征二次电子的发射系数，δ
e
为弹性散射电子的发射系数，δ
r
为反射电子的发射系数。
[0028]
进一步地，电子在电磁场中的运动轨迹求解方法为：
[0029]
直角坐标系中，电子的轨迹由决定；
[0030]
圆柱坐标系中，电子的轨迹由决定；
[0031]
其中，(x，y，z)为直角坐标系中电子在电磁场中的x、y、z向的位置，为圆柱坐标系中电子在电磁场中的径向、角向、轴向的位置，在电子所处位置上加一点代表了该方向的速度，加两点代表了该方向的加速度，e、b分别为该方向的电场和磁场，q为电子的单位电量，m为电子质量。
[0032]
进一步地，电子碰撞到边界后，更新电子进入电磁场的状态为：电子位置r
in
更新为电子碰撞的位置、能量e
in
更新为步骤6中出射电子的能量、角度(θ
in
,ф
in
)更新为步骤6中出射电子的角度，相位t
in
更新为电子碰撞到边界的时刻。
[0033]
一种计算机可读存储介质，所述的计算机可读存储介质存储有计算机程序，所述的计算机程序被处理器执行时实现所述一种多粒子多碰撞微放电阈值蒙特卡罗计算方法的步骤。
[0034]
本发明与现有技术相比的优点在于：
[0035]
(1)本发明克服了现有蒙特卡罗方法在计算过程中无论碰撞电子的二次电子发射系数为何值，均仅重新产生一个电子参与下一次碰撞的问题，本发明按照实际中二次电子发射的过程，可能出射实际个数的二次电子，其出射电子个数符合一定的概率分布，但多个相同状态的电子碰撞所产生的电子个数的统计值一定等于二次电子发射系数，提升了计算结果的准确性；
[0036]
(2)本发明以多个初始电子的随机状态作为起始状态，计算中又包含了电子碰撞部件表面后的二次电子发射过程中所有能量的电子，显著提升了计算结果的精度；
[0037]
(3)本发明在电磁场中推进电子运动轨迹的过程中，不需要按照粒子模拟的方法按照时间步长更新剖分网格内的粒子状态和网格节点上的场分布，因而避免了计算硬件开销大、费时长等问题。；
附图说明
[0038]
图1为本发明方法流程示意图；
[0039]
图2为本发明实施例微放电阈值计算与其他计算方法的对比结果。
具体实施方式
[0040]
为了更好的理解上述技术方案，下面通过附图以及具体实施例对本技术技术方案做详细的说明，应当理解本技术实施例以及实施例中的具体特征是对本技术技术方案的详细的说明，而不是对本技术技术方案的限定，在不冲突的情况下，本技术实施例以及实施例中的技术特征可以相互组合。
[0041]
以下结合说明书附图对本技术实施例所提供的一种多粒子多碰撞微放电阈值蒙特卡罗计算方法做进一步详细的说明，具体实现方式可以包括(如图1～2所示)：
[0042]
(1)按照初始电子能量、角度、相位所符合的概率分布，随机生成n个初始电子，n个初始电子状态为：位置r
i
、能量e
i
、角度(θ
i
,ф
i
)，相位t
i
；
[0043]
(2)按照输入频率f，输入功率w，计算微波部件内部的电磁场；
[0044]
(3)更新电子进入电磁场的状态，电子状态为：位置r
in
、能量e
in
、角度(θ
in
,ф
in
)，相位t
in
，对于初始电子两者状态相同；
[0045]
(4)采用四阶龙格
‑
库塔法计算时间步长δt内的电子运动轨迹；
[0046]
(5)判断时间步长δt内的电子所处位置是否到达部件边界，如果没有到达，按照(4)不断推进，如果到达边界，记录到达边界处电子碰撞的状态：位置r
p
、能量e
p
、角度(θ
p
,ф
p
)，相位t
p
；
[0047]
(6)采用碰撞时刻的电子能量e
p
和角度(θ
p
,ф
p
)，根据微波部件表面材料特性，基于furman模型计算出二次电子发射系数(δ
ts
、δ
e
、δ
r
)、二次电子电子的各能谱分布(f
1,e
、f
1,r
、f
n,ts
)及产生n个电子的概率
[0048][0049]
(6a)产生0～1内的随机数rand1，比较rnd1与p
n
，判断出实际产生的电子个数n；
[0050]
(6b)若n＝0，不产生二次电子；
[0051]
(7)若n＝1，根据x(1)＝f
1,e
f
1,r
f
1,ts
的概率分布计算出射电子的能量e
s1
，其中，f
1,e
为产生1个弹性散射电子的能谱，f
1,r
为产生1个反射电子的能谱，f
1,ts
为产生1个本征二次电子的能谱；；
[0052]
若n≥2，根据x(n)＝f
n,ts
的概率分布计算出射电子的能量e
s1
、e
s2
…
e
sn
；
[0053]
(8)判断是否达到设定的截止条件，如果达到判断是否发生微放电，如果没有达到，返回步骤(3)。
[0054]
具体的，在本技术实施例所提供的方案中，包括如下步骤：平行平板，1mm间隙，表面采用furman
‑
ag材料
[0055]
(1)产生n＝1000个初始电子，1000个初始电子能量e0满足分布
σ＝3ev,μ＝5ev；初始电子角度满足分布f(θ0)
‑
sin2θ0；初始电子方位角在[0,2π]随机分布，初始电子相位t0在[0,t
wave
]随机分布，则出射电子速度为
[0056]
(2)计算x＝1mm间隙的平行平板电磁场，频率依次设定为f＝1.06ghz、2.14ghz、3.71ghz、5.33ghz、7.97ghz，其内部电磁场由下式确定e＝[
‑
vsin(ωt
‑
βz),0,0]
[0057]
(3)电子进入电磁场的状态初始时刻
[0058]
(4)电子在电磁场中运动轨迹的推进由四阶龙格
‑
库塔法求解，
[0059][0060]
设定时间步长δt＝t
wave
/50，在dt范围内按照上式求解电子轨迹，第一次求解时[v
x v
y v
z
]为步骤(1)求解的初始电子初速度；
[0061]
(5)每在δt时间步长内求解一次电子轨迹就判断电子的x向轨迹是否满足x>1mm或者x<0mm，如果不满足，进入步骤(4)，[v
x v
y v
z
]设定为步骤(4)上一次求解结果；
[0062]
(6)电子此时与器件表面发生碰撞，碰撞位置x
p
，碰撞时刻t
p
，碰撞次数collsion_num＝1(后续每碰撞一次，collsion_num 1)，碰撞速度通过步骤(4)求解的结果得到[v
px v
py v
pz
]，则相对于材料表面入射电子能量入射电子角度根据furman模型计算出出射的电子个数、二次电子发射系数δ1及能量，假设出射了2个电子，对应能量为[e
s1
,e
s2
]，角度
[0063]
(7)判断collsion_num<10000，如果没有则返回步骤(3)，电子进入电磁场状态分别为和
[0064]
判断是否发生微放电的方法包括如下步骤：
[0065]
若以循环次数collsion_num作为判断微放电阈值方法，则当collsion_num>10000
时，计算δ
equl
＝(δ1 δ2
…
δ
collsion_num
)/collsion_num；其中，δ1、δ2、
…
、δ
collsion_num
分别为第i次碰撞的二次电子发射系数平均值；如果δ
equl
>1，则该计算条件下发生微放电；如果δ
equl
<1，则该计算条件下不发生微放电；如果δ
equl
＝1，该电压为微放电阈值电压。
[0066]
进一步，二次电子发射系数δ＝δ
ts
δ
e
δ
r
，δ
ts
为本征二次电子的发射系数，δ
e
为弹性散射电子的发射系数，δ
r
为反射电子的发射系数。
[0067]
若以随时间演化的电子数目变化趋势作为判断微放电阈值方法，则电子数目随时间而减少，则该计算条件下不发生微放电；如电子数目随时间而增多，则该计算条件下发生微放电；如果电子数目随时间不变，则该电压为微放电阈值。
[0068]
进一步，电子在电磁场中的运动轨迹求解方法为：
[0069]
直角坐标系中，电子的轨迹由决定；
[0070]
圆柱坐标系中，电子的轨迹由决定；
[0071]
其中，(x，y，z)为直角坐标系中电子在电磁场中的x、y、z向的位置，为圆柱坐标系中电子在电磁场中的径向、角向、轴向的位置，在电子所处位置上加一点代表了该方向的速度，加两点代表了该方向的加速度，e、b分别为该方向的电场和磁场，q为电子的单位电量，m为电子质量。
[0072]
可选的，在一种可能实现的方式中，电子碰撞到边界后，更新电子进入电磁场的状态为：电子位置r
in
更新为电子碰撞的位置、能量e
in
更新为步骤6中出射电子的能量、角度(θ
in
,ф
in
)更新为步骤6中出射电子的角度，相位t
in
更新为电子碰撞到边界的时刻。
[0073]
基于蒙特卡罗方法的传统mc、单粒子
‑
多碰撞(mc1)、多粒子
‑
单碰撞(mc2)、多粒子
‑
多碰撞(mc3)的计算结果来自文献(林舒,闫杨娇,李永东,刘纯亮2014物理学报63 147902)，cst软件的计算结果来自文献(王洪广,翟永贵,李记肖,李韵,王瑞,王新波,崔万照,李永东2016物理学报65 237901)，本方法计算结果如图2所示，计算效率如表1所示。
[0074]
表1平板传输线微放电阈值计算效率比较
[0075][0076]
本技术提供一种计算机可读存储介质，所述计算机可读存储介质存储有计算机指令，当所述计算机指令在计算机上运行时，使得计算机执行图1所述的方法。
[0077]
本领域内的技术人员应明白，本技术的实施例可提供为方法、系统、或计算机程序产品。因此，本技术可采用完全硬件实施例、完全软件实施例、或结合软件和硬件方面的实施例的形式。而且，本技术可采用在一个或多个其中包含有计算机可用程序代码的计算机可用存储介质(包括但不限于磁盘存储器和光学存储器等)上实施的计算机程序产品的形式。
[0078]
本技术是参照根据本技术实施例的方法、设备(系统)、和计算机程序产品的流程图和/或方框图来描述的。应理解可由计算机程序指令实现流程图和/或方框图中的每一流程和/或方框、以及流程图和/或方框图中的流程和/或方框的结合。可提供这些计算机程序指令到通用计算机、专用计算机、嵌入式处理机或其他可编程数据处理设备的处理器以产生一个机器，使得通过计算机或其他可编程数据处理设备的处理器执行的指令产生用于实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能的装置。
[0079]
这些计算机程序指令也可存储在能引导计算机或其他可编程数据处理设备以特定方式工作的计算机可读存储器中，使得存储在该计算机可读存储器中的指令产生包括指令装置的制造品，该指令装置实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能。
[0080]
这些计算机程序指令也可装载到计算机或其他可编程数据处理设备上，使得在计算机或其他可编程设备上执行一系列操作步骤以产生计算机实现的处理，从而在计算机或其他可编程设备上执行的指令提供用于实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能的步骤。
[0081]
显然，本领域的技术人员可以对本技术进行各种改动和变型而不脱离本技术的精神和范围。这样，倘若本技术的这些修改和变型属于本技术权利要求及其等同技术的范围之内，则本技术也意图包含这些改动和变型在内。
[0082]
本发明说明书中未作详细描述的内容属本领域技术人员的公知技术。