針-板空氣間隙流注放電起始過程的三維PIC/MCC仿真研究?

李晗蔚 孫安邦 張幸 姚聰偉 常正實 張冠軍

(西安交通大學,電力設備電氣絕緣國家重點實驗室,西安 710049)

1 引 言

流注放電(streamer discharge)多呈現出明亮、絲狀通道結構,通常在施加電場大于介質本身擊穿電場的區域形成,且在其頭部存在強烈的電場增強效應,具有向低于擊穿電場的未電離區域傳播的特性[1?3].作為高電壓領域氣體放電的一種常見初始形式,流注放電在自然界中大氣放電、特高壓輸電線路長空間間隙放電以及工業應用中都扮演著非常重要的角色.

流注放電理論于20世紀40年代初由Raether,Loeb和Meek提出[4,5],長期以來,為了解流注放電的物理機理,國內外學者進行了大量的實驗與數值仿真研究.在實驗研究方面,最具有代表性的為荷蘭埃因霍溫理工大學應用物理系,采用棱鏡法等對流注光電特性進行了大量實驗研究,包括其三維結構、傳播速度、流注尺寸等[6,7].然而,由于流注放電的產生與發展過程非常快,在大氣壓下一般為納秒量級,且實驗主要集中在毫米至厘米尺度放電間隙內進行,因此實驗手段的精確測量具有較大難度,當前流注物理特性的研究仍主要以數值仿真為主.

流注放電的數值仿真模型主要分為三類[8]:流體模型、粒子模型和流體粒子混合模型.流體模型從宏觀角度描述等離子體的放電過程,計算效率較高.如,Kulikovsky[9]首次將光電離加入到流注放電的二維軸對稱漂移-擴散近似流體模型中,研究了兩平板電極空氣間隙的流注放電特性,得到流注傳播行為主要受頭部電荷層寬度影響的結論;Luque和Ebert[10]引入自適應可變網格技術,采用三維流體模型研究了空氣中多個流注在壓強、氮氧混合比等因素影響下的運動行為,發現相鄰流注互相作用的機制主要是電荷的排斥和非局域光電離造成的相互吸引.但是,流體模型的使用范圍具有很大局限性,也很難描述放電過程的隨機行為(比如分叉等).基于第一性原理建立的粒子模型從微觀角度直接跟蹤放電過程中帶電粒子的運動過程,進而統計得到其宏觀參數的變化.如,西安交通大學李永東教授課題組采用粒子模型研究了大氣中光電離對流注形成和傳播的影響,獲得了逃逸電子產生過程的物理圖像以及流注的分叉結構[11];我們與Teunissen博士和Ebert教授合作,開發了一套流注放電三維粒子仿真模型,研究了均勻電場下自然界中光電離、背景電離(自由電子在強電場作用下從O?2解吸附)對流注形成的影響[12?14];粒子模型的主要缺點在于計算量相對較大.結合粒子和流體模型各自的優點,流體-粒子混合模型在電場梯度較小區域應用流體模型,在流注外圍及其他電場梯度大的區域應用粒子模型[15],然而,由于存在粒子流體的過渡區域,常規的并行計算以及可變網格等新技術較難實現,總體計算量依然比較大.

縱觀空氣中流注放電的數值研究現狀,宏觀流體模型仍是流注數值仿真的常用手段,且大多局限于二維仿真,僅適用于對單個流注的傳播特性進行研究.本文側重于對針-板電極結構中流注放電開展三維粒子模擬研究,研究價值在于:首先,針-板電極結構產生不均勻電場,電場在間隙中呈不均勻發散狀,影響流注演化特性與分叉結構;第二,針型電極附近產生高強度電場區域,具有較高的電子數密度,流注放電易于從針電極處產生;第三,實際放電大多是在不均勻電場結構中產生,而現有數值研究大多局限于二維仿真或均勻電場的仿真,亟需對以針-板電極結構類的非均勻電場下流注放電過程開展三維流注放電研究,反映真實放電形貌和結構特征.

2 仿真物理模型介紹

2.1 三維PIC/MCC數值模型

