陳俊霖，徐浩軍，魏小龍，陳增輝，呂晗陽

空軍工程大學 航空航天工程學院，西安 710038

等離子體隱身技術是一種新概念的隱身技術，具有吸波頻帶寬、隱形效果好、不改變目標的外形、使用時間長、可以通過開關迅速地讓等離子體產(chǎn)生和消失等優(yōu)點[1-3]。在飛行器中天線罩是一種強散射源，電磁波可透過天線罩照射在雷達天線的復雜結構上形成很強的后向散射，感性耦合等離子體(ICP)可以在較低的射頻功率和放電氣壓下獲得穩(wěn)定均勻、大面積高密度的等離子體，且裝置結構簡單，參數(shù)易于調(diào)節(jié)[4-6]，因此研究ICP在天線罩隱身技術上的應用具有較大的實際意義。

國內(nèi)外相關學者在等離子體隱身效果的計算上做了大量工作，文獻[9]采用Wenzel-Kramers-Brillouin (WKB)方法和分層近似研究了電磁波在非均勻等離子體中的傳輸特性，文獻[10]采用Z變換時域有限差分 (ZT-FDTD) 方法研究了不同電子密度和碰撞頻率的等離子體覆蓋圓柱目標的雷達散射截面(Radar Cross Section，RCS)，文獻[11]采用分段線性電流密度卷積時域有限差分(FDTD)方法分別計算了電磁波在磁化和非磁化、均勻和非均勻等離子體中的散射特性，并分析了在假設等離子體分布模型的覆蓋下三維目標的RCS變化，文獻[12]計算了薄層等離子體的衰減率，獲得了大氣等離子體在中高碰撞情況下的衰減系數(shù)。

前述的研究對等離子體參數(shù)空間分布普遍采用假設分布，如線性分布、二次分布、指數(shù)分布，直接利用假設分布構建WKB或時域有限差分模型計算等離子體的電磁散射特性，沒有較好地解決等離子體放電與電磁散射特性計算的耦合問題，不能較好地反映實際分布的等離子體源對覆蓋目標電磁散射特性的影響。

針對前述工作存在的不足，本文設計了一種石英夾層ICP天線罩模型，采用有限元和ZT-FDTD聯(lián)合仿真的方法，建立二維等離子流體模型，得到了穩(wěn)態(tài)下ICP中與電磁散射相關的電子密度空間分布，在此基礎上建立ZT-FDTD模型，利用自編程序ZT-FDTD計算了不同放電條件下石英夾層ICP隱身天線罩的電磁散射特性，利用微波干涉法實驗及XFDTD對軟件計算方法進行了驗證，該方法可為等離子體隱身實驗的放電參數(shù)選擇提供依據(jù)。

1 物理模型

本文研究的石英艙夾層ICP隱身天線罩模型如圖1所示。天線罩外形是簡單的旋轉空心橢球體，外側長半軸為10 cm，短半軸為8 cm，內(nèi)側長半軸為8 cm，短半軸為6 cm，中間夾層為等離子體產(chǎn)生區(qū)；在底部的石英窗下安裝平面型線圈天線，該天線采用空心銅管繞制而成，外徑為8 mm，內(nèi)徑為6 mm，匝數(shù)為一圈，直徑為15 cm；為保證放電天線的長時間穩(wěn)定工作，射頻線圈連接自動循環(huán)水冷系統(tǒng)；射頻電源為RSG-1000型，工作頻率為13.56 MHz，額定輸出功率為1 000 W，在功率源和射頻線圈之間連接有自動射頻阻抗匹配器以調(diào)節(jié)反射功率使其達到最小。在石英窗下側對稱安裝2個直徑為15 mm的石英接口，用于連接工質氣體和真空泵。

工質氣體設置為氬氣，腔室內(nèi)氣壓為2 Pa和20 Pa，放電功率為400～700 W，氣溫為常溫293 K，根據(jù)天線罩模型設置放電特性的計算區(qū)域(圖1)。

2 計算方法

2.1 ICP流體模型

在低溫等離子體中，電子和電子能量的輸運可以用一對擴散漂移方程來描述[13]，電子連續(xù)性方程為

(1)

