大氣壓甲烷針-板放電等離子體中粒子密度和反應路徑的數(shù)值模擬?

趙曰峰1) 王超1)2) 王偉宗3) 李莉4) 孫昊2) 邵濤2) 潘杰1)

1)(山東師范大學物理與電子科學學院,濟南 250014)

2)(中國科學院電工研究所,北京 100190)

3)(比利時安特衛(wèi)普大學化學學院,安特衛(wèi)普 2610,比利時)

4)(國家電網(wǎng)濟南供電公司,濟南 250012)

(2017年10月10日收到;2018年2月11日收到修改稿)

1 引 言

甲烷作為天然氣(氣田氣、油田氣與頁巖氣等)、沼氣和礦井氣的主要成分,不僅是制造氫氣、炭黑、甲醛、乙炔和一氧化碳等物質(zhì)的基本原料,而且是民用與工業(yè)加熱的重要燃料,在供熱取暖、交通運輸、燃氣發(fā)電、合成氣制取和農(nóng)作物施肥等領域得到了廣泛的關注與應用[1,2].隨著我國原油資源的日益緊缺,采用重油加氫技術處理高硫、高殘?zhí)俊⒏呓饘俚牧淤|(zhì)渣油,將原油中占比較高的重油轉(zhuǎn)化為輕質(zhì)油,從而優(yōu)化能源結構并降低污染物排放,已經(jīng)成為國民經(jīng)濟和社會發(fā)展中迫切需要解決的關鍵科技問題[3].傳統(tǒng)的重油加氫技術由于在高溫高壓條件下使用臨氫環(huán)境和催化劑,因而存在能耗高、輕油收率低和催化劑易結焦失活等缺點.千伏級高壓電源驅(qū)動的針-板放電能夠在寬廣的氣壓范圍內(nèi)產(chǎn)生高化學活性與非平衡性的甲烷放電等離子體,適宜開展放電等離子體甲烷轉(zhuǎn)化的科學研究和工業(yè)應用[4?12].由甲烷針-板放電與重油加氫反應耦合形成的甲烷轉(zhuǎn)化重油加氫技術,可以實現(xiàn)重油高效加氫并增產(chǎn)高附加值低碳烯烴,具有廣闊的應用前景和重要的研究意義.

比利時安特衛(wèi)普大學Bogaerts課題組是較早對放電等離子體甲烷轉(zhuǎn)化進行系統(tǒng)研究的科研團隊.2013年,該課題組整合實驗診斷數(shù)據(jù)與數(shù)值模擬結果,分析甲烷氮氣混合氣體介質(zhì)阻擋放電(dielectric barrier discharge,DBD)等離子體中氮氣濃度對甲烷轉(zhuǎn)化制氫的產(chǎn)額、機理和動力學過程的影響[13].2015年,該課題組與荷蘭埃因霍溫大學van Dijk等[3]開展合作,探討甲烷-氧氣和甲烷-二氧化碳混合氣體DBD等離子體中甲烷轉(zhuǎn)化產(chǎn)生合成氣以及產(chǎn)生高分子碳氫化合物和含氧化合物的主導路徑.2017年,美國德雷克塞爾大學Fridman等[1]使用納秒脈沖DBD和大氣壓輝光放電(atmospheric pressure glow discharge,APGD)等離子體活化甲烷,將氣態(tài)甲烷分子直接并入液態(tài)碳氫化合物燃料,從而實現(xiàn)了非平衡放電等離子體作用下的甲烷液化.上述研究表明,非平衡放電等離子體技術為甲烷轉(zhuǎn)化提供了新的發(fā)展方向,但在已經(jīng)公開發(fā)表的放電等離子體甲烷轉(zhuǎn)化工作中大多采用DBD和APGD的放電形式,鮮有針對大氣壓甲烷針-板放電等離子體的實驗診斷和數(shù)值模擬研究報道.

