復(fù)雜媒質(zhì)目標(biāo)電磁散射問題的有限元分析

彭文峰 宛 汀 郁美艷

（南京理工大學(xué)，江蘇 南京210094）

引 言

隨著微波和光學(xué)工程的迅速變革，以及人工合成媒質(zhì)技術(shù)的不斷發(fā)展，復(fù)雜電磁媒質(zhì)的研究受到越來越多的重視.選擇和運用恰當(dāng)?shù)碾姶艌鰡栴}的數(shù)值計算方法對復(fù)雜電磁媒質(zhì)進行研究有著重要的理論意義和工程應(yīng)用價值.在一些傳統(tǒng)的數(shù)值計算方法中，使用積分方程法如矩量法的算法，往往導(dǎo)致復(fù)雜的格林函數(shù)，有時候格林函數(shù)甚至是不可得的.使用時域有限差分法（FDTD）時，因為需要計算額外的時間導(dǎo)數(shù)而變得不穩(wěn)定.用于分析平面分層結(jié)構(gòu)的譜域法與直線法也需要重新計算格林函數(shù).雖然利用場分解特性［1］，可以將手征材料等效成兩個沒有耦合的簡單材料進行計算，但是這種方法卻不能應(yīng)用到更具有通用形式的一般雙各向異性媒質(zhì)中.眾所周知，有限元方法［2］是一種通用性很強的數(shù)值算法，它能靈活地適應(yīng)物體幾何結(jié)構(gòu)和材料的變化.與FDTD相比，有限元法更適用于分析任意不均勻復(fù)雜結(jié)構(gòu).而與矩量法相比，使用有限元算法計算雙各向異性媒質(zhì)的優(yōu)點在于推導(dǎo)雙各向異性的泛函公式相對比較簡單，有較強的復(fù)雜媒質(zhì)的處理能力.而且，有限元方法生成的線性方程組的系數(shù)陣為稀疏矩陣，節(jié)省了內(nèi)存.因此，本文研究將有限元方法發(fā)展應(yīng)用到含有復(fù)雜媒質(zhì)結(jié)構(gòu)的電磁散射問題中，包括各向異性媒質(zhì)、雙各向同/異性媒質(zhì)，并驗證其準(zhǔn)確性和通用性.

本文建立了具有一般通用形式的線性媒質(zhì)——雙各向異性媒質(zhì)［3］的有限元泛函公式模型，采用完全匹配層［4］來截斷散射問題開放區(qū)域，并推導(dǎo)了有限元矩陣方程的具體表達式.然后，運用該方法對四種典型的復(fù)雜媒質(zhì)結(jié)構(gòu)散射問題進行了數(shù)值計算，并將計算結(jié)果和相關(guān)文獻進行比較.結(jié)果證明，有限元方法能很好地解決含有各向異性、雙各向同性以及雙各向異性的各種復(fù)雜媒質(zhì)目標(biāo)電磁散射的數(shù)值仿真問題.本文的研究工作具有較強的通用性，能夠?qū)Χ喾N復(fù)雜媒質(zhì)目標(biāo)進行分析，在研究雷達目標(biāo)隱身和反隱身技術(shù)、復(fù)雜天線系統(tǒng)設(shè)計、現(xiàn)代電子系統(tǒng)的電磁兼容性分析等領(lǐng)域均能有效地發(fā)揮作用.

1 基本理論

對于復(fù)雜媒質(zhì)，其本構(gòu)關(guān)系為

圖1給出了一個典型的采用有限元方法分析電磁散射問題的示意圖.我們采用完全匹配層（PML）作為區(qū)域截斷邊界條件.在PML區(qū)域內(nèi)，有效磁導(dǎo)率和介電常數(shù)為對角矩陣的形式.對于z方向的PML，強加如下約束條件

圖1 三維目標(biāo)的散射

假定波的傳播方向為z，當(dāng)bc＝1，a＝b時，可以得到0反射.為了使入射波充分的衰減，a、b、c應(yīng)取作復(fù)數(shù).令a＝b＝1/c＝α-jβ，α決定了此介質(zhì)中的波長，β決定了波的衰減程度.則PML內(nèi)的相對磁導(dǎo)率與電導(dǎo)率可以表示為

PML參量的取值依賴于波的傳播方向和PML層本身的位置（PML層可以設(shè)在面、邊或者角落）.針對二維和三維問題時，我們選擇不同的ε＝r和＝μr應(yīng)用于不同的方向.棱邊和體角落含有特殊的張量，例如放于棱邊，則對于z方向上的一條棱邊，應(yīng)選擇r和r的值為r＝r＝［μr］x×［μr］y；放在體角上，則為三個方向上的各張量的乘積.

