标准有限元方法模拟大气压流注等离子体中数值振荡与伪扩散问题研究与改进

气体物理

２０２１年３月

ＰＨＹＳＩＣＳＯＦＧＡＳＥＳ

Ｖｏｌ．Ｍａｒ．６　Ｎｏ．２０２１

２　　ＤＯＩ １０．１９５２７／ｊ．ｃｎｋｉ．２０９６⁃１６４２．０８７１

标准有限元方法模拟大气压流注等离子体中

数值振荡与伪扩散问题研究与改进

朱益飞１ ２（１．空军工程大学航空工程学院等离子体动力学重点实验室 　陈贤聪１

２．西安交通大学机械工程学院 陕西西安７１００４９）

陕西西安７１００３８

ＳｏｌｕｔｉｏｎｓＡｔｍｏｓｐｈｅｒｉｃｆｏｒＮｕｍｅｒｉｃａｌＰｒｅｓｓｕｒｅＩｎｓｔａｂｉｌｉｔｙＺＨＵＹｉ⁃ｆｅｉＰｌａｓｍａａｎｄＰｓｅｕｄｏ１ ２ 　ＢａｓｅｄＣＨＥＮｏｎＤｉｆｆｕｓｉｏｎＸｉａｎ⁃ｃｏｎｇＦｉｎｉｔｅ１

ＥｌｅｍｅｎｔＰｒｏｂｌｅｍｓＭｅｔｈｏｄ

ｉｎＭｏｄｅｌｉｎｇｏｆ

２．Ｓｃｈｏｏｌ（１．ＳｃｉｅｎｃｅｏｆＭｅｃｈａｎｉｃａｌａｎｄＴｅｃｈｎｏｌｏｇｙＥｎｇｉｎｅｅｒｉｎｇｏｆＰｌａｓｍａ Ｘｉ′ａｎＤｙｎａｍｉｃｓＪｉａｏｔｏｎｇＵｎｉｖｅｒｉｓｔｙＬａｂ Ｘｉ′ａｎ Ｘｉ′ａｎ７１００３８ ７１００４９ ＣｈｉｎａＣｈｉｎａ

）

摘　要 基于商用软件的等离子体数值模拟方案大多是基于标准Ｇａｌｅｒｋｉｎ格式的有限元方法 在计算大气压流注放电时 存在数值振荡和伪扩散的问题 计算结果可靠性不高．针对针⁃板体流注放电结构 构建了修正的ＣＯＭＳＯＬ基于弱形式的人工稳定项Ｍｕｌｔｉｐｈｙｓｉｃｓ等离子体模型．利用修正模型复现了经典的雪崩．修正模型与ＣＯＭＳＯＬＭｕｌｔｉｐｈｙｓｉｃｓ⁃流注转捩算例内置等离子体模块的区别主要在于增加了 并与经典纳秒秒冲过电压针板放电实验结果进行了对比验证．结果表明 引入修正后的商用软件模型能够克服数值振荡和伪扩散缺陷 确保等离子体数值计算的准确性和可复现性．关键词 流注 等离子体 针板放电

气体放电 有限元数值计算

　中图分类号Ａｂｓｔｒａｃｔ Ｔｈｅ ＴＬ６１ｓｔａｎｄａｒｄ＋２．３

Ｇａｌｅｒｋｉｎｆｉｎｉｔｅ　文献标志码ｅｌｅｍｅｎｔ　ｍｅｔｈｏｄ Ａ

（ＦＥＭ）ｉｓｗｉｄｅｌｙｕｓｅｄｉｎｃｏｍｍｅｒｃｉａｌｐｌａｓｍａｍｏｄｅｌｉｎｇｓｏｆｔｗａｒｅ．

ＴｈｅＣＯＭＳＯＬｎｕｍｅｒｉｃａｌｔｒｏｄｕｃｔｉｏｎＭｕｌｔｉｐｈｙｓｉｃｓｉｎｓｔａｂｉｌｉｔｙｏｆａｒｔｉｆｉｃｉａｌｓｔａｂｉｌｉｚａｔｉｏｎ．ｐｌａｓｍａａｎｄｐｓｅｕｄｏｍｏｄｅｌＴｈｅｗａｓｄｉｆｆｕｓｉｏｎｐｒｏｐｏｓｅｄｃｏｎｓｔｒｕｃｔｅｄｈａｖｅｌｏｎｇｍｅｔｈｏｄｆｏｒｔｈｅｂｅｅｎｗａｓｐｉｎ⁃ｔｏ⁃ｐｌａｎｅａｐｒｏｂｌｅｍｖａｌｉｄａｔｅｄｂｙｄｉｓｃｈａｒｇｅｒｅｓｔｒｉｃｔｉｎｇａｃｌａｓｓｉｃａｌｃａｓｅ．ｔｈｅａｃｃｕｒａｃｙａｖａｌａｎｃｈｅ⁃ｔｏ⁃ｓｔｒｅａｍｅｒＴｈｅｍａｉｎｏｆｉｍｐｒｏｖｅｍｅｎｔｔｈｅｍｏｄｅｌ．Ａｂｅｎｃｈｍａｒｋｗａｓｃｏｒｒｅｃｔｅｄｔｈｅａｎｄｉｎ⁃ｂｙｅｘｐｅｒｉｍｅｎｔａｌｒｅｓｕｌｔｓｏｆａｎｏｖｅｒｖｏｌｔａｇｅｄｉｓｃｈａｒｇｅ．Ｔｈｅｒｅｓｕｌｔｓｓｈｏｗｔｈａｔ ｔｈｅｃｏｒｒｅｃｔｅｄＦＥＭｍｏｄｅｌｃａｎｅｆｆｅｃｔｉｖｅｌｙｏｖｅｒ⁃ｓｉｍｕｌａｔｉｏｎｃｏｍｅｔｈｅｐｒｏｂｌｅｍｏｆｐｌａｓｍａｏｆｄｉｓｃｈａｒｇｅ．

ｎｕｍｅｒｉｃａｌｉｎｓｔａｂｉｌｉｔｙａｎｄｐｓｅｕｄｏｄｉｆｆｕｓｉｏｎ ａｎｄｅｎｓｕｒｅｔｈｅａｃｃｕｒａｃｙａｎｄｒｅｐｅａｔａｂｉｌｉｔｙｏｆｔｈｅｎｕｍｅｒｉｃａｌＫｅｙｗｏｒｄｓ ｓｔｒｅａｍｅｒ ｐｌａｓｍａ ｐｉｎ⁃ｔｏ⁃ｐｌａｎｅｄｉｓｃｈａｒｇｅ ｇａｓｄｉｓｃｈａｒｇｅ ＦＥＭ

引　言

助燃、薄膜工艺 低温等离子体的应用领域正在迅大气压低温等离子体由于能够可控、节能、高速扩展［１⁃８］效地产生活性粒子 近年来越来越受到工业和科学在大气压下．

产生等离子体的基本方法之一是

界重视．从臭氧产生、聚合物处理、激光制造、污在由金属制造的针板结构的气隙间施加高电压．高染控制、生物除菌、医药治疗 到流动控制、点火

电压在针尖形成局部强电场 驱动电子加速并引发

收稿日期 ２０２０⁃０８⁃１２ 修回日期 ２０２０⁃０８⁃２６

基金项目 国家自然科学基金（５１９０７２０４ ５１７９０５１１ ９１９４１１０５ ９１９４１３０１）

第一作者简介 朱益飞（１９８９⁃）　男 博士 助理研究员 主要研究方向为纳秒脉冲等离子体与数值计算．Ｅ⁃ｍａｉｌ ｓｙｊｃｏｏｌｒ＠

１６３．ｃｏｍ　　　引用格式　　　 朱益飞１⁃８．　　　 陈贤聪　　　．标准有限元方法模拟大气压流注等离子体中数值振荡与伪扩散问题研究与改进　　　　　　　　　　　　　　　　　　　　　　　　　　　［Ｊ］．　气体物理　　　 　２０２１ 　　６（２） 　　　　　　　Ｃｉｔａｔｉｏｎ 　　　　ｆｉｎｉｔｅＺｈｕ　Ｙ　ｅｌｅｍｅｎｔＦ 　Ｃｈｅｎ　ｍｅｔｈｏｄ［Ｊ］．　Ｘ　Ｃ．　Ｓｏｌｕｔｉｏｎｓ　Ｐｈｙｓｉｃｓ　　ｆｏｒ　ｏｆｎｕｍｅｒｉｃａｌＧａｓｅｓ 　　　２０２１ ｉｎｓｔａｂｉｌｉｔｙ　　６（２） 　ａｎｄ　１⁃８．　ｐｓｅｕｄｏ　　ｄｉｆｆｕｓｉｏｎ　　　ｐｒｏｂｌｅｍｓ　　　ｉｎｍｏｄｅｌｉｎｇ　　　ｏｆ　ａｔｍｏｓｐｈｅｒｉｃ　　　　ｐｒｅｓｓｕｒｅ　　　ｐｌａｓｍａ　　ｂａｓｅｄ　　　　ｏｎ　.com.cn. All Rights Reserved.２

