薛 創(chuàng)，寧 成，彭先覺

（北京應(yīng)用物理與計算數(shù)學(xué)研究所，北京 100094）

電爆炸絲是加載強電流脈沖的金屬絲在短時間內(nèi)產(chǎn)生大量焦耳熱從而發(fā)生爆炸的實驗技術(shù)[1-4]，它具有測量數(shù)據(jù)可靠、重復(fù)性好、操作簡便等優(yōu)點，在Z箍縮等離子體物理、爆炸力學(xué)、物性參數(shù)研究等方面有著重要應(yīng)用。其中，在爆炸與沖擊應(yīng)用方面，Antonov 等[5]利用40 根銅絲組成的球形絲陣在水中產(chǎn)生了2 TPa 的高壓力，Tang 等[6]采用電爆炸銅絲方法代替?zhèn)鹘y(tǒng)的起爆技術(shù)改善了柱面爆轟波的品質(zhì)。國內(nèi)早期關(guān)于電爆炸絲和電爆炸導(dǎo)體的實驗研究是測量電壓、電流等電參數(shù)[1-2]；隨著診斷技術(shù)的進步，發(fā)展為測量絲邊界運動軌跡、環(huán)境介質(zhì)中的沖擊波軌跡[3-4,7]，近年來進一步發(fā)展為測量絲等離子體的溫度和密度分布[8-10]，取得了大量研究成果。

在零維理論建模方面，Tucker[11]建立了比作用量模型，采用該模型可以計算金屬絲的電阻率變化和能量沉積，至今仍在實驗和理論研究中發(fā)揮作用[12-13]，但基于初始半徑和橫截面積計算電流密度，原則上僅適用于爆炸前的過程，需要合適的模型描述絲的膨脹過程，獲得隨時間變化的半徑。Bloomberg 等[14]建立了一個動力學(xué)模型來描述Z箍縮絲陣等離子體的形成過程，假設(shè)單絲等離子體的溫度、密度和壓力均勻分布，膨脹是自相似運動過程，則動能方程和內(nèi)能方程可簡化為兩個常微分方程，采用修正的理想氣體狀態(tài)方程和Spitzer 電阻率作為封閉條件，能夠獲得半徑的時間變化。郭軍等[13]建立的零維模型將金屬絲電爆炸過程分為膨脹前和膨脹后兩個階段，分別采用比作用量模型和動力學(xué)模型描述，其中膨脹后的動力學(xué)過程由動量方程與總能量方程描述，但能量方程形式與Bloomberg 的不同，內(nèi)能及總能量轉(zhuǎn)換過程中熱壓的貢獻不同。晁攸闖等[15]建立了一個描述水中電爆炸絲的零維模型，采用了數(shù)據(jù)庫形式的電導(dǎo)率，為了合理利用理想氣體狀態(tài)方程，采用了比作用量模型計算爆炸時刻的等離子體狀態(tài)。Rososhek等[16]建立了一個描述水中電爆炸絲膨脹及沖擊波傳輸?shù)牧憔S模型，但該模型依賴于實驗測量的電功率和一個調(diào)節(jié)能量沉積和做功份額的參數(shù)作為輸入條件，也采用了理想氣體狀態(tài)方程。

在一維理論和數(shù)值模擬方面，通常采用磁流體建模[17-18]。Chung 等[18]基于考慮液氣相變的實際氣體物態(tài)方程和半經(jīng)驗的電導(dǎo)率公式，給出了水中銅絲電爆炸過程的冷啟動物理建模，獲得了與實驗數(shù)據(jù)相符合的數(shù)值模擬結(jié)果。本文基于含人工黏性的差分方法研制一維拉格朗日程序SWEET (single wire electrically exploding test)，采用實際氣體狀態(tài)方程QEOS[18-19]模型和修正的Lee-More[20-23]電導(dǎo)率模型作為封閉條件，進行金屬絲電爆炸過程的冷啟動數(shù)值模擬。

為了方便理解和分析水中金屬絲電爆炸動力學(xué)過程和能量轉(zhuǎn)化過程，本文在一維數(shù)值模擬工作基礎(chǔ)上，依據(jù)Bloomberg 的基本假設(shè)發(fā)展了零維模型，采用實際氣體狀態(tài)方程和電導(dǎo)率作為封閉條件，采用Rososhek 的激波模型作為邊界條件，進行從電容器放電到金屬絲電爆炸的冷啟動及連續(xù)過程計算：首先，確定零維模型的控制方程和定解條件；其次，將計算結(jié)果與一維模型以及實驗結(jié)果進行比較，驗證零維模型基本假設(shè)的合理性。

1 物理模型

