ARTICLE · INTELLIGENCE

战地情报 · 详情页

来自尧图项目组的一线实战观察与深度解析

基于COMSOL的棒-棒电极流注放电仿真与电子密度演化分析

基于COMSOL的棒-棒电极流注放电仿真与电子密度演化分析

流注放电,这四个字放在一起,搞高电压、搞脉冲功率、搞雷电物理的人看到都会多看两眼。工程里遇到的绝缘击穿、电网闪络,物理上最锋利的阶段恰恰就是流注放电这一段——它不像汤森放电那样缓慢暗涨,也不像火花击穿那样瞬间势不可挡,而是以每秒钟几十万米的速度在气体间隙里撕开一条高电离的通道。为了把这段演化过程看得清清楚楚,“棒-棒电极 + 空气流注放电 + 电子密度演化”这个组合就成了非常经典的研究对象。我这次就用 COMSOL 做了一版棒-棒电极流注放电仿真,写完确实有点心得想说。

先说明一下,这篇内容适合谁看:如果你要用仿真复现放电演化过程、需要观察电子密度怎么从种子电子一路涨到 (10^{19})、(10^{20}/\text{m}^3) 量级,或者你在做高压设备绝缘配合时想评估一段空气间隙的自由穿与阻断条件,那这篇文章的建模思路、参数选型和踩坑记录,应该能帮你节省一到两周的调试时间。已经熟到能自己写等离子体程序的人,可以重点看后面的数值鲁棒性处理部分;纯新手,也可以跟着完整做一遍,比一头扎进手册里有效率得多。

1. 棒-棒电极流注放电:从物理图像到仿真建模

1.1 流注放电是什么、解决什么问题

放电这个词在气体里其实更像一场“接力赛”。外电场把种子电子加速,电子和中性分子碰撞电离出新的电子和正离子,这个雪崩过程不断放大。早期阶段电流很小,雪崩也比较稀疏,这叫汤森放电;但电场一旦高到某个程度,雪崩头部积累的空间电荷开始反过来扭曲电场,局部场增强到临界水平,电离就出现了“时空自组织”——高速传播的高电离通道从电极尖端长出,这就是流注。

流注的核心机制说起来不复杂,就是空间电荷场的正反馈。正流注从阳极往阴极爬,头部是一团净正电荷,这些正电荷在它前方加强了朝向阳极方向的电场,让前方雪崩疯狂生长;通道内部因为带电位差被空间电荷屏蔽,场强反而很低。看到这你会发现,流注物理的本质其实就是电子连续性方程、电场演化方程和化学反应源项之间的强耦合,纯实验很难把每个量都测出来,数值仿真就成了观察内部过程最顺手的工具。

流注放电研究在工程上是有实打实的意义的。输电线路空气间隙的击穿路径预测、气体开关的绝缘恢复时间优化、飞机或风机叶片遭遇雷击时的先导起始判据,还有脉冲功率装置里纳秒级放电的时空演化,背后都离不开流注模型。棒-棒电极又是其中最“刁钻”的一种构型——两根针对着打,电场极度不均匀,正负流注都能从各自电极尖端发起,你不仅能看到单支流注的传播,还能观察两支流注相向而行、谁先贯穿间隙的竞争过程,信息量比针-板构型大不少。

1.2 为什么用 COMSOL、为什么选棒-棒电极

说句实话,流注放电仿真用纯编程也能做,经典的做法是自编有限差分或 PIC 程序,再跑一个玻尔兹曼求解器算输运系数。但这条路门槛不低,光是把流体方程组和泊松方程在非结构网格上稳定耦合起来,就够磨一阵子了。COMSOL 在这里的优势是,等离子体模块把漂移-扩散-反应方程、静电方程、边界条件和自适应网格都放在同一个框架下,你可以把精力集中在物理建模和参数调试上,而不是纠结矩阵组装和 MPI 并行通信。

