等離子射流近場特性的數(shù)值模擬*

張 琦，余永剛，鄭四發(fā)

(1.清華大學(xué)汽車研究院(相城)，江蘇蘇州 215134；2.南京理工大學(xué)能源與動力工程學(xué)院，江蘇南京 210094)

1 引 言

等離子體是由帶電的正、負(fù)粒子組成的集合體，其中正、負(fù)電荷電量相等，宏觀上呈電中性，為物質(zhì)存在的第四態(tài)。隨著等離子體物理技術(shù)的不斷發(fā)展，它的應(yīng)用范圍也越來越廣泛，包括軍事、航空航天、冶金、噴涂、金屬切割等領(lǐng)域。

等離子射流的自身特性及其與工質(zhì)相互作用特性是等離子射流點火技術(shù)中的基礎(chǔ)科學(xué)問題，國內(nèi)外學(xué)者對此進行了大量研究。國外等離子射流研究主要集中在放電能量較高的情況，如Kim等人[1]在3.1 kJ左右放電能量下對電熱毛細(xì)管等離子射流在大氣中某時刻的膨脹特性做了實驗研究，發(fā)現(xiàn)等離子體溫度和電子密度在經(jīng)過馬赫盤處出現(xiàn)跳躍。Kohel等人[2]實驗研究了3 kJ放電能量條件下毛細(xì)管等離子射流在大氣中自由膨脹的序列過程特性，認(rèn)為脈沖等離子射流為欠膨脹等離子射流，射流中有弓形激波和高亮的馬赫盤形成。Taylor[3]研究了約30 kJ放電能量條件下電爆炸絲電爆炸的細(xì)節(jié)過程，以及等離子射流在大氣中的自由膨脹過程及其溫度、速度等參量的分布情況。Chang等人[4]研究了在3 kJ放電能量條件下，放電脈沖長度對等離子射流在大氣中擴展以及沖擊平板特性的影響，并認(rèn)為放電脈沖長度可能是進行優(yōu)化點火的一個關(guān)鍵參數(shù)。Das等人[5]研究了放電能量為0.6 kJ、不同材料的毛細(xì)管和電爆炸絲條件下，等離子射流沖擊平板過程中的輻射傳熱特性。國內(nèi)劉東堯等人[6]研究了毛細(xì)管等離子射流壓力隨時間變化的規(guī)律，認(rèn)為等離子體及其射流壓力幅值取決于放電能量，射流產(chǎn)生的脈沖推力具有沖擊震蕩效應(yīng)。趙文華等人[7-8]對電弧等離子射流核脈動的實驗研究發(fā)現(xiàn)，射流并不存在一個處于穩(wěn)定狀態(tài)的核心區(qū)域，他們認(rèn)為核心區(qū)是脈動性的，而這是由電弧電壓脈動及電弧分流現(xiàn)象共同造成的。余永剛等人[9-12]研究了燃?xì)馍淞饕约暗入x子射流在漸擴型及圓柱型藥室中的擴展特性。張琦等人[13]采用高速錄像記錄了小放電能量(100 J以下)條件下等離子射流在大氣中的擴展過程，發(fā)現(xiàn)等離子射流徑向尺寸出現(xiàn)先增大后減小的現(xiàn)象。本研究在前人工作的基礎(chǔ)上，針對文獻(xiàn)[13]中的某一實驗工況進行數(shù)值模擬，分析小放電能量條件下等離子射流近場的壓力、速度、溫度以及形態(tài)演化特性，并分析文獻(xiàn)[13]中等離子射流徑向尺寸出現(xiàn)先增大后減小現(xiàn)象的原因。

2 物理模型

由等離子射流在大氣中擴展的特點，對物理過程做如下假設(shè):(1) 等離子射流在大氣中的膨脹擴展過程可以看作二維軸對稱非穩(wěn)態(tài)過程，采用RNGk-ε湍流模型模擬流場中的湍流效應(yīng)。(2) 小放電能量條件下，低溫等離子射流溫度約為8 000 K左右，電離度較小，將等離子體混合物近似作為理想可壓氣體。(3) 忽略等離子體與大氣中其它成分的化學(xué)反應(yīng)。(4) 忽略等離子體的電磁作用力、質(zhì)量力、體積力等次要因素的影響。(5) 等離子體處于局部熱力學(xué)平衡，光學(xué)厚度較薄，忽略輻射效應(yīng)。