由復(fù)雜媒質(zhì)的本構(gòu)關(guān)系以及麥克斯韋方程組推導(dǎo)出雙各向異性媒質(zhì)的電場波動Helmholtz方程為

對于復(fù)雜媒質(zhì)，有限元泛函可以寫為［2］

由于式（6）的工作變量是總場，在散射問題中有總場E＝Esc＋Einc，其中Esc為散射場，Einc為入射場，因此式（6）也可以寫為

又因Einc為已知場，在變分公式求偏導(dǎo)時等于零，因此式（7）中只含Einc的項可以去掉，應(yīng)用第一矢量格林定理，可得

式中：Vsc表示散射體的體積；Ssc表示散射體的表面；Einc為激勵平面波.

在獲得整個分析區(qū)域的泛函之后，接下來要進行的是區(qū)域離散的工作.選用靈活的四面體單元對整個區(qū)域進行離散，生成網(wǎng)格之后，根據(jù)選定的切向連續(xù)的矢量基函數(shù)（9）來建立有限元方程.

式中：i＝1，2，…，6，i表示第i條邊；i1i2表示第i條邊兩個端點的編號；Li為體積坐標(biāo).

選擇相同的基函數(shù)和加權(quán)函數(shù)Ni，可以得到復(fù)雜媒質(zhì)的有限元散射公式為

根據(jù)伽遼金方法，對于其單元每一棱邊元的殘數(shù)加權(quán)積分為0，結(jié)合式（10），得到整個e單元的矩陣形式為

式中：

其中，

式中，

由方程（11）解得散射場Esc（近場）之后，在遠處r的散射場Esc（r）（遠場）就可以使用等效原理求解，即

式中：S是包圍散射體的任意閉合曲面；Er為等效的無窮遠處入射的平面波，其表達式為

其中α為極化角，當(dāng)α＝0時入射波為θ極化，α＝π/2時入射波為φ 極化.kinc是傳播矢量，θinc、φinc為激勵平面波的入射角，有

式中：Ni為基函數(shù)；ai為有限元方程已求出的基函數(shù)的系數(shù).

將式（17）代入式（15）得到

求出了遠處r的散射場，就可以用雷達散射截面（RCS）表示物體的散射特性，雷達散射截面的定義為

由于入射場是平面波，因此有｜Einc（r）｜＝1，則式（19）可以寫為

2 數(shù)值算例

下面計算幾個典型的例子來證明本文方法的正確性和有效性.本文均采用PML作為區(qū)域截斷邊界條件，為保證PML良好的截斷效果，算例中PML的厚度均取為0.25λ0，距離散射體為0.3λ0，λ0為自由空間波長，參數(shù)值取為α＝β＝1.5.

第一個例子分析的是電尺寸為k0a＝0.2π的鐵氧體球.直角坐標(biāo)系的原點位于球心，設(shè)外加偏置磁場在方向得到，即H0＝H0，入射的平面波沿z軸方向，其極化方向沿正x軸方向.各向異性鐵氧體的介質(zhì)參數(shù)為ε0（真空中的介電常數(shù)），其相應(yīng)的導(dǎo)磁率為［5-6］，可表示為

本例中，鐵氧體的參數(shù)為μ1＝5μ0，μ2＝j(luò)μ0，μ3＝7μ0.圖2為本文計算結(jié)果與文獻［7］的比較，可以看出曲線有細小差別，其原因是本程序采用PML作為邊界截斷條件，一般情況下PML并不百分之百吸收，并且-35dB已經(jīng)是一個很小的值了，所以該誤差在允許的范圍之內(nèi).

第二個例子是電尺寸為k0a＝0.5的等離子球，入射的平面波沿z軸正方向入射，其極化方向沿正x軸方向.各向異性等離子體的導(dǎo)磁率為μ0（真空中的導(dǎo)磁率），電容率可表示為

圖3為本文結(jié)果與文獻［8］的比較.可以看出，本程序計算的結(jié)果和文獻結(jié)果在-40dB之上還是非常吻合的.由于本程序采用PML作為邊界截斷條件，一般來講，PML是不能達到百分百吸收的，另外，-40dB是可以忽略的十分小的誤差.

第三個例子是尺寸為k0a＝1.5的手征介質(zhì)球.為了與文獻對比，手征的本構(gòu)關(guān)系可以寫成如下形式

根據(jù)文獻［3］，手征參數(shù)設(shè)為εDBF＝4ε0、μDBF＝μ0、ε＝εDBF/（1-k2r）、μ＝μDBF/（1-k2r）.入射的平面波沿z軸正方向，其極化方向沿正x軸方向.圖4中的曲線對比再次驗證了本程序的正確性.