气体物理２０２１年　第６卷

雪崩⁃流注放电 形成一条细长的放电通道 并在数十个纳秒内击穿整个气隙［９⁃１０］数值计算是测试理论假设．

、验证实验结论和预

测放电参数与行为的有力工具．雪崩⁃流注放电是一个多尺度系统 包含了诸多过程 包括化学反应、组分输运、电磁场传播以及流体耦合等 高精度、高效率地模拟这个系统是一个巨大的挑战．

首先 流注的动力学过程是高度非线性的．由于在流注头部以及与介质间隙 电场与空间电荷分离层直接耦合 对应区域的计算需要极高的空间精度．另外 电荷分离层通常具有弧度 因而纯一维１．５模型不足以准确描述流注的传播．文献［１１］使用了

程结合实现了径向方向流注的模拟维模型 通过将二维Ｐｏｉｓｓｏｎ方程和一维输运方 但是该方法无法模拟负流注和表面流注．文献［１２⁃１３］引入了ＳＧ方法模拟二维轴对称几何下的火花和辉光放电．文献ＳＧ［１４⁃１５］板放电方法献分析了正流注和负流注的动力学过程．文献［１８］［１６⁃１７］进一步对算法进行改良的改良 用于精确模拟针． 用于计⁃算流注与介质接触问题．其次 流注传播是一个多空间尺度问题．大气压下 流注头部的电荷分离层

厚度仅为几十个电子自由程大小（几十微米）．与此同时 流注的传播距离却可达数厘米．因此 流注模型须解析跨越了４个量级的两种空间尺度．最常用的解决方案是使用非均匀网格 并在等离子体区域加密．文献［１９⁃２０］等采用了最新的结构化自适应网格 使得计算效率得到了显著提高．流注也是一个多时间尺度问题．流注的传播时间一般为几十到数百纳秒 而大气压下捕捉流注输运过程所需的最小时间尺度为１０－１４献［２１⁃２３］使用了半隐式推进ｓ 跨越了 文献７～８［２４］个数量级使用异步．文时间推进的方法克服了这些困难．再次 流注具有强非局域特性 使得并行加速变得十分困难．这种非局域性主要是由于电场和光电离的全局特性．在有介质的情况下 边界特性的变化会导致局部电场的陡增和间断．为了克服非局域性问题 一种解决方案是使用局域能量近似 将化学反应速率和输运特性表达为局域平均电子能量的函数．另一种方法是引入对电离反应速率系数的修正项．更加直接的方法是采用混合模拟 考虑快速电子（逃逸电子）的形成．文献［２５］将粒子模型用于模拟流注头部 将流体模型用于模拟流注主体．这种“电子束⁃主体区”方法综合了电子ＭｏｎｔｅＣａｒｌｏ模型和流体模型 得到了大量研究［２６⁃２８］注分叉现象的研究工作中．在一些专注于逃逸电子 纯粒子⁃ＭｏｎｔｅＣａｒｌｏ、流模拟也得到了应用［２９⁃３２］ＣＯＭＳＯＬ近年来 很多课．

题组采用子体放电中有一些独特的优势Ｍｕｌｔｉｐｈｙｓｉｃｓ开展流注模拟商用 （１）．具有友善的用其在模拟等离有限元软件

户界面 便于快速入门并获得一些初始结果 （２）有限元方法使其可以适应几乎任何复杂几何结构（３）量 自带的等离子体模块能够大大减少编程工作 （４）自适应网格特性使其能够在一定程度上克服有限元商用软件的通病———网格数量太小 （５）自带的后处理系统为计算结果分析提供了有力的工具 （６）拥有通用性的偏微分方程模块 可以较为自由地组织模型．

尽管ＣＯＭＳＯＬＭｕｌｔｉｐｈｙｓｉｃｓ具有以上优势 但是也存在很多非常显著且难以解决的缺陷 导致其无法作为研究在大气压下通常呈流注形态的纳秒脉冲放电的理想工具 （１）ＣＯＭＳＯＬＭｕｌｔｉｐｈｙｓｉｃｓ的等离子体模块采用了标准Ｇａｌｅｒｋｉｎ方法 该方法无法捕捉输运占优的流注头部区域的强间断面 导致计算结果严重失真 呈伪扩散形体 （２）商用软件的改“黑盒子 因而在模拟遇到问题时只能凭借猜测处理”模式使得无法对程序运行模式进行任何修

几乎没有办法调试 （３）全隐式的数值格式使得时间积分对计算域任何一点扰动都十分敏感 无关紧要处的一个奇异点就会导致整个模型计算效率下降 收敛极其困难 （４）与基于质量守恒的有限体积法不同 有限元方法本质是基于能量守恒原理的 不具有保证质量守恒的特性 在等离子体计算中不可避免会出现负浓度情况 （５）ＣＯＭＳＯＬＭｕｌｔｉｐｈｙｓｉｃｓ直接求解大规模矩阵而不充分利用并行计算机结构 使其存在计算速度缓慢、内存消耗巨大的问题．

ＣＯＭＳＯＬ本文介绍一种以弱形式偏微分方程形式在以克服等离子体模拟中出现的数值振荡和伪扩散Ｍｕｌｔｉｐｈｙｓｉｃｓ软件中引入人工稳定项 用问题的解决方案．文章分３个部分 第１部分以弱形式的方式推导人工稳定项的输入形式 第２部分给出基准验证和实验验证 最后给出结论．

１　弱形式与人工稳定项

分析引言所述的ＣＯＭＳＯＬＭｕｌｔｉｐｈｙｓｉｃｓ缺陷 发现缺陷（３）和（４）是当等离子体与介质、金属接触时存在的最主要问题 而对于针⁃板型的体放电

.com.cn. All Rights Reserved.第２期朱益飞 等 标准有限元方法模拟大气压流注等离子体中数值振荡与伪扩散问题研究与改进

３

而言 这些缺陷不会对计算结果产生太大影响．缺陷（５）无法克服 而缺陷（１）和（２）可以通过引入弱形式偏微分方程 添加人工稳定项来解决．本节将对在ＣＯＭＳＯＬ中构建弱形式ＰＤＥ引入流线扩散人工稳定性ＳＵＰＧ）（ｓｔｒｅａｍｌｉｎｅｕｐｗｉｎｄＰｅｔｒｏｖＧａｌｅｒｋｉｎｍｅｔｈｏｄ 首先进行阐述 利用．

ＣＯＭＳＯＬＭｕｌｔｉｐｈｙｓｉｃｓ建立一套最

简单的等离子体方程体系 Ｐｏｉｓｓｏｎ方程、Ｈｅｌｍｈｏｌｔｚ方程和输运方程．其中 Ｐｏｉｓｓｏｎ方程和Ｈｅｌｍｈｏｌｔｚ方程是单纯的椭圆形方程 可以直接使用偏微分方程模块以“黑盒子”的形式求解．输运方程则用弱形式偏微分方程表达．

在纳秒放电时间尺度下 可以假设离子和中性粒子不运动 仅求解不包含输运项的反应方程即可

∂∂ｎｔ

ｉ

－Δ（ＤｉΔｎｉ）＝Ｓｉ＋Ｓｐｈ

式中 ｎｉ为离子密度 Ｄｉ为离子扩散率 Ｓ Ｓｉ为化学反应主导的离子产生率．

ｐｈ为光电离主导的离子产生率电子的输运方程为∂∂ｎｔ

ｅ

－Δ（ＤｅΔｎｅ＋μｅｎｅΔΦ）＝Ｓｅ＋Ｓｐｈ

式中 ｎ迁移率 ｅ为电子密度 ＤΦ为电势 Ｓｅ为电子扩散率 μｅ为电子率 Ｓｅ为化学反应主导的电子产生迁移率和化学反应电子产生率均可通过开源电子ｐｈ为光电离主导的电子产生率．电子扩散率、