電子能量密度連續(xù)性方程為

(2)

表1 模型涉及的化學反應Table 1 Chemical reaction of model

注：a為依賴于電子能量的分布函數(shù)。

重粒子的質量守恒采用混合平均方法進行計算，對于粒子k，滿足：

(3)

式中：jk為擴散流矢量；Rk為粒子k的產(chǎn)生率；u此刻代表重粒子的平流流速；ρ為重粒子的混合密度；ωk為第k種粒子的質量分數(shù)。

腔體中的矢量磁勢分布由頻域的安培定律求解得到，即

(4)

式中：σ為電導率，σ=nee2/[me/(jω+vm)];me為電子質量，vm為碰撞頻率;w為電磁角頻率；ε0為真空中的介電常數(shù)；εr為相對介電常數(shù)；A為矢量磁勢；μ0為真空中的磁導率；μr為相對磁導率；Je為線圈電流。

由法拉第定律，該磁場會在腔體中感應出電場：E=-jωA。

因此等離子體輸入功率為

Pind=(1/2)real(E·J)

(5)

式中：J=σEe為腔體內(nèi)的感應電流。電子和電子能量的邊界條件設置不考慮電子的二次發(fā)射：

-n·Γe=(1/4)·veth·ne

(6)

-n·Γε=(5/3)·(veth·nε)/4

(7)

式中：n為邊界表面的法線；veth為電子熱運動速度。

由式(1)～式(7)可知該流體模型包含一組高度非線性且相互耦合的偏微分方程，用解析方法無法進行求解，本文利用軟件COMSOL-Multiphysics中的等離子體模塊對該流體模型進行求解。

COMSOL在求解過程中采用的是有限元方法。有限元法是基于變分原理和加權余量法的一種差分方法，其基本思想是將模型的求解域分割為有限多個互不重疊且在頂角處相互聯(lián)接的小單元，這些小單元又稱為基元，基元的角點稱為節(jié)點，在每個基元內(nèi)選擇恰當?shù)牟逯岛瘮?shù)，將偏微分方程組中的所求變量分解為各節(jié)點上插值函數(shù)的線性組合，從而將非線性偏微分方程組近似為矩陣方程并進行求解。

網(wǎng)格是進行有限元計算的關鍵因素，對計算結果精度有著很大的影響。由于等離子體放電過程涉及的物化反應主要在放電腔室區(qū)域內(nèi)進行，同時由于感應耦合放電具有趨膚效應等特點，在腔室壁面附近反應較為強烈。因此，本文采用了不同方式對各區(qū)域進行網(wǎng)格生成，如圖2所示。對放電腔室區(qū)域進行加密處理生成高密度網(wǎng)格，并采用了邊界層網(wǎng)格，其他區(qū)域網(wǎng)格較為稀疏，以減小計算量。

2.2 ZT-FDTD模型

電磁波在等離子體中的傳播過程采用Maxwell方程表示為

(8)

(9)

(10)

式中：H為磁場矢量；Jd為極化電流密度wp為等離子特徑頻率。

本文采用ZT-FDTD方法作為等離子體與電磁波相互作用的計算分析方法。Z變換法可直接將電磁場方程在Z域內(nèi)表征為離散差分方程，避免了積分項的引入，易于計算機編程且具有較高的精度[16]。等離子體在Z域上的本構關系為

D(z)=E(z)+S(z)

(11)

式中:D(z)為電位移矢量在Z域上的表達式，S(z)為輔助參量，表達式為

S(z)=z-1S(z)[1-exp(-vmΔt)]-

Δt[1-exp(-vmΔt)]

(12)

其中:含z-1和z-2的項分別代表輔助參量S(z)的前兩個時間步的取值。ZT-FDTD對等離子體中電磁波的各場量的遞推方法如式(13)～式(16)所示：

En=Dn-Sn

(13)

Sn=Sn-1[1-exp(-vmΔt)]-Sn-2·

(14)

(15)

(16)

該模型計算區(qū)域分為3個功能區(qū)：總場區(qū)、散射場區(qū)和完全匹配層(UPML)，如圖3所示。在等離子體層下方有一塊圓形銅板模擬天線，為降低建模難度，忽略了用于約束等離子體的石英結構，入射電磁波為覆蓋主要雷達波頻段的高斯脈沖，入射方向為天線罩軸線方向，極化方向為θθ極化。

