1 前言

随着微电子技术不断迭代升级，电子器件向着集成化、微小型化发展。随着晶体管密度的逐渐提高，其热流密度也在不断上升，传统的散热手段无法完全满足要求，于是，人们开始关注通过电晕放电产生离子风强化传热的方法。电晕放电现象是在曲率半径很小的尖锐电极中施加一定的电压时，气体会被电离，从而产生放电的现象^[1]。而当被电离的正负离子由于电场力的作用开始移动时，则会产生气体的流动，也就是所谓的“离子风现象”^[2-3]。离子风的应用范围很广 ^[4]。首先，利用离子风的流体性质，可以进行强化传热^[5]，其在电子设备与元器件中得到了广泛应用^[6-7]，具有震动小、噪声低的特点；此外，离子风在航天推进^[8-9]、回收微小颗粒^[10-12]中也得到了应用。电晕放电过程是一种存在多个化学反应的过程，实验观测难度较大，而数值模拟研究可以更好地获得电场内部的详细情况与离子风的流场分布。众多学者对不同的电极结构，采用了不同的模拟模型与模拟方法进行数值模拟研究。Bouazza等^[13]采用半解析的方法对线—板型离子风装置进行模拟，并分析了电压、导线半径对风速的影响。Monrolin等^[14]采用多尺度模型进行建模，并且将求解域分为内部域与离子漂移域两个域，分别计算。Chen等^[15]对筒电极结构的电晕放电以及离子风现象进行了模拟，发现在离尖端较近的电离区和离尖端较远的离子漂移区，流速和体积力方向相反，计算得到的平均电流、脉冲频率和流速与实验结果非常接近。Qu等^[16]设计了T型多导线离子风泵，并将实验与模拟2种手段相结合，分析了其对发光二极管(light-emitting diode，LED)芯片的散热效率，结果发现导线在18 kV的固定电压下放电，LED外壳温度可以降到74 ℃，低于最高额定温度85 ℃。Guan等^[17]研究了针-环系统中的离子风现象，采用无因次参数来表示电场力与惯性力的比值，模拟得到阴极电流和出口速度分布与实验数据的误差在5% 以内，最大速度与电晕电压成正比。Wang等^[18]采用模拟的方法对离子风对单个以及多个热源强化散热进行分析，综合考虑对局部热点的冲击效应、电极之间的阻隔效应以及电极的入口效应3种散热效应。结果表明，对单电极而言，局部接地的效果更佳；而对多电极而言，强化电极对结构的散热效应更强。Ladan等^[19]采用微极性流体模型对离子风引起的自然对流换热进行了数值模拟，并与常用的k-ε模型进行了对比，结果发现当流体处于过渡流时，微极性流体模型可以较好地模拟离子风的状态，但相关的材料参数的选取则高度依赖于外加的电压与流体的瑞利数。此外，Ma等^[20]对环型离子风装置进行了实验分析，强化传热效果可以提升200%；Feng等^[21]比较了环型离子风装置与翅片型离子风装置的散热情况，建立了Nu-Re关系的经验公式；Wang等^[22]分析了离子风的产生过程与强化传热的原理，总结了近年来离子风强化传热的进展，并且对其发展趋势进行了分析。

虽然前人对离子风强化传热已开展较多研究，但大多数针对某一固定结构的离子风装置，在研究中大多也只改变了操作参数，对结构参数的分析不够系统全面。由于离子风的发生可以由不同的结构产生，大多数文献并未说明其设计依据，尤其是在数值模拟领域，针对不同结构离子风装置的对比研究仍然不足。为此，本研究针对不同构型的离子风发生装置，讨论不同结构参数与操作参数下离子风对热源的强制对流传热的性能，为优化离子风电极形状，选择传热效率更好、更经济的电极参数，对离子风在强化传热方面的应用提供理论依据。

2 计算模型 2.1 数学模型

由于电晕放电过程中的化学反应众多，机理十分复杂，在确保模拟精确性的情况下，为了节约计算时间，对计算区域的参数做出如下假设：

(1) 只考虑空气中一种电荷的移动，忽略电化学反应；

(2) 空气视为不可压缩流体；

(3) 忽略黏性耗散引起的热效应；

(4) 忽略电晕放电过程中的焦耳热。

根据上述假设，电势的泊松方程表达式为

$ {\nabla ^{\text{2}}}V{\text{ = }} - \frac{\rho }{\varepsilon } $

(1)

式中：∇为微分算子；V为电势，V；ρ为空间电荷密度，C⋅m⁻³；ε为介电常数，F⋅m⁻¹。

忽略由热扩散产生的电荷移动，电流的连续性方程可以表示为

$ \nabla \boldsymbol{J}{\text{ = 0}} $

(2)

$ \boldsymbol{J}{\text{ = }}\mu \rho \boldsymbol{E}{\text{ + }}\rho \boldsymbol{u} $

(3)