能量分布函数计算软件ＢＯＬＳＩＧ＋计算得到．为了便于比对 本文中的反应速率系数和输运系数均采用了文献［１６］中给定的公式计算

０＝Δ（ＤｅΔｎｅ＋μｅｎｅΔΦ）－∂∂ｎｔ

ｅ＋Ｓｉ＋Ｓｐｈ

以上方程乘以“试函数”ｕ（ｘ ｙ ｚ） 并对有限单元

进行体积分 可以得到

∭Ｖ

（

Δ·（ＤｅΔｎｅ））ｕｄＶ＋

∭ΔΔΦ））ｕｄＶ－

Ｖ

（

·（ｎｅμｅ　

è􀆟ｔøｄＶ＋＝０

对上式前两项进行分部积分∭

æç􀆟ｎｅö÷ｕＶ

∭Ｖ

（Ｓ

ｉ

＋Ｓｐｈ）ｕｄＶ 并简化可得

∯Ｓ

（Ｄｅ

Δｎｅｕ）·ｎｄＳ－

∭Ｄ

Ｖ

ｅ

Δｎｅ·ΔｕｄＶ＋

∭Ｖ

（Δｎｅμｅ·（ΔΦ））ｕｄＶ＋∯Ｖ

ｅ

μｅ·Φ）Ｓ

ｅｅ

Δ）·ｎｄＳ－

∭（ｎΔΔｕｄＶ＋

∭

∭（Ｖ

ΔｎμΦｕ（ｎｅμｅ·（ΔΦ））ｕｄＶ－ΔΦ·Δｕ）ｄＶ＋

Ｖ

（ｎｅμｅ由此 可以得到能够直接输入∭æ－∭Ｖ

（ｎｅ

Δμｅ·ΔΦ）ｕｄＶ－

ç

􀆟ｎｅè􀆟ｔ＋ＳＶ

ｉ＋Ｓｐｈö

÷ø

ｕｄＶ＝０ＣＯＭＳＯＬＭｕｌｔｉｐｈｙｓｉｃｓ的弱形式方程－∭

∭Ｄ

Ｖ

ｅ

Δｎｅ·ΔｕｄＶ－

　æ

∭Ｖ

（ｎ

ｅ

μｅΔΦ·Δｕ）ｄＶ＋

ç

ｎＶ

è－􀆟　

􀆟ｎｔｅ＋Ｓｉ＋Ｓｐｈö÷ø

ｕｄＶ＋∯∯Ｓ

（Ｄ

ｅ

Δｅｕ）·ｎｄＳ＋

Ｓ

（ｎ

ｅ

μｅΔΦｕ）·ｎｄＳ＝０

（１）

式中 第１行是等离子体计算域的弱形式偏微分方程 第２行代表自然满足的边界条件．式（１）可以写成更明晰的形式

∭

（（－Ｄｅｎｅ－ｎｅμｅΦ）·ｕ＋（－

＋ＳＶ

Δ　Ｓｐｈ）ｕ）ｄＶ＋

∯ΔΔ􀆟ｉ＋

Ｓ

（Ｄ

ｅ

Δｎｅ＋ｎｅμｅΔΦ）ｕ􀆟ｎ·ｔ

ｅｎｄＳ＝０

方程（１）所描述的漂移⁃扩散过程 如果求解的区域被扩散过程主导Ｇａｌｅｒｋｉｎ 那么传统的ＣＯＭＳＯＬ标准算区域是漂移主导差分方法可以获得稳定可靠的解 那么方程（１）则表现出双曲型．如果计方程特性 传统的数值方法变得非常不稳定 计算结果要么会在粗糙的网格下剧烈振荡 要么在精细的网格下形成伪扩散解．在前述方程的基础上 利用ＣＯＭＳＯＬ直接复现文献［１６］中的轴对称针⁃板放电基准算例 采用相同的方程（连续方程＋Ｐｏｉｓｓｏｎ方程） 分别使用５０μｍ和１５μｍ两种尺寸网格开展计算 结果如图１（ａ）～（ｄ）所示．

图１（ａ）和（ｃ）对比了使用ＣＯＭＳＯＬＭｕｌｔｉｐｈｙｓｉｃｓ在两种网格条件下计算的流注轴向电子密度分布（线蓝线）．在粗网格下）和基准算例中的轴向流注通道电子密度 ＣＯＭＳＯＬＭｕｌｔｉｐｈｙｓｉｃｓ计算所得流（红注能够复现基准算例流注传播速度 如图１（ｂ）所示 但是整个放电通道内出现了图１（ａ）中剧烈的非物理数值振荡 计算结果无法分析理解．如果仅仅加密网格 如图１（ｃ）所示 计算轴向电子密度不再产生数值振荡 但是流注传播距离远远落后于基准算

.com.cn. All Rights Reserved.４

气体物理２０２１年　第６卷

例 图１（ｄ）表明 ＣＯＭＳＯＬＭｕｌｔｉｐｈｙｓｉｃｓ标准有限元不可能正确计算Ｍｕｌｔｉｐｈｙｓｉｃｓ内置的等离子体模块进行测试大气压流注．本文使用方法无法捕捉到基准算例中的电子雪崩⁃流注转捩 计算结果呈扩散流注状．这两个案例清楚地表明 使 ＣＯＭＳＯＬ了同样的错误结果．

也得到

用ＣＯＭＳＯＬＭｕｌｔｉｐｈｙｓｉｃｓ默认的数值格式和计算方法

（ａ）Ｅｌｅｃｔｒｏｎｄｅｎｓｉｔｙｉｎｃｏａｒｓｅｍｅｓｈｅｓ

（ｂ）Ｅ／Ｎｉｎｃｏａｒｓｅｍｅｓｈｅｓ

（ｃ）Ｅｌｅｃｔｒｏｎｄｅｎｓｉｔｙｉｎｆｉｎｅｍｅｓｈｅｓ

（ｄ）Ｅ／Ｎｉｎｆｉｎｅｍｅｓｈｅｓ

图１　使用ＣＯＭＳＯＬＴａｂｌｅ１　Ｒｅｓｕｌｔｓ标准Ｇａｌｅｒｋｉｎｃｏｍｐａｒｉｓｏｎ格式计算针板放电并与ｂｅｔｗｅｅｎｓｔａｎｄａｒｄＣＯＭＳＯＬＫｕｌｉｋｏｖｓｋｙｍｏｄｅｌ流注基准案例ａｎｄＲｅｆ［１６］

［１６］对比结果

　　这种数值不稳定性 源于计算网格与算法不匹严苛的网格条件在ＣＯＭＳＯＬ下基本不可能实现（除配 有限元中的标准Ｇａｌｅｒｋｉｎ方法与有限体积／有非计算域非常小）．为此 可以通过引入人工稳定限差分中的标准ＳＧ方法一样 在高气压条件下 性［３３⁃３７］提升计算精度．最广泛使用的人工稳定性方计算等离子体流注的准确性会显著下降．假设电子法之一 就是流线扩散稳定ＳＵＰＧ方法［３３⁃３４］扩散系数符合Ｅｉｎｓｔｅｉｎ关系 那么为了保证计算结弱形式方程 在ＣＯＭＳＯＬＭｕｌｔｉｐｈｙｓｉｃｓ中对各个有．通过果的准确性 须保证以下关系

［１７］

Δ限单元添加人工稳定项

２ＥＴｋｈｋｅｋ

≪１（２）

Ｐ（ｕｈ）τｅＲ式中 ｋ表示第ｋ个网格节点 ΔＥ为电场强度梯∑ｎｅｌｅ＝

１

式中 各项的含义如下

∫Ωｅｌ

（ｎｅｈ）ｄΩ（３）

度 ｈ为网格尺寸 Ｔｅ为电子温度．上式表明 两个相邻节点的电势差须远小于电子温度．对于高气压Ｐ（ｕｈ∂情况 若要满足该条件 则需要极高的网格密度 ）＝μｅ∂Φｘ∂∂ｕｘ＋μｅ

∂∂Φｙ∂∂ｕｙ

这使得高性能计算几乎不可能．方程（２）所规定的

τｅ＝

１２æ

çｈｈöèαｘμｅＥ

＋αｙμ÷

ｅＥø.com.cn. All Rights Reserved.第２期朱益飞 等 标准有限元方法模拟大气压流注等离子体中数值振荡与伪扩散问题研究与改进

５

Ｒ（ｎｅｈ）＝

∂Γ∂ｘ

ｅｘ＋

∂Γ∂ｙ

ｅｙ＋

∂∂ｎｔ

ｅｈ－Ｓｉ－Ｓｐｈ