3 數(shù)學(xué)模型

根據(jù)第2節(jié)的物理模型，建立如下數(shù)學(xué)模型。

(1) 連續(xù)性方程

?α2/?t+v·α2=Sα2/ρ2

(1)

式中：v為速度矢量;t為時間;模型中有兩個體積相，分別為大氣相和等離子體相，其中α1、ρ1為大氣相的體積分?jǐn)?shù)和密度，α2、ρ2是等離子體相的體積分?jǐn)?shù)和密度；Sα2為源項，由于不考慮化學(xué)反應(yīng)，該項為零。

(2) 動量方程

?(ρv)/?t+·(ρvv)=-p+·[μ(v+vT)]+ρg

(2)

(3) 能量方程

?(ρE)/?t+·[v(ρE+p)=·(keffT)+Sh

(3)

(4) 狀態(tài)方程

p=ρRT

(4)

式中：R為氣體常數(shù)。

(5) 湍流基本方程：包括湍流動能方程和湍流耗散率方程

式中：k為湍動能；ε為湍流耗散率；vi、vj是速度矢量分量，xi、xj是位移矢量分量；Gk代表由于平均速度梯度而產(chǎn)生的湍動能;Gb代表由浮力引起的湍動能;YM表示可壓縮湍流中波動擴張對整體耗散率的貢獻(xiàn);σk、σε為湍流普朗特數(shù)，σk=1.0，σε=1.3;湍流粘度μt=Cμk2/ε，其中Cμ為常數(shù)，Cμ=0.08；C1ε、C2ε、C3ε是湍流模型方程中的常數(shù)，C1ε=1.44，C2ε=1.92。

(6) 計算域及初邊界條件

由于流場為軸對稱結(jié)構(gòu)，取一半流場為計算域。為了減少計算網(wǎng)格，根據(jù)射流擴展特性取梯形計算域,網(wǎng)格均為結(jié)構(gòu)網(wǎng)格，入口處網(wǎng)格尺寸為0.33 mm，如圖1所示。壓力插值格式采用PRESTO!，密度、動量、湍動能、湍流耗散率及能量方程離散格式均采用二階迎風(fēng)格式，體積分?jǐn)?shù)采用Geo-reconstruct格式，時間離散格式采用一階隱式，離散方程的壓力和速度耦合采用SIMPLE方法。

開始時，計算區(qū)域內(nèi)還沒有等離子射流(α2=0)，因此初始參數(shù)即為環(huán)境參數(shù)。計算區(qū)域中，出口邊界參數(shù)可認(rèn)為與初始環(huán)境參數(shù)相同；入口條件選用壓力入口條件，即p=p(t),T=T(t)。

圖1 計算域及網(wǎng)格Fig.1 Computational domain and mesh grid

4 模擬結(jié)果及分析

圖2 不同網(wǎng)格尺寸下等離子射流軸向擴展位移Fig.2 Axial displacement of plasma jet under different grid sizes

采用文獻(xiàn)[13]中放電電壓為3 000 V、噴孔直徑為4 mm的實驗工況作為模擬工況。運用Fluent軟件模擬等離子射流在大氣中的擴展過程。

首先進行網(wǎng)格無關(guān)性驗證，分別對原網(wǎng)格(入口處網(wǎng)格尺寸為0.33 mm)及加密1.5倍后的網(wǎng)格模型進行計算，得到如圖2所示的等離子射流軸向擴展位移的對比曲線，可以看出，網(wǎng)格加密后對計算結(jié)果影響較小，說明原網(wǎng)格符合要求。

4.1 等離子射流擴展界面特性

圖3為噴嘴近場處等離子射流擴展界面的演化過程序列圖,可以看出，等離子射流從噴嘴噴出后，沿徑向(r)和軸向(x)擴展。0.097 ms時射流中出現(xiàn)卷吸。0.097～0.213 ms時，射流中心處擴展出現(xiàn)微小倒退，而中心射流外側(cè)依然向前擴展，形成分叉，射流中心出現(xiàn)空化現(xiàn)象。0.273 ms后，分叉射流逐漸向中心合并。射流擴展過程中頭部破碎，摻混較為劇烈。從圖3中還可以看出，在射流分叉出現(xiàn)過程中，射流的徑向尺寸先逐漸變大，分叉合并后，在中心形成一條主射流，前期形成的分叉流在主射流的作用下被拉伸，并向主射流靠攏，使射流整體徑向尺寸逐漸變小。這種變化趨勢與文獻(xiàn)[13]中實驗結(jié)果定性吻合。