計算中Yee元胞邊長為3 mm，時間步長設為電磁波通過半個網(wǎng)格空間步長的時間。在等離子體區(qū)域迭代計算時，從流體模型計算得到的ωp的采樣點數(shù)量是等離子體區(qū)域網(wǎng)格數(shù)量的4倍，每一個元胞中有9個采樣點。介質界面處電磁參數(shù)的選取及完全匹配層的設置均參考文獻[17-18]。

由上述ZT-FDTD方法計算得到近場后利用惠更斯遠-近場變換公式可得到等離子體天線罩在某一方向的散射電場(設為Es(t))。經(jīng)過傅里葉變換，可得其對應的頻譜Es(f)為

(17)

同樣，經(jīng)過傅里葉變換也可以得到入射波Ei(t)對應的頻譜Ei(f)為

(18)

根據(jù)RCS的定義，可得到目標在某一方向的RCS隨頻率的變化為

(19)

式中:r為目標到雷達的矩離。

為驗證該算法及開發(fā)程序的正確性和有效性，分別用自編的ZT-FDTD程序和XFDTD軟件計算無等離子體及等離子體均勻分布且電子密度為1×1017m-3，氣壓為2 Pa時天線罩的后向RCS，結果數(shù)據(jù)對比如圖4所示。

可以看出，在有/無等離子體的條件下，采用ZT-FDTD得到的計算結果與XFDTD軟件結果都能較好地吻合，在某些頻段存在一定的差異，這主要是因為網(wǎng)格的剖分密度不同造成的。

3 計算結果分析

3.1 ICP電子密度空間分布

運用COMSOL軟件中的等離子體模塊對上述流體模型進行運算，圖5(a)和圖5(b)為在放電功率PD為400 W，氣壓為2 Pa和20 Pa時由流體模型得到的在放電處于穩(wěn)態(tài)時的電子密度分布,橫坐標R為柱坐標的徑向坐標?？梢钥闯?，氣壓為2 Pa時其夾層區(qū)域電子密度沿腔體分布較均勻，這是因為此時電子的平均自由程(約為2.2 cm)和放電的夾層區(qū)域的尺寸差不多，此時雖然電源能量主要注入到射頻電磁場趨膚層內(nèi)，但由于電離過程的變化尺度和電子平均自由程相近使得電子能擴散至整個腔體，從而使夾層內(nèi)電子密度的空間分布仍然比較均勻[19]。而在氣壓為20 Pa時電子密度明顯提高，同時沿著腔體分布具有較大梯度，出現(xiàn)這種現(xiàn)象的原因是氣壓為20 Pa時電子平均自由程相比于2 Pa時下降了10倍，遠小于夾層區(qū)域尺寸，碰撞電離主要被限制在射頻電磁場的趨膚層內(nèi)，導致ne在腔體夾層區(qū)域中的分布出現(xiàn)較大梯度，同時氣壓升高導致碰撞頻率有效提高，使電子密度有效增強[20]。

為驗證流體模型的準確性，基于圖1的模型開展放電實驗，運用微波干涉法對等離子體電子密度進行診斷，其實驗系統(tǒng)裝置如圖7所示，該系統(tǒng)由矢量網(wǎng)絡分析儀和喇叭天線組成，在石英腔底面放置鋁板，由矢量網(wǎng)絡分析儀記錄微波信號通過等離子體后的相移Δφ，在鋁板下方放置吸波材料。

微波探測路徑上等離子體的平均電子密度為[21]

(20)

圖8為不同放電條件下實驗結果與計算結果的對比。可以看出實驗結果與仿真結果符合較好，微波干涉法診斷結果略低于流體模型結果，這主要是由于在流體模型中沒有考慮功率損失、容性分量等對的影響，同時也可以看出隨著功率提高，模型的誤差減小，說明在低功率情況下，容性分量在總耦合功率中占比較高。

3.2 等離子體天線罩隱身效果分析