式中 ｎ有限元素中的电子密度ｅｌ是有限元素的总数 ｎ ｕｅｈ是特征尺寸为ｈ的ΓΔｈ是对应的试函数Δ．Γ在ｘ和ｅｙｘ和ｅｙ是电子通量Γ＝（－Ｄｅｎｅ－ｎｅμｅΦ）方向上的分量．αｘ与αｙ分别是ｘ和ｙ方向上与Ｐｅｃｌｅｔ数Ｐｅ相关的参数

αｘ＝ｃｏｔｈ（Ｐｅｘ）－æ１ö

μç

÷ Ｐｅｅ

∂èＰｅｘøｘ＝１２Ｄ∂Φ

ｈｅ

ｘ　

αｙ＝ｃｏｔｈ（Ｐｅｙ）－æ

çèＰｅ１ö÷ Ｐｅｙø

ｙ＝１μ２ｅ∂Ｄ∂Φｅｙｈ

将方程（３）以弱形式方程的方式代入ＣＯＭＳＯＬ中电子密度输运方程弱形式 即可得到稳定的解．针对该改进的ＣＯＭＳＯＬ模型 具体的验证工作将在下节详细介绍．

２　基准验证与实验验证

针对改进的ＣＯＭＳＯＬＭｕｌｔｉｐｈｙｓｉｃｓ流注放电模型 基于经典计算结果［１６ ３８］和实验结果［３９］开展验证．２．１　计算基准验证

针板放电是一种经典和基础的等离子体放电结构 国外前期开展了大量研究 主要采用有限体积方法编程求解连续方程、Ｐｏｉｓｓｏｎ方程与光电离．其中 俄罗斯科学院Ｋｕｌｉｋｏｖｓｋｙ等对针板放电开展了数值计算研究 并构建了一系列便于验证的数值模拟基准算例．文献［１６］计算了１ｃｍ间隙条件下 大气压针板放电构型在１３ｋＶ恒定电压下的击穿过程 考虑了电子碰撞电离反应、三体和双体吸附反应、电子⁃离子复合和离子⁃离子复合反应 成功复现了雪［１６］崩⁃流注转捩现象和流注传播特性复现 的参数被国际上广泛认可为大气压低温等离子体数、方程完整详实 得到了众多课题组的．文献值计算的重要基准算例之一．

本文第１部分对比了使用标准有限元方法计算的结果与文献［１６］的结果 在图１中展示了数值振荡和伪扩散问题．本节使用修正后的ＣＯＭＳＯＬＭｕｌ⁃运参数和反应体系再次开展计算ｔｉｐｈｙｓｉｃｓ模型 采用与文献［１６］相同的控制方程 结果如图２所示、输

．　同时刻轴向电子密度分布和电子密度空间分布云　图２（ａ）和（ｃ）分别为文献［１６］计算得到的不图．图２（ｂ）和（ｄ）为采用修正ＣＯＭＳＯＬＭｕｌｔｉｐｈｙｓｉｃｓ模型计算结果．由图可见 添加人工稳定性的修正模型计算结果与基准算例几乎没有区别 电场强度、电子密度与传播过程均实现了准确的复现．注意到使用修正ＣＯＭＳＯＬＭｕｌｔｉｐｈｙｓｉｃｓ模型计算的流注直径比基准算例稍细 引起这一区别的原因可能是光电离计算方式的不同 文献［１６］采用了简化的“积分环”方式 仅计算流注头部强电场附近的光电离源项 并将全局积分简化为对流注头部附近的环形体积进行积分．本文采用了近年来使用较为广泛的Ｈｅｌｍｈｏｌｔｚ三方程光电离模型［４０］准确性已经得到了广泛验证［４１］ 本文不再赘述 该模型的１．图

（１９中基准算例的红线与图型计算结果中已经没有了图ｎｓ时刻）对应．对比图２（ａ）（ｂ）１ １２蓝线出现的数值振可见中左数第 修正后的模３条线荡和伪扩散的问题．

（ａ）Ａｘｉａｌｅｌｅｃｔｒｏｎｄｅｓｎｉｔｙｆｒｏｍｂｅｎｃｈｍａｒｋ

（ｂ）ＡｘｉａｌｅｌｅｃｔｒｏｎｄｅｓｎｉｔｙｆｒｏｍｃｏｒｒｅｃｔｅｄＣＯＭＳＯＬｍｏｄｅｌ

.com.cn. All Rights Reserved.６

气体物理２０２１年　第６卷

（ｃ）Ｐｌｏｔｏｆｅｌｅｃｔｒｏｎｄｅｎｓｉｔｙｃｏｎｔｏｕｒｓｆｒｏｍｂｅｎｃｈｍａｒｋ

（ｄ）ＥｌｅｃｔｒｏｎｄｅｎｓｉｔｙｃｏｎｔｏｕｒｓｆｒｏｍｃｏｒｒｅｃｔｅｄＣＯＭＳＯＬｍｏｄｅｌ

图２　修正的ＣＯＭＳＯＬ流注模型复现Ｋｕｌｉｋｏｖｓｋｙ

Ｆｉｇ．２　ＣＯＭＳＯＬＲｅｓｕｌｔｓ基准流注算例ｃｏｍｐａｒｉｓｏｎ［１６］

ｍｏｄｅｌａｎｄｂｅｔｗｅｅｎＲｅｆ［１６］

ｃｏｒｒｅｃｔｅｄ

需要注意的是 本节进行的计算基准验证的核心目的是数值计算角度确保 在相同的输运参数、几何结构和边界条件下 求解相同的控制方程能够得到完全相同的解 从而确保数值计算研究的一致性．文献［１６］的结果作为数值计算的良好基准 却并未得到实验的充分验证．在克服标准有限元方法引发的数值振荡、伪扩散问题的基础上 使用修正的模型与实验结果进行了进一步比对验证．２．２　实验结果验证

利用修正后的ＣＯＭＳＯＬＭｕｌｔｉｐｈｙｓｉｃｓ等离子体

模型［３９］ 利用相增强高速相机针对一例过电压实验进行直接比对［３９］ＩＣＣＤ研究了大气压下．文献１．６

用下的光学特性ｃｍ间隙针板结构在纳秒脉冲过电压．在实验（见图３（ａ））和后续的计（８６ｋＶ）作算中［４２］均发现 如果电压上升时间明显短于电场屏蔽效应的特征时间尺度 局域场就会高于空气电离阈值 引起独特的倒锥形放电结构．

（ａ）Ｍｅａｓｕｒｅｄｄｉｓｃｈａｒｇｅｍｏｒｐｈｏｌｏｇｙ

（ｂ）Ｃａｌｃｕｌａｔｅｄ（ａｆｔｅｒｄｉｓｔｒｉｂｕｔｉｏｎｄｉｒｅｃｔＡｂｅｌｏｆｔｒａｎｓｆｏｒｍａｔｉｏｎ）ｓｅｃｏｎｄａｒｙｐｏｓｉｔｉｖｅｓｙｓｔｅｍ

图３　Ｆｉｇ．修正的３　ＯｖｅｒｖｏｌｔａｇｅＣＯＭＳＯＬｃｏｒｒｅｃｔｅｄｄｉｓｃｈａｒｇｅ流注模型复现过电压放电实验ＣＯＭＳＯＬｅｘｐｅｒｉｍｅｎｔｍｏｄｅｌ

ｂａｓｅｄｏｎ

在与实验完全相同的条件下 本文试图利用标准有限元模型计算该倒锥形放电结构 模型因为未知（黑盒子）原因不收敛 未能完成计算．本文仅展示了利用修正后的ＣＯＭＳＯＬＭｕｌｔｉｐｈｙｓｉｃｓ等离子体模型的计算结果．

在实验中 光学辐射的主要来源是氮气分子的第二正带系Ｎ分密度 为了与实验测量结果直接对照２（Ｃ３Πｕ） 为此利用模型计算了该组 对计算结果进行了Ａｂｅｌ变换以考虑放电区域的厚度．将变换后的计算结果在０～１之间进行归一化 如图３（ｂ）所示．实验和计算均捕捉到了阳极和阴极区域的局域强辐射１．２ｃｍ 实验结果与计算结果对比吻合良好 以及在大约Ｚ＝１ｃｍ处的最大放电直径

．

３　结论

针对商用有限元软件模拟流注放电时存在的

.com.cn. All Rights Reserved.第２期朱益飞 等 标准有限元方法模拟大气压流注等离子体中数值振荡与伪扩散问题研究与改进

７

