高级搜索

电弧等离子体非平衡输运数值模拟和关键参数主动调控研究进展

,  ,  , 

downloadPDF

上一篇

下一篇

曾云海, 张子明, 方川, 李和平. 电弧等离子体非平衡输运数值模拟和关键参数主动调控研究进展[J]. 物理学报, 2026, 75(9): 090501. doi: 10.7498/aps.75.20251759
引用本文: 曾云海, 张子明, 方川, 李和平. 电弧等离子体非平衡输运数值模拟和关键参数主动调控研究进展[J]. 物理学报, 2026, 75(9): 090501. doi: 10.7498/aps.75.20251759
shu
Yunhai ZENG, Ziming ZHANG, Chuan FANG, Heping LI. Research progress on numerical simulations of non-equilibrium transports in arc plasmas and active controls of key plasma parameters[J]. Acta Physica Sinica, 2026, 75(9): 090501. doi: 10.7498/aps.75.20251759
Citation: Yunhai ZENG, Ziming ZHANG, Chuan FANG, Heping LI. Research progress on numerical simulations of non-equilibrium transports in arc plasmas and active controls of key plasma parameters[J]. Acta Physica Sinica, 2026, 75(9): 090501. doi: 10.7498/aps.75.20251759
shu

电弧等离子体非平衡输运数值模拟和关键参数主动调控研究进展

  • 清华大学工程物理系, 北京　100084

    • 通讯作者: E-mail: liheping@tsinghua.edu.cn.

  • 中图分类号: 05.60.-k, 05.70.Ln, 52.50.Nr, 94.20.Fg