棒-棒电极在 COMSOL 里建起来也非常省事:二维轴对称就够了。两根细棒尖端相对,间隙长度取几个毫米到十几个毫米即可,计算域形状仍然简单,但电场不均匀性、流注起始条件、传播规律这些关键物理全都能保留下来。相比针-板电极,棒-棒构型在高电压工程里对应的是更真实的“间隙放电”场景,研究正负流注相互作用时也更有代表性。这也是我最终选这个构型的原因。

提示:做流注仿真,物理内核的优先级永远高于软件操作。网格再密、求解器再花哨,如果空间电荷效应和电离-附着平衡没设对,结果也只是“画得很漂亮的不对”。后面我会反复提到这一点。

2. 模型搭建的关键思路与方案选型

2.1 几何与边界条件设计的考量

流注放电仿真几何参数直接影响起始电压和电场增强倍数,不能随便拍脑袋。我这一版的等效结构如下:

参数取值说明
棒半径0.2 mm尖端电场增强与曲率半径成反比,太细则网格压力大
间隙长度10 mm常压空气下典型流注传播长度,时间尺度好分辨
棒端形状半球头曲率连续,避免锐边界带来的数值振荡
计算域半径25 mm至少为间隙的2倍以上,削弱外部边界对电场的约束
气压1 atm默认常压,迁移率、扩散系数按 (N=2.5\times 10^{25}/\text{m}^3) 换算
施加电压20 kV 正极性对于10 mm间隙,起始场强需要十几 kV,可在此附近调

尖端电场增强这里有个经验值:对于半球头细棒,最大电场大约为

[ E_{max}\approx\frac{V}{r\cdot \ln(2d/r)} ]

取 (r=0.2\text{ mm})、(d=10\text{ mm}),(E_{max}) 就能到每毫米几千伏的水平,足以驱动常压空气的流注起始。如果你把半径改成 0.05 mm,需要的时间步长会小一个量级,收敛压力立刻上来,建议从 0.2 mm 起步。

边界条件上我会分成三类:阳极棒和阴极棒表面是高压与接地,周向远端设置为电绝缘(或开放边界),对称轴用轴对称边界。等离子体壁面这里有一条容易被忽略——二次电子发射系数,常压放电里它对流注传播的影响不如低压明显,但建议还是给个 0.01~0.05 的默认值,不要设成 0,否则会造成数值上电子通量的不连续。

2.2 等离子体化学机制:简化空气模型的取舍

全免空气化学模型可以涉及几十个组分、上百个反应,但你做流注仿真,特别是想稳定观察电子密度的演化,一开始就上全套反应集基本是给自己找麻烦。最值得先跑通的,是下面的“三组分+三反应”简化模型:

  • 组分:电子 (e)、正离子 (p)、负离子 (n)(ODEs输出这三个变量)
  • 电离反应:(e + M \rightarrow 2e + p)
  • 附着反应:(e + M \rightarrow n)
  • 复合反应:(e + p \rightarrow M)、(n + p \rightarrow M)

注意这里的“M”代表空气分子,而不是某个具体组分。电离速率系数和附着速率系数要写成约化电场 (E/N) 的函数,可以通过查询公开的等离子体化学数据集,或者在本地用玻尔兹曼求解器预先算好一把表,再以插值函数形式灌进 COMSOL。

为什么不做全反应集?原因有三个:

  1. 大量增益反应会让方程组刚性陡增,时间步长被迫缩小,流注传播还没看到,已经算了几个星期。
  2. 流注物理的主要量级由电离-附着差值和迁移率决定,附加反应在常压空气流注尺度上主要起微调作用。
  3. 简化模型能让你更容易看清楚每一步物理,后续要加光电离、激励态或更精细的化学,再逐步扩展也不迟。

电子输运参数里,最核心的是电子迁移率 (\mu_e) 和扩散系数 (D_e)。常压空气里约化场强超过 45 Td 左右,电离速率就开始超过附着速率,电子密度进入净增长区间。流注头部附近的约化场强通常在 100~200 Td,迁移率大约在 300~500 cm²/(V·s) 量级。直接用常数也可以,但仿真精度会打折;我建议的折中方案是,用插值表给出 (\mu_e(E/N)),扩散系数通过爱因斯坦关系 (D_e = \mu_e T_e) 推算,电子温度也做一个约化场强的弱函数。