数值问题 推导了基于流线扩散人工稳定技术的修正项 并在ＣＯＭＳＯＬＭｕｌｔｉｐｈｙｓｉｃｓ软件中进行了测试．计算发现 采用标准有限元方法计算的针板流注结果存在显著的数值振荡和伪扩散流注问题 而引入修正后的模型克服了这些数值问题 能够成功复现１ｃｍ气隙、１３ｋＶ电压条件下大气压空气流注基准算例．进一步计算了１．６ｃｍ气隙、８６ｋＶ过电压大气压空气放电形态 并与实验结果进行了对比验证 证明了模型本身的可靠性．计算与验证结果表明 以弱形式偏微分方程的方式将人工稳定项引入基于有限元（标准Ｇａｌｅｒｋｉｎ方法）的商用软件中 能够有效克服放电计算中的数值振荡和流注伪扩散问题 确保等离子体计算准确性和可重复性．

致谢５１７９０５１１ 　本文受到国家自然科学基金９１９４１１０５ ９１９４１３０１）支持．

（５１９０７２０４ 参考文献（Ｒｅｆｅｒｅｎｃｅｓ）

［１］　ＥｌｉａｓｓｏｎｐｌａｓｍａｃｈｅｍｉｃａｌＢ

Ｋｏｇｅｌｓｃｈａｔｚｐｒｏｃｅｓｓｉｎｇ［Ｊ］．Ｕ．ＮｏｎｅｑｕｉｌｉｂｒｉｕｍＩＥＥＥＴｒａｎｓａｃｔｉｏｎｓｖｏｌｕｍｅ

ＰｌａｓｍａＳｃｉｅｎｃｅ １９９１ １０６３⁃１０７７．

ｏｎ［２］ＢｉｅｄｅｒｍａｎＨ．Ｐｌａｓｍａ１９（６）

ｐｏｌｙｍｅｒｆｉｌｍｓ［Ｍ］．Ｓｉｎｇａｐｏｒｅ ［３］ＷｏｒｌｄＧｕｓｅｖＳｃｉｅｎｔｉｆｉｃ ２００４．

［４］

ＮＡＳＡＶＫｉｍＳＴＩＥ ／ＲｅｃｏｎＫａｒａｂｕｔｏｖＴｅｃｈ．ＡＲｅｐ．Ａ．Ａ Ｌａｓｅｒｏｐｔｏａｃｏｕｓｔｉｃｓ［Ｒ］．ｃｏｎｔｒｏｌ ＨＨ．ａＮｏｎｔｈｅｒｍａｌｈｉｓｔｏｒｉｃａｌｒｅｖｉｅｗ ｐｌａｓｍａｃｕｒｒｅｎｔｐｒｏｃｅｓｓｉｎｇ１９９１ ９３．

ｉｓｓｕｅｓ ｆｏｒａｉｒ⁃ｐｏｌｌｕｔｉｏｎａｎｄｆｕｔｕｒｅ（２） ｐｒｏｓｐｅｃｔｓ［Ｊ］．９１⁃１１０．

ＰｌａｓｍａＰｒｏｃｅｓｓｅｓａｎｄＰｏｌｙｍｅｒｓ ２００４ １

［５］ＰｉｌｌａｔｕｒｂｕｌｅｎｔＧ Ｇａｌｌｅｙｐｌａｓｍａ［Ｊ］．ｐｒｅｍｉｘｅｄＤ ｆｌａｍｅＬａｃｏｓｔｅＩＥＥＥｕｓｉｎｇＤＡ Ｔｒａｎｓａｃｔｉｏｎｓａｎａｎｏｓｅｃｏｎｄｅｔａｌ．ＳｔａｂｉｌｉｚａｔｉｏｎｏｎＰｌａｓｍａｒｅｐｅｔｉｔｉｖｅｌｙｏｆａｅｎｃｅ ｐｕｌｓｅｄ２００６ ３４（６） ２４７１⁃２４７７．

Ｓｃｉ⁃

［６］ＭａｓｓｉｎｅｓｃｈｅｍｉｓｔｒｙＦ ｉｎａＳéｇｕｒｇｌｏｗＰ ｄｉｅｌｅｃｔｒｉｃＧｈｅｒａｒｄｉａｎｄｂａｒｒｉｅｒＮ ｍｏｄｅｌｌｉｎｇ［Ｊ］．ｄｉｓｃｈａｒｇｅｅｔａｌ．ＰｈｙｓｉｃｓａｔＳｕｒｆａｃｅ

ａｔｍｏｓ⁃ａｎｄ［７］ａｎｄｐｈｅｒｉｃＦｒｉｄｍａｎＣｏａｔｉｎｇｓｐｒｅｓｓｕｒｅ Ｔｅｃｈｎｏｌｏｇｙ ｄｉａｇｎｏｓｔｉｃｓ２００３ １７４⁃１７５ ８⁃１４．ｐｈｅｒｉｃ［８］

ｐｌｉｅｄｐｒｅｓｓｕｒｅＡ Ｃｈｉｒｏｋｏｖｄｉｓｃｈａｒｇｅｓ［Ｊ］．Ａ ＧｕｔｓｏｌＪｏｕｒｎａｌＡ．Ｎｏｎ⁃ｔｈｅｒｍａｌｏｆａｔｍｏｓ⁃ＳｔａｒｉｋｏｖｓｋｉｉＰｈｙｓｉｃｓ ＳＤＢＤｐｌａｓｍａＡＹ ２００５ ａｃｔｕａｔｏｒＮｉｋｉｐｅｌｏｖ３８（２） ｗｉｔｈＡＲ１⁃Ｒ２４．

ＰｈｙｓｉｃｓＤ Ａｐ⁃ｎａｎｏｓｅｃｏｎｄＡ Ｎｕｄｎｏｖａｐｕｌｓｅ⁃ｐｅｒｉｏｄｉｃＭＭ ｅｔａｌ．ｄｉｓｃｈａｒｇｅ［Ｊ］．［９］２００９ 邵涛１８（３） ＰｌａｓｍａＳｏｕｒｃｅｓＳｃｉｅｎｃｅａｎｄＴｅｃｈｎｏｌｏｇｙ

离子体应用 章程 ［Ｊ］．王瑞雪０３４０１５．

高电压技术 等．大气压脉冲气体放电与等 ２０１６ ４２（３） ６８５⁃７０５．

ＳｈａｏｐｒｅｓｓｕｒｅＴ ＺｈａｎｇＣ ＷａｎｇＲＸ ｅｔａｌ．Ａｔｍｏｓｐｈｅｒｉｃ⁃ｔｉｏｎ［Ｊ］．ｐｕｌｓｅｄ７０５（ｉｎＣｈｉｅｓｅ）．

ＨｉｇｈＶｏｌｔａｇｅｇａｓｄｉｓｃｈａｒｇｅＥｎｇｉｎｅｅｒｉｎｇ ａｎｄｐｕｌｓｅｄ２０１６ ｐｌａｓｍａ４２（３） ａｐｐｌｉｃａ⁃６８５⁃［１０］李和平研究进 展于达仁综述［Ｊ］． 孙文廷高电 压等．技大气压放电等离子体

术 ２０１６ ４２（１２） ３６９７⁃３７２７．

ＬｉｍｏｓｐｈｅｒｉｃＨＰ Ｙｕｎｅｅｒｉｎｇ ２０１６ ｄｉｓｃｈａｒｇｅＤＲ Ｓｕｎ４２（１２） ｐｌａｓｍａｓＷＴ ｅｔ３６９７⁃３７２７（ｉｎ［Ｊ］．ａｌ．Ｓｔａｔｅ⁃ｏｆ⁃ｔｈｅ⁃ａｒｔＨｉｇｈＣｈｉｅｓｅ）．ＶｏｌｔａｇｅｏｆＥｎｇｉ⁃ａｔ⁃［１１］ＭｏｒｒｏｗＲ ＬｏｗｋｅＪＪ．Ｓｔｒｅａｍｅｒｐｒｏｐａｇａｔｉｏｎｉｎａｉｒ［Ｊ］．［１２］Ｂｏｅｕｆ６１４⁃６２７．

ＪｏｕｒｎａｌｏｆＰｈｙｓｉｃｓＤ ＡｐｐｌｉｅｄＰｈｙｓｉｃｓ １９９７ ３０（４） ｃｏｍｐｕｔｅｒＪＰ Ｐｉｔｃｈｆｏｒｄ［１３］Ｓｃｉｅｎｃｅ ＦｉａｌａＡ １９９１ ｓｉｍｕｌａｔｉｏｎ［Ｊ］．ＬＣ．Ｐｉｔｃｈｆｏｒｄ１９（２） ＬＣ Ｂｏｅｕｆ２８６⁃２９６．