圖3 噴嘴近場等離子-液體兩相界面演化過程Fig.3 Evolution of the plasma-liquid interface at the nozzle near field

圖4 等離子射流軸向擴展位移隨時間的變化曲線Fig.4 Axial displacement of plasma jet varying with time

由等離子射流界面演化序列圖可以處理得到等離子射流前端面的軸向擴展位移。圖4為等離子射流前端面軸向擴展位移隨時間變化的模擬曲線與文獻(xiàn)[13]中實測值的比較，二者定性吻合。

4.2 壓力分布特性

圖5為初期等離子射流近場處壓力分布序列圖。等離子射流前端與大氣相互作用時，兩相流體相互壓縮，產(chǎn)生壓力波向前推進。隨著射流的擴展，壓力沿軸向逐漸出現(xiàn)周期性脈動。射流分叉時，壓力也出現(xiàn)分叉特性，如0.153和0.213 ms時刻。

圖5 初期等離子射流場中壓力分布序列圖Fig.5 Pressure distribution of plasma jet at early time

圖6 不同時刻x軸上的壓力分布Fig.6 Pressure on the x-axis at different times

下面以0.022、0.097、0.153及0.633 ms時刻為例，取軸線上的壓力進行分析，如圖6所示。射流從噴嘴噴出時，壓力迅速下降至低于環(huán)境大氣壓力，其后又突然上升，并在該處出現(xiàn)激波。0.022 ms時，軸線上的壓力相對較大，激波位置離噴嘴較近。隨著射流的擴展，壓力變化逐漸減小，且沿軸向逐漸出現(xiàn)脈動性，如0.633 ms時壓力沿軸向呈現(xiàn)脈動性；激波位置也向下游運動，并基本固定在x=15 mm附近。各時刻激波上游的壓力均相等。

再看壓力沿徑向的分布特性，以x=1和10 mm處的徑向截面為例，分析不同時刻壓力在所取截面上的分布特性，如圖7所示?？傮w上看，壓力沿徑向迅速遞減，并逐漸趨于環(huán)境壓力。x=1 mm處的徑向截面上，初始(如0.012 ms時)壓力相對較大,0.052 ms后壓力較小，且在壓力突降段基本不隨時間變化，而在壓力突降段之外，壓力隨時間有較小的變化。x=10 mm處徑向截面處于射流場中的低壓區(qū)，壓力沿徑向先出現(xiàn)一段較小的衰減，再迅速遞增，最后趨于環(huán)境壓力；0.052 ms時，壓力相對較??；0.153 ms后，壓力衰減段和迅速上升段基本不隨時間變化。

結(jié)合圖6和圖7可以得出，隨著射流的擴展，流場中的壓力隨時間變化幅度逐漸變小，但具有波動性。

圖7 不同時刻x=1和10 mm處壓力的徑向分布Fig.7 Radial distribution of pressure at x=1 and 10 mm at different times

4.3 速度分布特性

圖8為初期等離子射流近場處的軸向速度(vx)分布序列圖。噴嘴前方逐漸形成一個馬赫盤，馬赫盤上游速度較大，下游速度較小。射流經(jīng)過馬赫盤時速度突然降低，逐漸形成分叉空心射流。分叉空心射流合并后，在馬赫盤下游一直存在一個錐形低速區(qū)。圖8中流速為負(fù)值的區(qū)域為回流區(qū)。

圖8 等離子射流場軸向速度分布Fig.8 Axial velocity distribution of plasma jet

圖9為初期等離子射流近場處的徑向速度(vr)分布序列圖。由圖9可知，軸線上的徑向速度較小，靠近噴嘴的近場處射流徑向擴展速度較大，并且擴展速度沿徑向先增大后減小。0.069 ms以前，噴嘴近場vr較大的區(qū)域逐漸增大;0.069～0.153 ms過程中，vr較大區(qū)域逐漸減小;0.153 ms后，該區(qū)域基本不變。

圖9 等離子射流場徑向速度分布Fig.9 Radial velocity distribution of plasma jet

圖10 不同時刻軸線上的馬赫數(shù)分布Fig.10 Ma on the x axis at different times