圖9給出放電氣壓為2 Pa，不同放電功率PD下和無等離子體時天線罩的后向RCS對比。加上ICP后，在低頻段(0～0.7 GHz)，RCS變化很小，這是由于此時電磁波波長遠大于天線罩尺寸，天線罩處在電磁散射的瑞利區(qū)，RCS主要由天線罩體積決定；而在2 GHz以上較寬的頻率范圍內(nèi)加上ICP后可使天線罩的RCS明顯下降。

由于氣壓在2 Pa時碰撞頻率處于106Hz量級，而對應的頻率處于1 GHz量級：ωp?vm，碰撞衰減較弱，因此在2～20 GHz頻段內(nèi)，對RCS的平均衰減較低，但在接近等離子體頻率的局部窄帶衰減較高，衰減的峰值隨著功率的增大向高頻方向移動，這是因為當入射電磁波的頻率接近時，會產(chǎn)生較強的共振衰減作用，衰減峰值區(qū)出現(xiàn)在4.2～5.1 GHz頻帶內(nèi)，功率增加衰減峰值不一定增加。

圖10為放電氣壓為20 Pa，不同放電功率PD下和無等離子體時天線罩的后向RCS對比，相比于氣壓為2 Pa，此時RCS衰減超過10 dBsm的頻帶更寬，同時衰減效果也有所增加，這是由于此時電子密度分布不均勻使得等離子體振蕩頻率在空間有較大梯度，因此在更寬的頻段范圍均能對電磁波造成共振衰減，同時，等離子體碰撞頻率提高了10倍，使得碰撞衰減增強，從而有效增強衰減效果。主衰減頻帶隨著功率的升高向高頻方向移動，同時曲線的波動特性加強，即出現(xiàn)“振鈴”現(xiàn)象[22]，這是由于ωp空間梯度的擴大加劇了散射波的繁雜性，導致RCS的波動特性加強。

4 結 論

1) 感性耦合等離子體可產(chǎn)生較高電子密度的等離子體，能有效實現(xiàn)RCS的減縮。

2) 氣壓為2 Pa時，碰撞衰減較弱，等離子體密度分布較均勻，衰減帶寬集中在共振頻率附近，功率增加會使衰減峰值向高頻方向移動。

3) 氣壓為20 Pa時，碰撞衰減增強，且等離子體密度分布有較大梯度，衰減帶寬有效增加，RCS曲線的波動特性加強。

參 考 文 獻

[1] SINGH A, DESTLER W W, CATRAVAS P. Experimental study of interaction of microwave with a nonmagnetized pulsed-plasma column[J]. Journal of Applied Physics, 1992, 72(5): 1707-1719.

[2] DAVID H L. Plasma stealth[J]. New Scientist, 2000, 168(2264): 60.

[3] 白希堯, 張芝濤, 楊波, 等. 用于飛行器的強電離放電非平衡等離子體隱身方法研究[J]. 航空學報, 2004, 25(1): 52-54.

BAI X Y, ZHANG Z T, YANG B, et al. Study on the method of non-Equiblium plasma stealth by using strong ionization discharge[J]. Acta Aeronautica et Astronautica Sinica, 2004, 25(1): 52-54 (in Chinese).

[4] 趙文鋒, 楊洲, 王衛(wèi)星, 等. 基于CFDRC 的感應耦合等離子體離子數(shù)密度空間分布仿真[J]. 高電壓技術, 2014,40(1): 206-211.

ZHAO W F, YANG Z, WANG W X, et al. Simulation on the spatial distribution of ICP ion number density based on CFDRC[J]. High Voltage Engineering, 2014, 40(1): 206-211 (in Chinese).

[5] NAKAGAWA H, MORISHITA S, NODA S, et al. Characterization of 100 MHz inductively coupled plasma ICP by comparison with 13.56 MHz ICP[J]. Journal of Vacuum Science & Technology A, 1999, 17(4): 1514-1519.

[6] MAHONEY L J, WENDT A E, BARRIROS E, et al. Electron-density and energy distributions in a planar inductively coupled discharge[J]. Journal of Applied Physics, 1994, 76: 2041-2047.

[7] 常雨, 陳偉芳, 孫明波, 等. 等離子體渦電磁散射特性及隱身性能[J]. 航空學報, 2008, 29(2): 304-308.