ＩＥＥＥＰｓｅｕｄｏｓｐａｒｋＴｒａｎｓａｃｔｉｏｎｓｄｉｓｃｈａｒｇｅｓｏｎＰｌａｓｍａｖｉａ

ｂｒｉｄＪＰ．Ｔｗｏ⁃ｄｉｍｅｎｓｉｏｎａｌ ｈｙ⁃

［１４］ＲｅｖｉｅｗｍｏｄｅｌＢａｂａｅｖａＥ ｏｆＮ１９９４ ｌｏｗ⁃ｐｒｅｓｓｕｒｅＹ Ｎａｉｄｉｓ４９（６） ＧＶ．５６０７⁃５６２２．

ｇｌｏｗｄｉｓｃｈａｒｇｅｓ［Ｊ］．ＰｈｙｓｉｃａｌｐｏｓｉｔｉｖｅＴｗｏ⁃ｄｉｍｅｎｓｉｏｎａｌｍｏｄｅｌｌｉｎｇｏｆ

ｉｎａｉｒ［Ｊ］．ｓｔｒｅａｍｅｒＪｏｕｒｎａｌｄｙｎａｍｉｃｓｏｆＰｈｙｓｉｃｓｉｎＤ ｎｏｎ⁃ｕｎｉｆｏｒｍＡｐｐｌｉｅｄＰｈｙｓｉｃｓ ｅｌｅｃｔｒｉｃ１９９６ ｆｉｅｌｄｓ［１５］２９（９） Ｂａｂａｅｖａ２４２３⁃２４３１．

ａｔｉｖｅＮＹ ＮａｉｄｉｓＩＥＥＥｓｔｒｅａｍｅｒｓＴｒａｎｓａｃｔｉｏｎｓｉｎａｉｒＧＶ．Ｄｙｎａｍｉｃｓｏｆｐｏｓｉｔｉｖｅａｎｄｎｅｇ⁃

ｏｎｉｎＰｌａｓｍａｗｅａｋｕｎｉｆｏｒｍＳｃｉｅｎｃｅ ｅｌｅｃｔｒｉｃ１９９７ ｆｉｅｌｄｓ［Ｊ］．２５（２） ［１６］３７５⁃３７９．

Ｋｕｌｉｋｏｖｓｋｙｍｏｖｉｎｇａｖａｌａｎｃｈｅ⁃ｔｏ⁃ｓｔｒｅａｍｅｒＡＡ．Ｐｏｓｉｔｉｖｅｓｔｒｅａｍｅｒｔｒａｎｓｉｔｉｏｎ［Ｊ］．ｉｎａｗｅａｋｆｉｅｌｄＰｈｙｓｉｃａｌｉｎａｉｒ

［１７］ＲｅｖｉｅｗＡＥ １９９８ ５７（６） ７０６６⁃７０７４．

ｇｏｒｉｔｈｍＫｕｌｉｋｏｖｓｋｙｏｆｅｌｅｃｔｒｏｎＡＡ．ＡｍｏｒｅａｃｃｕｒａｔｅＳｃｈａｒｆｅｔｔｅｒ⁃Ｇｕｍｍｅｌａｌ⁃

ｃｈａｒｇｅｓｉｍｕｌａｔｉｏｎ［Ｊ］．ｔｒａｎｓｐｏｒｔＪｏｕｒｎａｌｆｏｒｓｅｍｉｃｏｎｄｕｃｔｏｒｏｆＣｏｍｐｕｔａｔｉｏｎａｌａｎｄＰｈｙｓｉｃｓ

ｇａｓｄｉｓ⁃［１８］１９９５ Ｐｅｃｈｅｒｅａｕ１１９（１） １４９⁃１５５．

ｒｅｉｇｎｉｔｉｏｎａｔｍｏｓｐｈｅｒｉｃｏｆＦ ａＪáｎｓｋｙＪ ＢｏｕｒｄｏｎＡ．Ｓｉｍｕｌａｔｉｏｎｏｆｔｈｅ

ｐｒｅｓｓｕｒｅ［Ｊ］．ｄｉｓｃｈａｒｇｅｂｅｈｉｎｄＰｌａｓｍａａｄｉｅｌｅｃｔｒｉｃｌａｙｅｒｉｎａｉｒａｔＴｅｃｈｎｏｌｏｇｙ ２０１２ ２１（５） Ｐ Ｃａｐｅｉｌｌèｒｅ０５５０１１．

ＳｏｕｒｃｅｓＳｃｉｅｎｃｅａｎｄ［１９］ＰａｎｃｈｅｓｈｎｙｉＳ ＳéｇｕｒＪ ｅｔｐａｒａｌｌｅｌａｌ．Ｎｕｍｅｒｉｃａｌ

ａｄａｐｔｉｖｅ

［２０］２００８ ｍｅｓｈｓｉｍｕｌａｔｉｏｎｒｅｆｉｎｅｍｅｎｔ［Ｊ］．ｏｆｆｉｌａｍｅｎｔａｒｙＵｎｆｅｒ２２７（１３） ６５７４⁃６５９０．ＪｏｕｒｎａｌｄｉｓｃｈａｒｇｅｓｏｆＣｏｍｐｕｔａｔｉｏｎａｌｗｉｔｈＰｈｙｓｉｃｓ ｒｅｆｉｎｅｍｅｎｔ［Ｊ］．

ｃｈａｒｇｅＴ ｓｉｍｕｌａｔｉｏｎｓＢｏｅｕｆＪＰ Ｃｏｍｐｕｔｅｒ

ｕｓｉｎｇＲｏｇｉｅｒａｓｙｎｃｈｒｏｎｏｕｓＦ ｅｔａｌ．Ｍｕｌｔｉ⁃ｓｃａｌｅａｄａｐｔｉｖｅｇａｓｍｅｓｈ

ｄｉｓ⁃

［２１］２０１０ Ｖｅｎｔｚｅｋ１８１（２） ２４７⁃２５８．

Ｐｈｙｓｉｃｓ

Ｃｏｍｍｕｎｉｃａｔｉｏｎｓ

ｓｉｏｎａｌｍｏｄｅｌｉｎｇＰＬ ＨｏｅｋｓｔｒａｏｆｈｉｇｈＲＪ ｐｌａｓｍａＫｕｓｈｎｅｒｄｅｎｓｉｔｙＭＪ．Ｔｗｏ⁃ｄｉｍｅｎ⁃

ｉｎｄｕｃｔｉｖｅｌｙ

.com.cn. All Rights Reserved.８

气体物理

２０２１年　第６卷ｃｏｕｐｌｅｄＶａｃｕｕｍｓｏｕｒｃｅｓＮａｎｏｍｅｔｅｒＳｃｉｅｎｃｅｆｏｒＳｔｒｕｃｔｕｒｅｓ＆Ｔｅｃｈｎｏｌｏｇｙｍａｔｅｒｉａｌｓｐｒｏｃｅｓｓｉｎｇ［Ｊ］．Ｐｒｏｃｅｓｓｉｎｇ Ｂ Ｍｅａｓｕｒｅｍｅｎｔ ＭｉｃｒｏｅｌｅｃｔｒｏｎｉｃｓＪｏｕｒｎａｌａｎｄａｎｄｏｆＰｈ⁃［２２］ｅｎｏｍｅｎａ ＨａｇｅｌａａｒｇａｓｄｉｓｃｈａｒｇｅｓＧ１９９４ Ｊ Ｋｒｏｅｓｅｎ１２（１） ｂｙｉｍｐｌｉｃｉｔＧＭ．４６１⁃４７７．

ｔｒｅａｔｍｅｎｔＳｐｅｅｄｉｎｇｏｆｕｐｔｈｅｆｌｕｉｄｅｌｅｃｔｒｏｎｍｏｄｅｌｓｅｎｅｒｇｙ

ｆｏｒ

［２３］２０００ ｓｏｕｒｃｅＬｉｎ１５９（１） ｔｅｒｍ［Ｊ］．１⁃１２．

ＪｏｕｒｎａｌｏｆＣｏｍｐｕｔａｔｉｏｎａｌＰｈｙｓｉｃｓ