電子溫度和粒子密度的時空演化特性與反應路徑貢獻是非平衡等離子體及其應用研究方向的長期學術熱點[14?21].對大氣壓甲烷針-板放電等離子體而言,由于相應的原位實時診斷方法尚不完善并且所適用的粒子種類有限,數(shù)值模擬成為探索上述問題的必要研究方法.天津大學Wang等[22]使用零維全局模型,研究等離子體能量和放電氣隙間距對大氣壓DBD微等離子體反應器中甲烷轉(zhuǎn)化產(chǎn)生高分子碳氫化合物的影響.美國德克薩斯大學奧斯汀分校Raja和Levco[23]結合全局模型與流體模型,研究在正極性和負極性觸發(fā)電壓作用下,甲烷流注放電中擊穿機理、等離子體組分和活性粒子生成的共性與區(qū)別.然而,現(xiàn)有的數(shù)值模擬通常采用零維全局模型或一維流體模型[24,25],缺乏與實驗診斷結果相符合的二維模擬以及對大氣壓甲烷針-板放電等離子體中活性粒子的產(chǎn)生、輸運、時空分布和反應路徑等問題的全面分析與合理闡釋.

本文對負直流高壓電源激勵下的大氣壓甲烷針-板放電等離子體建立二維流體模型并進行實驗驗證,開展數(shù)值模擬研究,呈現(xiàn)電場強度、電子溫度和粒子密度的空間與軸向分布,計算產(chǎn)生帶電粒子和中性粒子的反應貢獻,探討生成各種粒子的關鍵路徑,為調(diào)控大氣壓甲烷針-板放電等離子體的產(chǎn)物選擇性、提高甲烷轉(zhuǎn)化重油加氫技術的能量利用效率,以及優(yōu)化放電等離子體甲烷轉(zhuǎn)化與重油加氫反應的耦合時間和作用方式等問題提供依據(jù).

2 模型描述

大氣壓甲烷針-板放電等離子體二維流體模型的幾何結構如圖1所示,其中直徑為2 mm的圓柱金屬針(材料為鎢)電極位于平行金屬板(材料為銅)電極上方,針電極的長度為3 mm,針尖的曲率半徑為1 mm,針電極尖端距離板電極表面的距離為2 mm,電壓振幅為5 kV的負直流高壓電源經(jīng)過阻值為10 k?的限流電阻接針電極,板電極接地.求解區(qū)域大小為10 mm×5 mm.

圖1 大氣壓甲烷針-板放電等離子體二維流體模型的幾何結構Fig.1.Geometry structure of two-dimensional fluid model of methane needle-plane discharge plasma at atmospheric pressure.

工作氣體甲烷在求解區(qū)域內(nèi)設為均勻分布,氣體壓強設為1 atm(1 atm=101325 Pa),氣體溫度設為300 K.除基態(tài)甲烷分子CH4以外,數(shù)值模型中包含的21種帶電和中性粒子有e,CH,CH2,CH3,C2H2,C2H4,C2H5,C2H6,C3H8,H和H2,粒子種類的選擇參照文獻[3,13].甲烷放電等離子體中的等離子體物理化學反應及其反應速率系數(shù)如表1所列,其中電子碰撞反應R1—R16的反應速率系數(shù)由碰撞截面計算,R1,R2,R3,R8和R13的碰撞截面取自Bolsig+軟件,R4—R7,R9—R12和R14—R16的碰撞截面取自LXcat數(shù)據(jù)庫,其他反應的反應速率系數(shù)取自參考文獻[26—29].

表1 甲烷放電等離子體中的等離子體物理化學反應及其反應速率系數(shù)Table 1.Plasma physical and chemical reactions and the corresponding reaction rate coefficients in methane discharge plasma.

流體模型中的基本控制方程主要包括連續(xù)性方程、泊松方程和電子能量平衡方程[14,23].甲烷放電等離子體中各種粒子的動態(tài)行為滿足如下連續(xù)性方程:

式中nk為第k種粒子的粒子數(shù)密度;Γ為粒子通量;S為粒子源項.

在連續(xù)性方程中,帶電粒子通量為Γk=μknkE?Dk?nk;中性粒子通量僅考慮擴散項,表示為Γk=?Dk?nk.其中E為電場強度,μ和D分別為粒子的漂移率和擴散系數(shù).

電場強度E由泊松方程描述如下:

式中φ為電位;ε0為真空介電常數(shù);εr為氣體的相對介電常數(shù);e為電子電荷;qk為第k種粒子的電量.