下面分析馬赫數(shù)(Ma)的時空分布特性。圖10為不同時刻軸線上的馬赫數(shù)分布。結(jié)合圖8可知,馬赫數(shù)突降處為一正激波。激波上游，馬赫數(shù)沿軸向迅速遞增。激波下游，初始階段馬赫數(shù)基本上沿軸向遞減；后期,如0.633 ms時，馬赫數(shù)沿軸向先脈動性遞增，再脈動性遞減。0.022 ms時激波位置離噴嘴較近，0.097 ms時激波位于x=14 mm處。激波上游馬赫數(shù)基本不隨時間變化。

再看馬赫數(shù)沿徑向的分布特性。圖11為0.633 ms 時軸線上不同位置處馬赫數(shù)的徑向分布，其中x=1 mm位于噴嘴近場，x=50 mm位于射流中部，x=110 mm位于射流頭部?？傮w上看，離噴嘴越近，初始馬赫數(shù)越大，沿徑向的衰減也越多。圖12 為x=10 mm處不同時刻馬赫數(shù)沿徑向的分布情況。0.052 ms時馬赫數(shù)較大，0.153 ms后，馬赫數(shù)隨著時間的增加而逐漸變小。

圖11 軸線上不同位置處馬赫數(shù)的徑向分布Fig.11 Radial distribution of Ma at different axial positions

圖12 不同時刻馬赫數(shù)沿徑向的分布Fig.12 Radial distribution of Ma at different times

4.4 溫度分布特性

圖13為初期噴嘴近場處等離子射流場中的溫度分布序列圖，總體上看,流場溫度在空間上具有非單調(diào)性。0.022 ms時，等離子射流剛開始進入大氣，與大氣相互壓縮時，形成了一個環(huán)形高溫區(qū)。隨著射流的擴展，噴嘴前方逐漸出現(xiàn)馬赫盤，射流經(jīng)過馬赫盤時受到強烈的壓縮，溫度迅速升高，如圖13中0.052 ms后的圖像所示。0.273 ms后處于射流分叉的合并過程，合并區(qū)域形成一個較大的高溫區(qū)。

圖13 等離子射流場中溫度分布Fig.13 Temperature distribution of plasma jet

下面以0.022、0.097、0.153及0.633 ms時刻為例，取軸線上的溫度進行量化分析，如圖14所示。經(jīng)過激波時，射流溫度迅速上升。0.022～0.153 ms，軸線上等離子射流與大氣的界面距離激波較近，導(dǎo)致射流溫度到達(dá)最高峰后，又迅速下降。0.633 ms時，兩相界面距離激波較遠(yuǎn)，射流溫度達(dá)到高峰后，沿軸線逐漸非單調(diào)下降，在兩相界面處迅速下降。根據(jù)不同時刻的溫度曲線，在激波下游，溫度隨時間呈非單調(diào)變化。

再看溫度沿徑向的分布特性。圖15(a)給出了0.633 ms時軸線上不同位置處溫度的徑向分布,可以看出，溫度沿徑向遞減；且從x=1 mm到x=110 mm過程中，溫度沿徑向的衰減率趨于平緩。圖15(b)為x=10 mm處，不同時刻等離子射流沿徑向的溫度分布圖。x=10 mm處位于正激波的上游，結(jié)合圖13可知，溫度沿徑向先逐漸升高，再迅速下降；0.153 ms時刻，在溫度下降過程中出現(xiàn)非單調(diào)性，是由空氣卷吸引起的。

圖14 不同時刻軸線上溫度分布Fig.14 Temperature of plasma jet on x-axis at different times

圖15 溫度沿徑向的分布Fig.15 Temperature distribution of plasma jet along radial direction

圖16(a)為Das等人[5]用高速CCD攝像機拍攝的0.12 ms時的等離子射流圖像，其中放電能量為0.6 kJ，噴嘴直徑為3.2 mm。圖16(b)為本研究模擬的0.097 ms時噴嘴近場的溫度分布圖，模擬實驗條件為：放電能量85 J，噴嘴直徑4 mm。相較于文獻(xiàn)[5]的實驗條件，模擬條件的放電能量較低，噴嘴直徑較大。對比圖16(a)和圖16(b)可以看出，雖然工況不同，但模擬結(jié)果與實驗照片中的射流結(jié)構(gòu)類似。

5 結(jié) 論