2.3 数值求解的参数选择与时间尺度

流注放电是典型的多尺度问题:电离反应的特征时间在皮秒量级,流注跨越 10 mm 间隙只需要几十纳秒,但整个过程的电磁瞬态可能又需要微秒级才有完整响应。COMSOL 求解这类问题,我会直接选时域瞬态,并打开自适应时间步进。

关键参数里,我建议时间步长下限设在 (10^{-14}) 秒量级。原因很简单:流注头部的高场区域宽度只有几微米,电子以 (10^5\sim10^6) m/s 的速度穿过一个网格单元,所花时间就是 (10^{-10}\sim10^{-11}) 秒等分到网格尺度后,大概到 (10^{-14}\sim10^{-12}) 秒。如果不限制最小步长,求解器很可能会跳过快速电离过程,导致“电子密度被负值污染”甚至直接发散。

空间网格上,流注通道内的网格尺度不低于 5~20 μm,流注头部前方的高梯度区建议通过自适应加密主动压到 1~5 μm。棒尖端的曲率半径处也建议局部加密,否则起始电场计算误差会导致起始时间大幅偏离实验值。

3. 核心实现:物理场设置与求解配置

3.1 全局定义与几何参数

打开 COMSOL,先建立二维轴对称模型。全局定义里把上面表格中的参数全部写成参数节点,这样后期调电压、调间隙时可以不动几何,只改参数表,跑批量工况非常方便。

几何本身可以这样做:以对称轴为 (r=0),在轴向上绘制两根细棒。棒体用矩形命名“Anode rod”和“Cathode rod”,棒尖用圆弧段连接形成半球头。两棒尖端相对距离设为 10 mm。外围计算域用一个大矩形或梯形包围,边界距离尖端至少 20 mm 以上。整个几何可以一次创建完成,不涉及CAD导入,更不容易出现拓扑错误。

这一阶段的重点检查项有两个。第一,轴向上两棒是否严格对称——棒-棒构型没有对称轴方向的对称性要求,但几何两侧的曲率半径需要一致,否则电场分布会偏。第二,半导体接头或尖锐角——任何小于毫米尺寸的奇点都会让局部电场在网格细化时无限升高,物理上不存在,数值上导致伪击穿。我的做法是把尖端倒角做成均匀圆角,宁可损失一点点“锋利”,也要保住电场分布的可信度。

3.2 等离子体模块的方程与边界条件

物理场层面,我用“等离子体-瞬态”接口配合“静电”接口。核心变量包括电子密度 (n_e)、电子能量密度 (n_\varepsilon)、正离子密度 (n_p)、负离子密度 (n_n),以及电位 (V)。方程组就是经典的漂移-扩散-反应形式:

[ \frac{\partial n_e}{\partial t} + \nabla \cdot \left(-D_e\nabla n_e - \mu_e \mathbf{E} n_e\right) = S_e ]

[ \frac{\partial n_\varepsilon}{\partial t} + \nabla \cdot \left(-D_\varepsilon\nabla n_\varepsilon - \mu_\varepsilon \mathbf{E} n_\varepsilon\right) = S_\varepsilon ]

[ \frac{\partial n_p}{\partial t} + \nabla \cdot \left(-D_p\nabla n_p + \mu_p \mathbf{E} n_p\right) = S_p ]

负离子方程形式类似。电场通过“静电”接口求解:

[ -\nabla\cdot(\varepsilon_r\varepsilon_0\nabla V) = e(n_p - n_n - n_e) ]

边界条件我总结如下:

  • 阳极棒:电位 (V = V_0),电子密度取“二次发射边界”,二次发射系数 ( \gamma=0.01)。
  • 阴极棒:接地 (V=0),电子通量允许二次发射。
  • 计算域外边界:默认绝缘即 (\partial V/\partial n=0),而粒子通量也设置为“出流”。
  • 对称轴:轴对称,所有量的径向导数为0。