本文采用我們與Ebert教授課題組合作開發的三維PIC/MCC模型對針-板電極放電起始過程開展數值研究.該模型前期已用于研究均勻電場空氣中背景電離對流注形成的影響[12?14].PIC/MCC數值模擬的流程如圖1所示.

為了減少計算時間,本模型中采用的模擬環境和基本假設如下:

1)模擬氣體為干燥空氣,主要為氮氧混合氣體組成;

2)由于離子的運動速度遠小于電子,在所關注的較短的放電時間內 (大氣壓下ns時間尺度),離子假設為靜止;

3)整個模擬過程處于弱電離狀態,只考慮電子-中性粒子之間的碰撞,忽略電子-電子碰撞以及電子-離子之間的庫侖碰撞,中性氣體假設為背景,參與電子之間的碰撞行為,而本身的密度不做改變;電子與中性氣體的碰撞界面由Lxcat數據庫獲得[16];

4)忽略電荷運動產生的自洽磁場,通過求解Poisson方程更新模擬區域的電場變化.

圖1 PIC/MCC模擬流程圖Fig.1.Flow chart of PIC/MCC simulation.

2.2 電子運動與碰撞

2.2.1 電子運動

電子的運動遵循牛頓-洛倫茲定律,不考慮磁場作用,運動方程為:

式中,me為電子質量,q為電子電荷,v和x分別為電子速度及其位置矢量,E為電場.

采用Verlet算法描述電子運動,具體計算方法為:

式中,vk,xk為電子運動k個時間步長后的速度和位置;ai為k個時間步長后的電子加速度.

2.2.2 粒子碰撞

模型中,考慮粒子之間的彈性碰撞和非彈性碰撞(例如電子激發、電離、吸附等),主要反應方程如下:

彈性碰撞

電離碰撞

激發碰撞

電子吸附

電子在不同碰撞截面積σi與背景氣體分子發生碰撞時,其碰撞頻率可以表示為

式中:N為背景氣體分子密度;v為電子相對速度.

在Δt時間內,電子與中性氣體碰撞概率可表示為

2.3 光電離(photoionization)

流注從正極性針型電極起始和發展過程中,其傳播方向和外加電場方向相同,與電子漂移方向相反,這意味著外部電子需要從其流注前方反向運動到流注通道內部.在空氣中,光電離是產生外部自由電子的主要機制:空氣中的N2受激發之后,處于激發態的N2分子退激發為基態N2分子時釋放出的光子與O2分子發生碰撞電離,從而釋放出自由電子[8].雖然光電離機制被廣泛采納并且應用已久,然而其內部具體的物理過程還不甚清楚.本模型采用Zheleznyak等提出的光電離模型,計算流注放電過程中的光電離過程:

定義湮滅因子ηq:

式中:p為壓力;pq為湮滅壓力,pq=30 mbar.

當p=pq時,一半激發態的分子會由于與其他分子的碰撞而退激發,并不會對外釋放出光電子.在整個模型中,平均每單位光電離會產生η個光電子,定義η=ηE·ηq,ηE是一個與電場強度有關的因子[18].

2.4 粒子數值模擬加速技術

三維粒子模型的計算量極大,因此采取并行計算的方法,加速電勢及電場的計算速度.同時,也為了在計算過程中保持計算準確度,使用了自適應可變網格(adaptive mesh re fi nement,AMR)、粒子可變權重和并行計算技術.

2.4.1 自適應可變網格

根據電子的密度梯度和電場大小對計算區域的網格大小進行調整.當區域內電子數密度超過設定閾值時,區域內網格將自動加密處理,新網格大小為原網格的一半;電子數密度低于閾值時,網格將自動粗化處理.圖2是多重自適應網格示意圖.在電場強度高、粒子數密度梯度較大的區域內,網格劃分的較細;在電場較低、粒子數密度梯度較低的區域內,網格劃分的相對較粗.

圖2 多重自適應網格Fig.2.Multi-scope adaptive grid.

不同尺寸網格數量、最大和最小網格尺寸決定著模擬的計算精度和計算時間.本文計算中,允許網格細化次數為5,在經過5次細化后,達到最小網格,不再繼續進行細化.網格進行細化時的最小電子數密度為1×1011m3;每隔2.0×10?11s判斷網格是否需要細化.