式中：J为电流密度，A⋅m⁻²；μ为离子迁移率，m²⋅V⁻¹⋅s⁻¹；E为电场强度，V⋅m⁻¹；u为空气流速，m⋅s⁻¹。

空气的流动采用不可压缩流体质量守恒方程与纳维—斯托克斯方程进行描述，并加入电晕放电产生的电场力，方程表达式为

$ \nabla \boldsymbol{u}{\text{ = 0}} $

(4)

$ (\boldsymbol{u} \cdot \nabla )\boldsymbol{u} = - \frac{{\nabla \boldsymbol{p}}}{{{\rho _{\text{air}}}}} + \nu {\nabla ^2}\boldsymbol{u} + \frac{{Q \boldsymbol{E}}}{{{\rho _{\text{air}}}}} $

(5)

式中：p为压差，Pa；ρ_air为空气的密度，kg⋅m⁻³；ν为空气的运动黏度，m²⋅s⁻¹；Q为电荷的带电量，C。

传热的能量方程由下式表达：

$ \nabla (\boldsymbol{u}T) = \nabla (\frac{\lambda }{{{\rho _{\text{air}}}{c_p}}}\nabla T) $

(6)

式中：T为温度，K；λ为空气的导热系数，W⋅m⁻¹⋅K⁻¹；c_p为空气的比定压热容，J⋅kg⁻¹⋅K⁻¹。

2.2 几何模型与边界条件

本研究建立的离子风模型包含3种不同的几何构型，如图 1所示。3个模型的下部均为1个长方体热源，热源上方放置一定宽度的铜散热片用来辅助散热，而模型最上方则为电晕放电的发射极。不同的是，图 1(a)中的散热片直接接地，充当电晕放电的收集电极；图 1(b)添加了2个互相平行的板电极充当收集电极；图 1(c)则采用4根互相平行的电晕线充当收集电极。为了保证散热面积相同，图 1(b)与(c)的散热片也得以保留。

由于本次模拟涉及电场、湍流流场以及温度场的耦合计算，求解需要相当大的计算资源，且本次模拟模型存在高度对称性，为了便于求解且不失准确性，将图 1的三维模型在对称面上进行纵切，从而简化为图 2的二维模型，同时设定外部的较大空间作为计算域。模型的几何尺寸见表 1。表中，L_s为计算域的长度，H_s为计算域的宽度，L_h为热源长度，H_h为热源宽度，L_c为散热片长度，H_c为散热片宽度，r_e为放电电极半径，H_e为放电电极高度，单位均为mm。

边界条件的划分如下，对于静电场与电荷场，求解域的边界条件采用第二类边界条件，电荷与电压梯度为0；起晕电极与接地电极采用第一类边界条件，起晕电极电压的初始值为输入电压。由于无法准确计算起晕电极的初始电荷值，所以采用Kaptzov假设来进行计算。Kapzov假设的步骤为先假定电场中的电荷密度，通过迭代计算的方式得到电场强度，再与Peek公式的电场值进行比较，对电荷密度进行修正。起晕电极的临界电场强度则由Peek定律求得^[2]：

$ {E_{\text{0}}}{\text{ = 3}} \times {\text{1}}{{\text{0}}^{\text{6}}}\delta (1 + \frac{{0.308}}{{\sqrt {\delta r} }}) $

(7)

式中：E₀为临界电场强度，V⋅m⁻¹；δ为空气相对密度，无量纲参数；r为电晕线半径，cm。

流场的求解域边界的设置为开放边界条件，初始速度为0，而电极表面则设置为无滑移边界条件。

2.3 网格划分与模拟结果验证

图 3为离子风模型的网格划分示意图。主要采用三角形网格；对于离子风的发射极与接收极，则采用相对密集的网格；此外，在电晕极与收集极壁面添加了若干层边界层网格，保证计算的准确性；而对于外部计算域的精度要求不高，因此网格尺寸较大，以提高计算效率。

根据网格划分的疏密采用7套不同的网格进行网格无关性验证，验证参数为热源的温度，计算结果如图 4所示。从图 4可知，热源的温度逐渐降低，且当网格数量为28 000个左右时，热源温度趋于不变。为了兼顾计算的准确度与计算速度，采用28 000个网格进行后续的计算。

本研究以文献[24]中的实验为考核算例，对实验进行数值模拟，将模拟得到的平均离子风风速与文献中的实验值对比，结果如图 5所示。对比结果显示，本研究建立的离子风模型在5~7 kV准确性非常高，误差小于5%；最大误差出现在8 kV位置处，最大偏差为11.7%，基本满足工程应用中的误差范围。

3 计算结果分析 3.1 模拟结果