在一維柱對稱、單溫、單速度、單流體和磁擴散等近似下，等離子體的質(zhì)量方程、運動方程、內(nèi)能方程和磁場演化方程分別為：

式中：ρ 為質(zhì)量密度，v為徑向速度，t為時間，r為徑向位置，p為壓力，j為電流密度，B為磁感應(yīng)強度，ei為比內(nèi)能，η 為電阻率，時間導(dǎo)數(shù)都是拉格朗日隨體導(dǎo)數(shù)。

根據(jù)Bloomberg 的假設(shè)條件：密度沿徑向均勻分布，除邊界點外，空間偏導(dǎo)數(shù)為零，式(1)的左端密度不含空間位置信息，對其從0 到r積分，得到速度沿徑向的分布：

即速度沿徑向呈線性分布，如果絲等離子體邊界半徑為rp，邊界速度為vp，則速度分布為：

式(5)和式(6)聯(lián)立消去速度后，對時間積分的結(jié)果滿足質(zhì)量守恒條件。

對運動方程（式(2)）左右兩邊乘以2v就能得到動能方程，再對動能方程和內(nèi)能方程兩邊求空間體積積分，消去空間導(dǎo)數(shù)，分別得到零維模型的動能和內(nèi)能方程：

式中：m= ρSh為絲質(zhì)量，S=πrp2為絲等離子體柱的橫截面積，h為絲長度；eˉk=v2p/4 ，為單位質(zhì)量的平均動能；pˉ 和ēi分別為平均壓力和平均比內(nèi)能，在壓力和內(nèi)能均勻分布的假設(shè)下，可以去掉上標，即p=pˉ ，ei=ēi；p0為絲等離子體邊界處的壓力，等于水中的激波波后壓力；Rload為絲的電阻，在本文的研究范圍內(nèi)電流隨時間的上升速率較小，趨膚深度較大，絲半徑較小，電流密度均勻，因此Rload= ηh/S；I為電流強度，pB=μ0I2/(8πS)為磁壓力，μ0為真空磁導(dǎo)率，在本文單位制下μ0= 4π。

動能方程（式(7)）可以簡化為等離子體邊界的運動方程：

一維模型方程組（式(1)～(4)）及零維模型方程組（式(7)和式(8)）的求解需要狀態(tài)方程和電導(dǎo)率參數(shù)作為封閉條件，以密度、溫度作為基本熱力學(xué)變量，有：p=p(ρ,T)，ei=e(ρ,T)，η = 1/σ(ρ,T)，其中σ 為電導(dǎo)率。本文采用的銅介質(zhì)狀態(tài)方程和電導(dǎo)率參數(shù)隨溫度、密度的變化關(guān)系見圖1。

圖1 銅介質(zhì)物性參數(shù)以及絲等離子體的軌跡Fig. 1 Physical properties parameters of copper and trajectory of wire plasma

圖1 還給出了典型電爆炸絲等離子體的密度和溫度的狀態(tài)變化，可以看出從固態(tài)變?yōu)榈入x子體，密度變化了約4 個量級，溫度變化了約3 個量級。圖2 比較了電導(dǎo)率參數(shù)的插值公式和實驗結(jié)果。

圖2 銅介質(zhì)的電導(dǎo)率參數(shù)與實驗結(jié)果[22]比較Fig. 2 Comparison of the electrical conductivity for copper with experiment results[22]

計算磁壓力和焦耳熱還需要已知電流I，可以由實驗測量得到，本文根據(jù)等效電路模型計算得到。通常，脈沖功率裝置對金屬絲負載放電的等效電路如圖3 所示。

圖3 電爆炸金屬絲的等效電路模型Fig. 3 Circuit model for the electrically exploding wire

脈沖功率裝置的等效電容為C0，初始充電電壓為U0，放電過程中隨時間變化的電壓為Ubank，等效電阻為R0，等效電感為L0。金屬絲負載的等效電感[18]為Lload= 2h[ln(2h/rp)-0.75]，負載電感和電阻都隨著電爆炸過程而變化，通常在電爆炸過程中的電感變化較小，電阻的變化很大。裝置電壓Ubank、負載電壓Uload和電流I之間滿足電路方程：

水介質(zhì)環(huán)境下，等離子體的膨脹受到環(huán)境介質(zhì)的影響。假設(shè)邊界壓力p0與壓縮后的水介質(zhì)壓力pw平衡，如圖4 所示，根據(jù)波陣面的激波關(guān)系式和水介質(zhì)的質(zhì)量守恒條件可以獲得這個邊界壓力。

圖4 絲等離子體邊界及水中激波間斷示意圖Fig. 4 Schematic profile of the wire plasma interface and shock wave in the water