2.4.2 粒子可變權重

放電過程中,計算域內粒子數目多,隨放電的發展,粒子數目也急劇增長.為了控制模擬速度并減小計算量,對計算粒子權重進行可變調節,如圖3所示,每個計算粒子可以分解為兩個計算粒子,或與其他計算粒子合并成為一個計算粒子.

圖3 計算粒子的分裂與合并Fig.3.The splitting and merging of calculated particles.

詳細的權重調整方法以及權重變化后粒子的速度和位置分布處理等,可參見文獻[19,20].

3 結果與討論

本文主要探討電壓幅值、氣體組分和電極尖端形狀尺寸三個參數對流注放電的影響,著重關注放電過程中形成的初始電離云、流注分叉結構和電場強度差異等.

3.1 仿真初始條件

3.1.1 計算域和電極

本模型中,計算域為圖4(a)所示的(5×5×5)mm3立方體區域,電極尖端正對平面為地電極.

圖4 計算域示意圖及電極尖端放大圖Fig.4.The schematic of the computational domain and the enlarged electrode tip.

計算域的選擇主要出于以下考慮:

1)在立方體計算域中,電極尖端處于立方體正中心處,與各平面距離相同,流注向空間各方向發展概率均等;若采用長方體區域,電極尖端與各平面距離不等,導致空間中電場強度分布存在較大差別,可能也會影響流注的發展趨勢;

2)在一個相對較小的區域內,較低電壓即可實現流注放電.

針型電極由球體、圓柱與圓錐模擬而成,如圖4(b)所示,以在電極表面放置點電荷的形式描述外施加電壓.

3.1.2 初始參數設置

模型中,模擬區域內背景電離密度為1012m?3,等離子體在整個計算域空間內隨機分布.為加快流注放電仿真速度,在電極尖端正上方0.05 mm處人為放置一定數量的電子和正離子(宏觀為中性),帶電粒子數密度呈高斯分布,半徑為50μm,包含3×106個電子并呈現為高斯分布.

模擬中,氣體壓強為1 bar,溫度293 K,使用N2和O2的混合氣體模擬空氣,比例為N2:O2=80%:20%.

3.2 電壓幅值對流注放電的影響

選擇電壓參數分別為:3.0,3.5,4.0和4.5 kV,分別在此四個電壓幅值下進行模擬.

3.2.1 初始電離云

流注放電開始階段,在針尖附近觀察到形成球狀類似于電暈的初始電離云(inception cloud),這與之前實驗中觀測到的現象是一致的[6].如圖5所示,電離云大小隨施加電壓的增大而增大.

圖5 不同電壓幅值下初始電離云圖像 (a)3.0 kV,1.0 ns;(b)3.5 kV,0.8 ns;(c)4.0 kV,0.6 ns;(d)4.5 kV,0.4 nsFig.5.Photos of inception clouds under different voltage amplitudes:(a)3.0 kV,1.0 ns;(b)3.5 kV,0.8 ns;(c)4.0 kV,0.6 ns;(d)4.5 kV,0.4 ns.

假定初始電離云是一個半徑為R的球體,其表面近似等電位U0,電極間隙的擊穿場強Ec,則其最大半徑可近似估算為R0=V0/Ec.因此,電壓越高,初始電離云半徑越大,當電離云發展到最大半徑R0位置時,單個流注通道開始出現并向前發展,在此類低場強區域,流注的傳播主要依靠頭部空間電荷產生的較強自洽電場(可達到3Ec—5Ec)的作用.

圖6給出了不同電壓下產生初始電離云時的電場強度分布.所有圖像均是x-z軸平面,在y軸方向取計算區域的中心位置截斷取得的,即y=2.5 mm處.后文中所有電場強度圖像均是按照此方法得到.

可以發現:各電壓條件下,針型電極表面存在一個較薄的高場強區域,即為上述所描述的電場增強效應;在人為放置種子電子的影響下,電壓較低(3.0 kV及3.5 kV)時,種子電子處整體電場強度較高,而當電壓升高時(4.0 kV及4.5 kV),僅其頭部成為高電場強度區域.因此,高斯種子電子的放置會改變電極尖端附近的電場強度分布.