这里有一个很容易踩的坑:如果你把正离子或负离子的壁面通量设成“零通量”,它们会在电极表面积累成一层数值边界层,导致电极附近的电场和粒子密度分布失真。我建议壁面统一用“通量条件”让离子自然吸附,电极本体不参与化学反应。

3.3 网格、求解器与时间步进策略

网格划分上,这个模型不是均匀网格能解决的。流注头部的高电场区域很小,但动力学剧烈,所以我采用“分层加密”策略:

  • 尖端周围和间隙轴线附近:固定网格尺寸 5 μm;
  • 计算域外围:粗到 1~2 mm;
  • 中间过渡区:用最大单元增长率控制在 1.3 以下,避免网格尺寸跳变引起数值反射。

COMSOL 的等离子体接口里可以开“自适应网格”,但我不建议一开始就开。原因是你还没把物理场跑出正当形态时,自适应判断很容易被数值噪声误导,把网格全堆在无意义的区域。先把固定网格跑通一轮,确认流注传播特征出现,再开自适应提高分辨率,这样每一步的变量都在你控制范围内。

时间步进上,选择 BDF 二阶格式,最大阶数限制为 2 更稳定。初始时间步设为 (10^{-13}) 秒,最大时间步不要超过 1 ns。求解器容差我一般设成相对 (10^{-4}) 绝对 (10^{-8}) 到 (10^{-6}) 之间与各变量尺度匹配。COMSOL 里不同变量量级差异极大,电子密度在低场区可能只有 (10^{10}),流注通道里却是 (10^{20}),绝对容差必须按变量分别设置,不能一个值通吃所有方程。

4. 放电演化过程的观察与分析

4.1 电子密度的时空演化:从种子电子到流注通道

初始值不能设成零。真实气体里总有宇宙射线和背景辐射维持着一个微弱的本底电离,常压空气中电子密度大概在 (10^{10}\sim10^{13}/\text{m}^3) 量级。我在模型里赋了均匀的初始电子密度 (10^{12}/\text{m}^3),同时为了确保阳极附近第一代雪崩有足够的“启动者”,在棒尖附近额外叠了一个局部高斯分布的种子电子团,中心电子密度 (10^{15}/\text{m}^3),半径 0.3 mm。

不加这个局部种子,仿真也能跑,但流注起始时间会明显延后,而且起始位置受网格非对称性影响,可能从棒尖侧壁而不是最尖端发起,物理上就不对了。

放电开始后的演化大致分三段:

  • 第一阶段(0~2 ns):棒尖高场区内种子电子原地雪崩,电子密度从 (10^{12}) 快速涨到 (10^{16}) 以上,由于电离源项强烈依赖约化场强,雪崩区集中在强场几十微米范围内;
  • 第二阶段(2~20 ns):空间电荷场开始主导。正离子堆积在雪崩头部,头部前方场增强,电离区开始离开电极表面,形成向前移动的高电离波前,这就是流注头;
  • 第三阶段(20 ns 以后):流注通道建立。通道内部电子密度达到 (10^{19}\sim10^{20}/\text{m}^3),电场被屏蔽到远低于击穿阈值,电子密度曲线出现典型的“头部陡降、通道平缓”结构。

建议后处理时把电子密度显示切换为对数坐标。线性坐标下头部和通道内密度的巨大差异会让图像糊成一团,对数坐标才能清楚看到通道延伸和头部推进位置。我通常用“Surface plot”加“log10”表达式,色标范围 (10^{10}\sim10^{20}),这样每一帧的视觉对比度都很舒服。

4.2 电场增强与空间电荷效应

流注放电的灵魂不只是电子密度,而是电场如何被空间电荷改造。仿真里最能说明问题的是沿间隙轴线的电场分布。