在隨激波波陣面運動的參考系中，水中激波的間斷關(guān)系為：

式中：ρw0= 1 g/cm3為水的波前密度，ρw為水的波后密度；pw0為水的波前壓力，pw為水的波后壓力；Sw為激波波陣面的速度，vw為波后水的速度。波前狀態(tài)已知，式(11)包含4 個未知量，還需要補充兩個方程。第一個是水的狀態(tài)方程，即壓力和密度的關(guān)系為：

式中新引入了兩個常數(shù)：參數(shù)n= 7.15，代表壓縮性，數(shù)值越大，表示越難壓縮；Cw0= 1 478 m/s，為水中的聲速。第二個方程是水介質(zhì)波后與絲等離子體界面位置之間的質(zhì)量守恒關(guān)系：

式中：rs為激波位置，r0為絲的初始半徑。聯(lián)立式(11)和式(12)消去波后壓力與水的速度，得到激波速度與密度之間的關(guān)系式：

式中：δ 為水的壓縮比。已知等離子體位置rp，就能根據(jù)式(13)和式(14)確定激波位置和速度；而等離子體邊界位置rp只需要知道式(12)的邊界壓力，求解控制方程（式(8)和式(9)）就能得到。

2 模擬應(yīng)用

選取Chung 等[18]的水下銅絲電爆炸實驗作為零維模型的應(yīng)用，等效電路參數(shù)為：電容C0= 5.22 μF，裝置電阻R0= 65.45 mΩ，裝置電感L0= 3 050 nH，充電電壓U0= 13 kV；銅絲參數(shù)為：長度h= 4 cm，初始半徑r0= 25, 50, 100 μm。首先給出50 μm 半徑下的基本動力學(xué)過程以及能量轉(zhuǎn)化關(guān)系，對比一維磁流體程序的模擬結(jié)果和實驗結(jié)果，然后改變絲參數(shù)，考察放電模式的變化。

絲等離子體的電壓及電流波形是實驗的直接測量量，見圖5。對比時間波形、峰值、峰值時刻等特征，本文的計算結(jié)果與實驗結(jié)果符合較好，而Tucker 模型[11]計算的電流在上升沿至爆炸時刻之前與本文中結(jié)果符合，爆炸時刻之后出現(xiàn)明顯偏離。

圖5 電壓和電流的時間演化曲線Fig. 5 Temporal evolution of voltages and currents

絲等離子體的膨脹軌跡及水中激波軌跡同樣是實驗關(guān)注的測量量，實驗結(jié)果與零維模型結(jié)果的比較見圖6。圖7 在一維壓力分布上疊加了絲邊界運動軌跡、激波軌跡的實驗數(shù)據(jù)和部分絲質(zhì)點流線?？傮w上零維模型和一維模型均給出了與實驗符合的動力學(xué)結(jié)果。由于零維模型始終以激波波后壓力作為邊界條件，但實際上水中的激波波后壓力并非均勻，后期該邊界壓力偏大（如圖7 所示），這導(dǎo)致絲等離子體的膨脹速度及激波速度偏低。

圖6 等離子體邊界和激波陣面的時間演化曲線Fig. 6 Temporal evolution of wire radius and shock front

圖7 絲等離子體和水介質(zhì)壓力分布隨時間的演化Fig. 7 Pressure histories of wire plasma and water

密度和溫度是零維模型的基本熱力學(xué)狀態(tài)量，也是壓力、電阻率等其他物理量的基礎(chǔ)，假設(shè)絲等離子體密度均勻，則根據(jù)實驗測量的邊界軌跡也能獲得體積和密度，圖8 給出了它們的時間演化曲線。對比了零維模型、一維模型及實驗的密度結(jié)果，發(fā)現(xiàn)三者在一定范圍內(nèi)較符合；一維程序采用了數(shù)十個質(zhì)點描述絲等離子體，它們的熱力學(xué)狀態(tài)在早期較長的時間內(nèi)基本重合，圖7 也能看出這些網(wǎng)格在膨脹過程中均勻分布，說明在該實驗應(yīng)用中零維模型的基本假設(shè)是合理的。另外，前期零維模型的溫度比一維模型的更高，后期更低，說明二者的能量轉(zhuǎn)化過程存在一定的差異。

圖8 等離子體密度和溫度的時間演化曲線Fig. 8 Temporal evolution of density and temperature of the wire plasma

圖9 給出了早期焦耳熱（qJ）的沉積、熱壓力做功（pdV）與絲等離子體內(nèi)能和動能的關(guān)系，與能量守恒關(guān)系較為吻合，與本文中式(7)和式(8)的描述完全一致。