圖6 初始電離云電場強度截面圖 (a)3.0 kV,1.0 ns;(b)3.5 kV,0.8 ns;(c)4.0 kV,0.6 ns;(d)4.5 kV,0.4 nsFig.6.Cross Section view of inception clouds’electrical field strength:(a)3.0 kV,1.0 ns;(b)3.5 kV,0.8 ns;(c)4.0 kV,0.6 ns;(d)4.5 kV,0.4 ns.

3.2.2 流注形態結構

圖7 不同電壓幅值下流注放電電子數密度演化過程 (a)3.0 kV;(b)3.5 kV;(c)4.0 kV;(d)4.5 kVFig.7.Electron density evolution of streamer discharge under different voltage amplitudes:(a)3.0 kV;(b)3.5 kV;(c)4.0 kV;(d)4.5 kV.

圖7是流注放電中電子數密度隨時間演化過程的模擬結果,為了便于比較,所有流注圖像的觀測角度均相同.因不同電壓下的流注放電持續時間不同,圖中進行比較的時間選取略有不同.

常用者為頭孢哌酮、頭孢哌酮舒巴坦、頭孢克肟，菌素對腸桿菌科細菌等革蘭陰性桿菌具有強大抗菌作用，頭孢他啶和頭孢哌酮除腸桿菌科細菌外對銅綠假單胞菌亦具高度抗菌活性。

可以看出,在流注的發展過程中,易產生分叉形成類似樹枝狀結構.不同電壓幅值作用下,流注發展形態有很大差別.電壓越高,流注分叉時間越早,其分支結構越多.

電壓為3.0 kV時,流注形態與其他電壓下流注形態有較大差別.第一,如3.2.1節所述,幾乎未形成初始電離云;第二,其他電壓工況下流注分叉,形成多流注通道;而此工況下形成單流注通道,直到5.4 ns時流注開始分叉、產生多流注通道的跡象.此外,在此電壓工況下,流注通道內電子密度比其他更高電壓下略高,可能是由于流注發展過程中通道相對集中所導致的.隨放電發展,現有流注通道將會形成更多的流注通道.

3.3 氣體組分對流注放電的影響

根據文獻,N2-O2氣體混合物中,N2與O2的比例直接影響流注放電的結構特征.本節在外加4.0 kV電壓下,確定另外兩種工況下N2和O2的比例為:0:100和99.9:0.1.前一種情況等效為純凈O2,后一種情況等效為N2中摻雜少量O2.

Nijdam等[21]在實驗中觀察了不同O2濃度下流注放電的發展過程,發現在100%O2中幾乎看不到初始電離云的存在;N2體積分數越大,流注分叉越多,且分支越細.

3.3.1 初始電離云

在本文的仿真條件下,流注放電在不同氣體組分下形成的初始電離云形貌如圖8所示:初始電離云都具有較高的電子數密度,隨著O2濃度提高,電離云半徑增大,且單獨流注通道形成時間提早.根據流注理論,光電離是空氣中產生二次電子崩的主要原因,而二次電子崩是流注向前發展的重要因素;在本模型中同時考慮了電子的吸附與解吸附過程.因此,O2所占體積越高,光電離和電子解吸附過程越容易提供自由電子,利于正流注的向前發展.

圖8 不同氣體組分影響下的初始電離云圖像 (a)0.1%O2,0.8 ns;(b)20%O2,0.6 ns;(c)100%O2,0.4 nsFig.8.Comparison of inception clouds with different gas components:(a)0.1%O2,0.8 ns;(b)20%O2,0.6 ns;(c)100%O2,0.4 ns.

3.3.2 流注形態結構

三個工況下流注的發展形態如圖9所示.

圖9 不同氣體組分流注發展圖像 (a)0.1%O2;(b)20%O2;(c)100%O2Fig.9.Photos of the development of streamer discharge with different gas components:(a)0.1%O2;(b)20%O2;(c)100%O2.