ａＫＭ ＨｕｎｇＣＴ ＨｗａｎｇＦＮ ｅｔａｌ．ｍｏｄｅｌｉｎｇｐａｒａｌｌｅｌｃｏｄｅｓｅｍｉ⁃ｉｍｐｌｉｃｉｔｕｓｉｎｇｆｉｎｉｔｅ⁃ｖｏｌｕｍｅｔｗｏ⁃ｄｉｍｅｎｓｉｏｎａｌＤｅｖｅｌｏｐｍｅｎｔｏｆ

［２４］ＰｈｙｓｉｃｓＵｎｆｅｒＴ Ｃｏｍｍｕｎｉｃａｔｉｏｎｓ ＢｏｅｕｆＪＰ Ｒｏｇｉｅｒ２０１２ Ｆ １８３（６） ｍｅｔｈｏｄ［Ｊ］．ｐｌａｓｍａｅｔａｌ．Ａｎ１２２５⁃１２３６．

Ｃｏｍｐｕｔｅｒｆｌｕｉｄｓｃｈｅｍｅｗｉｔｈｌｏｃａｌａｓｙｎｃｈｒｏｎｏｕｓ

ｐｒｏｂｌｅｍｓ ［２５］ＣｏｍｐｕｔａｔｉｏｎａｌＡｐｐｌｉｃａｔｉｏｎｔｉｍｅＬｉＣ ＴｅｕｎｉｓｓｅｎＰｈｙｓｉｃｓ ｔｏｓｔｅｐｐｉｎｇＪ Ｎｏｏｌ２００７ ｇａｓｄｉｓｃｈａｒｇｅｓ［Ｊ］．ｆｏｒｍｕｌｔｉ⁃ｓｃａｌｅＭ ｅｔ２２７（２） ａｌ．Ａｃｏｍｐａｒｉｓｏｎ８９８⁃９１８．Ｊｏｕｒｎａｌｔｒａｎｓｐｏｒｔｏｆｐａｒｔｉｃｌｅ

ｏｆ３Ｄ

２０１２ ｓｔｒｅａｍｅｒｓ［Ｊ］．ｆｌｕｉｄａｎｄｈｙｂｒｉｄｓｉｍｕｌａｔｉｏｎｓｆｏｒｎｅｇａｔｉｖｅ

［２６］Ｃｈａｎｒｉｏｎ２１（５） ＰｌａｓｍａｌｅｃｔｒｏｎｓｆｒｏｍＯ Ｂｏｎａｖｅｎｔｕｒａ０５５０１９．

ＳｏｕｒｃｅｓＳｃｉｅｎｃｅａｎｄＴｅｃｈｎｏｌｏｇｙ

ａ′ｂｅａｍ⁃ｂｕｌｋ′Ｚ ｍｏｄｅｌＣｉｎａｒＤ ｏｆｓｔｒｅａｍｅｒ ｅｔａｌ．Ｒｕｎａｗａｙｅ⁃

９（５） ｔｉｏｎｔｏＴＧＦｓ［Ｊ］．０５５００３．

ＥｎｖｉｒｏｎｍｅｎｔａｌＲｅｓｅａｒｃｈＬｅｔｔｅｒｓ ａｐｐｌｉｃａ⁃

２０１４

［２７］ＢａｂａｅｖａｈｙｂｒｉｄｏｆｐｏｌａｒｉｔｙｍｏｄｅｌｉｎｇＮＹ ＴｅｒｅｓｈｏｎｏｋａｎｄｓｅｃｏｎｄａｒｙｏｆｎａｎｏｓｅｃｏｎｄＤＶ ｅｌｅｃｔｒｏｎｓｓｕｒｆａｃｅＮａｉｄｉｓｅｍｉｓｓｉｏｎ［Ｊ］．ｄｉｓｃｈａｒｇｅｓ ＧＶ．Ｆｌｕｉｄｅｆｆｅｃｔａｎｄ

ＳｏｕｒｃｅｓＳｃｉｅｎｃｅａｎｄＴｅｃｈｎｏｌｏｇｙ Ａ Ａｇａｒｗａｌ２０１６ Ａ ｅｔ２５（４） ａｌ．Ｔｈｒｅｅ⁃ｄｉｍｅｎ⁃

０４４００８．Ｐｌａｓｍａ［２８］ｓｉｏｎａｌＲａｕｆＳ ＰｌａｓｍａｍｏｄｅｌＢａｌａｋｒｉｓｈｎａＳｏｕｒｃｅｓｏｆｅｌｅｃｔｒｏｎＳｃｉｅｎｃｅａｎｄｂｅａｍＴｅｃｈｎｏｌｏｇｙ ｇｅｎｅｒａｔｅｄｐｌａｓｍａ［Ｊ］．（６） ０６５００６．

２０１７ ２６［２９］ＣｈａｎｒｉｏｎｏｆｓｔｒｅａｍｅｒＯ ＮｅｕｂｅｒｔＴ．ＡＰＩＣ⁃ＭＣＣｃｏｄｅｆｏｒｓｉｍｕｌａｔｉｏｎ［３０］ＬｅｖｋｏｔｉｏｎａｌＰｈｙｓｉｃｓ ｐｒｏｐａｇａｔｉｏｎＤ Ｙａｔｏｍ２００８ Ｓ Ｖｅｋｓｅｌｍａｎ２２７（１５） ｉｎａｉｒ［Ｊ］．Ｖ ７２２２⁃７２４５．

ＪｏｕｒｎａｌｏｆＣｏｍｐｕｔａ⁃ｕｌａｔｉｏｎｓｏｆｒｕｎａｗａｙｅｌｅｃｔｒｏｎｇｅｎｅｒａｔｉｏｎｅｔａｌ．Ｎｕｍｅｒｉｃａｌｉｎｐｒｅｓｓｕｒｉｚｅｄ

ｓｉｍ⁃

［３１］０１３３０３．

ｇａｓｅｓ［Ｊ］．ＪｏｕｒｎａｌｏｆＡｐｐｌｉｅｄＰｈｙｓｉｃｓ ２０１２ １１１（１） ＭａｃＬａｃｈｌａｎＣＳ ＰｏｔｔｓＨＥ Ｄｉｖｅｒｓｕｂ⁃ｎｓＤＡ．ｓｔｒｅａｍｅｒＳｉｍｕｌａｔｉｏｎｆｏｒｍａ⁃ｏｆ

［３２］２２（１） ｔｉｏｎ［Ｊ］．ｔｒａｎｓｉｅｎｔＴｅｕｎｉｓｓｅｎ０１５０２５．

ＰｌａｓｍａｅｎｅｒｇｙＳｏｕｒｃｅｓｄｉｓｔｒｉｂｕｔｉｏｎｓＳｃｉｅｎｃｅｉｎａｎｄＴｅｃｈｎｏｌｏｇｙ ２０１３

ｃｈａｒｇｅｉｎｃｅｐｔｉｏｎＪ ＥｂｅｒｔａｒｏｕｎｄＵ．３ＤａｓｈａｒｐＰＩＣ⁃ＭＣＣａｎｏｄｅｓｉｍｕｌａｔｉｏｎｓｉｎｎｉｔｒｏｇｅｎ／ｏｆｏｘｙｇｅｎ

ｄｉｓ⁃

［３３］２０１６ ｍｉｘｔｕｒｅｓ［Ｊ］．Ｃｏｄｉｎａ２５（４） ＰｌａｓｍａＲ Ｏñａｔｅ０４４００５．

ＳｏｕｒｃｅｓＳｃｉｅｎｃｅａｎｄＴｅｃｈｎｏｌｏｇｙ

Ｅ ＣｅｒｖｅｒａＭ．Ｔｈｅｉｎｔｒｉｎｓｉｃｔｉｍｅｆｏｒｔｈｅ

ｓｔｒｅａｍｌｉｎｅｑｕａｄｒａｔｉｃＭｅｃｈａｎｉｃｓｅｌｅｍｅｎｔｓ［Ｊ］．ｕｐｗｉｎｄ／Ｐｅｔｒｏｖ⁃Ｇａｌｅｒｋｉｎ

ａｎｄＥｎｇｉｎｅｅｒｉｎｇ Ｃｏｍｐｕｔｅｒ１９９２ Ｍｅｔｈｏｄｓｆｏｒｍｕｌａｔｉｏｎ

９４（２） ２３９⁃２６２．

ｉｎＡｐｐｌｉｅｄｕｓｉｎｇ