CHANG Y, CHEN W F , SUN M B, et al. Scattering and stealth of plasma vortex[J]. Acta Aeronautica et Astronautica Sinica, 2008, 29(2): 304-308 (in Chinese).

[8] 陳楠. 基于FDTD等離子體天線隱身及輻射性能研究[D]. 衡陽： 南華大學， 2010： 31-32.

CHEN N. Research on plasma antenna stealth and radiation performance based on FDTD method[D]. Hengyang: University of South China, 2010： 31-32 (in Chinese).

[9] BRUSKI L G, MASE A, TAMANO T, et al. Application of one-dimensional Wentzel-Kramers-Brillouin approximation in microwave reflectometry of plasma density profiles[J]. Review of Scientific Instruments, 1998, 65(5): 2184-2185.

[10] 劉少斌, 莫錦軍, 袁乃昌. 非磁化等離子體密度與目標雷達隱身的關系[J]. 電波科學學報, 2003, 18(1): 57-61.

LIU S B, MO J J, YUAN N C. Research on the relation between the unmagnetized plasma density and the stealth of target[J]. Chinese Journal of Radio Science, 2003, 18(1): 57-61 (in Chinese).

[11] SULLIVAN D M. Frequency-dependent FDTD method using Z transforms[J]. IEEE Transactions on Antennas Propagation, 1992, 40(10): 1223-1230.

[12] LIU M, HU X, JIANG Z. Attenuation of wave in a thin plasma layer by finite-difference time-domain analysis[J]. Journal of Applied Physics, 2007, 101(5): 1661.

[13] 張文茹. 氬氣放電的流體力學模擬及其COMSOL軟件的驗證[D]. 大連： 大連理工大學， 2013: 8-9.

ZHANG W R . The fluid simulation of argon discharge and its verification with COMSOL software[D]. Dalian: Dalian University of Technology, 2013: 8-9 (in Chinese).

[14] ANGEL O B, CORNELIA B. Fast and reliable simulations of argon inductively coupled plasma using COMSOL[J]. Vacuum, 2015, 116: 65-72.

[15] MORAVEJ M. Properties of an atmospheric pressure radio-frequency argon and nitrogen plasma[J]. Plasma Sources Science & Technology, 2006, 15(2): 204-210.

[16] 楊利霞, 許紅蕾, 孫棟, 等. 雙各向異性色散介質電磁波傳播Z-時域有限差分分析[J]. 電波科學學報, 2015, 30(3): 423-428.

YANG L X, XU H L, SUN D, et al. Electromagnetic scattering by bianisotropic dispersive media by using Z-FDTD method[J]. Chinese Journal of Radio Science, 2015, 30(3): 423-428 (in Chinese).

[17] ELSHERBENI A Z, DEMIR V. The finite-difference time-domain method for electromagnetics with MATLAB simulation[M]. North Carolina: SciTech publishing Inc., 2006: 31-40.

[18] 葛德彪, 閆玉波. 電磁波時域有限差分方法[M]. 3版. 西安: 西安電子科技大學出版社, 2011: 83-96.

GE D B, YAN Y B. Finite-difference time-domain method for electromagnetic waves[M]. 3rd ed. Xi’an: Xidian University Press, 2011: 83-96 (in Chinese).

[19] BRAINTHWAITE N St J. Introduction to gas discharge[J]. Plasma Sources Science and Technology, 2000, 9(4): 517-527.

[20] LICHTENBERG A J, LIEBERMAN M A. Principles of plasma discharges and materials processing[M]. New York: John Wiley & Sons, Inc., 2005: 304-305.

[21] 杜寅昌, 曹金祥, 汪建, 等. 射頻電感耦合夾層等離子體中的模式轉換[J]. 物理學報. 2012, 6(19): 337-342.

DU Y C, CAO J X, WANG J, et al. Mode transition of inductively coupled plasma in interlayer chamber[J]. Acta Physica Sinica, 2012, 6(19): 337-342 (in Chinese).

[22] LIU S B, YUAN N C, MO J J. A novel FDTD formulation for dispersive media[J]. IEEE Microwave and Wireless Compoents Letters, 2003, 13(5): 187-189.