圖10 不同氣體組分流注放電電場強度截面圖 (a)0.1%O2;(b)20%O2;(c)100%O2Fig.10.Sectional views of electrical field strength of streamer discharge with different gas components:(a)0.1%O2;(b)20%O2;(c)100%O2.

可以看出,在0.1%O2濃度工況中,流注的電子數密度明顯大于高O2濃度工況下流注的電子數密度;在3.2 ns內,形成的分支較少,前期近似為單流注向板電極發展,在2.2 ns時可見明顯流注分支結構;流注表面整體較粗糙,類似于存在許多細小凸起.隨O2濃度的增高,20%O2工況下,流注分叉數目,尤其是向兩側分叉數目逐漸變多;100%O2工況下,除形成向上的主流注通道外,還向四周形成了細小的流注放電通道.

圖10給出了三種工況下流注發展時的場強變化圖.流注放電初期,尖端電極表面存在高場強區域薄層;隨著流注發展,薄層消失;此后,流注頭部電場強度較高.初期流注可能向空間各個方向發展,后期由于流注通道之間的相互抑制,部分流注由于電場強度過低或少量自由電子注入而熄滅,最終產生隨機形態的分叉結構.

在0.1%O2工況時,流注頭部始終存在大片高場強區域,如圖10(a)所示,由于光電離項的影響,該工況下流注發展較慢,流注通道沒有完全形成;隨著流注繼續向前發展,會在通道內部產生大量電荷,電場強度相應降低.

3.4 電極形狀尺寸對流注放電的影響

上文模擬設置中,電極尖端球體半徑為50μm,本節考慮了另外兩種工況為尖端球體半徑5μm和250μm,如圖11所示.

圖11 不同尺寸電極尖端示意圖 (a)5μm;(b)250μmFig.11. Photos of different electrode tip sizes:(a)5μm;(b)250μm.

3.4.1 初始電離云

三種工況下,初始電離云圖像及其電場強度截面圖分別如圖12和圖13所示.

圖12 不同電極尖端尺寸工況下初始電離云圖像 (a)5μm,0.4 ns;(b)50μm,0.6 ns;(c)250μm,1.2 nsFig.12.Photos of inception clouds with different electrode tip sizes:(a)5μm,0.4 ns;(b)50μm,0.6 ns;(c)250μm,1.2 ns.

圖13 不同電極尖端尺寸形成的初始電離云電場強度截面圖 (a)5μm,0.4 ns;(b)50μm,0.6 ns;(c)250μm,1.2 nsFig.13.Sectional views of electrical field strength of inception clouds with different electrode tip sizes:(a)5 μm,0.4 ns;(b)50μm,0.6 ns;(c)250μm,1.2 ns.

可見,尖端半徑為250μm的工況下,幾乎沒有形成初始電離云,且電極尖端處電場強度較低.在5μm和50μm電極尺寸尖端均有密度較高的電場強度區域.在尖端半徑較大時,低電場強度會導致碰撞電離減少,放電不足以形成高電子數密度聚集的初始電離云.

3.4.2 流注形態結構

尖端電極尺寸影響下流注的發展形貌如圖14所示.

圖14 不同電極尖端流注放電發展圖 (a)5μm;(b)50μm;(c)250μmFig.14.Photos of streamer discharge with different electrode tip sizes:(a)5 μm;(b)50 μm;(c)250 μm.

對于電極尖端半徑為5μm的工況,與原先50μm時流注形態對比,流注后期發展無明顯區別.這是由于流注分叉且獨立通道形成后,空間電荷產生的自洽電場成為影響流注發展的主要因素,因此后期流注發展與電極尖端形狀尺寸關系不大.尖端半徑為250μm的工況下,流注放電形態與前兩種工況出現了明顯區別,未形成明顯的初始電離云,且隨放電時間發展,基本呈現單流注形式,并在5 ns左右產生分叉.可以看出,針型電極尺寸可顯著影響流注發展過程分叉的時間和結構.圖15給出了5 ns后,250μm半徑流注發展的圖像.

綜合圖15和圖14(c)可知,由于電極尺寸較大,電極上方始終存在一片高電子數密度區域;然而隨時間發展,流注后期仍會產生分叉結構,與其他工況下的流注發展無異.這同樣是由于后期流注的發展主要取決于空間電荷的自洽電場,與電極附近的電場強度高低關系不大.