初始阶段,即使没有流注,棒尖附近电场也远高于平均值,呈指数衰减;流注形成后,在流注头部位置会出现一个幅值远高于背部电场的电场尖峰,典型峰值可达到背部的 5~10 倍。这个尖峰代表的就是净正电荷积聚区域,它的存在驱动着头部前方不断产生新雪崩,维持流注自持传播。

而通道内部则相反,空间电荷屏蔽让电场掉到接近甚至低于临界电离场强的水平。你从电场曲线图上能很直观地看到“一个高峰,后面一片洼地”的形态——高峰是流注头,洼地是已贯穿的导电通道。

这里有个物理细节需要特别注意:正流注和负流注的空间电荷结构不同。正流注头部是正离子堆积,电场增强发生在头部前方;负流注头部则是电子层和正离子分离,电场结构更复杂,传播也更不稳定。从实验上看,负流注更容易分支;我的仿真是正极性棒-棒结构,模拟结果是以单个主通道向前推进为主。如果你想研究分支,需要额外引入光电离项和随机扰动,光靠流体方程加高斯背景电离,很难自发出现树状分支。

4.3 流注传播速度、通道半径与诊断

做完一组仿真后,要从原始数据里提炼准物理量,不能只停留在“视频好看”层面。我最常用的三个诊断量:

第一是流注传播速度。以最大电子密度位置(头部)为追踪对象,沿着对称轴提取每一时刻头部位置,然后做位置-时间曲线的数值微分。放电中期速度大约在 (2\sim5\times10^5) m/s,这个量级和文献中的常压空气正流注速度基本吻合。如果仿出来的速度超过 (10^7) m/s,多半是电离系数给得太大或者电压超出了真实条件;如果低于 (10^4) m/s,那可能流注还没真正形成,只是在爬行。

第二是通道半径。以某一时刻的电子密度分布为基准,取头部后方 1~2 mm 处的径向剖面,用半高峰宽(FWHM)来定义通道半径。因为这个量对可视化参数不敏感,适合做定量比较。常压空气流注通道半径一般在 0.1~0.5 mm 量级,如果你的剖面画出来半径超过 2 mm,基本可以断定网格太粗或者扩散系数被高估了。

第三是电流演化。电极上的传导电流可以用 (I = \int_S (e \mu_e \mathbf{E} n_e) dS) 来提取,看它随时间是否出现“阶梯式”上升,那是新流注分支或通道贯穿的典型信号。电流曲线的形状在工程上还直接对应着击穿危险的诊断窗口。

5. 实操记录:从求解失败到稳定收敛

5.1 常见错误特征与排查

仿真流注的过程,说“十次有八次不能一次跑通”一点也不夸张。我在这一版模型上调试时遇到的典型问题,整理成下表:

现象可能原因解决方案
计算一开始就发散初始电子密度太低或网格太粗,泊松方程无法耦合把背景电子密度提到 (10^{12}/\text{m}^3),加密棒尖附近网格
电子密度出现负值对流占优时的数值振荡,常见于头部高梯度区对电子方程采用一阶迎风格式,或把头部网格尺寸降一个量级
流注传播速度异常偏慢电离速率系数被低估,或有效场强被边界削弱复核 (E/N) 插值表,确认是否达到了45 Td以上的净增长区间
流注头部在网格间“跳跃”时间步太大导致头部一帧跑过多个网格固定最大时间步为 0.1 ns 以下,再配合自适应加密
两个棒尖同时出现伪放电计算域边界离电极太近,边界电场不真实扩大外部计算域或改用开放边界条件
数值结果在某个时间点突然爆炸复合项或附着项在极端密度下数值溢出给源项加密度低限保护,比如ne>=0,同时调整绝对容差

最容易被新手忽略的是第一行的“初始电子密度”。有人为了追求“干净的物理”,把初始电子密度设成 0 或 (10^5/\text{m}^3),结果求解器在第一步就报奇异,因为电势方程右侧没有空间电荷,但左侧却有很强的拉普拉斯项,方程变得病态。真实空气里常温常压下总有过量种子电子,设得太低反而与物理不符。

5.2 收敛性调试的具体技巧