電子溫度Te由電子能量平衡方程描述如下:

式中KB為玻爾茲曼常數(shù);下標e表示電子;ψe為電子能量通量;me為電子質(zhì)量;ne為電子數(shù)密度;M為氣體的平均分子質(zhì)量;vel為彈性碰撞頻率;Tg為氣體溫度;ri和?εi分別為第i次反應的速率和由非彈性碰撞引起的能量損失.

在針電極和板電極處,電子通量Γe和電子能量通量ψe滿足如下邊界條件:

式中n為向外單位向量;re為反射系數(shù);νe,th為電子熱速度;γk為二次電子發(fā)射系數(shù);Γk為第k種粒子的粒子通量;Γt為熱發(fā)射通量;?k為第k種粒子的熱離子能量;?為平均熱離子能量.

使用Comsol 5.0軟件進行網(wǎng)格剖分和數(shù)值運算,在求解區(qū)域內(nèi)使用軸對稱結構和不均勻三角形網(wǎng)格,剖分結果如圖2所示,其中生成的網(wǎng)格數(shù)量約為2.3×105,總自由度個數(shù)約為3.7×106.此外,初始電子密度設為1.0×1016m?3,初始正離子密度總和設為1.0×1016m?3,中性粒子密度的初始值設為1.0×1011m?3,電子溫度的初始值設為3 eV.需要指出的是,對于大氣壓針-板放電等離子體,上述帶電粒子密度和電子溫度的初始值設置很高.

圖2 求解區(qū)域內(nèi)使用軸對稱結構和不均勻三角形網(wǎng)格的剖分結果Fig.2. Subdivision results in calculation domain using axial symmetry structure and asymmetry triangle mesh.

為了驗證建立的二維流體模型的正確性,使用Stark展寬法估算大氣壓甲烷針-板放電等離子體中的電子密度[30].需要指出的是,對于大氣壓針-板放電,Stark展寬法估算等離子體中的電子密度存在一定誤差.圖3為在二維流體模型中使用的放電條件下,開展實驗測量得到的針電極下方區(qū)域的Hα譜線,測量時光纖垂直于放電通道放置,與放電等離子體反應器外邊沿的距離為2 cm,將圖中半高全寬區(qū)域的Stark展寬代入可得相應位置的電子密度數(shù)量級為1021m?3.

圖3 Hα譜線Fig.3.Hαline.

圖4 給出了外施電源的作用時間為6 ns時刻電子密度的空間分布的數(shù)值模擬結果,表明大氣壓甲烷針-板放電等離子體具有接近“等離子體子彈”的空間分布特性.從圖中可以看出,在針電極尖端的下方區(qū)域內(nèi),電子密度最大值超過1.0×1020m?3,表明數(shù)值模擬結果與實驗診斷數(shù)據(jù)基本符合,因此本文呈現(xiàn)的大氣壓甲烷針-板放電等離子體中電場強度、電子溫度和粒子密度的空間與軸向分布以及反應路徑貢獻均為6 ns時刻的數(shù)值模擬結果.

圖4 電子密度的空間分布Fig.4.Spatial distribution of electron densities.

3 數(shù)值模擬結果與討論

圖5 (a)電場強度和(b)電子溫度的空間分布Fig.5.Spatial distributions of(a)electric intensity and(b)electron temperature.

圖5 為大氣壓甲烷針-板放電等離子體中的電場強度和電子溫度的空間分布,圖5(a)中的箭頭表示電場強度矢量的方向.如圖所示,在半球形針電極頭部的外側(cè)邊緣(與針電極外側(cè)的距離小于0.1 mm)區(qū)域,電場強度和電子溫度均達到最大值,分別高于5×106V/m和4.5 eV.已有研究指出,大氣壓甲烷DBD等離子體中電子密度和電子溫度最大值分別約為5.5×1019m?3和3 eV[13],與之相比,大氣壓甲烷針-板放電等離子體有更高的電子密度和電子溫度最大值.此外,在由針電極尖端指向板電極中心位置的軸向方向,電場強度和電子溫度的變化趨勢均為在與針電極尖端的距離約為0.1—0.6 mm的區(qū)域內(nèi)先減小后增大、在與針電極尖端的距離大于0.6 mm的區(qū)域內(nèi)再逐漸減小.上述電場強度和電子溫度的數(shù)值模擬結果表明,在負直流高壓電源的激勵和空間電荷的作用下,大氣壓甲烷針-板放電的等離子體鞘層不明顯,在靠近陰極的區(qū)域內(nèi)能夠形成陰極暗區(qū)和負輝區(qū),在遠離陰極的區(qū)域內(nèi)則主要由正柱區(qū)構成[31].