圖9 能量及其源項隨時間的演化曲線Fig. 9 Temporal evolution of Joule heating, work, internal energy and kinetic energy of the wire plasma

上述結(jié)果展示了電爆炸絲等離子體的基本能量轉(zhuǎn)化過程。焦耳熱（qJ）提供了內(nèi)能增加和溫度上升的能量源，膨脹過程的做功（pdV,V為體積）會損失等離子體內(nèi)能（Ei），使得溫度下降；膨脹過程中的熱壓力克服邊界水壓力做功（pwdV）及磁壓力做功（pBdV）是動能的源項，對于動能（Ek）與內(nèi)能之和即總能量而言，焦耳熱是正的源項，磁壓和邊界壓力做功是負的源項。從數(shù)值上看，磁壓力做功和絲等離子體的動能相對較小，焦耳熱主要轉(zhuǎn)化為內(nèi)能和對水介質(zhì)做功。

進一步考察本文模型與其他模型中采用的理想氣體壓力以及Tucker 電阻率[11]的對比，見圖10。膨脹早期絲密度較大，考慮電離度修正的理想氣體壓力[13]（p=RρT(1+Z)/A，R為氣體常數(shù)，Z為電離度，A=63.5 為銅的原子序數(shù)）將得到較大的壓力；在爆炸時刻之后，根據(jù)實際氣體狀態(tài)方程得到的壓力近似與電離度為1 的修正理想氣體壓力相等。在熔化時刻前，本文零維模型采用的電阻率與Tucker 模型[11]的數(shù)值完全相同，但在氣化階段產(chǎn)生了明顯分離；爆炸時刻之后的Tucker 電阻率[11]較小，與此對應(yīng)，圖5的電流在第一個峰之后下降不明顯，從而與實驗測量的電流波形不相符。

圖10 壓力及電阻率的時間演化曲線Fig. 10 Temporal evolution of pressure and resistivity

最后，給出絲半徑變化對放電模式的影響，見圖11。在固定充電電壓和絲長度的情況下，增大絲半徑會減小絲電阻和電感，從而減小電壓峰值；同時，由于絲質(zhì)量增大，升高相同的溫度需要更多的焦耳熱，因此電壓峰值時刻推遲；電流的第一個峰值時刻也隨之推遲，且第一個電流峰值更大。爆炸發(fā)生之后，電阻率在一定時間內(nèi)維持較大的數(shù)值，絲電阻較大，使得電流減小，當(dāng)絲半徑小于某個閾值時，電流出現(xiàn)“停滯”，直到絲膨脹后截面積增大，絲電阻減小，電流再次增大。另外，當(dāng)絲半徑足夠大時，由于脈沖功率源的能量有限，絲的溫度不會升高到發(fā)生液氣相變的溫度，此時電流波形將呈現(xiàn)短路電流的狀態(tài)。

圖11 絲半徑改變下的電流及電壓波形Fig. 11 Current and voltage waveform for wires with different radii

3 結(jié) 語

本文發(fā)展了一個描述水中金屬絲電爆炸過程的零維模型，包含了動能和內(nèi)能的演化方程、實際氣體狀態(tài)方程、修正的李-莫爾電導(dǎo)率參數(shù)，以及水中激波關(guān)系式等，實現(xiàn)了以實驗室溫度和密度為初始狀態(tài)的冷啟動和連續(xù)過程模擬計算。當(dāng)采用實際氣體狀態(tài)方程時，初始壓力遠小于理想氣體壓力，早期不會發(fā)生非物理膨脹；當(dāng)采用本文的電導(dǎo)率模型時，早期電阻率與Tucker 比作用量模型的結(jié)果符合較好，因此不需要將電爆炸絲過程分為熔化前和熔化后兩個階段建模。

對比分析了基于一維磁流體模型的SWEET 程序模擬結(jié)果，認為在電流上升階段，水環(huán)境中質(zhì)量密度、溫度和電流密度的均勻分布假設(shè)較為合理，采用零維模型描述絲等離子體對于校驗物性參數(shù)更方便。然而，當(dāng)電流下降后，絲等離子體的壓力可能小于水中激波的波后壓力，將會產(chǎn)生稀疏波追趕激波的圖像，零維建模的基本假設(shè)就不再滿足了。盡管如此，在本文的參數(shù)范圍內(nèi)，零維模型的等離子體邊界軌跡及激波軌跡與實驗結(jié)果相符合；能量轉(zhuǎn)化關(guān)系清楚，守恒性好；改變絲直徑，獲得了不同的放電模式，同樣與實驗規(guī)律相符合，證明了該零維模型的合理性。