圖15 尖端球體半徑為250μm的電極5—8 ns流注圖像Fig.15.Photos of streamer discharge from 5 ns to 8 ns with 250μm radius of electrode tip.

4 結 論

本文針對針-板電極結構空氣中流注放電現象,采用三維粒子仿真模型,確定計算域和相關計算參數,調整電壓幅值、氣體組分和電極尖端形狀尺寸參數并進行相關仿真,分析各不同工況下的初始電離云及后續流注發展形態結構,得出以下結論.

1)放電初始階段,電極尖端附近形成近似為球體的初始電離云,具有很高的電子數密度.不同工況下的初始電離云半徑不同,其發展為多個單獨流注通道的時間也不同.后期流注的形態發展方向具有很大的隨機性.

2)隨電壓幅值增大,電離云的半徑越大,電壓工況為4.0 kV和4.5 kV時流注半徑無明顯差異;流注放電所形成的樹枝狀分叉結構也隨著電壓的增高變多.

3)N2-O2氣體混合物的比例顯著影響流注放電形貌.隨N2-O2混合物中O2濃度的增大,初始電離云的半徑逐步增大,分叉結構和數目明顯變化.

4)改變電極尖端形狀尺寸對放電初期影響較大,而在放電后期,電極尖端的改變對流注發展無明顯影響,這主要是因為流注后期發展過程中自洽電場起主導作用.

[1]Nijdam S 2011Experimental Investigations on the Physics of Streamers(Eindhoven:Technische Universiteit Eindhoven)pp2–4

[2]Nijdam S,van Veldhuizen E,Peter B 2012Plasma Chemistry&Catalysis in Gases&Liquids(Germany:WILEY-ICH)pp1–44

[3]Zhang C,Tarasenko V F,Shao T,Beloplotov D V,Lomaev M I,Wang R,Sorokin D A,Yan P 2015Phys.Plasmas22 033511

[4]Raether H 1939Z.Phys.112 464

[5]Loeb L B,Meek J M 1940J.Appl.Phys.11 438

[6]Briels T M P,van Veldhuizen E M,Ebert U 2008J.Phys.D:Appl.Phys.41 234008

[7]Nijdam S,Moerman J S,Briels T M P,van Veldhuizen E M,Ebert U 2008Appl.Phys.Lett.92 101502

[8]PengQ J2012Ph.D.Dissertation(Chongqing:Chongqing University)(in Chinese)[彭慶軍 2012博士學位論文(重慶:重慶大學)]

[9]Kulikovsky A A 1997J.Phys.D:Appl.Phys.30 441

[10]Luque A,Ebert U 2014New J.Phys.16 013039

[11]Li Y D,Wang R P,Zhang Q G,Zhou Y,Wang H G,Liu C L 2011IEEE Trans.Plasma Sci.39 2226

[12]Sun A B,Teunissen J,Ebert U 2013Geophys.Res.Lett.40 2417

[13]Sun A B,Teunissen J,Ebert U 2014J.Phys.D:Appl.Phys.47 445205

[14]Teunissen J,Sun A B,Ebert U 2014J.Phys.D:Appl.Phys.47 365203

[15]Li C,Ebert U,Hundsdorfer W 2010J.Comput.Phys.229 200

[16]Lxcat,http://lxcat net laplace univ-tlse fr[2017-3-22]

[17]Sun A B,Li H W,Xu P,Zhang G J 2017Acta Phys.Sin.66 195101(in Chinese)[孫安邦,李晗蔚,許鵬,張冠軍2017物理學報66 195101]

[18]Teunissen J,Ebert U 2016Plasma Sources Sci.Technol.25 044005

[19]Sun A B,Becker M M,Loffhagen D 2016Comput.Phys.Commun.206 35

[20]Teunissen J,Ebert U 2014J.Comput.Phys.259 318

[21]Nijdam S,van de Wetering F M J H,Blanc R,van Veldhuizen E M,Ebert U 2010J.Phys.D:Appl.Phys.43 145204