圖6 帶電粒子密度的空間分布Fig.6.Spatial distributions of charged particle densities,including(a)

圖7 帶電粒子密度的軸向分布Fig.7.Axial distributions of charged particle densities,including

圖8和圖9分別為數(shù)值模擬結果中密度較高的6種中性粒子(CH,CH2,CH3,C2H4,C2H5和H)密度的空間分布和軸向分布.從圖中可以看出,甲基CH3和氫原子H是大氣壓甲烷針-板放電等離子體中粒子密度最高的兩種中性粒子,其粒子密度最大值達到1020m?3數(shù)量級.并且,CH3和H的粒子密度空間分布和軸向分布幾乎相同.與甲烷-二氧化碳混合氣體DBD等離子體的數(shù)值模擬結果相比,甲烷針-板放電等離子體中的CH3密度更高,與H密度較為接近[3].亞甲基CH2、乙烯C2H4和乙基C2H5是另外3種粒子密度較高的中性粒子,其粒子密度分布在靠近陰極的區(qū)域內(nèi)明顯不同,最大值分別處于1019,1018和1017m?3數(shù)量級,而在正柱區(qū)內(nèi)則呈現(xiàn)均勻且相像的分布規(guī)律,粒子密度分別約為5.6×1015,7.2×1015和3.7×1015m?3.在與針電極尖端的距離約為0.1 mm處,次甲基CH具有與C2H5接近的粒子密度最大值(約為2.1×1017m?3),但在正柱區(qū)內(nèi)CH的粒子密度降低至1012m?3數(shù)量級.此外,CH粒子的軸向分布與其他中性粒子有顯著區(qū)別,而與帶電粒子更為相似,是本文模擬結果中僅有的在與針電極尖端距離約為0.1—0.6 mm的區(qū)域內(nèi)粒子密度先減小后增大的中性粒子.

圖8 中性粒子(a)CH,(b)CH2,(c)CH3,(d)C2H4,(e)C2H5和(f)H密度的空間分布Fig.8.Spatial distributions of neutral particle densities,including(a)CH,(b)CH2,(c)CH3,(d)C2H4,(e)C2H5and(f)H.

圖9 中性粒子CH,CH2,CH3,C2H4,C2H5和H密度的軸向分布Fig.9.Axial distributions of neutral particle densities,including CH,CH2,CH3,C2H4,C2H5and H.

圖10 為大氣壓甲烷針-板放電等離子體中帶電粒子的反應路徑,其中實線表示生成箭頭所指向的帶電粒子的惟一反應或主要反應,虛線表示生成箭頭所指向的帶電粒子的次要反應,紅色和黑色的線段分別表示是否有工作氣體甲烷參與的反應,線段起點和線段中注明的粒子表示參與反應的粒子,線段終點的粒子表示該反應生成的粒子,線段中注明的百分比表示該反應產(chǎn)額占生成該粒子的所有反應總產(chǎn)額的百分比,括號內(nèi)的序號為該反應在表1中的反應序號.如圖10所示,電子與基態(tài)甲烷分子CH4間的電子碰撞電離反應R5和R4是生成的主要路徑.分別與CH4發(fā)生分子碰撞解離反應R18和R17生成大部分的電子與氫氣分子H2間的電子碰撞電離反應R10生成發(fā)生分子碰撞解離反應R20生成分別與CH4和C2H4分子發(fā)生反應R21和R23生成少量的與乙烷分子C2H6間的反應R19也生成少量的由于電子碰撞電離反應R5和R4是生成的主導反應,而電場強度和電子溫度為在與針電極尖端距離約為0.1—0.6 mm的區(qū)域內(nèi)先減小后增大,因此和的粒子密度在該區(qū)域內(nèi)也呈現(xiàn)先減小后增大的變化規(guī)律.