Research progress on numerical simulations of non-equilibrium transports in arc plasmas and active controls of key plasma parameters

  • Department of Engineering Physics, Tsinghua University, Beijing 100084, China

    • Corresponding author: E-mail: liheping@tsinghua.edu.cn.

  • MSC: 05.60.-k, 05.70.Ln, 52.50.Nr, 94.20.Fg

  • 摘要: 电弧等离子体由于其具有高能量密度和丰富化学活性粒子的特点, 在能源环保、航空航天、先进材料和国防等诸多领域具有广泛的应用. 而等离子体非平衡特性则是影响应用范围和效果的关键因素之一. 本文以自由燃烧弧作为研究电弧等离子体非平衡协同输运的典型模式体系, 系统地总结了直流电弧等离子体非平衡输运物理-数学模型的建立以及非平衡数值模拟工作的研究进展, 讨论了电弧等离子体体系非平衡质能协同输运机制及其对等离子体特性的影响规律, 总结了目前通过发生器结构设计、外加电磁场位型设计及工作参数调节等手段对电弧等离子体关键参数进行主动调控的方法; 并基于此, 简要讨论了推动电弧等离子体基础研究和工业应用需要解决的若干关键科学问题.
    Abstract: Arc plasma, characterized by its high energy density and abundant chemically reactive species, has wide applications in numerous fields including energy and environmental protection, aeronautics and astronautics, advanced materials processing, and national defense. The non-equilibrium feature of arc plasmas is one of the key factors influencing the scope and effectiveness of plasma applications. Taking a direct-current (DC) free-burning arc as a model system for studying non-equilibrium synergistic transports in arc plasmas, this paper provides a comprehensive review of the state-of-the-art of numerical modeling of transport mechanisms and active modulation methods of key parameters of arc plasmas. First, the development of physical-mathematical models and the related research progress of numerical simulations for non-equilibrium transports in DC arc plasmas are summarized, highlighting a shift in research focus from local thermodynamic equilibrium modeling for the arc column region to considering thermal, chemical, or even electrical non-equilibrium effects within the cold wall boundary layer and/or the plasma-cold gas interaction region, and further extending to multi-region-coupled modeling with consideration of the arc column, electrode boundary layer, solid electrode, and external circuit. Consequently, based on systematic discussions of the synergistic mass-momentum-energy transport mechanism and its influences on the characteristics of arc plasma systems, a summary of current active regulation methods for the key parameters of arc plasmas is presented, e.g., with the aid of geometrical design of plasma generators, variation of external electromagnetic field configurations, and adjustment of operating parameters. Finally, the key scientific issues for promoting fundamental research and industrial applications of arc plasmas are discussed briefly, e.g., development of multi-region, multi-phase coupled numerical models; high-efficiency numerical methods with high spatiotemporal resolutions; high-efficiency chemical reaction pathway screening methods; reliable databases for chemical reaction kinetics and related thermodynamic and transport properties of plasmas; and even creation of a new research paradigm for analyzing complicated mass-energy synergistic transport mechanisms and for developing novel arc plasma sources with desired parameters facing various application demands.
  • 加载中

  • 图 1  电弧等离子体的(a)非平衡模式体系及(b)主要应用领域示意图

    Figure 1.  Schematics of (a) non-equilibrium model system and (b) major application fields of arc plasmas.

    下载: 全尺寸图片 幻灯片

    图 2  电弧等离子体数值模拟不同模型示意图[38]

    Figure 2.  Schematic of different numerical simulation models for arc plasmas[38].

    下载: 全尺寸图片 幻灯片

    图 3  氩自由燃烧电弧典型放电工况下的二维温度场和流场分布(工作气压p = 1 atm, 弧电流I = 200 A, 电极间距L = 10 mm)[51]

    Figure 3.  Two-dimensional distributions of the temperature and flow fields for a free-burning argon arc (Operating pressure p = 1 atm, arc current I = 200 A, electrode gap spacing L = 10 mm)[51].

    下载: 全尺寸图片 幻灯片

    图 4  (a)采用2T-LCE模型计算得到的氩自由燃烧电弧等离子体中电子温度和重粒子温度二维空间分布及与实验测量结果的比较[17]; (b)当气体温度低于5000 K时采用提高电子与重粒子间碰撞频率方法计算得到的等离子体弧柱区电子与重粒子温度分布图[95]. 此处, 图中的数字1, 5, 7, ···, 21给出了箭头所指的等值线的温度值; 例如, 数字“5”表示所指的电子或重粒子温度为5×103 K

    Figure 4.  (a) Calculated 2-D distributions of the electron and heavy-particle temperatures, as well as comparison with the experimental data[17]; (b) the calculated electron and heavy-particle temperature distributions with an artificially increased collision frequency between electrons and heavy particles when the gas temperature is lower than 5000 K in a free-burning argon arc[95] based on the 2T-LCE model. Here, each number such as 1, 5, 7, ···, 21 refers to the temperature value of the corresponding contour; for example, “5” means that the electron temperature or the heavy particle temperature of the corresponding contours is 5×103 K.

    下载: 全尺寸图片 幻灯片

    图 5  采用(a) 2T-LCE模型和(b) 2T-nLCE模型计算得到的氩电弧等离子体气体温度(Th)、电子温度(Te)及相应的热非平衡度(θ)沿电弧径向的变化规律. 工作参数: 弧电流I = 100 A, 气压p = 1 atm, 电极间距L = 8 mm[26]

    Figure 5.  Calculated radial distributions of the heavy-particle temperature (Th) and electron temperature (Te), as well as the corresponding thermal non-equilibrium parameter (θ) for argon arc plasmas using (a) 2T-LCE model and (b) 2T-nLCE model, respectively. Operating parameters: arc current I = 100 A, pressure p = 1 atm, electrode gap spacing L = 8 mm[26].

    下载: 全尺寸图片 幻灯片

    图 6  考虑与不考虑激发态氩原子对氩电弧射流推力器内电流和电离区空间分布的影响[100]. 工作参数: 氩气流量Q = 115 mg/s, 弧电流I = 8 A

    Figure 6.  Comparison of the predicted enclosed current contours and ionization regions with and without considering the excited argon atoms within the arcjet nozzle for argon discharges. Operating parameters: argon flowrate Q = 115 mg/s, arc current I = 8 A[100].

    下载: 全尺寸图片 幻灯片

    图 7  基于多区域耦合模型的(a)电弧等离子体体系示意图和(b)数值模拟流程图[38,120]

    Figure 7.  (a) Schematics of the arc plasma system and (b) flowchart for numerical modeling based on the multi-region-coupled non-equilibrium fluid model[38,120].

    下载: 全尺寸图片 幻灯片

    图 8  不同磁化方向下(a)母材表面处电弧温度场计算结果, (b)焊接电弧形貌图及(c)母材焊缝形貌图[125,126]

    Figure 8.  (a) Modeling results of arc temperature distributions on the base metal surface, (b) arc morphology of welding, and (c) weld morphology of base metal under different magnetization directions[125,126].

    下载: 全尺寸图片 幻灯片

    图 9  (a)不同气体旋流注入角以及(b)不同外加磁场强度条件下的阳极弧根运动速度、阳极表面最高温度[128]

    Figure 9.  Anode arc attachment velocity and maximum anode surface temperature under (a) different gas swirl injection angles and (b) different external magnetic fields[128].

    下载: 全尺寸图片 幻灯片

    图 10  (a)电子-原子和电子-离子碰撞频率沿电弧几何轴线的变化规律[38]; (b)壁面驻点处气体温度和热流密度随真空腔压强的变化规律[36]

    Figure 10.  (a) Axial distributions of the electron-atom and electron-ion collision frequencies[38]; (b) variations of the gas temperature and heat flux density at the wall stagnation point with pressures[36].

    下载: 全尺寸图片 幻灯片

    图 11  (a)采用环形阳极结构的电弧等离子体发生器示意图及(b)反向冷气体注入条件下的典型氩电弧放电照片, 计算得到的平板型阳极和环形阳极发生器对应的(c)热非平衡度(θ)二维分布图, 以及(d)电子温度、重粒子温度、电子数密度沿电弧几何轴线的分布图[122]

    Figure 11.  (a) Schematic of the annular anode arc plasma generator, and (b) typical discharge image of an argon arc plasma with the counter injection of a cold argon gas, (c) two-dimensional distributions of the thermal non-equilibrium degree (θ), and (d) the axial profiles of the electron and heavy-particle temperatures and electron number density along the geometrical axis of the arc plasma generator using the solid and annular anodes[122].

    下载: 全尺寸图片 幻灯片

    图 12  利用外加电磁场调控等离子体形貌和粒子通量　(a)采用S型外加磁场分离等离子体与固体颗粒[131]; (b)采用磁阱结构改变电弧等离子体形貌从而调节基板处的离子电流分布[130,132]

    Figure 12.  Control of morphologies and particle fluxes of arc plasmas with the aid of external electromagnetic fields: (a) Separation of macroparticles and cathodic arc plasma using an open S-shaped magnetic filter[131]; (b) control of plasma arc jet shape and regulating ion current distribution to substrate by using a magnetic trap[130,132].

    下载: 全尺寸图片 幻灯片

    图 13  复合电弧(a)发生装置示意图及(b)平均弧宽和峰值压强随自由弧电流的变化规律[133]

    Figure 13.  (a) Schematic of a hybrid arcing torch, and (b) variations of averaged arc width and arc pressure peak value with the outer free arc current[133].

    下载: 全尺寸图片 幻灯片

    图 14  电弧等离子体体积的调控: 采用聚焦气体的(a)电弧等离子体发生器结构示意图及(b)对应的电弧形貌和焊缝形貌照片[138]

    Figure 14.  Volume control of arc plasmas: (a) Schematic of a gas focusing plasma arc welding torch; (b) the images of arcs and their corresponding weld formations[138].

    下载: 全尺寸图片 幻灯片
    下载: 全尺寸图片 幻灯片

    表 1  不同状态下等离子体控制方程的完整形式

    Table 1.  Governing equations corresponding to different states of plasmas

    状态 参数
    条件    		 质量守恒方程 动量守恒方程 能量守恒方程 电流连续
    性方程    		 可能的
    适用区域
    2T-nLCE-
    nQCN    		 Te ≠ Th
    Rne ≠ 0
    ne ≠ np    		 $ \begin{aligned} &\dfrac{\partial {\rho }_{i}}{\partial t}+\nabla \cdot \left({\rho }_{i}{\boldsymbol{v}}_{0}\right)\\ =\;&-\nabla \cdot \left({\rho }_{i}{\boldsymbol{V}}_{i}\right) \\& +\displaystyle\sum \limits_{k=1}^{r}{\xi }_{ki}{G}_{ki} \end{aligned} $ $ \begin{aligned} &\dfrac{\partial }{\partial t}\left(\rho {\boldsymbol{v}}_{0}\right)+\nabla \cdot \left(\rho {\boldsymbol{v}}_{0}{\boldsymbol{v}}_{0}\right)\\ =\;&-\nabla p+ \nabla \cdot {\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\leftrightarrow$}}\over{\tau}}}+{\boldsymbol{j}}_{\text{cond}}\times \boldsymbol{B}\\ &+{\rho }^{\text{c}}\left(\boldsymbol{E}+{\boldsymbol{v}}_{0}\times \boldsymbol{B}\right) \end{aligned} $ $\left\{\begin{aligned} & 电子: \dfrac{\partial \left({\rho }_{\text{e} }h_{\text{e} }^{*}\right)}{\partial t}+\nabla \cdot \left({\rho }_{\text{e} }h_{\text{e} }^{*}{\boldsymbol{v} }_{0}\right)= \\& \nabla \cdot \left({\lambda }_{\text{e} }\nabla {T}_{\text{h} }\right) +\nabla \cdot \left(\lambda _{\text{e} }^{\theta }{T}_{\text{e} }\nabla \ln \theta \right)-\nabla \cdot \left({\rho }_{\text{e} }h_{\text{e} }^{*}{\boldsymbol{V} }_{\text{e} }\right) +\dfrac{\text{d}{p}_{\text{e} } }{\text{d}t} \\ & -({\rho }_{\text{e} }{\boldsymbol{V} }_{\text{e} })\cdot \dfrac{\text{d}{\boldsymbol{v} }_{0} }{\text{d}t}+ {\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\leftrightarrow$}}\over{\tau}}} _{\text{e} }\colon \nabla {\boldsymbol{v} }_{0}+{n}_{\text{e} }\left\langle {\boldsymbol{X} }_{\text{e} }{\boldsymbol{C} }_{\text{e} }\right\rangle +Q_{\text{eh} }^{\text{el} }+Q_{\text{e} }^{\text{inel} } \\& 重粒子: \dfrac{\partial \left({\rho }_{\text{h} }h_{\text{h} }^{*}\right)}{\partial t}+\nabla \cdot \left({\rho }_{\text{h} }h_{\text{h} }^{*}{\boldsymbol{v} }_{0}\right) = \\& \nabla \cdot \left({\lambda }_{\text{h} }\nabla {T}_{\text{h} }\right) +\displaystyle\sum \limits_{i=2}^{N}\nabla \cdot \left(\lambda _{i}^{\theta }{T}_{i}\nabla \ln \theta \right) -\displaystyle\sum \limits_{i=2}^{N}\nabla \cdot \left({\rho }_{i}h_{i}^{*}{\boldsymbol{V} }_{ {i} }\right) \\ & +\displaystyle\sum \limits_{i=2}^{N}\dfrac{\text{d}{p}_{i} }{\text{d}t} -\displaystyle\sum \limits_{i=2}^{N}({\rho }_{ {i} }{\boldsymbol{V} }_{i})\cdot \dfrac{\text{d}{\boldsymbol{v} }_{0} }{\text{d}t}+\displaystyle\sum \limits_{i=2}^{N} {\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\leftrightarrow$}}\over{\tau}}}_{i}\colon \nabla {\boldsymbol{v} }_{0} \\ & +\displaystyle\sum \limits_{i=2}^{N}{n}_{i}\langle {\boldsymbol{X} }_{i}{\boldsymbol{C} }_{i}\rangle -Q_{\text{eh} }^{\text{el} }+Q_{\text{h} }^{\text{inel} }  \end{aligned}\right.$ (12d) 鞘层区
    2T-nLCE-
    QCN    		 Te ≠ Th
    Rne ≠ 0
    ne = np    		 $ \begin{aligned} & \dfrac{\partial {\rho }_{i}}{\partial t}+\nabla \cdot ({\rho }_{i}{\boldsymbol{v}}_{0}) \\ =\; &-\nabla \cdot ({\rho }_{i}{\boldsymbol{V}}_{i} ) \\& +\displaystyle\sum \limits_{k=1}^{r}{\xi }_{ki}{G}_{ki} \end{aligned} $ $ \begin{aligned} &\dfrac{\partial }{\partial t}\left(\rho {\boldsymbol{v}}_{0}\right)+\nabla \cdot \left(\rho {\boldsymbol{v}}_{0}{\boldsymbol{v}}_{0}\right)\\=\; &-\nabla p+\nabla \cdot {\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\leftrightarrow$}}\over{\tau}}} +{\boldsymbol{j}}_{\text{cond}}\times \boldsymbol{B} \end{aligned} $ $ \left\{\begin{aligned} &电子: \dfrac{\partial \left({\rho }_{\text{e}}h_{\text{e}}^{*}\right)}{\partial t}+\nabla \cdot \left({\rho }_{\text{e}}h_{\text{e}}^{*}{\boldsymbol{v}}_{0}\right) = \\& \nabla \cdot \left({\lambda }_{\text{e}}\nabla {T}_{\text{h}}\right)+ \nabla \cdot \left(\lambda _{\text{e}}^{\theta }{T}_{\text{e}}\nabla \ln \theta \right)-\nabla \cdot \left({\rho }_{\text{e}}h_{\text{e}}^{*}{\boldsymbol{V}}_{\text{e}}\right) +\dfrac{\text{d}{p}_{\text{e}}}{\text{d}t} \\ & -({\rho }_{\text{e}}{\boldsymbol{V}}_{\text{e}})\cdot \dfrac{\text{d}{\boldsymbol{v}}_{0}}{\text{d}t}+ {\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\leftrightarrow$}}\over{\tau}}} _{\text{e}}\colon \nabla {\boldsymbol{v}}_{0} + {n}_{\text{e}}\left\langle {\boldsymbol{X}}_{\text{e}}{\boldsymbol{C}}_{\text{e}}\right\rangle +Q_{\text{eh}}^{\text{el}}+Q_{\text{e}}^{\text{inel}} \\&重粒子: \dfrac{\partial \left({\rho }_{\text{h}}h_{\text{h}}^{*}\right)}{\partial t}+\nabla \cdot \left({\rho }_{\text{h}}h_{\text{h}}^{*}{\boldsymbol{v}}_{0}\right) = \\& \nabla \cdot \left({\lambda }_{\text{h}}\nabla {T}_{\text{h}}\right) +\displaystyle\sum \limits_{i=2}^{N}\nabla \cdot \left(\lambda _{i}^{\theta }{T}_{i}\nabla \ln \theta \right)-\displaystyle\sum \limits_{i=2}^{N}\nabla \cdot \left({\rho }_{i}h_{i}^{*}{\boldsymbol{V}}_{{i}}\right) \\ &+\displaystyle\sum \limits_{i=2}^{N}\dfrac{\text{d}{p}_{i}}{\text{d}t} -\displaystyle\sum \limits_{i=2}^{N}({\rho }_{{i}}{\boldsymbol V}_{i})\cdot \dfrac{\text{d}{\boldsymbol{v}}_{0}}{\text{d}t} +\displaystyle\sum \limits_{i=2}^{N} {\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\leftrightarrow$}}\over{\tau}}}_{i}\colon \nabla {\boldsymbol{v}}_{0} \\ & +\displaystyle\sum \limits_{i=2}^{N}{n}_{i}\langle {\boldsymbol{X}}_{i}{\boldsymbol{C}}_{i}\rangle -Q_{\text{eh}}^{\text{el}}+Q_{\text{h}}^{\text{inel}}  \end{aligned}\right. $ $ \nabla \cdot (\sigma \nabla {\varphi }^{\text{ext}})=0 $ 电离区、
    非平衡
    过渡区、
    弧柱区
    2T-LCE-
    QCN    		 Te ≠ Th
    Rne = 0
    ne = np    		 $ \dfrac{\partial \rho }{\partial t}+\nabla \cdot \left(\rho {\boldsymbol{v}}_{0}\right)=0 $ $ \begin{aligned} &\dfrac{\partial }{\partial t}\left(\rho {\boldsymbol{v}}_{0}\right)+\nabla \cdot \left(\rho {\boldsymbol{v}}_{0}{\boldsymbol{v}}_{0}\right)\\ =\;&-\nabla p+\nabla \cdot {\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\leftrightarrow$}}\over{\tau}}} +{\boldsymbol{j}}_{\text{cond}}\times \boldsymbol{B} \end{aligned} $ $ \left\{\begin{aligned} &电子: \dfrac{\partial \left({\rho }_{\text{e}}h_{\text{e}}^{*}\right)}{\partial t}+\nabla \cdot \left({\rho }_{\text{e}}h_{\text{e}}^{*}{\boldsymbol{v}}_{0}\right) = \\& \nabla \cdot \left({\lambda }_{\text{e}}\nabla {T}_{\text{h}}\right) +\nabla \cdot \left(\lambda _{\text{e}}^{\theta }{T}_{\text{e}}\nabla \ln \theta \right) -\nabla \cdot \left({\rho }_{\text{e}}h_{\text{e}}^{*}{\boldsymbol{V}}_{\text{e}}\right) \\ & +\dfrac{\text{d}{p}_{\text{e}}}{\text{d}t}-({\rho }_{\text{e}}{\boldsymbol{V}}_{\text{e}})\cdot \dfrac{\text{d}{\boldsymbol{v}}_{0}}{\text{d}t}+ {\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\leftrightarrow$}}\over{\tau}}}_{\text{e}}\colon \nabla {\boldsymbol{v}}_{0} +{n}_{\text{e}}\left\langle {\boldsymbol{X}}_{\text{e}}{\boldsymbol{C}}_{\text{e}}\right\rangle \\ & +Q_{\text{eh}}^{\text{el}}-{U}_{\text{rad}} \\ & 重粒子: \dfrac{\partial \left({\rho }_{\text{h}}h_{\text{h}}^{*}\right)}{\partial t}+\nabla \cdot \left({\rho }_{\text{h}}h_{\text{h}}^{*}{\boldsymbol{v}}_{0}\right) = \\& \nabla \cdot \left({\lambda }_{\text{h}}\nabla {T}_{\text{h}}\right) +\displaystyle\sum \limits_{i=2}^{N}\nabla \cdot \left(\lambda _{i}^{\theta }{T}_{i}\nabla \ln \theta \right)-\displaystyle\sum \limits_{i=2}^{N}\nabla \cdot \left({\rho }_{i}h_{i}^{*}{\boldsymbol{V}}_{i}\right) \\ & +\displaystyle\sum \limits_{i=2}^{N}\dfrac{\text{d}{p}_{i}}{\text{d}t} -\displaystyle\sum \limits_{i=2}^{N}({\rho }_{i}{\boldsymbol{V}}_{i})\cdot \dfrac{\text{d}{\boldsymbol{v}}_{0}}{\text{d}t}+\displaystyle\sum \limits_{i=2}^{N} {\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\leftrightarrow$}}\over{\tau}}}_{i}\colon \nabla {\boldsymbol{v}}_{0} \\ & +\displaystyle\sum \limits_{i=2}^{N}{n}_{i}\langle {\boldsymbol{X}}_{i}{\boldsymbol{C}}_{i}\rangle -Q_{\text{eh}}^{\text{el}}  \end{aligned}\right. $ $ \nabla \cdot (\sigma \nabla {\varphi }^{\text{ext}})=0 $ 非平衡
    过渡区、
    弧柱区
    LTE-LCE-
    QCN    		 Te = Th
    Rne = 0
    ne = np    		 $ \dfrac{\partial \rho }{\partial t}+\nabla \cdot \left(\rho {\boldsymbol{v}}_{0}\right)=0 $ $ \begin{aligned} &\dfrac{\partial }{\partial t}\left(\rho {\boldsymbol{v}}_{0}\right)+\nabla \cdot \left(\rho {\boldsymbol{v}}_{0}{\boldsymbol{v}}_{0}\right) \\ =\;& -\nabla p+\nabla \cdot {\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\leftrightarrow$}}\over{\tau}}} +{\boldsymbol{j}}_{\text{cond}}\times \boldsymbol{B} \end{aligned} $ $ \begin{aligned} &\dfrac{\partial \left(\rho {h}^{*}\right)}{\partial t}+\nabla \cdot \left(\rho {h}^{*}{\boldsymbol{v}}_{0}\right)\\ =\;&\nabla \cdot (\lambda \nabla T)+\dfrac{5}{2}\dfrac{{k}_{\mathrm{B}}}{e}\cdot {\boldsymbol{j}}_{\text{cond}}\cdot \nabla T\\ &-\displaystyle\sum \limits_{i=1}^{N}({\rho }_{i}{\boldsymbol{V}}_{i})\cdot \dfrac{\text{d}{\boldsymbol{v}}_{0}}{\text{d}t}+\displaystyle\sum \limits_{i=1}^{N}\dfrac{\text{d}{p}_{i}}{\text{d}t}\\ &+\displaystyle\sum \limits_{i=1}^{N} {\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\leftrightarrow$}}\over{\tau}}}_{i}\colon \nabla {\boldsymbol{v}}_{0}+\displaystyle\sum \limits_{i=1}^{N}{n}_{i}\langle {\boldsymbol{X}}_{i}{\boldsymbol{C}}_{{i}}\rangle -{U}_{\text{rad}} \end{aligned} $ $ \nabla \cdot (\sigma \nabla {\varphi }^{\text{ext}})=0 $ 弧柱区
    注: Rne为电子净产生率; $n_{\rm p} = \displaystyle\sum\nolimits_{j=2}^NZ_jn_j$; Gki为反应速率的绝对值; ξki =1, –1和0分别代表在反应k中, 粒子i生成、消亡和不变; $ Q_{\text{eh}}^{\text{el}} $为电子-重粒子弹性碰撞能量交换项; $ Q_{\text{e}}^{\text{inel}} $和$ Q_{\text{h}}^{\text{inel}} $分别为非弹性碰撞过程中电子和重粒子子体系获得的能量, 可以由唯象理论获得, 具体表达式可参考文献[39]; Urad为等离子体单位体积的辐射热损失; $ \rho {h}^{*}=\displaystyle\sum\nolimits _{i=1}^{N}{\rho }_{i}h_{i}^{*} $为等离子体比焓; λ为等离子体热导率.
    下载: 导出CSV
  • [1] Crookes W 1879 Am. J. Sci. 3 241 doi: 10.2475/ajs.s3-18.106.241
    [2] 王晓钢 2014 等离子体物理基础 (北京: 北京大学出版社) 第3, 4页 Wang X G 2014 Fundamentals of Plasma Physics (Beijing: Peking University Press) pp3, 4
    [3] Nehra V, Kumar A, Dwivedi H K 2008 Int. J. Eng. 2 53
    [4] Samal S 2017 J. Clean Prod. 142 3131 doi: 10.1016/j.jclepro.2016.10.154
    [5] Gururani P, Bhatnagar P, Bisht B, Kumar V, Joshi N C, Tomar M S, Pathak B 2021 Environ. Sci. Pollut. Res. Int. 28 65062 doi: 10.1007/s11356-021-16741-x
    [6] Kalra C S, Gutsol A F, Fridman A A 2005 IEEE Trans. Plasma Sci. 33 32 doi: 10.1109/TPS.2004.842321
    [7] Zhang H, Li L, Li X, Wang W, Yan J, Tu X 2018 J. CO2 Util. 27 472 doi: 10.1016/j.jcou.2018.08.020
    [8] George A, Shen B, Craven M, Wang Y, Kang D, Wu C, Tu X 2021 Renewable Sustainable Energy Rev. 135 109702 doi: 10.1016/j.rser.2020.109702
    [9] Sahoo A, Tripathy S 2021 Mater. Today Proc. 41 363 doi: 10.1016/j.matpr.2020.09.562
    [10] Colombo V, Concetti A, Ghedini E, Dallavalle S, Vancini M 2009 Plasma Sources Sci. Technol. 18 023001 doi: 10.1088/0963-0252/18/2/023001
    [11] Fauchais P 2004 J. Phys. D: Appl. Phys. 37 R86 doi: 10.1088/0022-3727/37/9/R02
    [12] Laux C O, Spence T G, Kruger C H, Zare R N 2003 Plasma Sources Sci. Technol. 12 125 doi: 10.1088/0963-0252/12/2/301
    [13] Jiang L, Chen Y, Mao C, Han J, Chen A, Ye J 2025 Plasma Sci. Technol. 27 015503 doi: 10.1088/2058-6272/ad92f8
    [14] Lee J, Shin S, Jang I, Bae S, Jeong S 2025 Plasma Sci. Technol. 27 035502 doi: 10.1088/2058-6272/ad9bfd
    [15] Yan J, Liu K, Ma J, Li Y, Li K, Wei H 2025 Plasma Sci. Technol. 27 035503 doi: 10.1088/2058-6272/ad9e90
    [16] Ye J, Li Y J, Zhang Z, Wang X, Tao K, Zeng Q, Chen L, Qian D, Zhang S, Yang L, Ma X 2024 Plasma Sci. Technol. 26 095501 doi: 10.1088/2058-6272/ad5067
    [17] Hsu K C, Pfender E 1983 J. Appl. Phys. 54 4359 doi: 10.1063/1.332672
    [18] Chen D M, Hsu K C, Pfender E 1981 Plasma Chem. Plasma Process. 1 295 doi: 10.1007/BF00568837
    [19] Trelles J P, Heberlein J V R, Pfender E 2007 J. Phys. D: Appl. Phys. 40 5937 doi: 10.1088/0022-3727/40/19/024
    [20] Zhang Z, Wang C, Sun Q, Xia W 2022 Plasma Sci. Technol. 24 055402 doi: 10.1088/2058-6272/ac55fe
    [21] Hsu K C, Pfender E 1983 J. Appl. Phys. 54 3818 doi: 10.1063/1.332606
    [22] Heberlein J, Mentel J, Pfender E 2010 J. Phys. D: Appl. Phys. 43 023001 doi: 10.1088/0022-3727/43/2/023001
    [23] Cao M, Proulx P, Boulos M I, Mostaghimi J 1994 J. Appl. Phys. 76 7757 doi: 10.1063/1.357953
    [24] Zhou W, Guo H, Jiang W, Li H P, Li Z Y, Lapenta G 2016 Plasma Sources Sci. Technol. 25 05LT01 doi: 10.1088/0963-0252/25/5/05LT01
    [25] Maecker H 1955 Z. Physik 141 198 doi: 10.1007/BF01327300
    [26] Guo H, Zhang X N, Chen J, Li H P, Ostrikov K 2018 Sci. Rep. 8 4783 doi: 10.1038/s41598-018-22911-8
    [27] Sakura N, Yoshida M, Tanaka M, Takayuki W 2019 J. Phys. D: Appl. Phys. 52 404002 doi: 10.1088/1361-6463/ab3139
    [28] Murphy A B, Tanaka M, Tashiro S, Sato T, Lowke J J 2009 J. Phys. D: Appl. Phys. 42 115205 doi: 10.1088/0022-3727/42/11/115205
    [29] Villarreal-Medina R, Murphy A B, Méndez P F, Ramírez-Argáez M A 2023 Plasma Chem. Plasma Process. 43 787 doi: 10.1007/s11090-023-10328-9
    [30] Kim Y H, Gunasegaram D R, Cleary P W, Murphy A B 2025 J. Phys. D: Appl. Phys. 58 143001 doi: 10.1088/1361-6463/adb3b0
    [31] Le D K, Tashiro S, Trinh Q N, Suga T, Sawamura N, Fukuda K, Sasakura S, Alvarez-Rocha J E, Mendez P F, Murphy A B, Bui V H, Tanaka M 2025 J. Manuf. Processes 139 105 doi: 10.1016/j.jmapro.2025.02.021
    [32] Trinh N Q, Le K D, Tashiro S, Suga T, Sawamura S, Fukuda K, Murphy A B, Bui, H V, Tanaka M 2024 J. Manuf. Processes 124 590 doi: 10.1016/j.jmapro.2024.06.047
    [33] Song I, Jeong G H, Kim S K, Kim Y H, Murphy A B, Park T K, Kim D, Park H, Cho D W 2024 Processes 12 1757 doi: 10.3390/pr12081757
    [34] Goebel D M, Katz I 2008 Fundamentals of Electric Propulsion: Ion and Hall Thrusters (New Jersey: John Wiley & Sons.
    [35] Zha B L, Su Q D, Xie F, Shi Y A, Wang J J 2018 IOP Conf. Ser. : Earth Environ. Sci. 186 012016 doi: 10.1088/1755-1315/186/4/012016
    [36] 方川, 张子明, 汪耀庭, 罗岚月, 曾实, 李志辉, 李和平 2023 力学学报 55 2818 doi: 10.6052/0459-1879-23-337 Fang C, Zhang Z M, Wang Y T, Luo L Y, Zeng S, Li Z H, Li H P 2023 Chin J. Theor. Appl. Mech. 55 2818 doi: 10.6052/0459-1879-23-337
    [37] Mandilas C, Daskalos E, Karagiannakis G, Konstandopoulos A G 2013 Mater. Sci. Eng. , B 178 22 doi: 10.1016/j.mseb.2012.10.004
    [38] Fang C, Chen J, Li J, Zhang Z M, Guo H, Li Z H, Zeng S, Li H P 2022 Plasma Sources Sci. Technol. 31 015015 doi: 10.1088/1361-6595/ac2c8d
    [39] Li H P, Zhang X N, Xia W D 2013 Phys. Plasmas 20 033509 doi: 10.1063/1.4794970
    [40] 张子明, 赵恒欣, 吕自坚, 李志辉, 曾实, 李和平 2025 高电压技术 51 5696 doi: 10.13336/j.1003-6520.hve.20241624 Zhang Z M, Zhao H X, Lyu Z J, Li Z H, Zeng S, Li H P 2025 High Volt. Eng. 51 5696 doi: 10.13336/j.1003-6520.hve.20241624
    [41] Hirschfelder J O, Curtiss C F, Bird R B 1964 Molecular Theory of Gases and Liquids (New York: John Wiley & Sons.
    [42] Devoto R S 1966 Phys. Fluids 9 1230 doi: 10.1063/1.1761825
    [43] Devoto R S 1967 Phys. Fluids 10 2105 doi: 10.1063/1.1762005
    [44] Devoto R S 1968 J. Plasma Phys. 2 617 doi: 10.1017/S0022377800004086
    [45] Devoto R S 1973 Phys. Fluids 16 616 doi: 10.1063/1.1694396
    [46] Devoto R S 1965 Ph. D. Dissertation (Michaigan: Stanford Univeristy
    [47] Bonnefoi C, Aubreton J, Mexmain J M 1985 Z. Naturfors. , A: Phys. Sci. 40 885 doi: 10.1515/zna-1985-0904
    [48] Colombo V, Ghedini E, Sanibondi P 2008 Prog. Nucl. Energy 50 921 doi: 10.1016/j.pnucene.2008.06.002
    [49] Rat V, André P, Aubreton J, Elchinger M F, Fauchais P, Lefort A 2001 Phys. Rev. E 64 026409 doi: 10.1103/PhysRevE.64.026409
    [50] Zhang X N, Li H P, Murphy A B, Xia W D 2013 Phys. Plasmas 20 033508 doi: 10.1063/1.4794969
    [51] Hsu K C, Etemadi K, Pfender E 1983 J. Appl. Phys. 54 1293 doi: 10.1063/1.332195
    [52] Kovitya P, Lowke J J 1985 J. Phys. D: Appl. Phys. 18 53 doi: 10.1088/0022-3727/18/1/009
    [53] Wang H, Chen X, Pan W, Murphy A B, Geng J, Jia S 2010 Plasma Sci. Technol. 12 692 doi: 10.1088/1009-0630/12/6/11
    [54] Lowke J J, Morrow R, Haidar J Murphy A B 1997 IEEE Trans. Plasma Sci. 25 925 doi: 10.1109/27.649597
    [55] Murphy A B, Arundell C J 1994 Plasma Chem. Plasma Process. 14 451 doi: 10.1007/BF01570207
    [56] Boulos M I, Fauchais P L, Pfender E (eds) 2023 Handbook of Thermal Plasmas (Cham: Springer
    [57] Murphy A B 2012 Chem. Phys. 398 64 doi: 10.1016/j.chemphys.2011.06.017
    [58] Murphy A B 1995 Plasma Chem. Plasma Process. 15 279 doi: 10.1007/BF01459700
    [59] Capitelli M, Celiberto R, Gorse C, Giordano D 1995 Plasma Chem. Plasma Process. 16 S267 doi: 10.1007/BF01512639
    [60] Wang C, Wu Y, Chen Z, Yang F, Feng Y, Rong M, Zhang H 2016 Plasma Sci. Technol. 18 732 doi: 10.1088/1009-0630/18/7/06
    [61] Chen S Q, Wang H X 2012 Chin. Phys. Lett. 29 025202 doi: 10.1088/0256-307X/29/2/025202
    [62] Cressault Y, Murphy A B, Teulet Ph, Gleizes A, Schnick M 2013 J. Phys. D: Appl. Phys. 46 415207 doi: 10.1088/0022-3727/46/41/415207
    [63] Wang X, Zhong L, Cressault Y, Gleizes A, Rong M 2014 J. Phys. D: Appl. Phys. 47 495201 doi: 10.1088/0022-3727/47/49/495201
    [64] Hsu K C, Pfender E 1984 Plasma Chem. Plasma Process. 4 219 doi: 10.1007/BF00566842
    [65] Wang H, Chen X, Pan W 2009 Plasma Sci. Technol. 11 163 doi: 10.1088/1009-0630/11/2/07
    [66] Siddanathi S L, Westerberg L G, Åkerstedt H O, Wiinikka H, Sepman A 2025 Appl. Sci. 15 8367 doi: 10.3390/app15158367
    [67] Nishiyama H, Sato T, Takamura K 2003 ISIJ Int. 43 950 doi: 10.2355/isijinternational.43.950
    [68] Wu G Q, Li H P, Bao C Y, Chen X 2008 Prog. Comput. Fluid Dyn. 8 424 doi: 10.1504/PCFD.2008.021318
    [69] Long N P, Katada Y, Tanaka Y, Uesugi Y, Yamaguchi Y 2012 J. Phys. D: Appl. Phys. 45 435203 doi: 10.1088/0022-3727/45/43/435203
    [70] Xiang J, Park H, Tanaka K, Shigeta M, Tanaka M, Murphy A B 2020 J. Phys. D: Appl. Phys. 53 044004 doi: 10.1088/1361-6463/ab51f3
    [71] Wang X, Tashiro S, Tanaka M, Murphy A B 2024 Plasma Chem. Plasma Process. 44 25 doi: 10.1007/s11090-023-10359-2
    [72] Gonzalez J J, Gleizes A, Proulx P, Boulos M 1993 J. Appl. Phys. 74 3065 doi: 10.1063/1.354624
    [73] Li H P, Chen X. 2001 J. Phys. D: Appl. Phys. 34 L99 doi: 10.1088/0022-3727/34/17/102
    [74] Li H P, Chen X 2001 Thin Solid Films 390 175 doi: 10.1016/S0040-6090(01)00915-4
    [75] Li H P, Chen X 2002 Plasma Chem. Plasma Process. 22 27 doi: 10.1023/A:1012988430995
    [76] Li H P, Chen X 2002 Chin. Phys. 11 44 doi: 10.1088/1009-1963/11/1/310
    [77] Li H P, Heberlein J, Pfender E 2005 IEEE Trans. Plasma Sci. 33 402 doi: 10.1109/TPS.2005.844955
    [78] Chen X, Li H P 2003 Surf. Coat. Technol. 171 124 doi: 10.1016/S0257-8972(03)00252-4
    [79] Li H P, Pfender E, Chen X 2003 J. Phys. D: Appl. Phys. 36 1084 doi: 10.1088/0022-3727/36/9/306
    [80] Benilov M S, Benilova L G, Li H P, Wu G Q 2012 J. Phys. D: Appl. Phys. 45 355201 doi: 10.1088/0022-3727/45/35/355201
    [81] Huang R, Fukanuma H, Uesugi Y, Tanaka Y 2012 J. Therm. Spray Technol. 21 636 doi: 10.1007/s11666-011-9710-5
    [82] Huang R, Fukanuma H, Uesugi Y, Tanaka Y 2013 J. Therm. Spray Technol. 22 183 doi: 10.1007/s11666-012-9860-0
    [83] Rat V, André P, Aubreton J, Elchinger M F, Fauchais P, Vacher D 2002 J. Phys. D: Appl. Phys. 35 981 doi: 10.1088/0022-3727/35/10/306
    [84] Zhang X N, Li H P, Murphy A B, Xia W D 2015 Plasma Sources Sci. Technol. 24 035011 doi: 10.1088/0963-0252/24/3/035011
    [85] Ghorui S, Heberlein J V R, Pfender E 2007 Plasma Chem. Plasma Process. 27 267 doi: 10.1007/s11090-007-9053-7
    [86] 王海兴, 孙素蓉, 陈士强 2012 物理学报 61 195203 doi: 10.7498/aps.61.195203 Wang H X, Sun S R, Chen S Q 2012 Acta Phys. Sin. 61 195203 doi: 10.7498/aps.61.195203
    [87] Wang W Z, Rong M Z, Wu Y, Spencer J W, Yan J D, Mei D H 2012 Phys. Plasmas 19 083506 doi: 10.1063/1.4739778
    [88] 潘子晗, 陈仙辉, 王城, 夏维东 2021 物理学报 70 085201 doi: 10.7498/aps.70.20202040 Pan Z H, Chen X H, Wang C, Xia W D 2021 Acta Phys. Sin. 70 085201 doi: 10.7498/aps.70.20202040
    [89] Ghorui S, Heberlein J V R, Pfender E 2008 Plasma Chem. Plasma Process. 28 553 doi: 10.1007/s11090-008-9141-3
    [90] Colombo V, Ghedini E, Sanibondi P 2011 Plasma Sources Sci. Technol. 20 035003 doi: 10.1088/0963-0252/20/3/035003
    [91] Wang G, Zhang B, Cao M, Deng J, Shi Y, Li B, Li X, Murphy A B 2025 J. Phys. D: Appl. Phys. 58 165502 doi: 10.1088/1361-6463/adbcbd
    [92] Wu Y, Chen Z, Yang F, Cressault Y, Murphy A B, Guo A, Liu Z, Rong M, Sun H 2015 J. Phys. D: Appl. Phys. 48 415205 doi: 10.1088/0022-3727/48/41/415205
    [93] Bhigamudre V G, Trelles J P 2019 J. Phys. D: Appl. Phys. 52 015205 doi: 10.1088/1361-6463/aae643
    [94] Charrada K, Zissis G, Aubes M 1996 J. Phys. D: Appl. Phys. 29 2432 doi: 10.1088/0022-3727/29/9/030
    [95] Haidar J 1999 J. Phys. D: Appl. Phys. 32 263 doi: 10.1088/0022-3727/32/3/014
    [96] Ghorui S, Heberlein J V R, Pfender E 2007 J. Phys. D: Appl. Phys. 40 1966 doi: 10.1088/0022-3727/40/7/020
    [97] Baeva M, Benilov M S, Almeida N A, Uhrlandt D 2016 J. Phys. D: Appl. Phys. 49 245205 doi: 10.1088/0022-3727/49/24/245205
    [98] Baeva M, Kozakov R, Gorchakov S, Uhrlandt D 2012 Plasma Sources Sci. Technol. 21 055027 doi: 10.1088/0963-0252/21/5/055027
    [99] Wang H X, He Q S, Murphy A B, Zhu T, Wei F Z 2017 Plasma Chem. Plasma Process. 37 877 doi: 10.1007/s11090-017-9807-9
    [100] Wei F Z, Wang H X, Murphy A B, Sun W P, Liu Y 2013 J. Phys. D: Appl. Phys. 46 505205 doi: 10.1088/0022-3727/46/50/505205
    [101] Wang H X, Sun W P, Sun S R, Murphy A B, Ju Y 2014 Plasma Chem. Plasma Process. 34 559 doi: 10.1007/s11090-013-9501-5
    [102] Baeva M, Uhrlandt D 2013 J. Phys. D: Appl. Phys. 46 325202 doi: 10.1088/0022-3727/46/32/325202
    [103] Golubovskii Y, Kalanov D, Gortschakow S, Baeva M, Uhrlandt D 2016 J. Phys. D: Appl. Phys. 49 475202 doi: 10.1088/0022-3727/49/47/475202
    [104] Golubovskii Y, Kalanov D, Baeva M, Gorchakov S, Uhrlandt D 2015 J. Phys. D: Appl. Phys. 48 225203 doi: 10.1088/0022-3727/48/22/225203
    [105] Wang H X, Sun S R, Sun W P 2015 Plasma Chem. Plasma Process. 35 543 doi: 10.1007/s11090-015-9610-4
    [106] 李和平, 吴贵清, 李鹏, 包成玉 2013 高电压技术 39 1549 doi: 10.3969/j.issn.1003-6520.2013.07.001 Li H P, Wu G Q, Li P, Bao C Y 2013 High Volt. Eng. 39 1549 doi: 10.3969/j.issn.1003-6520.2013.07.001
    [107] Trelles J P 2014 Plasma Sources Sci. Technol. 23 054002 doi: 10.1088/0963-0252/23/5/054002
    [108] Pekker L, Hussary N 2014 J. Phys. D: Appl. Phys. 47 445202 doi: 10.1088/0022-3727/47/44/445202
    [109] Pekker L, Hussary N 2015 Phys. Plasmas 22 083510 doi: 10.1063/1.4928882
    [110] Pekker L, Murphy A B 2016 J. Phys. D: Appl. Phys. 49 375202 doi: 10.1088/0022-3727/49/37/375202
    [111] Franklin R N 2003 J. Phys. D: Appl. Phys. 36 R309 doi: 10.1088/0022-3727/36/22/R01
    [112] Hershkowitz N 2005 Phys. Plasmas 12 055502 doi: 10.1063/1.1887189
    [113] Benilov M S 2008 J. Phys. D: Appl. Phys. 41 144001 doi: 10.1088/0022-3727/41/14/144001
    [114] Benilov M S 2009 Plasma Sources Sci. Technol. 18 014005 doi: 10.1088/0963-0252/18/1/014005
    [115] Robertson S 2013 Plasma Phys. Controlled Fusion 55 093001 doi: 10.1088/0741-3335/55/9/093001
    [116] Li H P, Benilov M S 2007 J. Phys. D: Appl. Phys. 40 2010 doi: 10.1088/0022-3727/40/7/024
    [117] Chen T, Wang C, Liao M R, Xia W D 2016 J. Phys. D: Appl. Phys. 49 085202 doi: 10.1088/0022-3727/49/8/085202
    [118] Chen T, Wang C, Zhang X N, Zhang H, Xia W D 2017 Plasma Sources Sci. Technol. 26 025002 doi: 10.1088/1361-6595/26/2/025002
    [119] Sun Q, Wang C, Chen T, Xia W D 2017 J. Phys. D: Appl. Phys. 50 425202 doi: 10.1088/1361-6463/aa882b
    [120] Zhang Z M, Fang C, Wang Y T, Luo L Y, Li H P 2024 Plasma Sci. Technol. 26 115402 doi: 10.1088/2058-6272/ad6705
    [121] Li H P, Ostrikov K, Sun W 2018 Phys. Rep. 770-772 1 doi: 10.1016/j.physrep.2018.08.002
    [122] Fang C, Zhang Z M, Wang Y T, Luo L Y, Li Z H, Zeng S, Li H P 2023 J. Phys. D: Appl. Phys. 56 11LT01 doi: 10.1088/1361-6463/acbc88
    [123] Sun L, Zhou W, Jiang W, Li H P, Li Z Y 2020 J. Phys. D: Appl. Phys. 53 455201 doi: 10.1088/1361-6463/ab9fc2
    [124] Sun L, Sun X P, Guo B Y, Zhou W, Li Z Y 2022 Phys. Fluids 34 067120 doi: 10.1063/5.0096325
    [125] Nomura K, Morisaki K, Hirata Y 2009 Weld. World 53 R181 doi: 10.1007/BF03266730
    [126] Nomura K, Ogino Y, Hirata Y 2011 Weld. Int. 26 759 doi: 10.1080/09507116.2011.592691
    [127] Zhukovskii R, Chazelas C, Vardelle A, Rat V, Distler B 2020 J. Therm. Spray Technol. 29 894 doi: 10.1007/s11666-020-01052-3
    [128] Zhukovskii R, Chazelas C, Vardelle A, Rat V 2020 J. Therm. Spray Technol. 29 3 doi: 10.1007/s11666-019-00969-8
    [129] Zhukovskii R, Chazelas C, Rat V, Vardelle A, Molz R 2022 J. Therm. Spray Technol. 31 28 doi: 10.1007/s11666-021-01253-4
    [130] Baranov O, Bazaka K, Kersten H, Keidar M, Cvelbar U, Xu S, Levchenko I 2017 Appl. Phys. Rev. 4 041302 doi: 10.1063/1.5007869
    [131] André A 2002 IEEE Trans. Plasma Sci. 30 108 doi: 10.1109/TPS.2002.1003947
    [132] Baranov O, Fang J, Keidar M, Lu X, Cvelbar U, Ostrikov K 2014 IEEE Trans. Plasma Sci. 42 2464 doi: 10.1109/TPS.2014.2323263
    [133] Qiu J Y, Liu Z M 2021 J. Manuf. Processs 64 766 doi: 10.1016/j.jmapro.2021.02.028
    [134] Kato K, Iizuka S, Sato N 1994 Appl. Phys. Lett. 65 816 doi: 10.1063/1.112240
    [135] Yao Y, Yatsuda K, Watanabe T, Matsuura T, Yano T 2009 Plasma Chem. Plasma Process. 29 333 doi: 10.1007/s11090-009-9182-2
    [136] Matsuura T, Taniguchi K, Watanabe T 2007 Thin Solid Films 515 4240 doi: 10.1016/j.tsf.2006.02.086
    [137] Liu Y, Tanaka M, Choi S, Watanabe T 2014 Int. J. Appl. Glass Sci. 5 443 doi: 10.1111/ijag.12081
    [138] Li T Q, Yang X M, Chen L, Zhang Y, Lei Y C, Yan J C 2020 Sci. Technol. Weld. Joining 25 329 doi: 10.1080/13621718.2019.1702284
    [139] Li J, Fang C, Chen J, Li H P, Wang Z B, Nie Q Y, Guo H, Zhao X, Zhao L X, Zhang Y, Ostrikov K 2022 Plasma Sources Sci. Technol. 31 055015 doi: 10.1088/1361-6595/ac6d0c
  • 加载中
图( 15) 表( 1)
计量
  • 文章访问数:  746
  • HTML全文浏览数:  746
  • PDF下载数:  8
  • 施引文献:  0
出版历程
  • 收稿日期:  2025-12-25
  • 刊出日期:  2026-05-05

电弧等离子体非平衡输运数值模拟和关键参数主动调控研究进展

  • 收稿日期:  2025-12-25

摘要: 电弧等离子体由于其具有高能量密度和丰富化学活性粒子的特点, 在能源环保、航空航天、先进材料和国防等诸多领域具有广泛的应用. 而等离子体非平衡特性则是影响应用范围和效果的关键因素之一. 本文以自由燃烧弧作为研究电弧等离子体非平衡协同输运的典型模式体系, 系统地总结了直流电弧等离子体非平衡输运物理-数学模型的建立以及非平衡数值模拟工作的研究进展, 讨论了电弧等离子体体系非平衡质能协同输运机制及其对等离子体特性的影响规律, 总结了目前通过发生器结构设计、外加电磁场位型设计及工作参数调节等手段对电弧等离子体关键参数进行主动调控的方法; 并基于此, 简要讨论了推动电弧等离子体基础研究和工业应用需要解决的若干关键科学问题.

English Abstract

    全文HTML

    1.   引　言

    • 等离子体是一种包含有足够数量的非束缚态电荷的多粒子体系, 具有带电性和多体性的独特性质; 带电性使得等离子体中粒子间存在长程相互作用, 而多体性则体现为多自由度特征; 二者的结合导致等离子体中运动模式的集体效应和多时空尺度效应, 成为不同于固、液、气三态的物质第四态[1,2]. 根据气体温度的不同, 等离子体可分为高温等离子体和低温等离子体. 其中, 高温等离子体的气体温度可达108 K量级, 如聚变等离子体; 而低温等离子体又可进一步分为热等离子体和冷等离子体[3-5], 也有学者将其分为热等离子体、暖等离子体和冷等离子体[6-8]. 其中, 热等离子体由于其重粒子温度(Th)或气体温度(Tg)较高, 且与电子温度(Te)接近, 即Th = Tg ≈ Te ~104 K量级, 因此, 其核心区往往处于局域热力学平衡(local thermodynamic equilibrium, LTE)状态和局域化学平衡(local chemical equilibrium, LCE)状态, 而且由于高气压条件下粒子间的频繁碰撞, 粒子速度/能量分布可认为满足Maxwell分布. 也正是由于处于LTE-LCE状态的热等离子体所具有的高焓值或高能量密度的特性, 使得其在焊接、切割、喷涂、光谱检测等领域得到了广泛应用[9-16]. 而另一方面, 就整个热等离子体体系而言, 由于核心区等离子体与冷壁面或环境冷气体间存在复杂的相互作用, 从而使得靠近冷壁面或与环境冷气体交界的这些区域的等离子体往往显著偏离平衡状态[17-20]. 对于本文所讨论的电弧等离子体, 通常可采用直流(direct current, DC)或交流(alternating current, AC)电源驱动产生, 其典型的放电形式包括直流自由燃烧电弧、转移电弧或非转移电弧等[4].

      以直流自由燃烧电弧为例, 如图1(a)所示, 整个等离子体发生系统包括了外电路和放电单元两个部分, 而等离子体放电单元又可划分为五个区域, 即阴极和阳极固体区、对应的两个电极边界层区域, 以及位于二者之间的弧柱区. 其中, 电极边界层由一个厚度在10–6 m量级的空间电荷区(即鞘层)和一个厚度在10–3—10–4 m量级的电离区组成[21,22]; 弧柱区的空间尺度则通常在mm到cm量级, 甚至更大[23]. 已有研究表明, 即使在弧电流较大(如I > 100 A)、弧柱中心区域等离子体满足LTE-LCE的条件下, 一方面, 在阴极和阳极边界层区域, 由于固体冷壁面的作用, 等离子体往往会显著偏离局域热力学平衡(non-local thermodynamic equilibrium, nLTE)和/或偏离局域化学平衡(non-local chemical equilibrium, nLCE)状态, 甚至在鞘层区由于复杂的等离子体-壁面相互作用而导致带电粒子数密度不再满足准电中性(quasi-charge neutrality, QCN)条件、粒子速度/能量分布亦不再满足各向同性Maxwell分布[24]; 另一方面, 由于Maecker效应产生的引射作用使得环境冷气体被大量卷吸进入弧柱中心区[25], 于是在弧柱中心“热”平衡区和周围环境气体的“冷”平衡区之间也会产生处于nLTE-nLCE状态的非平衡过渡区[26]. 由此可见, 在如图1(a)所示的直流自由燃烧电弧等离子体中同时存在着热、电、化学非平衡以及粒子非Maxwell速度/能量分布四种非平衡效应, 而且在非平衡过渡区的电子温度和重粒子温度同时存在着显著的空间温度梯度[26], 这是明显不同于冷等离子体条件下较小的气体温度空间梯度的突出特点. 因此, 自由燃烧电弧作为一种典型的高气压等离子体非平衡输运模式体系, 不仅可用于研究热等离子体中复杂的非平衡协同输运机制, 而且对于其他类型等离子体源(如前述的冷等离子体和暖等离子体)的非平衡输运机制和特性的研究亦具有重要的参考价值.

      从实际应用来看, 电弧等离子体在诸多领域有着广泛的应用, 如微纳米材料制备、高危废弃物绿色处置、增材制造、临近空间飞行器鞘套环境地面模拟、空间电推进等, 如图1(b)所示. 随着电弧等离子体应用领域的不断拓展, 对于等离子体关键参数的调控需求亦不断提升. 例如, 在早期的等离子体焊接、切割、喷涂等材料加工领域, 其关键参数是电弧等离子体的气体温度(Th或Tg)或能量密度[27-33]; 在电推进领域则主要利用气体放电后加速喷出的等离子体羽流来提供维持飞行器运动所需的动力, 其关键参数是等离子体的质量平均速度($ {\boldsymbol{v}}_{0} $)[34]; 在高速飞行器再入大气的等离子体鞘套环境地面模拟中, 由于存在着等离子体与飞行器壁面间强烈的热-力-化学耦合效应, 涉及了气体温度(Th或Tg)、电子温度(Te)和密度(ne)、化学活性粒子浓度和壁面气动剪切力等多个等离子体关键参数[35,36]; 而在微纳米材料制备过程中, 既需要大体积、高温、高化学活性的等离子体区以实现材料的快速成核与团聚, 又需要气体温度陡降的非平衡区以避免颗粒尺寸的过度增长, 提高粒径均匀度, 相应地, 等离子体关键参数包括了气体温度(Th或Tg)、电子温度(Te)和密度(ne)、等离子体非平衡区体积和非平衡度的空间梯度等[37].

      由此可见, 建立高气压等离子体非平衡输运的模式体系, 从理论上深入研究电弧等离子体中质量-动量-能量非平衡协同输运机制, 并以实际应用需求为牵引, 探索电弧等离子体关键参数(如温度、活性粒子浓度、非平衡区尺度等)主动调控方法和技术, 对于丰富低温等离子体物理的基础理论, 推动现有等离子体技术的优化升级, 甚至催生一系列战略性、颠覆性的新兴应用, 具有重要的科学意义和工程应用价值. 因此, 本文将主要从数值模拟的角度综述直流电弧等离子体非平衡协同输运和关键参数主动调控的研究进展. 本文第2部分主要讨论直流电弧等离子体非平衡输运物理-数学模型的研究进展, 第3部分和第4部分主要讨论有关直流电弧等离子体非平衡输运和等离子体特性主动调控的研究进展, 第5部分针对上述两方面的研究进展简要讨论本领域未来发展需要解决的关键科学问题, 第6部分给出本综述的主要结论. 此处需要指出的是, 限于篇幅, 本文无意成为详尽的综述, 但我们尽可能基于直流电弧等离子体数值模拟方面的研究, 展现从机制解析到关键参数主动调控的研究思路的转变和新的研究进展.

    2.   直流电弧等离子体非平衡输运物理-数学模型研究进展

    • 通常来讲, 一套完整描述直流电弧等离子体体系中质量-动量-能量非平衡协同输运过程的物理-数学模型主要包括基本控制方程、等离子体输运性质参数和相应的所有待求解变量的边界条件. 随着人们对等离子体中传热、流动、电磁、化学反应耦合输运过程理解的不断深入, 数值模拟的研究对象从等离子体弧柱区扩展到包括固体电极区、电极边界层区, 乃至外电路; 相应地, 直流电弧等离子体物理-数学模型亦经历了从单温度平衡态(LTE-LCE)流体模型到双温度化学平衡(two temperature-LCE, 2T-LCE)、双温度化学非平衡或热力学-化学非平衡(2T-nLCE/nLTE-nLCE)流体模型, 再到考虑电极鞘层和外电路能量注入的非平衡耦合模型[38]. 我们认为, 用于电弧等离子体非平衡协同输运过程数值模拟的物理-数学模型应该满足如下三方面的要求.

      1)自洽性: 即描述等离子体中传热、流动、化学反应及电磁耦合过程的基本控制方程与所涉及的输运性质参数(如黏性系数、扩散系数、热导率、电导率等)在定义上应保持一致;

      2)一致性: 即在局域热力学和化学非平衡(nLTE-nLCE)条件下的基本控制方程和输运性质参数在满足局域化学平衡(2T-LCE)和局域热力学-化学平衡(LTE-LCE)条件下, 应能回归到相应状态下的基本控制方程表达形式和输运性质参数计算方法.

      3)守恒性: 即尽管由于电子质量远小于重粒子质量, 从而使得电子与重粒子间的碰撞过程对于重粒子子体系的动量和能量输运的改变非常小, 但仍须满足电子与重粒子间的质量输运(即扩散通量)守恒.

      基于上述三方面的考虑, 描述准电中性(QCN)等离子体区域非平衡协同输运的基本控制方程和输运性质参数通常从Boltzmann方程((1)式)出发进行理论推导(在本文以下的讨论中, 仅以单原子气体放电产生的等离子体为例进行讨论[39], 该思路也可以推广到考虑分子转动和振动等非平衡效应的情形[40]):

      式中,

      其中$ {\boldsymbol{c}}_{i} $和fi分别表示粒子i的相空间速度和速度分布函数; t为时间; mi和$ {\boldsymbol{X}}_{i} $分别表示粒子i的质量和受到的外力; g为粒子间的相对运动速率; σij为粒子i与j碰撞过程的微分散射截面, i或j = 1或e代表电子, i或j = 2, 3, ···, N代表重粒子, N为体系中所包含的粒子种类总数, Ω为立体角. 上标“′”表示碰撞后的对应物理量.

      等离子体输运性质参数的计算通常基于Chapman-Enskog方法, 通过对粒子速度分布函数进行微扰展开, 将扰动函数表示为包含温度梯度、压力梯度、速度梯度和浓度梯度的扩散驱动力的线性组合, 并通过Sonine多项式对扰动函数线性展开, 最终通过求解括号积分和碰撞积分得到输运性质参数的计算公式[41]. 在等离子体条件下, 根据电子子体系与重粒子子体系是否解耦, 可以将输运性质参数的计算方法分为完全解耦方法、完全耦合方法和半解耦方法. Devoto[42-46]最先基于等离子体中电子质量远远小于重粒子质量这一物理事实, 提出了完全解耦的计算方法, 即将电子子体系和重粒子子体系间的碰撞过程完全解耦, 从而极大地简化了对等离子体输运性质参数的求解, 并开展了大量有关平衡态等离子体和双温度等离子体输运参数的计算. 随后, Bonnefoi等[47]修正了Devoto所给出的双温度等离子体条件下有关扩散驱动力($ {\boldsymbol{d}}_{i} $)的定义, 使其满足$ \displaystyle\sum {\boldsymbol{d}}_{i}=0 $这一条件, 进一步完善了该完全解耦方法. 尽管基于完全解耦方法得到的大部分等离子体输运性质参数, 如热导率、电导率和黏性系数等与完全耦合方法的计算结果吻合良好[43,48], 但在Rat等[49]的研究中发现, 基于完全解耦方法给出的粒子扩散通量不满足质量通量守恒, 特别是随着等离子体温度的升高、电离度的增大, 粒子扩散通量偏离质量通量守恒的程度也愈加显著, 作者认为产生这一问题的主要原因在于完全解耦方法忽略了电子-重粒子间由于碰撞所产生的输运过程. 基于此, Rat等[49]认为即使在等离子体条件下, 依然应该采用与常温气体相同的完全耦合计算方法, 即在Boltzmann方程等号右端的碰撞项中完整保留电子与重粒子间的碰撞耦合项. 虽然采用完全耦合方法给出的等离子体中粒子扩散通量的守恒性得到了满足, 但计算过程复杂、耗时; 特别是Colombo等[48]通过计算分析再次表明, 除了组分扩散系数外, 由完全解耦方法和完全耦合方法计算得到的等离子体热导率、电导率和黏性系数的差别微乎其微.

      为了兼顾等离子体条件下组分扩散通量在物理上的守恒性和等离子体输运性质参数在计算上的复杂性, Zhang和Li等[50]在考虑电子质量远小于重粒子质量这一物理事实的同时, 考虑了电子和重粒子这两个子体系间的相互作用过程, 提出了针对碰撞项的半解耦处理方法, 即根据Chapman-Enskog微扰理论, 给出了如下一阶近似的Boltzmann方程:

      其中$\phi_i $为粒子i的扰动函数, 上标“0”代表速度分布函数的零阶近似解.

      一方面, 在双温度条件下, 等离子体中电子(i = 1)和重粒子(i ≥ 2)分别服从以电子温度(Te)和重粒子温度(Th)为特征温度的Maxwell速度分布, 即

      其中, nj和mj分别为粒子j (j =1, 2, ···, N)的数密度和质量; $ {\boldsymbol{C}}_{i} {\equiv} {\boldsymbol{c}}_{i}-{\boldsymbol{v}}_{0} $为粒子i的特定速度, $ {\boldsymbol{v}}_{0} $为等离子体的质量平均速度; kB = 1.381×10–23 J/K为Boltzmann常数.

      于是, 由$m_1/m_j \ll $1这一物理事实, 可简化得到

      另一方面, 对于(2)式等号右边的一阶扰动函数, 同样存在:

      这意味着在电子和重粒子碰撞前后, 重粒子一阶扰动函数的变化与电子相比可忽略不计. 于是, 对于等离子体中的电子子体系和重粒子子体系, (2)式可进一步简化为

      此处需要强调的是, 相较于Devoto[42-46]提出的完全解耦方法, Zhang等[50]提出的半解耦方法一定程度上在简化输运性质参数计算复杂性的同时, 包含了电子与重粒子两个子体系间的质量输运过程, 即通过保留(5b)式中等号右边的第一个碰撞积分项$ \displaystyle\iint f_{i}^{(0)}f_{1}^{(0)}\left({\phi }'_{1}-{\phi }_{1}\right)g{\sigma }_{i1}\text{d}{\boldsymbol{\varOmega }}\text{d}{\boldsymbol{c}}_{1} $使得整个等离子体体系中粒子扩散通量的守恒性得到了满足.

      于是, 根据输运通量矢量的定义, 粒子i的扩散速度可写为

      其中Ti为粒子i的温度($ {T}_{1}={T}_{\text{e}}, {T}_{i}{}_{=2, 3, \cdots , N}={T}_{\text{h}} $); Dij和$ D_i^T $分别为双温度等离子体的多组分寻常扩散系数和热扩散系数; $D_{ij}^{{\theta }} $和$ D_{i}^{{\theta }^{*}} $为相应的非平衡扩散系数; ω1和ωj分别定义为ω1 = x1p(1–x1)/D2和$ {\omega }_{j} = -{x}_{j}{x}_{1}p/D^2\;(j=2,3,\cdots,N) $; $ {\boldsymbol{d}}_{i} $为粒子i的扩散驱动力, 仍采用与Rat等[49]相同的表达式, 即

      其中ρi = mini为粒子i的质量密度, ρ $ =\displaystyle\sum\nolimits _{i=1}^{N}{\rho }_{i} $和nt $ =\displaystyle\sum\nolimits _{i=1}^{N}{n}_{i} $分别为等离子体总的质量密度和总数密度; θ =Te/Th为表征等离子体偏离局域热力学平衡程度的无量纲参数(即热非平衡度); xi = ni/nt为粒子i的摩尔分数; D = 1 + x1(θ – 1); p = $ \displaystyle\sum\nolimits_{i=1}^{N}{p}_{i}=\displaystyle\sum\nolimits _{i=1}^{N}\left({n}_{i}{k}_{\text{B}}{T}_{i}\right) $为等离子体的总压强, pi = nikBTi为粒子i的分压强. (6)式和(7)式在保证扩散驱动力守恒($ \displaystyle\sum\nolimits _{i=1}^{N}{\boldsymbol{d}}_{i}=0 $)的情况下考虑了等离子体中热非平衡度(θ)空间不均匀性对体系非平衡质能输运过程的影响[50].

      等离子体的黏性系数可由胁强张量计算得到, 即

      其中μ为等离子体的黏性系数, ${\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\leftrightarrow$}}\over{\tau}}}$为黏性应力张量, ${\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\leftrightarrow$}}\over{I}}} $为单位张量, ${\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\leftrightarrow$}}\over{S}}} $为应变率张量. 在等离子体条件下, 由于电子质量远小于重粒子质量, 故忽略了电子对黏性系数的贡献, 而仅考虑重粒子的贡献.

      在等离子体条件下, 热流密度矢量的表达式为

      其中λi和$ \lambda _{i}^{\theta } $分别为双温度等离子体条件下的平动热导率和相应的非平衡热导率; $ \varepsilon _{i}^{\text{r}} $和$ \varepsilon _{{i}}^{\text{int}} $则分别为粒子i的反应能和内能; $ \lambda _{ij}^{\text{D}} $为与扩散驱动力相关的导热率; Ejk为一矩阵的元素, 其具体表达式可参见Zhang等[50]的论文; 等号右边的前三项表征了由于等离子体中组分扩散引起的粒子焓输运过程, 后三项则分别表征由于重粒子温度的空间梯度、热非平衡度的空间梯度以及扩散驱动力引起的热流密度.

      在外加电场的作用下, 等离子体中的传导电流密度可写为

      其中e = 1.602×10–19 C为基本电荷量, Zi为粒子i净电荷量, σ和$ {\boldsymbol{E}}^{\text{ext}} $分别为等离子体电导率和外加电场强度. 在等离子体条件下, 同样由于电子质量远小于重粒子质量, 仅考虑了自由电子对传导电流的贡献.

      与上述等离子体输运性质参数的定义和计算方法相自洽, 描述等离子体非平衡协同输运过程的基本控制方程可做如下的推导. 若用ψi = mi, $ {m}_{i}{\boldsymbol{C}}_{i} $, $ {m}_{i}C_{i}^{2}/2+\varepsilon _{i}^{\text{r}}+\varepsilon _{{i}}^{\text{int}} $分别表示粒子i的质量、动量和能量, 将Boltzmann方程(1)两边同时乘以ψi, 并对其在整个速度空间积分, 即

      则可分别得到如下形式的粒子i的质量、动量和能量守恒方程:

      其中$ {\boldsymbol{v}}_{i} $为粒子i的运动速度, 且满足$ {\boldsymbol{v}}_{i}={\boldsymbol{V}}_{i}+{\boldsymbol{v}}_{0} $; $ {\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\leftrightarrow$}}\over{P}}}_{i} $, $ {\boldsymbol{q}}_{i} $, ${\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\leftrightarrow$}}\over{\tau}}}_{i} $分别为粒子i对应的胁强张量、热流密度、黏性应力张量; $ h_{i}^{*} = {h}_{i}/{m}_{i} = (5{k}_{\mathrm{B}}{T}_{i}/2+\varepsilon _{i}^{\text{r}}+\varepsilon _{{i}}^{\mathrm{int}})/m_i $为粒子i的比焓; $ {\varDelta }_{i}^{\text{el}}\left[{\psi }_{i}\right] $和$ {\varDelta }_{{i}}^{\text{inel}}\left[{\psi }_{i}\right] $分别为弹性碰撞和非弹性碰撞过程中粒子属性的变化; $ \left\langle {\psi }_{i}\right\rangle = \dfrac1{n_i} \displaystyle \int {\psi }_{i}{f}_{i}\text{d}{\boldsymbol{c}}_{i} $表示对参数ψi取平均值; $ \dfrac{\text{d}}{\text{d}t} = \dfrac{\partial }{\partial t} + {\boldsymbol{v}}_{0}\cdot \nabla $; (12)式中各项的物理含义可参见Li等[39]的论文.

      当存在外加电磁场时, 数值模拟还须考虑电磁场的时空演化, 即需要与上述方程(12a)—(12c)联立求解如下方程:

      其中φ和φext分别为总电势和外加电场的电势, 且φ = φint + φext, 此处φint为内电场的电势; $\boldsymbol E= -\nabla \varphi$为总电场; ρc为净电荷密度; ε0为真空介电常数; $ {\boldsymbol{j}}_{\text{tot}} $为总电流密度; $ {\boldsymbol{j}}_{\text{disp}} $为位移电流密度.

      当放电电流较大时, 需考虑电流流经等离子体产生的自感应磁场($ \boldsymbol{B} $)对等离子体特性的影响, 即在动量守恒方程((12b)式)右端的外力作用项中需考虑磁场产生的洛伦兹力. 在通常情况下, 可引入磁矢势$ \boldsymbol{A} $(即$ \boldsymbol{B}=\nabla \times \boldsymbol{A} $). 于是, 从Maxwell方程组出发, 可以得到如下形式的磁矢势方程:

      其中μ0为真空磁导率.

      方程(12)为等离子体处于完全非平衡(即热力学非平衡、化学非平衡和非电中性, nLTE-nLCE-nQCN或2T-nLCE-nQCN)条件下的控制方程的一般形式. 有关等离子体体系处于不同条件(如是否存在外加电磁场、是否满足电中性条件、是否具有宏观流动等)下的方程具体形式的讨论, 有兴趣的读者可参见Li等[39]的论文. 表1列出了当外力仅为电场力和磁场力时, 对应于不同平衡状态下的控制方程的具体形式. 针对表1中的方程, 此处做如下简要讨论:

      1)控制方程与输运性质参数的自洽性: 一方面, 等离子体的基本控制方程只有在给定对应的输运通量(如扩散速度、胁强张量、热流密度等)及输运性质参数表达式时, 才能满足计算上的封闭性; 而另一方面, 在推导控制方程和对应的输运性质参数时, 只有基于相同的基本假定, 才能保证数值模拟所采用的物理-数学模型的自洽性. 由上述讨论可以看到, 基本控制方程(12)及其中所涉及的输运通量和输运性质参数均基于Chapman-Enskog方法在相同的假设条件下推导得到, 守恒方程中的质量通量、动量通量和能量通量与对应的定义输运性质参数的三个通量完全一致, 从而保证了同一等离子体状态下控制方程和输运参数的自洽性. 另外, 值得注意的是, 热流密度表达式(9)中等号右边第六项相较于其他项可以忽略不计[42]; 前三项之和表征了粒子焓输运, 在LTE状态下, 这三项之和通常用来定义反应热导率, 但在双温度等离子体条件下, 选取电子温度还是重粒子温度来定义反应热导率仍有待商榷, 而采用(9)式则避免了定义双温度条件下的反应热导率. 类似地, 在描述非平衡态等离子体的能量守恒方程中, 选取温度为待求解变量, 也避免了双温度条件下定义定压比热时在温度选择上可能存在的困惑.

      2)不同状态下控制方程的一致性: 基于完全非平衡态(nLTE-nLCE-nQCN)等离子体条件下的控制方程(方程(12)), 在不同的假定下, 这套控制方程可逐步得到简化, 但不失一致性. 例如, 当满足准电中性(QCN)条件(如电极鞘层区以外的等离子体区)时, ρc = 0; 当满足局域化学平衡(LCE)条件(如大电流下的弧柱区等)时, $ {\varDelta }_{{i}}^{\text{inel}}\left[{m}_{i}\right]=0 $, 相应的等离子体组分则可由质量作用定律或自由能最小原理得到, 仅需求解总的质量守恒方程, 而无需求解每个粒子的质量守恒方程; 当等离子体处于完全平衡状态(LTE-LCE-QCN)时, 所有粒子的温度均相等(T = Ti, i = 1, 2, ···, N), 且有$ \displaystyle\sum \nolimits_{i=1}^{N}\nabla {\cdot} ({\rho }_{i}h_{i}^{*}{\boldsymbol{V}}_{i}) {=} -\nabla {\cdot} ({\lambda }_{\text{r}}\nabla T)-\frac{5}{2}\frac{{k}_{\text{B}}}{e}{\cdot} {\boldsymbol{j}}_{\text{cond}}{\cdot} \nabla T $, 于是, 将所有粒子的能量方程相加, 即可得到完全平衡态下的能量守恒方程, 相应地, 等离子体的热导率λ = λe + λh + λr (此处λe, λh, λr分别为电子热导率、重粒子热导率和反应热导率), 控制方程也回归到文献中普遍采用的完全平衡态下的方程形式, 且在光学薄假设下, $ Q_{\text{e}}^{\text{inel}}=-{U}_{\text{rad}} $[51].

      3)控制方程中扩散通量的守恒性: 相较于Devote提出的完全解耦方法会导致粒子扩散通量不满足质量通量守恒的要求, 半解耦方法考虑了电子与重粒子两个子体系间的质量输运过程, 从而在兼顾计算复杂性的同时, 保证了等离子体体系中粒子扩散通量的守恒性.

      4)控制方程在数值模拟中的简化形式: 表1列出了不同状态下等离子体控制方程的完整形式, 但从数值模拟的角度出发, 上述方程在满足一定条件时, 仍可进一步简化, 从而减小实际数值计算的复杂性、降低计算量. 例如, 对于低速流动的等离子体, $ \left({\rho }_{\text{e}}{\boldsymbol{V}}_{\text{e}}\right)\text{}\cdot \text{d}{\boldsymbol{v}}_{0}/\text{d}t $, $ \displaystyle\sum\nolimits _{i=2}^{N}\left({\rho }_{i}{\boldsymbol{V}}_{i}\right)\text{}\cdot \text{d}{\boldsymbol{v}}_{0}/\text{d}t $, $ \text{d}{p}_{\text{e}}/\text{d}t $, $ \displaystyle\sum\nolimits _{i=2}^{N}\text{d}{p}_{i}/\text{d}t $, ${\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\leftrightarrow$}}\over{\tau}}}_{\text{e}}\colon \nabla {\boldsymbol{v}}_{0} $, $ \displaystyle\sum\nolimits _{i=2}^{N} {\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\leftrightarrow$}}\over{\tau}}}_{i}\colon \nabla {\boldsymbol{v}}_{0} $均可忽略; 若不考虑自感应磁场诱导的电流对体系能量注入的贡献, 则

      其中, $ {\boldsymbol{j}}_{\text{e}} $和$ {\boldsymbol{j}}_{\text{ion}} $分别为电子电流和离子电流. 于是, 双温度化学平衡(2T-LCE)等离子体条件下的能量方程可简化为如下形式[39]:

    3.   电弧等离子体体系质量-动量-能量协同输运研究进展

    • 基于等离子体流体模型, 国内外学者开展了大量的电弧等离子体体系质量-动量-能量协同输运的数值模拟研究, Fang等[38]较为系统地总结了目前电弧等离子体数值模拟研究所采用的8类模型, 如图2所示. 其中, 模型1, 2, 5基于局域热力学平衡(LTE)或非平衡(2T-LCE, 2T-nLCE)假设, 通过在电极鞘层外边界人为给定边界条件的方法封闭方程组, 因而无法考虑电弧弧柱区与固体电极区以及电极鞘层间的耦合作用, 且准确边界条件的给定通常比较困难; 模型3, 4, 6通过对电极鞘层采取简化处理方式保证电流导通, 实现了电弧弧柱区与固体电极区之间的耦合处理, 但依然无法揭示电极鞘层对等离子体体系能量输运的影响; 为研究阴极鞘层对弧柱区等离子体行为的影响, 模型7在简化阳极鞘层的基础上引入了阴极鞘层模型, 从而可以研究阴极鞘层对阴极固体区及电弧弧柱区的影响; 而模型8则为考虑了阴、阳极固体区和相应的两个电极鞘层区, 以及弧柱区的完全非平衡耦合模型, 能够用于自洽分析等离子体体系中的质量-动量-能量非平衡协同输运机制.

      • 3.1.   平衡态电弧等离子体物理特性研究进展

    • 在等离子体焊接、切割、喷涂等应用领域, 电弧等离子体的研究主要集中于高温高能量密度的核心区, 而对电极边界层及电弧边缘过渡区的关注相对较少. 由于在等离子体弧柱区温度较高, 电子与重粒子间的频繁碰撞不仅使得这两个子体系间的能量交换达到近平衡态(Th ≈ Te), 而且化学反应亦接近于动态平衡; 也就是说, 在电弧弧柱区的等离子体通常满足LTE, LCE及QCN假设. 因此, 可采用表1中的平衡态控制方程, 并结合相应的边界条件开展弧柱区传热与流动特性的数值模拟研究(如在不考虑固体电极区条件下, 给定电极鞘层外边界的温度和速度边界条件[17,52]; 或在包含固体电极区条件下, 采取将靠近电极壁面的等离子体区第一个计算节点布置于电极鞘层外的近似处理方法[53,54]等). 相应地, 数值模拟主要采用图2中的模型1—4.

      在平衡态等离子体输运性质参数的计算方面, Devoto[42]研究了Chapman-Enskog方法求解等离子体输运参数的收敛性, 不同于普通气体只需要一阶Sonine多项式展开就已满足精确要求, 等离子体条件下黏性系数和电导率通常需要展开到二阶, 热导率和扩散系数则需要展开到三阶或更高阶. 以Devoto[42–46]的研究工作为基础, 美国的Pfender研究组、澳大利亚的Murphy研究组、意大利的Capitelli研究组、法国的Fauchais研究组以及中国的Rong和Wang等研究组在典型的电弧等离子体温度范围内计算了Ar[55,56], N2[55–57], O2[55], 空气[56,58,59,60], Li[61], Ar/Fe[62], Ar/Cu[62], Ar/Al[62], SF6/Cu[63]等多种平衡态等离子体的输运性质参数, 并建立了相应的数据表, 为平衡态条件下等离子体输运特性的数值模拟提供了重要的数据支撑.

      基于LTE-LCE假定, Hsu等[51]成功预测了自由燃烧电弧弧柱区的铃铛形气体温度分布和阴极上游冷气体的卷吸效应, 与实验测量结果吻合良好, 从而验证了平衡态流体模型在预测电弧弧柱区传热与流动特性的可靠性(如图3所示). 由于电弧等离子体的温度场和流场分布会直接影响等离子体材料加工的工艺效果, 因此, 基于平衡态流体模型, 国内外研究者就等离子体工作气压、放电电流、电极结构与间距, 以及工作气体化学成分、注入位置、流速等多种因素对电弧等离子体的物理特性(包括温度、壁面热流密度、粒子通量、电弧形貌与空间尺寸、三维效应等)的影响机制和规律开展了系统的数值模拟研究.

      1)工作气压和弧电流的影响: Hsu和Pfender[64]的研究表明, 电弧高温区会随着工作气压的升高而收缩, 从而导致焓通量的增大和阴极射流速度的降低. Wang等[65]进行了大气压和减压条件下非转移氩电弧等离子体发生器流动和传热特性的数值模拟和实验研究, 发现尽管1和0.1 atm (1 atm = 105 Pa)条件下的电弧流场和温度场分布整体相似, 但在0.1 atm条件下的流场很难达到充分发展, 而弧电流的升高会导致出口处电弧的最高温度和轴向速度均增大. Kovitya和Lowke[52]研究了不同弧电流下氩自由燃烧电弧动量和能量的输运特性, 结果表明, 随着弧电流的增大, 电弧中心区温度会逐步升高, 电弧体积扩展, 轴向电场增强; 而减小电流密度则会显著降低近阴极区域的等离子体温度和轴向速度, 但对远离阴极区等离子体的影响则很小; 等离子体温度和速度的计算结果与实验测量值吻合良好.

      2)电极结构与间距的影响: Siddanathi等[66]研究了锥形、渐缩形、平端、圆柱形阴极及阶梯式阳极结构对非转移弧等离子体发生器中传热与流动特性的影响, 结果表明, 为了实现恒定功率下发生器的稳定运行, 不同的阴极结构需要适配不同的弧电流水平, 且平端阴极表现出最高的热效率, 而阳极形状的影响则微乎其微. Nishiyama等[67]对阴极尖端角度的影响进行了数值模拟研究, 发现当阴极尖端角度为45°时电流分布最集中, 等离子体温度亦最高. Cao等[23]则对长达200 mm的氩转移电弧特性进行了数值模拟和实验研究, 结果表明, 电极间距主要影响阳极附近的温度场分布, 而对阴极附近的等离子体温度场分布则影响很小, 且长转移弧等离子体具有抛物线型的流动特征, 温度分布的数值模拟结果与实验测量结果吻合良好. Wu等[68]则研究了棒状阴、阳极间形成的直流电弧等离子体, 发现在发生器阴、阳极之间会形成轴向速度为零的滞止层, 且滞止点处的气体温度会随着弧电流的增大或电极间距的减小而升高.

      3)工作气体化学成分和流速的影响: Lowke等[54]研究了氩气和氩-氢混合气体条件下的电弧特性, 发现氢气的加入会导致气体温度最大值和弧电压的提高, 并使得等离子体高温区收缩. 此外, Villareal-Medina等[29]的研究表明, 尽管不同气体成分不会对大气压下焊接电弧中焦耳热、对流、辐射和导热等主要热传递方式产生影响, 但不同气体放电产生的等离子体的热导率、电导率和黏性系数等输运性质参数明显不同, 从而会对具体的能量输运特性产生较大的影响. Wang等[53]比较了氢气、氮气和氩气三种不同工作气体条件下的低功率电弧加热推进器中的等离子体流动和传热特性, 计算结果表明, 同其他两种气体相比, 氢气条件下的弧柱核心区等离子体温度最高, 且轴向速度最大, 但推进器出口处的温度则最低. Murphy等[28]研究发现, 焊接电弧中氦气、氢气和氮气等屏蔽气的引入会提高阳极处的热流和电流密度, 从而导致阳极及其附近区域的等离子体温度升高、工件的熔池深度增大. Long等[69]的研究结果则表明, 在氧气等离子体切割电弧中, 气体流量的增加会引起阴极附近电流密度与等离子体温度分布特性的显著变化, 其径向最大值的空间位置会随着气流的增大而向外侧移动; 同时, 气流的增强会对电弧等离子体产生更强的压缩效应, 导致电流密度峰值与电弧温度最大值上升, 加剧了阴极材料的蒸发, 从而导致阴极的总质量损失增大.

      4)电极金属蒸气的影响: Xiang等[70]研究了铁金属蒸气对氩电弧等离子体输运特性的影响, 结果表明, 在400 A的弧电流下, 虽然金属蒸气主要存在于阳极上方和阴极附近区域, 在弧柱区浓度则相对较低, 但仍会导致电弧温度和电流密度的显著降低, 从而造成熔池深度和体积的显著减小. Wang等[71]进一步探究了含氧量对氩电弧焊过程中铁蒸气的影响, 发现铁蒸气仅局限在靠近熔池的区域, 且会导致纯氩气电弧焊熔池表面的电流密度和热通量在轴线上减小并径向扩展, 而在加入0.5%的氧气后, 铁蒸气反而会造成轴线上的电流密度和热通量的增大, 并导致电弧的径向收缩. Gonzalez等[72]的研究则表明, 自由燃烧电弧中由于阳极蒸发产生的金属蒸气会显著地增大等离子体的电导率和辐射损失功率, 导致阳极附近等离子体温度降低, 进而又会反过来影响阳极材料的蒸发速率, 使得阳极附近的金属蒸气浓度降低, 这也意味着, 如果不考虑阳极材料的蒸发, 数值模拟可能会过高预测阳极附近的等离子体温度.

      5)等离子体发生器内传热与流动的三维效应: Chen和Li等针对直流非转移电弧等离子体发生器及其下游射流区的三维效应进行了较为系统的数值模拟研究[73-79]. 对于采用切向旋流进气的、具有完全轴对称几何结构的直流非转移电弧等离子体发生器, Li等通过给定阳极外壁面稳态传热的边界条件, 并应用Steenbeck弧压最小原理确定了阳极弧根的轴向和周向位置及电弧的三维形状, 自洽预测了发生器内部阳极弧根局部贴附的三维效应, 数值模拟结果与实验观察定性一致[76-79]. 而由于等离子体发生器出口电弧温度场和速度场的三维分布, 以及实际的等离子体喷涂等应用中载气的侧向注入, 在等离子体射流区同样会存在显著的三维效应, 从而影响随载气注入的颗粒的运动与加热历程, 进而影响涂层质量[74,75].

      此处需要说明的是, 尽管平衡态流体模型能够在整体上合理预测电弧等离子体体系中的传热与流动特性, 但由于在电极冷壁面附近以及高温等离子体与环境冷气体交界面附近存在着显著的非平衡效应, 从而使得数值模拟的计算结果与实验测量结果之间存在一定的差异, 例如, 采用平衡态流体模型计算得到的弧电压值通常会明显偏离实验测量值, 其中一个重要的原因即为平衡态数值模拟忽略了电极边界层区域的电位降, 计算得到的弧电压主要是电弧弧柱区的电位降[80]. 因此, Li等[79]指出, 为了准确地获得电极附近的等离子体参数, 需要在数值模拟中进一步考虑电极附近等离子体的非平衡效应和弧根运动导致的非稳态效应. Huang等[81]发展了三维非稳态流体模型, 研究了直流非转移弧等离子体发生器内部的传热与流动过程, 成功地预测了频率为11 kHz的阳极弧根波动, 但由于平衡态流体模型得到的电极附近等离子体温度低、电导率小, 从而导致计算得到的弧电压较实验测量值偏高. 于是, Huang等[82]进一步采用改进的平衡态流体模型计算得到标称电子温度, 从而在电极附近区域获得相对LTE模型更高且更加贴近实际的等离子体电导率, 进而大幅提高了电弧等离子体发生器内部等离子体温度场和速度场的数值模拟预测精度. 这也进一步说明了开展非平衡态电弧等离子体数值模拟的必要性.

      • 3.2.   非平衡态电弧等离子体物理特性研究进展

    • 随着微纳米材料制备与表面改性、临近空间飞行器鞘套环境地面模拟等新兴应用领域的快速发展, 电弧数值模拟的研究区域不只局限于弧柱区, 而是不断向非平衡的过渡区和电极边界层区域发展. 一方面, 由于等离子体与环境冷气体或电极冷壁面间复杂的相互作用, 这些区域中的等离子体往往不再处于LTE状态, 甚至不满足LCE、QCN假设, 而是表现出明显的非平衡特性. 另一方面, 电弧等离子体的核心参数, 如气体温度(Th或Tg)、电子温度(Te)和数密度(ne)、等离子体非平衡区体积和非平衡度的空间梯度等, 不仅与等离子体体系内的质量-动量-能量非平衡输运机制密切相关, 而且会直接影响等离子体在微纳米材料合成等诸多领域的应用效果. 因此, 在过去的几十年里, 国内外研究者对电弧等离子体的非平衡特性亦开展了大量的数值模拟研究工作. 针对不同区域等离子体的不同非平衡特性, 可采用表1中对应的非平衡态控制方程, 并结合相应的边界条件开展等离子体非平衡协同输运机制与特性的数值模拟研究, 其主要采用的数值模拟模型为图2中的模型5—8.

      与平衡态等离子体数值模拟类似, 非平衡态等离子体数值模拟同样需要完备的输运性质参数. 在双温度等离子体输运性质参数的计算方面, Devoto[42-46]首先将基于电子与重粒子子体系完全解耦的计算方法, 从LTE等离子体推广到了双温度等离子体. 之后, 美国的Pfender研究组、法国的Fauchais研究组、澳大利亚的Murphy研究组、意大利的Colombo研究组以及中国的荣命哲、王海兴、李和平、夏维东研究组等相继开展了Ar[83,84], O2[85], He[86], SF6[87], Ar/N2[88], N2/O2[89], C/O[90], C4F7N/CO2/O2[91]和SF6/Cu[92]等多种等离子体在双温度条件下的输运性质参数计算, 建立了大量数据表, 为双温度条件下非平衡态等离子体的数值模拟研究提供了基础数据.

      • 3.2.1.   等离子体弧柱区2T-LCE物理特性研究

    • 基于2T-LCE模型, Chen等[18]针对直流转移电弧约束型阳极区的传热与流动特性进行了二维数值模拟研究. 计算结果表明, 在温度超过12000 K的电弧等离子体核心区, 电子温度与重粒子温度基本相等; 而在电弧边缘和靠近阳极壁面附近的区域, 等离子体则显著偏离LTE状态, 电子温度与重粒子温度差异较大. 这一计算结果也在很大程度上表明, 双温度流体模型在精确预测电弧传热与流动特性方面的必要性. Hsu和Pfender[17]也在平衡态电弧特性研究的基础上, 进一步研究了大气压氩自由燃烧电弧的非平衡特性. 从图4(a)给出的电子和重粒子温度的二维空间分布数值模拟结果可以看到, 在等离子体与冷气体相互作用的电弧边缘区域, 电子温度明显高于重粒子温度, 等离子体显著偏离LTE状态. Trelles等[19]分别使用LTE-LCE和2T-LCE模型对直流非转移电弧等离子体发生器内部的三维、非定常传热与流动过程进行了数值模拟研究, 结果表明, 在等离子体与冷气流相互作用区域以及阳极壁面附近均存在着显著偏离平衡态的区域; 与LTE-LCE模型不同, 2T-LCE模型无需额外的再附着模型, 即可自洽地再现电弧再附着过程, 并自动限制总的弧电压与电弧长度的量级, 且2T-LCE模型的计算结果与实验测量结果的吻合度更高. Bhigamudre和Trelles[93]以雷诺数和焓数这两个无量纲参数为控制变量, 对有横向冷气体注入条件下的电弧特性进行了系统的研究, 揭示了在工作参数变化时, 由于冷气体与高温电弧间的强烈相互作用所导致的电弧形态、弧电压降、电弧功率、非平衡程度以及等离子体迎向来流方向的电弧前缘特征厚度等的变化规律. Charrada等[94]研究了不同弧电流和工作气压对汞电弧灯非平衡特性的影响, 结果表明, 电弧中心高温区总是较电弧边缘区更接近于热力学平衡态, 且随着弧电流或工作气压的降低, 等离子体会进一步偏离LTE状态. 而Haidar[95]在自由燃烧弧非平衡数值模拟的研究中发现, 由于在阴极尖端附近大量冷气体被卷吸进入电弧核心区, 从而导致该区域电子温度显著高于重粒子温度.

      尽管诸多数值模拟研究工作均揭示了不同参数条件下电弧等离子体在靠近电极冷壁面和环境冷气体区显著偏离LTE状态的特性, 但Haidar[95]的数值模拟研究工作亦指出, 如果仅采用2T-LCE模型, 而不考虑可能存在的化学非平衡效应, 沿电弧径向即使在物理空间已经处于冷气体区, 计算得到的电子温度仍然显著高于重粒子温度, 这在物理上显然是不合理的. 为此, 当气体温度低于5000 K时, Haidar通过人为增大电子与重粒子间的碰撞频率以提高两者间的能量交换, 从而使得在冷气体区电子温度尽可能接近重粒子温度, 如图4(b)所示. 类似的现象, 在有关直流非转移电弧双温度数值模拟的计算结果中也有报道[77]. 由此可见, 对于高温电弧等离子体与环境冷气体强烈相互作用的过渡区域, 化学非平衡效应可能对电弧特性具有显著的影响, 从而使得发展考虑局域化学非平衡效应的2T-nLCE(或nLTE-nLCE)模型变得十分必要.

      • 3.2.2.   等离子体弧柱区2T-nLCE物理特性研究

    • 为了能够自洽地预测沿电弧径向位于弧柱中心区和环境冷气体间非平衡过渡区的等离子体特性, Guo等[26]采用2T-nLCE模型开展了氩自由燃烧电弧非平衡特性的二维数值模拟研究, 特别是与相同条件下2T-LCE模型的计算结果进行了对比. 如图5(a)所示, 在弧电流为100 A的条件下, 若采用2T-LCE模型, 则计算得到的电子温度即使在冷气体区(如$ {T}_{\text{h}}\ll 5000\;\text{K} $@ r > 20.0 mm)也会远高于重粒子温度, 这在物理上显然是不合理的; 而采用2T-nLCE模型则可自洽地预测出等离子体弧柱区沿径向的非平衡特性变化规律, 即电弧中心“热”平衡区-非平衡过渡区-环境“冷”平衡区, 如图5(b)所示. Guo等[26]进一步分析了自由燃烧电弧中热力学与化学非平衡效应的协同作用, 即随着远离弧柱中心区、等离子体向外扩张的过程中, 一方面, 由于电流密度以及相应的焦耳热效应的降低, 在电子温度、重粒子温度及电子数密度均降低的同时, 电子与重粒子间的碰撞频率亦降低, 从而导致热力学非平衡度(θ)急剧增大, 等离子体显著偏离LTE状态; 而另一方面, 当Te降低到大约5000 K以下时, 电子与离子复合形成原子的反应速率相对于电离速率显著提高, 并在过渡区明显超过电离反应速率; 此时, 非弹性碰撞能量交换过程在控制子体系之间能量再分配过程中起着决定性的作用, 并与弹性碰撞能量交换过程一起, 共同驱动体系的Te很快趋向于Th; 相应地, 在电弧边界处, 冷气体与等离子体之间的强烈作用导致等离子体逐渐消失, 转变为平衡态冷气体, 而热非平衡度(θ)在这个过程中亦重新回到1. 而Fang等[38]则在此基础上进一步研究了Ar的三体复合反应对电弧非平衡特性的影响, 发现该反应在远离电弧几何轴线的低电子数密度区域主导非弹性碰撞过程, 显著地影响电子温度和数密度的空间分布, 以及热非平衡和化学非平衡区的宽度. 在Ghorui等[96]和Baeva等[97,98]的研究中亦指出, 等离子体密度梯度的变化对化学非平衡效应的影响非常显著, 在电弧收缩或扩张显著的区域(如电极附近、电弧边缘等), 等离子体偏离LCE的程度更加显著; 甚至即使在满足LTE状态的区域也可能会偏离LCE状态.

      另外, 等离子体体系中的化学反应动力学过程亦会对数值模拟所预测的电弧特性产生重要的影响. 例如, Wang等[99]采用2T-nLCE模型研究了电弧推力器中氢-氦混合气体放电时的粒子空间分布, 结果表明, 非平衡化学动力学过程导致的粒子分离主要发生在氮气和氢气解离和电离的区域, 且热扩散对电弧射流推进器内粒子扩散的影响比压力扩散大得多; Wei等[100]和Wang等[101]对高速等离子体电弧射流的数值模拟研究发现, 虽然激发态氩原子(Ar*)数密度远小于电子、氩离子等其他粒子的数密度, 但Ar*对整个等离子体体系的电离、复合反应过程以及电弧阳极弧根贴附模式起着重要作用, 如图6所示. Baeva和Uhrlandt[102], Golubovskii等[103,104]的研究亦指出, 数值模拟所采用的化学反应动力学过程会对计算得到的等离子体特性产生明显的影响, 如激发态粒子的加入会影响电弧的共振辐射. 由此可见, 在2T-nLCE数值模拟中, 除了描述非平衡态等离子体的控制方程和输运性质参数外, 化学反应动力学过程同样对数值模拟结果的准确性起着重要的作用[105].

      • 3.2.3.   考虑多区耦合效应的非平衡态等离子体物理特性研究

    • 本文3.1, 3.2.1和3.2.2节有关平衡态和非平衡态电弧等离子体的数值模拟研究, 主要关注电弧弧柱区的传热、流动和非平衡特性, 而忽略了弧柱区与固体电极区之间的相互作用. 一方面通过给定电极鞘层外缘的电势(φ)、电流密度(j)、温度(Te和Th)及粒子数密度(ni)等参数的边界条件, 使得数值模拟的难度和计算量相对较低; 而另一方面, 也导致边界条件的给定成为影响数值计算准确度的一个重要因素[106-111]. 于是, 拓展计算域成为解决这一问题的一个主要途径. 从图1可以看到, 电极鞘层区是联系电弧弧柱区和固体电极区的一个重要的物理区域. 尽管在空间尺度上, 电弧等离子体鞘层的特征尺度通常在10–6 m量级, 远小于mm甚至cm量级的电极间距, 但在鞘层区存在着温度、电势、电子和离子数密度等参数的巨大梯度, 成为影响电弧等离子体体系能量输运的一个关键区域[111-115]. 从目前发表的研究结果来看, 考虑电弧等离子体体系多区耦合效应的数值模拟主要包括如下几种方法:

      1)不考虑电极鞘层的弧柱区-固体电极区耦合模型: 即计算域包含电弧弧柱区和电极固体区, 但不包括电极鞘层区, 而仅采用一定的简化处理方法保证等离子体-电极交界面处的电流导通, 如本文3.1节提到的在平衡态流体模型中采取将靠近电极壁面的等离子体区第一个计算节点布置于电极鞘层区外的近似处理方法[53,54], 或采用标称电子温度计算等离子体电导率[82]的近似处理方法等.

      2)考虑阴极鞘层的弧柱区-阴极区耦合模型: 由于大气压条件下电弧等离子体的鞘层厚度通常小于粒子的平均自由程, 故可引入无碰撞鞘层模型[114], 并把鞘层看作一个用于在电弧弧柱区和固体电极区之间进行信息交换的、不计厚度的界面. 基于此, Li和Benilov等[80,116]提出了采用表面非线性加热模型的单向耦合模型, 即首先根据电极几何结构和弧电流, 采用表面非线性加热模型计算得到阴极固体区温度和电流密度分布、阴极鞘层电位降, 以及沿阴极鞘层-弧柱交界面处的电子温度、重粒子温度和电流密度分布, 然后以此作为弧柱区数值模拟的边界条件, 计算得到弧柱区等离子体参数的空间分布. 计算结果表明, 电弧等离子体中鞘层电压降贡献了整个弧电压的2/3甚至更多[80], 单向耦合模型计算得到的电弧伏安特性曲线与实验测量结果吻合良好; 从能量输运的角度看, 注入阴极鞘层的能量不仅用于加热阴极固体区, 而且用于加热电弧弧柱区, 这一能量传递机制是采用忽略鞘层的流体模型所无法揭示的. 在单向耦合模型(即仅考虑阴极鞘层对弧柱的影响)基础上, Chen等[117]提出了阴极区双向耦合模型(即考虑阴极鞘层与弧柱的相互作用), 通过在弧柱边缘计算热扰动层实现鞘层/预鞘层与弧柱的耦合, 研究表明, 热电子发射是阴极壁面的主要冷却机制; 阴极弧根贴附存在点模式(温度峰值显著)与扩散模式(温度分布均匀)两种形式. 进一步地, Chen及其合作者[118,119]评估了等离子体区对阴极加热的影响, 并对比了单向与双向耦合模型在计算结果上的差异, 发现双向耦合模型计算得到的弧电压值更接近于实验测量值.

      3)多区耦合模型: 基于无碰撞鞘层模型, Fang等[38]提出了包括等离子体弧柱区、2个电极鞘层区和2个固体电极区的多区域耦合模型; Zhang等[120]则在此基础上, 进一步耦合了外电路模块, 通过自洽求解包括质量、动量、能量守恒方程和电流连续性方程的非平衡流体模型以及电路模型, 研究了电弧等离子体体系中的非平衡输运机制、放电模式及其模式转换, 如图7所示, 其中变量的物理含义和表达式可参见文献[120]. 基于多区域耦合模型, 结合分析碰撞等离子体中非平衡协同输运的“能量树”概念[121], 数值计算结果表明: 1)多区耦合模型能够自洽地预测阴极和阳极鞘层区特性, 并在不同的弧电流下, 预测了放电从正常辉光、反常辉光到过渡模式, 再到电弧模式的转变; 总的鞘层电压降贡献了约70.6%的弧电压, 而阴极鞘层电压降则占总鞘层电压降的约84%, 这说明阴极鞘层在能量沉积中起关键作用; 而进一步的研究表明, 阴极鞘层区注入的功率密度(Pc)可作为不同放电模式转变的判据, 且对于氩气放电的辉光-电弧模式转变的临界窗口为Pc= (0.28—1.2) × 1012 W·m–3[38,120]. 2)对流换热过程会显著地影响等离子体体系的质量-动量-能量交换, 对流换热效应的强化会导致等离子体的非平衡度提高; 重粒子子体系的对流换热过程是其内部主要的能量传递过程及其向环境的能量耗散过程之一, 它的改变直接影响到重粒子温度和体系的热非平衡度的空间分布[38,120]. 3)工作压力降低会导致粒子平均自由程增大, 电子与重粒子间的弹性碰撞能量交换减弱, 同时会使得等离子体速度提高, 从而也会强化对流换热过程, 导致等离子体热非平衡度(θ)的增大和非平衡区域的扩展[38]. 上述有关电弧等离子体体系中非平衡质能协同输运机制的分析也为等离子体关键参数的主动调控提供了新的思路[36,122].

      尽管无碰撞鞘层模型能够在一定程度上反映鞘层行为对电弧特性的影响, 但却无法获得鞘层厚度以及鞘层区等离子体参数的空间分布, 亦无法完全揭示鞘层区与电极固体区以及弧柱区之间的复杂相互作用. 因此, 近几年来, 一些研究者开始尝试采用粒子模拟的方法深入地研究电弧等离子体的电极鞘层特性, 并针对大气压条件下的微放电开展了电弧特性的数值模拟研究. 例如, Zhou等[24]采用粒子模拟结合蒙特卡罗碰撞方法(particle-in-cell Monte Carlo collision, PIC-MCC)预测了直流电弧等离子体阴极鞘层的厚度和阴极鞘层区电子和离子的数密度、电场强度的空间分布, 以及空间平均的电子能量分布函数, 揭示了热阴极表面的热电子发射是维持电弧放电的主要电子发射过程, 这与Zhang等[120]给出的电弧放电模式下的阴极表面电流成分的计算结果是一致的. Sun等[123,124]基于该方法进一步开展了大气压条件下微电弧放电的数值模拟, 研究了热-场致发射驱动的微等离子体的击穿后特性, 发现热-场致发射驱动的微尺度放电可以在低至20 V的外加电压下产生稳定的微电弧放电等离子体; 与外加电压和电极间距相比, 阴极温度对热-场致发射的影响最为显著, 并导致电子能量分布函数呈现出多峰分布的特点[123]; 鞘层会显著地影响近阴极区的粒子输运以及弧柱向热阴极的能量输运, 焦耳热与阴极表面带电粒子沉积是近阴极区向阴极传递能量的主导机制; 相比自然对流, 阴极表面辐射热损失占主导地位[124].

    4.   电弧等离子体特性主动调控研究进展

    • 随着等离子体应用领域的不断拓展, 对等离子体关键参数(包括等离子体气体温度、电子温度、活性粒子数密度、粒子通量或能量通量、等离子体区体积等)的主动调控日益受到国内外研究者的广泛关注. 数值模拟为深入理解电弧等离子体中复杂的非平衡质能协同输运机制提供了重要的手段, 这不仅为针对已有电弧等离子体发生器中的质量-动量-能量非平衡协同输运机制的分析提供了一个强有力的工具, 而更为重要的是, 对数值模拟所揭示的质能协同输运机制的深入分析, 例如采用“能量树”概念[121]的等离子体非平衡协同输运机制的分析, 为主动调控电弧等离子体的结构参数、运行参数和/或环境参数, 从而创造全新的、满足特殊应用需求的等离子体源提供了重要的理论指导. 我们认为, 后者对于推动包括电弧等离子体在内的各类碰撞等离子体源的应用基础研究具有重要的意义. 在本节中, 我们着重讨论在深入分析平衡或非平衡数值模拟结果基础上的等离子体关键参数主动调控的研究进展, 同时, 也对仅从实验角度开展的一些有关等离子体关键参数主动调控的研究工作进行了简要的讨论.

      • 4.1.   基于数值模拟的关键参数主动调控研究进展

    • 在等离子体焊接中, 电弧等离子体形貌是影响焊缝质量和焊接效率的重要因素之一. Nomura等[125,126]基于局域热力学平衡(LTE)模型开展了外加磁场对焊接电弧等离子体形貌影响的数值模拟研究, 为实验工作提供了重要的理论指导. 如图8所示, 数值模拟结果表明, 无外加磁场时, 电弧关于发生器几何轴线呈对称的铃铛形、而在垂直于几何轴线的横截面上则为圆形; 当施加尖峰型磁场后, 因外加磁场与电弧自感应磁场在垂直方向叠加增强、水平方向抵消减弱, 从而导致电磁收缩力偏移, 电弧横截面转变为椭圆形. 基于此, Nomura等进一步给出了磁化方向对电弧横截面椭圆化效果的调控规律: 水平磁化对弧柱近电极侧椭圆化作用显著, 而垂直磁化则对近母材表面侧作用更强且椭圆化程度更大. 由于母材侧等离子体椭圆化会使得热量输入范围沿焊缝方向扩展, 可更好地防止电弧后方区域冷却过快的问题, 从而延缓熔池金属的凝固速度. 此外, 计算结果还表明, 电弧横截面椭圆化效果越强, 电弧作用在液态熔池表面的压力越低, 从而有利于抑制液态金属飞溅、降低驼峰缺陷的发生概率. 基于上述数值模拟结果, Nomura等确定了实验中永磁体的最优磁化方向为垂直磁化, 实验结果表明: 无永磁体时, 在100 cm/min焊接速度时会出现焊缝驼峰缺陷, 水平磁化使无驼峰缺陷焊接速度提升至140 cm/min, 而垂直磁化则将该速度上限进一步提高至180 cm/min, 这为在高速电弧焊接条件下实现高质量焊缝成形提供了重要的理论指导[126].

      在等离子体材料加工领域, 等离子体发生器的使用寿命、放电稳定性直接影响材料加工效果, 而阳极材料烧蚀是制约其性能提升的关键因素之一. Zhukovski等[127]针对该问题, 基于局域热力学平衡(LTE)假设建立了三维非稳态直流放电等离子电弧射流物理数学模型, 该模型在25º旋流进气、无外加磁场工况下的数值预测结果与实验测量结果吻合良好. 在此基础上, Zhukovski等研究了旋流气体注入与外加轴向磁场耦合作用对阳极弧根运动的调控机制[128], 结果发现(如图9所示): 在不破坏喷枪出口等离子体射流稳定性的前提下, 施加0.05—0.1 T的外加磁场或采用45º旋流气体注入, 可有效驱动弧根旋转, 实现阳极热负荷的均匀分散与阳极表面最高温度的降低, 从而有效地抑制了阳极材料烧蚀, 为延长实际应用中等离子体喷枪的使用寿命提供了可靠的理论指导.

      尽管基于局域热力学平衡(LTE)假设的数值模拟可为等离子体关键参数调控提供一定的理论参考, 但Zhukovski等[129]的研究亦表明, 相较于LTE模型, 2T-LCE非平衡模型对电弧等离子体发生器弧电压的预测结果更贴近实验测量结果, 且对阳极弧根贴附行为的描述也更符合实际过程. 因此, 开展非平衡数值模拟研究对于推动等离子体关键参数主动调控技术的发展具有重要的意义.

      在高速再入飞行器壁面热环境地面模拟研究领域, 如何精准地调控等离子体射流冲击壁面的热流密度对于更好地推动地面模拟实验的发展十分重要. 基于多区域耦合的非平衡自由燃烧电弧模型, Fang等[38]的研究表明: 提高工作气压可显著地提升电子-重粒子间的碰撞频率, 促使电子温度(Te)与重粒子温度(Th)趋于一致, 并提高电弧气体温度; 以70 A电流为例, 当工作气压由0.1 atm增至1 atm时, 电弧轴线上电子-原子与电子-离子碰撞最大频率分别从6×1015 s–1和2×1017 s–1提高至2×1017 s–1和2×1018 s–1(如图10(a)所示), 最高气体温度亦由15000 K升高至15600 K. 受此启发, 方川等[36]进一步提出了通过改变气压来调节等离子体射流冲击壁面热流密度的新方法. 非平衡数值模拟结果表明: 在其他条件保持不变时, 当工作气压从102 Pa提高至104 Pa时, 壁面处气体温度由约300 K升高至950 K, 热流密度由5 kW·m–2增加至25 kW·m–2, 如图10(b)所示, 从数值上验证了通过工作气压调控壁面热流密度的有效性, 从而可以根据飞行实验中测量得到的热流密度反演地面模拟实验中对应的真空腔工作气压, 进而实现对高速飞行器再入过程中壁面热环境的地面模拟. 基于该方法开展的典型飞行器部件烧蚀实验与NASA电弧风洞实验结果具有良好一致性, 进一步从实验上验证了该方法可作为高速飞行器再入过程壁面材料烧蚀特性的一种低成本预测试手段, 可用于材料性能风洞考核试验前的初步预评估试验[36].

      在微纳米材料加工领域, 高的材料制备速率不仅需要大体积的高气体温度区域来蒸发原料, 而且需要气体温度陡降区域来促进材料快速成核、团聚且避免颗粒尺寸过度增长, 同时还需要大体积的非平衡区域来保证等离子体较高的化学反应活性. Fang等[38]的非平衡数值模拟研究结果表明, 重粒子子体系的对流换热过程是等离子体非平衡过渡区的主要能量耗散通道, 强化重粒子的对流换热可显著地降低重粒子温度、加剧非平衡程度并扩大非平衡区体积. 基于此, Fang等提出了一种新的环形阳极结构设计, 通过从阳极中心孔反向注入冷气体的方式强化重粒子子体系的对流换热过程, 如图11(a)和图11(b)所示[122]. 相应地, 基于2T-nLCE模型的非平衡数值模拟结果亦表明: 相较于传统平板型阳极结构, 环形阳极结构设计使重粒子子体系的对流换热能量耗散占比由–3.9%(此处的负号“–”表征重粒子子体系获得能量)提升至+28.2%(此处的正号“+”表征重粒子子体系损失能量), 从而在发生器两电极间形成厚度仅0.7 mm、温度梯度达3.3×106 K/m的滞止层, 同时将非平衡区面积占比从61.6%扩展至92.2%, 如图11(c)所示; 在此区域内, 气体温度低至约2300 K, 而电子温度(约8000 K)与密度(约2.4×1020 m–3)仍保持较高的水平, 如图11(d)所示, 数值计算结果与实验测量结果符合良好, 这为产生兼具大体积与高活性的非平衡态热等离子体源提供了有效途径.

      • 4.2.   等离子体关键参数主动调控的实验研究进展

    • 除了上述基于数值模拟结果深入分析基础上的等离子体关键参数主动调控研究工作以外, 也有不少研究者主要基于实验手段进行了等离子体参数调控的研究, 并取得了若干成果. 在等离子体先进材料加工领域, 控制等离子体中的粒子通量和能量通量的大小及方向(如离子电流和电子电流的大小与方向)对所制备材料的物理化学特性至关重要[130]. 例如, Anders[131]通过弯曲的螺线管所产生的磁场从电弧等离子体中分离并去除固体颗粒, 用于制备黏附性好、致密、高质量的薄膜和涂层, 如图12(a)所示; Baranov等[130,132]则利用在等离子体源和基板之间放置磁场线圈产生磁阱来束缚电子, 并通过控制线圈的供电方式和线圈间距设计出不同结构的磁阱, 进而调节基板上的离子电流密度分布, 如图12(b)所示, 最终实现对等离子体射流形貌和材料在基板上的沉积区域的调控.

      等离子体的气体温度与电子温度亦会显著地影响等离子体的材料加工效果. Qiu和Liu[133]为了在提高电弧等离子体温度的同时保持电弧压力基本不变, 设计了自由电弧与约束电弧相结合的复合电弧放电装置. 研究表明, 在电弧焊接过程中, 阳极热流主要由自由电弧主导, 而电弧压力则主要由约束电弧主导, 从而避免了电弧吹力过大破坏熔池的稳定性, 导致熔池坍塌甚至烧穿, 在一定程度上实现了电弧直径与电弧压力的准独立调控, 如图13所示. 而Kato等[134]则使用一个施加负偏压的栅网将等离子体空间物理分隔为两个区域, 通过改变栅网的负偏压来控制下游区域的电子温度, 对于氩气或氦气放电, 可以使得电子温度降低约两个数量级.

      等离子体体积同样也是影响等离子体材料处理效果的一个关键参数. 合适体积的等离子体能够为材料处理提供充足的反应空间, 助力高效完成粉末球化、微纳米材料合成等工艺过程; 而通过对等离子体体积的精准压缩, 则可显著地提升等离子体的能量密度, 进而增强等离子体在焊接等应用场景中的穿透能力与焊接精度. 例如, Yao等[135]、Matsuura等[136]、Liu等[137]采用多相交流放电产生了大体积(直径约100 mm)、高气体温度(6000 K 以上)的等离子体区域, 用于粉末球化和微纳米材料合成. 在焊接领域, Li等[138]通过从阴极侧的小孔径向向内注入聚焦气体的方式约束电弧(如图14(a)所示), 在压缩电弧等离子体体积、提高等离子体能量密度的同时, 既提高了电弧焊接的精度, 也提高了等离子体的穿透性和焊接速度(如图14(b)所示).

      上述实验工作提供了通过外加磁场、等离子体发生器结构设计等方式实现等离子体通量、温度、体积等关键参数主动调控的新思路, 如果能够结合本文第3节所讨论的非平衡数值模拟方法, 通过对计算结果进行等离子体体系质能协同输运机制的深入分析, 或许有可能为优化面向应用的等离子体参数的物理设计提供新的思路.

    5.   关键科学问题讨论

    • 随着电弧等离子体应用领域的不断拓展, 面向不同应用领域的等离子体关键参数调控需求各不相同. 由于电弧等离子体中的非平衡质能输运过程具有显著的多物理场耦合、多相、多界面、多时空尺度效应, 因此, 从数值模拟的角度来看, 在未来的研究中有如下几方面的关键科学问题需要深入开展:

      1)由于电弧等离子体中非平衡输运的多界面和空间多尺度效应, 不仅需要研究放电单元内部的非平衡质能输运, 特别是电极鞘层区行为, 还需要考虑放电单元与外电路、高温等离子体与环境冷气体等不同相态间的耦合效应. 因此, 发展考虑多区域、多相耦合的数值模拟模型(如流体-粒子-解析耦合模型)和高空间分辨能力的数值计算方法, 是分析等离子体体系空间多尺度、多界面效应需要解决的一个关键科学问题. 这将有助于深入地理解和主动调控气体的放电模式(如电弧放电模式、辉光放电模式以及电弧-辉光过渡模式)和亚模式(如电弧放电模式下收缩和弥散两种不同的阳极弧根贴附模式、辉光放电模式下的正常和反常辉光放电等), 以及不同(亚)模式下的等离子体温度、密度和空间尺度等关键等离子体参数.

      2)由于电弧等离子体中非平衡输运的多物理场和时间多尺度效应, 不仅需要同时分辨具有较长特征时间的传热与流动过程和具有较短特征时间的化学反应过程, 甚至需要考虑驱动电源的频率特性, 而且需要针对等离子体工作气体不同的化学成分、工作气压等进行等离子体组分和化学反应动力学过程的筛选. 已有的研究表明, 即使对于同一等离子体体系, 当在数值模拟中考虑不同的等离子体组分和化学反应动力学过程时, 计算得到的等离子体参数时空分布可能会有明显的不同. 因此, 建立快速高效的等离子体化学反应动力学路径筛选方法, 以及可靠的化学反应动力学、热力学性质和输运性质参数数据库, 发展具有高时间分辨能力的数值计算方法, 是分析等离子体体系时间多尺度、多物理场效应需要解决的另一关键科学问题.

      3)由于等离子体体系中质量-动量-能量的协同输运特性, 如何系统、深入地分析电弧等离子体中的非平衡输运机制, 从而为发展需求牵引的等离子体关键参数主动调控提供有效的控制策略, 也是数值模拟研究需要解决的又一关键科学问题. 本团队[36,120,122,139]近几年的初步研究表明, 以等离子体体系的能量“注入-再分配-损失”为主线的“能量树”概念[121]或许能够为分析已有等离子体体系中的质能协同输运机制, 进而提出等离子体参数的主动调控策略具有一定的指导意义.

    6.   总结与展望

    • 电弧等离子体具有高温、高能量密度以及包含丰富的化学活性粒子等突出特点, 在环境、能源、航空航天、新材料、国防等诸多领域发挥着重要作用. 由于实际应用中的电弧等离子体体系中通常存在复杂的传热、流动、化学反应、电磁场等多物理场、多时空尺度、多相、多界面耦合效应, 因此, 以实际应用需求为牵引, 在发展先进实验诊断技术的同时, 发展能够描述电弧等离子体体系非平衡质能协同输运的数值模拟方法, 不仅有助于深入理解电弧等离子体体系的非平衡协同输运机制、丰富低温等离子体物理的基础理论, 而且有助于推动现有电弧等离子体技术的优化升级, 甚至有望催生一系列战略性、颠覆性新兴应用. 本文从数值模拟的角度对电弧等离子体的研究现状进行了调研分析, 主要得到如下结论:

      1)自由燃烧电弧中同时存在着显著的热力学、化学、电非平衡效应, 是研究高气压碰撞等离子体体系非平衡质能协同输运的典型模式体系. 本文系统地总结了直流电弧等离子体非平衡输运物理-数学模型的研究进展, 给出了等离子体处于完全非平衡(nLTE-nLCE-nQCN或2T-nLCE-nQCN)条件下的控制方程, 以及在2T-nLCE-QCN, 2T-LCE-QCN, LTE-LCE-QCN条件下的具体形式; 从自洽性、一致性和守恒性3个方面讨论了物理-数学模型应该满足的基本条件; 在保证计算正确性和准确性的同时, 从提高数值模拟效率的角度出发, 给出了上述控制方程在满足一定条件下的简化形式.

      2)在电弧等离子体特性研究方面, 重点总结了电弧等离子体体系质量-动量-能量非平衡协同输运机制研究的发展历程. 在早期, 对电弧等离子体的研究主要集中于高温高能量密度的核心区, 在数值模拟中主要基于平衡态流体模型, 研究了工作气压、放电电流、气体种类等多种因素对电弧等离子体特性的影响机制. 随着等离子体应用的快速发展, 对电弧等离子体的研究区域不断向非平衡的过渡区和电极边界层发展, 数值模拟主要基于2T-LCE和2T-nLCE非平衡流体模型以及考虑固体电极区和电极鞘层区的多区耦合模型, 着重揭示等离子体体系质量-动量-能量非平衡协同输运机制和特性.

      3)在电弧等离子体关键参数调控研究方面, 以应用需求为牵引, 以数值模拟为主线, 重点总结了通过等离子体发生器结构设计、外加电磁场位型设计以及工作参数调节等手段对等离子体形貌、温度、体积、粒子通量、热流密度等关键参数的主动调控策略.

      随着电弧等离子体应用领域的不断拓展, 对等离子体关键参数的调控需求也将不断提升. 从数值模拟的角度来看, 为了适应电弧等离子体应用领域不断拓展对等离子体关键参数进行主动、精准调控提出的迫切需求, 进一步深度耦合流体模型、粒子模型和解析模型, 发展考虑多区域、多相耦合的物理数学模型和可靠的等离子体化学反应动力学、热力学性质和输运性质参数数据库, 以及高效、高时空分辨能力的数值计算方法, 甚至借助于当前快速发展的人工智能技术, 形成系统分析和科学预测电弧等离子体体系中非平衡质能协同输运机制的理论方法, 对于推动电弧等离子体源的广泛应用、甚至催生新的应用具有重要的理论意义和工程应用价值.

    参考文献 (139)

目录

/

返回文章
返回