图 6为电极电压在30 kV时，直接接地结构离子风模型的电势、电荷和速度分布。从图中可以看到，与文献[25]中的电势、流速分布相类似，最大的电势分布都集中在离子风的发射极附近，且周围的梯度都相当大，说明电晕放电主要集中在发射极附近，电场线也均匀地由离子风的发射极指向收集极，且电势逐渐减小；空间电荷的密度由发射极指向收集极快速减小，而且除了从发射极到收集极的一小部分区域外，其余的部分几乎不存在电荷。离子风风速在中心位置呈最大值，但是吹向散热片时中心位置会有速度很小的区域(图 6(c)的底部中间部分)，这里的换热情况较差，最大风速为2.91 m⋅s⁻¹。

图 7为电极电压在30 kV时，线—板结构离子风模型的电势、电荷和速度分布。从图中可以看出，由于接地电极与发射电极距离很近，线—板结构离子风的电势分布更为紧凑，电荷的分布区域也明显小于直接接地结构的电极，但是电荷密度的最大值有所提高，且由于板的存在，离子风的风速也有提高，为5.04 m⋅s⁻¹。值得注意的是，还有一些风是从板的上部向外流出去的，这是由于板之间的压力较大而板边缘的风速较小的原因。

图 8为电极电压在30 kV时，线—网结构离子风模型的电势、电荷和速度分布。从图 8(a)和(b)中可以看出，线—网结构的发射极与接地极虽然接近，模型结构也与线—板型电极相近，但是由于接地极的面积较小，电势与电荷密度的分布依然与直接接地结构相近。电荷密度与线—板结构在同一数量级，离子风的最大速度处于直接接地与线—板两者之间，如图 8(c)所示，速度为4.23 m⋅s⁻¹。

3.2 电压的影响

本节模拟热流密度在6 000 W⋅m⁻¹、不同电压下3种构型的离子风的平均风速、功率以及热源温度的变化情况。电压范围为22~40 kV。模拟结果如图 9所示。结果表明，随着电压的升高，3种构型的离子风风速逐渐增高，传热系数逐渐升高，功率逐渐增高。这是因为电压的升高导致空间电荷密度的上升，增大了电场力与电流。线—板、线—网和直接接地电极的平均离子风风速依次降低，传热系数也逐渐降低，但是线—板构型的电功率明显大于直接接地电极以及线—网电极，并且随着电压的增大，功率的差距也会越来越大，但此时平均风速以及传热系数的差距并没有增大。所以当电压很大时，线—板电极由于耗能问题，可能不是装置构型的最优解，需要针对具体的应用情况进行具体分析。

3.3 发射电极高度的影响

本节模拟热流密度在6 000 W⋅m⁻¹、电晕电极电压为30 kV、不同发射电极高度下3种构型的离子风的平均风速、功率以及热源温度的变化情况。发射电极的高度为40~80 mm。模拟结果如图 10所示。从图 10(a)中可以看出，离子风的平均风速随发射电极高度的增加而降低，这是由于放电极与接地极之间的电场强度减小的缘故；从图 10(b)中可以看出，传热系数随发射电极高度增加而降低，这是由于离子风速减小的缘故；而在图 10(c)中，功率随发射电极高度增加而降低。说明减小电极高度可以有效增强换热效果，但是电极高度不能无限减小，否则当电场强度很大时，空气会被彻底击穿，造成安全隐患。

3.4 散热片宽度的影响

本节模拟热流密度在6 000 W⋅m⁻¹、电晕电极电压为30 kV、电极高度为40 mm、不同散热片宽度下3种构型离子风的平均风速、功率以及热源温度的变化情况。散热片的宽度范围为40~120 mm。结果如图 11所示，对于线—板式及线—网式电极，由于散热片与电晕放电的情况基本无关，所以二者的功率保持不变，而因为在离子风吹向散热片的中心地区存在流速很小的区域，导致由于散热片宽度的减小，散热片上的离子风平均风速降低，传热系数也降低。直接接地结构与上述2种结构有一些差别，因为散热片本身也是电晕放电的接地极，所以当散热片宽度减小时，功率也减小。因此当散热片的宽度比较小、散热需求比较低时，结构简单、功率低的直接接地式结构也不失为一种好选择。

4 结论

对3种不同构型的离子风模型进行模拟，获得了模型内部以及周围的电势、电荷和流场分布，并研究了放电电压、电极高度与散热片宽度对离子风强化传热的影响，得出以下结论：

(1) 在所模拟的参数下，在相同电压时，线—板型结构的平均离子风速度最高，线—网型结构次之，直接接地型结构平均离子风速最低。

(2) 当电晕极的电压上升时，3种离子风构型的风速均增加，热源温度下降，功率也会有所升高，线—板型的离子风构型增加的功率要明显大于其他2种构型。所以从耗能角度看，如果电压较小时，线—板构型的效果更好。

(3) 电极高度的减小也可有效强化离子风模型的换热性能，但注意电极高度不能太低，否则空气会被彻底击穿，造成安全隐患。

(4) 对于线—板型电极与线—网型电极，散热片面积的变化只会影响换热面积；而对于直接接地型电极，散热片还兼顾收集电极的功能，所以当设计时的散热片面积比较小时，直接接地电极的效果可能会更好。