圖10 帶電粒子的反應路徑Fig.10.Reaction pathway of charged particles.

圖11 碳氫化合物分子的反應路徑Fig.11.Reaction pathway of hydrocarbon molecules.

生成碳氫化合物分子的反應路徑如圖11所示,其中藍色和黑色的線段分別表示是否有電子參與的反應,其他的設置與圖10相同.圖11表明,電子與基態(tài)CH4分子間的電子碰撞分解反應R6,R7和R8分別是生成甲基分子CH3、亞甲基分子CH2和次甲基分子CH的主要路徑.CH2與CH4分子間的反應R30也生成CH3,該反應產(chǎn)額占CH3分子總產(chǎn)額的12.68%.CH與CH4分子發(fā)生反應R29生成乙基分子C2H5.CH3分別與CH3和C2H5分子發(fā)生反應R24和R27生成乙烷分子C2H6和丙烷分子C3H8.盡管有多條反應路徑可以產(chǎn)生乙烯分子C2H4,但CH2與CH4分子間的反應R31是生成C2H4的主導反應.電子與C2H4分子間發(fā)生的電子碰撞分解反應R14是生成乙炔分子C2H2的關鍵路徑.

產(chǎn)生氫原子H的反應路徑有R5,R6,R8,R9,R12,R20,R32和R33,其中電子與CH4間的電子碰撞分解反應R6的反應產(chǎn)額比其他反應的總產(chǎn)額高105數(shù)量級,是生成H的主導反應.反應R7,R11,R14,R18,R19,R21,R22,R23,R28,R31,R34和R35產(chǎn)生氫分子H2,其中電子與CH4間的電子碰撞分解反應R7和CH2與CH4發(fā)生的反應R31的產(chǎn)額各占H2總產(chǎn)額的52.15%和47.85%,是生成H2的主要路徑.

4 結 論

對由負直流高壓電源驅(qū)動的大氣壓甲烷針-板放電等離子體建立二維流體模型并進行實驗驗證,數(shù)值模擬研究其中的粒子密度和反應路徑.模擬結果表明,大氣壓甲烷針-板放電能夠形成陰極暗區(qū)、負輝區(qū)和正柱區(qū).粒子密度的最大值都超過1018m?3數(shù)量級.和密度在沿軸向達到最大值后逐漸減小.和的粒子密度軸向演化與電場強度和電子溫度的軸向演化接近.的粒子密度最大值出現(xiàn)在針電極尖端下方約0.1 mm處,在其他區(qū)域內(nèi)的粒子密度分布較為均勻.CH3和H的粒子密度最大值達到1020m?3數(shù)量級,且空間分布和軸向分布幾乎相同.CH2,C2H4和C2H5的粒子密度分布在靠近陰極的區(qū)域內(nèi)明顯不同而在正柱區(qū)內(nèi)較為相像.CH的粒子密度最大值與C2H5接近,但在正柱區(qū)內(nèi)粒子密度降低至1012m?3數(shù)量級.電子與CH4間的電子碰撞電離是生成和的主要路徑.分別與CH4發(fā)生分子碰撞解離生成電子與H2間的電子碰撞電離生成發(fā)生分子碰撞解離生成電子與CH4間的電子碰撞分解是生成CH3,CH2,CH和H的主導反應.CH與CH4反應生成C2H5.CH3與CH3和C2H5反應分別生成C2H6和C3H8.CH2與CH4和電子與C2H4發(fā)生的反應分別是生成乙烯C2H4和乙炔C2H2的關鍵路徑.電子與CH4間的電子碰撞分解反應和CH2與CH4發(fā)生的反應的產(chǎn)額各占H2總產(chǎn)額的52.15%和47.85%,是生成H2的主要路徑.

[1]Liu C,Chernets I,Ji H,Smith J,Rabinovich A,Dobrynin D,Fridman A 2017IEEE Trans.Plasma Sci.45 683

[2]Kang H,Lee D,Kim K,Jo S,Pyun S,Song Y,Yu S 2016Fuel Process.Technol.148 209