接下来聊几个调试技巧,都是我自己踩坑总结出来的。

第一个技巧是先跑短间隙,再跑目标间隙。我最初直接从 10 mm 间隙开始,结果流注传播到中段后由于网格和库朗数问题反复崩溃。换成 3 mm 间隙调试完成后,再拉长到 10 mm,整个流程顺了很多。短间隙里流注能在更短物理时间内贯穿,调试周期指数级缩短。

第二个技巧是用手动时间步扫描替换全自动步进。COMSOL 的自动时间步进在刚性问题下非常保守,步长掉到 (10^{-16}) 秒,进度条几乎不走。我的做法是手动设定一段“步进日志”:0~5 ns 用 (10^{-14}) 秒固定步长,5~20 ns 用 (10^{-13}) 秒,20 ns 以后逐步放宽。这么做反而比全自动策略更快更稳定。

第三个技巧是做一场“预热”仿真。真正研究流注之前,先只求解静电电场,把电势分布和最大场强算出来,确认棒尖最大电场在预期范围内。这一场不需要等离子体物理,但能快速排查几何和边界的问题。我碰到过棒尖倒角半径设成 0.05 mm 后最大电场高到 100 kV/mm 都不止,明显超出工程常识,就是预热仿真发现的。

第四个技巧和求解器相关:对每个变量单独设置绝对容差。COMSOL 默认会给你一个统一的绝对容差,比如 (10^{-6}),对电子密度来说这个值太小了。电子密度的典型值在 (10^{12}\sim10^{20}) 之间,统一容差会导致相对误差判断失真。我是把正负离子密度和电子密度的绝对容差设在 (10^{6}\sim10^{9}),电位容差设在 (10^{-3}) V,这样求解器不会把精力浪费在追平无意义的微小密度上,收敛速度有明显提升。

5.3 数据导出与可视化要点

数据后处理不要全堆在 COMSOL 界面里做,我做仿真的习惯是“原场可视化 + 外部数据导出”两层结构。

原场可视化重点看三张图:电子密度对数颜色的空间分布、沿轴的电场强度线图、电子密度沿轴的剖面图。这三张图基本能回答“流注走到哪了、头部场强多少、通道密度结构如何”三个最核心的问题。动画导出可以用 COMSOL 的“动画”功能,把时间序列导出成 avi 或图片序列,方便后续做展示或对比。

外部数据导出我一般用“一维截线”功能,沿对称轴提取 (n_e)、(E_z)、(V) 等量,导出成 CSV,再用数据处理脚本画传播速度曲线和通道半径曲线。这一步的好处是,你能把多组不同电压、不同间隙下得到的数据放在同一坐标系里对比,直接看出参数趋势,而不必反复在软件里切视图。

关于网格无关性验证,我给个务实建议:至少跑两套网格——一套“基础网格”峰值尺寸 10 μm,一套“加密网格”峰值尺寸 3 μm,比较流注传播速度和头部电子密度峰值。如果两者差异在 5% 以内,说明网格已经收敛;如果差异超过 20%,就继续加密,直到关键量稳定为止。这个验证得花时间,但这是让你的仿真结果能写进论文或者报告里的最低门槛。

另外,保存求解结果时勾选“所有时间步”。我吃过一次亏,默认设置只保存了部分求解步,结果后面的速度曲线恨不得只有三个点,根本画不成线,只能重跑一遍,浪费了大半天。

讲心里话,COMSOL 把流注放电这个强非线性、强耦合的问题封装到了“可以操作”的程度,但这不代表它是个傻瓜化软件。真正让结果变得可信的,还是你对空间电荷效应、电离阈值 45 Td、流注头部网格分辨率和刚性问题时间尺度这几个物理切面的理解。把简化三组分模型跑通、跑稳,你已经对放电演化过程有了很直观的感知;等哪天真要研究分支流注、光电离或者复杂流注-先导转换,再从这一版骨架往上加机制,就会顺畅得多。

RELATED READING

延伸阅读

更多一线实战笔记与深度复盘,助您持续精进