第四個例子是雙各向異性媒質(zhì)圓柱，底面半徑為a＝0.5λ0，高度為h＝0.2λ0，雙各向異性媒質(zhì)參數(shù)具有如下張量形式入射的平面波沿z軸正方向，其極化方向沿正x軸方向.Ω 媒質(zhì)參數(shù)為ε1＝2.0，ε2＝3.0，ε3＝2.0，μ1＝1.2，μ2＝1.2，μ3＝1.0，Ω＝0.0，0.5，1.0.圖5中與文獻［9］曲線的對比驗證了本程序計算雙各向異性媒質(zhì)的正確性.

3 結(jié) 論

本文介紹了分析復(fù)雜媒質(zhì)散射問題的有限元方法.有限元法適用于分析任意非均勻復(fù)雜結(jié)構(gòu)，有較強的復(fù)雜媒質(zhì)的處理能力，生成的線性方程組的系數(shù)陣為稀疏矩陣，節(jié)省了存儲量.算例的數(shù)值分析表明，該方法的計算結(jié)果與文獻結(jié)果吻合較好，從而證明了其正確性和有效性.此外，該方法適用范圍廣，不僅能夠用來計算各向異性媒質(zhì)的散射特性，對雙各向同性和雙各向異性媒質(zhì)同樣適用，是一種通用性較強的方法.

［1］AL-KANHAL M A，ARVAS E.Electromagnetic scattering from a chiral cylinder of arbitrary cross section［J］.IEEE Trans Antennas and Propagation，1996，44（7）：1041-1048.

［2］金建銘，著.電磁場有限元方法［M］.王建國，葛德彪，譯.西安：西安電子科技大學(xué)出版社，1998.

［3］周 平.脊加載橢圓波導(dǎo)傳輸特性的矢量有限元法分析［J］.電波科學(xué)學(xué)報，2009，24（6）：1164-1167.ZHOU Ping.Propagation characteristics analysis for ridged elliptical waveguide by vector finite-element method［J］.Chinese Journal of Radio Science，2009，24（6）：1164-1167.（in Chinese）

［4］孫向陽，聶在平，李愛勇，等.用高階疊層矢量有限元法計算隨鉆測井的三維電磁響應(yīng)［J］.電波科學(xué)學(xué)報，2009，24（2）：273-279.SUN Xiangyang，NIE Zaiping，LI Aiyong，et al.The modeling of logging-while-drilling tool’s three-dimensional electromagnetic response using the high order hierarchical vector finite element method［J］.Chinese Journal of Radio Science，2009，24（2）：273-279.（in Chinese）

［5］EDWARD K N，YUNG H Y，DING R S，et al.Finite-difference analysis of H-Plane waveguide Y-junction circulars［J］.Microwave and Optical Technology Letters，1999，20（6）：414-422.

［6］耿友林，吳信寶，官伯然.導(dǎo)體球涂覆各向異性鐵氧體介質(zhì)電磁散射的解析解［J］.電子與信息學(xué)報，2006，28（9）：1740-1743.GENG Youlin，WU Xinbao，GUAN Boran.The analytical solution to the electromagnetic scattering by an anisotropic ferrite-coated conducting sphere［J］.Journal of Electronics ＆Information Technology，2006，28（9）：1740-1743.（in Chinese）

［7］耿友林.球矢量波函數(shù)在各向異性介質(zhì)散射中的應(yīng)用［D］.西安：西安電子科技大學(xué)，2006.GENG Youlin.The application of the sphere vector wave function for the electro-magnetic scattering from the anisotropic media［D］.Xi’an：Xidian University，2006.（in Chinese）

［8］ZHANG Y，WEI X，LI E.Electromagnetic scattering from three-dimensional bi-anisotropic objects using hybrid finite element-boundary integral method［J］.J of Electromagn Waves and App，2004，18（11）：1549-1563

［9］OZGUN O，KUZUOGLU M.A non-iterative domain decomposition method for finite element analysis of 3D electromagnetic scattering problems［C］//IEEE Antennas and Propagation Society International Symposium.San Diego，July，5-11，2008：40-50.

［10］ZHU J，CHEN R S，F(xiàn)AN Z H，et al.Efficient preconditioner for the finite-element boundary integral analysis of electromagnetic scattering in the half space［J］.IET Microwaves Antennas ＆ Propagation，2010，4（3）：374-380.

［11］JACKSON S A，VOUVAKIS M N.Octree-based finite element method for electromagnetic scattering problems［C］//IEEE Antennas and Propagation Society International Symposium.Toronto，July，11-17 2010：1-4.

［12］ZHANG Ruigang，CHEN Huanzhen.An immersed finite element method for anisotropic flow models in porous medium［C］//International Conference on Information Science and Technology.Nanjing，March 26-28，2011：168-175.