［３４］ＣｏｄｉｎａｓｏｌｖｉｎｇＲ．ＣｏｍｐａｒｉｓｏｎｏｆｓｏｍｅｆｉｎｉｔｅｅｌｅｍｅｎｔＣｏｍｐｕｔｅｒｔｈｅＭｅｔｈｏｄｓｄｉｆｆｕｓｉｏｎ⁃ｃｏｎｖｅｃｔｉｏｎ⁃ｒｅａｃｔｉｏｎｉｎｅｑｕａｔｉｏｎ［Ｊ］．ｍｅｔｈｏｄｓｆｏｒ

ｒｉｎｇ １９９８ １５６（１／４） Ａｐｐｌｉｅｄ１８５⁃２１０．

ＭｅｃｈａｎｉｃｓａｎｄＥｎｇｉｎｅｅ⁃［３５］ＣｏｄｉｎａｓｙｓｔｅｍｓＲ．ＣｏｍｐｕｔｅｒｏｆＯｎｓｔａｂｉｌｉｚｅｄｆｉｎｉｔｅｅｌｅｍｅｎｔｍｅｔｈｏｄｓｆｏｒｌｉｎｅａｒ

Ｍｅｔｈｏｄｓｃｏｎｖｅｃｔｉｏｎ⁃ｄｉｆｆｕｓｉｏｎ⁃ｒｅａｃｔｉｏｎｉｎＡｐｐｌｉｅｄＭｅｃｈａｎｉｃｓｅｑｕａｔｉｏｎｓ［Ｊ］．ｒｉｎｇ ２０００ １８８（１／３） ６１⁃８２．ａｎｄＥｎｇｉｎｅｅ⁃［３６］ｐｅｎｄｅｎｔ

ＪｏｈｎＶ Ｓｃｈｍｅｙｅｒｃｏｎｖｅｃｔｉｏｎ⁃ｄｉｆｆｕｓｉｏｎ⁃ｒｅａｃｔｉｏｎ

Ｅ．Ｆｉｎｉｔｅｅｌｅｍｅｎｔｍｅｔｈｏｄｓｅｑｕａｔｉｏｎｓ

ｆｏｒｔｉｍｅ⁃ｄｅ⁃

ｓｍａｌｌｄｉｆｆｕｓｉｏｎ［Ｊ］．ＣｏｍｐｕｔｅｒＭｅｔｈｏｄｓｉｎＡｐｐｌｉｅｄｗｉｔｈ

［３７］ｃｈａｎｉｃｓＣｏｄｉｎａｄｉｆｆｕｓｉｏｎＲ．ａｎｄＥｎｇｉｎｅｅｒｉｎｇ ２００８ １９８（３／４） ４７５⁃４９４．

Ｍｅ⁃ｔｉｅｓＣｌａｖｅｒｏａｎｄｅｑｕａｔｉｏｎ ＦｉｎｉｔｅｅｌｅｍｅｎｔＣ ａｎｉｓｏｔｒｏｐｉｃｓｕｂｇｒｉｄ⁃ｓｃａｌｅａｐｐｒｏｘｉｍａｔｉｏｎＧｒａｃｉａＪ Ｌｉｓｂｏｎａｓｐａｃｅ⁃ｔｉｍｅｓｐａｃｅｓ ｏｆｔｈｅｃｏｎｖｅｃｔｉｏｎ⁃

ｄｉｓｃｒｅｔｉｚａｔｉｏｎｓ［Ａ］ｌｏｃａｌｉｎｓｔａｂｉｌｉ⁃／／ａｎｄＭｅｔｈｏｄｓ［Ｍ］．ＩｎｔｅｒｉｏｒＬａｙｅｒｓ

Ｂｅｒｌｉｎ Ｓｐｒｉｎｇｅｒ ＣｏｍｐｕｔａｔｉｏｎａｌＦ．ＢＡＩＬ２０１１ ８５⁃９７．ａｎｄ２０１０⁃Ｂｏｕｎｄａｒｙ

Ａｓｙｍｐｔｏｔｉｃ

［３８］Ｐａｎｃｈｅｓｈｎｙｉｄｉｆｆｅｒｅｎｔｐｍｅｎｔｏｆａｐｒｅｓｓｕｒｅｓ

ｃａｔｈｏｄｅ⁃ｄｉｒｅｃｔｅｄＳ ＮｕｄｎｏｖａＭ ｅｘｐｅｒｉｍｅｎｔｓｔｒｅａｍｅｒＳｔａｒｉｋｏｖｓｋｉｉａｎｄｄｉｓｃｈａｒｇｅＡ．ｃｏｍｐａｒｉｓｏｎｉｎＤｅｖｅｌｏ⁃

ａｉｒｗｉｔｈ

ａｔｄｉｒｅｃｔｎｕｍｅｒｉｃａｌｓｉｍｕｌａｔｉｏｎ［Ｊ］．ＰｈｙｓｉｃａｌＲｅｖｉｅｗＥ

［３９］Ｔａｒｄｉｖｅａｕ２００５ ７１（１） ０１６４０７．

ａｎｄｎａｎｏｓｅｃｏｎｄｄｉｆｆｕｓｅ⁃ｔｏ⁃ｆｉｌａｍｅｎｔａｒｙＰ ＭｏｒｅａｕＮ ＢｅｎｔａｌｅｂｓｃａｌｅｃｏｒｏｎａｔｒａｎｓｉｔｉｏｎＳ ｄｉｓｃｈａｒｇｅｉｎｅｔａｌ．ｕｎｄｅｒａｈｉｇｈＤｉｆｆｕｓｅｍｏｄｅ

ｈｉｇｈｐｒｅｓｓｕｒｅｖｏｌｔａｇｅ［４０］Ｂｏｕｒｄｏｎ（１７） ［Ｊ］．Ｊｏｕｒｎａｌ１７５２０２．

ｏｆＰｈｙｓｉｃｓＤ ＡｐｐｌｉｅｄＰｈｙｓｉｃｓ ２００９ ４２ｆｏｒｐｈｏｔｏｉｏｎｉｚａｔｉｏｎＡ ＰａｓｋｏｐｒｏｄｕｃｅｄＶＰ ＬｉｕｂｙＮＹ ｎｏｎ⁃ｔｈｅｒｍａｌｅｔａｌ．Ｅｆｆｉｃｉｅｎｔｇａｓｄｉｓｃｈａｒ⁃

ｍｏｄｅｌｓ

ｇｅｓｉｎａｉｒｂａｓｅｄｏｎｒａｄｉａｔｉｖｅｔｒａｎｓｆｅｒａｎｄｔｈｅＨｅｌｍｈｏｌｔｚ［４１］蔡新景２００７ ｅｑｕａｔｉｏｎｓ［Ｊ］．１６（３） Ｐｌａｓｍａ６５６⁃６７８．

ＳｏｕｒｃｅｓＳｃｉｅｎｃｅａｎｄＴｅｃｈｎｏｌｏｇｙ 流注放电过程光电离快速计算 王新新 邹晓兵 等．基于［Ｊ］．Ｈｅｌｍｈｏｌｔｚ中国电机工程学模型的

报 ２０１５ ３５（１） ２４０⁃２４６．

ＣａｉｐｈｏｔｏｉｏｎｉｚａｔｉｏｎＸＪ ＷａｎｇＸＸ ＺｏｕＸＢ ｅｔａｌ．Ｆａｓｔｃｏｍｐｕｔａｔｉｏｎｏｆ２４６（ｉｎｍｏｄｅｌ［Ｊ］．ＰｒｏｃｅｅｄｉｎｇｓｉｎｓｔｒｅａｍｅｒｏｆｔｈｅｄｉｓｃｈａｒｇｅｓＣＳＥＥ ２０１５ ｂａｓｅｄｏｎ３５（１） ＨｅｌｍｈｏｌｔｚＣｈｉｅｓｅ）．２４０⁃［４２］ＰｅｃｈｅｒｅａｕｃｏｎｉｃａｌｐｈｅｒｉｃｄｉｓｃｈａｒｇｅＦ ＬｅｓｔｒｕｃｔｕｒｅＤｅｌｌｉｏｕＰ ＪáｎｓｋｙＪ ｅｔａｌ．Ｌａｒｇｅ

ＴｒａｎｓａｃｔｉｏｎｓｐｒｅｓｓｕｒｅｏｎＰｌａｓｍａｉｎａｐｏｉｎｔ⁃ｔｏ⁃ｐｌａｎｅｏｆａｎａｉｒｄｉｓｃｈａｒｇｅａｔａｔｍｏｓ⁃２３４７．

Ｓｃｉｅｎｃｅ ２０１４ ｇｅｏｍｅｔｒｙ［Ｊ］．４２（１０） ２３４６⁃

ＩＥＥＥ.com.cn. All Rights Reserved.