本研究通過數(shù)值模擬，得到如下結(jié)論：(1) 等離子射流在大氣中擴展時，噴嘴近場逐漸有馬赫盤形成，射流出現(xiàn)分叉、聚合現(xiàn)象。射流頭部破碎，摻混較為劇烈。(2) 等離子射流近場壓力波系由球形波逐漸演化成高、低壓相間的周期性波系。(3) 正激波上游各物理參數(shù)不隨時間改變，下游參數(shù)隨時間變化。(4) 等離子射流前端面擴展位移的計算值與實測值吻合較好。

圖16 實驗結(jié)果[5]和模擬結(jié)果的比較Fig.16 Comparison of experimental results[5] and simulation results

[1] KIM J U，SUK H.Characterization of a high-density plasma produced by electrothermal capillary discharge [J].Appl Phys Lett，2002，80(3):368-370.

[2] KOHEL J M，SU L K,CLEMENTS N T，et al.Emission spectroscopic measurements and analysis of a pulsed plasma jet [J].IEEE Trans Magn，1999，35(1):201-204.

[3] TAYLOR M J.Measurement of the properties of plasma from ETC capillary plasma generators [J].IEEE Trans Magn，2001，37(1):194-198.

[4] CHANG L M,HOWARD S L.Influence of pulse length on electrothermal plasma jet impingement flow:ARL-TR-4348 [R].Chantilly,Virginia:American Systems Corp,2007.

[5] DAS M，THYNEL S T，LI J，et al.Material dependence of plasma radiation produced by a capillary discharge:AIAA 2006-3425 [R].San Fransico，California:9th AIAA/ASME Joint Thermophysics and Heat Transfer Conference,2006.

[6] 劉東堯，郭海波，余永剛，等.毛細(xì)管消融等離子體及其射流壓力特性研究 [J].南京理工大學(xué)學(xué)報(自然科學(xué)版)，2009(3):348-351.

LIU D Y，GUO H B，YU Y G，et al.Investigation on the pressure of ablation plasma and its jet [J].Journal of Nanjing University of Science and Technology (Natural Science)，2009(3):348-351.

[7] 趙文華，唐皇哉，田 闊，等.存在空氣卷吸時等離子體射流光譜診斷應(yīng)做的修正 [J].光譜學(xué)與光譜分析，2004，24(4):388-391.

ZHAO W H，TANG H Z，TIAN K，et al.The correction to spectroscopic diagnostics of plasma jet with air engulfment [J].Spectroscopy and Spectral Analysis，2004，24(4):388-391.

[8] 趙文華，劉 笛，田 闊.電弧等離子體射流核脈動及射流卷吸的實驗研究 [J].核聚變與等離子體物理，2002，22(1):1-5.

ZHAO W H，LIU D，TIAN K.Experimental study on the fluctuation and engulfment of arc plasma spraying jet [J].Nuclear Fusion and Plsama Physics，2002,22(1):1-5.

[9] 余永剛，昌學(xué)霞，周彥煌，等.整裝式含能液體燃燒推進實驗研究及數(shù)值模擬 [J].工程熱物理學(xué)報，2008，29(3):531-534.

YU Y G，CHANG X X，ZHOU Y H，et al.Experiment and numerical simulation of combustion and propulsion of bulk-loaded energetic liquid [J].Journal of Engineering Thermophysics,2008,29(3):531-534.

[10] YU Y G，CHANG X X，ZHAO N.Study of bulk-loaded liquid propellant combustion propulsion processes with stepped-wall combustion chamber [J].J Appl Mech，2011，78(5):051001-051008.

[11] YU Y G，YAN S H，ZHAO N.Influence of boundary shape on interaction process of plasma jet and liquid media [J].Int J Appl Electromagn Mech,2010，33:541-548.

[12] 余永剛，閆善恒，陸 欣，等.雙股高速燃?xì)馍淞髟谝后w中擴展及相互作用的實驗研究 [J].科學(xué)技術(shù)與工程，2008，8(1):43-46.

YU Y G，YAN S H，LU X，et al.Experimental study on expansion process and interaction of high speed twin combustion-gas jet in liquid [J].Science Thechnology and Engineering,2008,8(1):43-46.

[13] 張 琦，余永剛，陸 欣，等.等離子射流在大氣中擴展特性實驗研究 [J].彈道學(xué)報，2010,22(4):97-101.

ZHANG Q，YU Y G，LU X，et al.Experimental study on plasma jet propagation in air [J].Journal of Ballistics，2010,22(4):97-101.