[3]Bie C D,van Dijk J,Bogaerts A 2015J.Phys.Chem.C119 22331

[4]Xu Y,Zhang X,Yang C,Zhang Y,Yin Y 2016Plasma Sci.Technol.18 1012

[5]Wang C,Zhang Z,Cui H,Xia W,Xia W 2017Chin.Phys.B26 085207

[6]Liu J L,Park H W,Chung W J,Park D W 2016Plasma Chem.Plasma Proc.36 437

[7]Zhang Z B,Wu Y,Jia M,Song H M,Sun Z Z,Li Y H 2017Chin.Phys.B26 065204

[8]Wang W,Snoeckx R,Zhang X,Cha M S,Bogaerts A 2018J.Phys.Chem.CDOI:10.1021/acs.jpcc.7b10619

[9]Bi Z H,Hong Y,Lei G J,Wang S,Wang Y N,Liu D P 2017Chin.Phys.B26 075203

[10]Zhang D Z,Wang Y H,Wang D Z 2017Chin.Phys.B26 065206

[11]Shao T,Wang R X,Zhang C,Yan P 2018High Voltage3 14

[12]Gao Y,Zhang S,Liu F,Wang R X,Wang T L,Shao T 2017Trans.China Electrotechnical Soc.32 61(in Chinese)[高遠,張帥,劉峰,王瑞雪,汪鐵林,邵濤 2017電工技術學報32 61]

[13]Snoeckx R,Setareh M,Aerts R,Simon P,Maghari A,Bogaerts A 2013Int.J.Hydrogen Energy38 16098

[14]Pan J,Li L 2015J.Phys.D:Appl.Phys.48 055204

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

[16]Pan J,Li L,Wang Y,Xiu X,Wang C,Song Y 2016Plasma Sci.Technol.18 1081

[17]Niu Z T,Zhang C,Ma Y F,Wang R X,Chen G Y,Yan P,Shao T 2015Acta Phys.Sin.64 195204(in Chinese)[牛宗濤,章程,馬云飛,王瑞雪,陳根永,嚴萍,邵濤2015物理學報64 195204]

[18]Pan J,Li L,Chen B,Song Y,Zhao Y,Xiu X 2016Eur.Phys.J.D70 136

[19]Babaeva N Y,Zhang C,Qiu J,Hou X,Tarasenko V F,Shao T 2017Plasma Sources Sci.Technol.26 085008

[20]Yang D P,Li S Y,Jiang Y F,Chen A M,Jin M X 2017Acta Phys.Sin.66 115201(in Chinese)[楊大鵬,李蘇宇,姜遠飛,陳安民,金明星2017物理學報66 115201]

[21]Yang W B,Zhou J N,Li B C,Xing T W 2017Acta Phys.Sin.66 095201(in Chinese)[楊文斌,周江寧,李斌成,邢廷文2017物理學報66 095201]

[22]Wang B,Yan W,Ge W,Duan X 2013Chem.Eng.J.234 354

[23]Levko D,Raja L L 2017Plasma Sources Sci.Technol.26 035003

[24]Yin Z Q,Wang Y,Zhang P P,Zhang Q,Li X C 2016Chin.Phys.B25 125203

[25]Wang Q,Yu X L,Wang D Z 2017Chin.Phys.B26 035201

[26]Herrebout D,Bogaerts A,Yan M,Gijbels R,Goedheer W,Dekempeneer E 2001J.Appl.Phys.90 570

[27]Lefkowitz J K,Guo P,Rousso A,Ju Y 2015Phil.Trans.R.Soc.373 20140333

[28]Adamovich I V,Li T,Lempert W R 2015Phil.Trans.R.Soc.373 20140336

[29]Takana H,Nishiyama H 2014Plasma Sources Sci.Technol.23 034001

[30]Nikiforov A Y,Leys C,Gonzalez M A,Walsh J L 2015Plasma Sources Sci.Technol.24 034001

[31]Yao C W,Ma H C,Chang Z S,Li P,Mu H B,Zhang G J 2017Acta Phys.Sin.66 025203(in Chinese)[姚聰偉,馬恒馳,常正實,李平,穆海寶,張冠軍 2017物理學報 66 025203]