二維時(shí)域電場(chǎng)積分方程新平均穩(wěn)定法

摘 要:對(duì)二維理想金屬導(dǎo)體柱的時(shí)域電場(chǎng)積分方程(TDEFIE)后期振蕩問(wèn)題提出了一種新的時(shí)間平均方法，此方法用前一時(shí)刻和后一時(shí)刻的值來(lái)修正當(dāng)前時(shí)刻的電流值，消除了振蕩，而且相對(duì)以往的平均方法提高了MOT算法的穩(wěn)定性，對(duì)精度的影響很小。這里采用電場(chǎng)積分方程(EFIE)MOT算法的隱式格式(Implicit)，在TM，TE高斯平面波激勵(lì)情況下分別計(jì)算了兩個(gè)例子:圓柱和方柱。又采用一種新的時(shí)間基函數(shù)BLIFs來(lái)結(jié)合新的時(shí)間平均法，取得了很好的效果。經(jīng)過(guò)Matlab軟件的數(shù)值分析，可以看出此時(shí)間平均法更加穩(wěn)定和精確。

關(guān)鍵詞:時(shí)域電場(chǎng)積分方程;時(shí)間平均法;隱式算法;時(shí)間步進(jìn)算法

中圖分類號(hào):TN820.1文獻(xiàn)標(biāo)識(shí)碼:A

文章編號(hào):1004-373X(2010)04-174-03

A Novel Average Technique for the Stability of Two-demension TDEFIE

ZHANG Cheng，LUO Wei，WU Xianliang

(School of Electronic Science and Technology，Anhui University，Hefei，230039，China)

Abstract:A new time average method using time domain electric integral equation，which adds the front and back of the current electricity to eliminate oscillation，it is presented to calculate two-demension conducting cylinder for improving the stability more than other average approaches and the precision influenced by the method can be almost ignored.An implicit method using EFIE is adopted to calculate two examples apiece inspired by TM and TE gauss plane wave:circle cylinder and square cylinder.Except that，a new temporal basis funcion BLIFs used in the novel average technique acquire better impression.The results of Matlab simulation prove the stability and precision of the approach.

Keywords:time domain electric integral equation;time average method;implicit method;marching-on-in-time

0 引 言

TDIE是四種數(shù)值算法中發(fā)展較為滯后的一種，但由于其在電磁散射和輻射中無(wú)可比擬的優(yōu)勢(shì)，得到越來(lái)越多的關(guān)注。影響它發(fā)展的原因主要是后期不穩(wěn)定現(xiàn)象。許多學(xué)者經(jīng)過(guò)潛心研究對(duì)其做出了各種假設(shè)，其中具有代表性的有:A.G.Tihuis認(rèn)為在建模時(shí)對(duì)時(shí)間和空間進(jìn)行了一系列近似帶來(lái)誤差，從而在時(shí)間步進(jìn)過(guò)程中積累，最終導(dǎo)致發(fā)散;B.P.Rynne認(rèn)為目標(biāo)本身存在某些內(nèi)諧振頻率。以上只是許多猜想中的兩個(gè)假設(shè)，人們基于以上假設(shè)對(duì)TDIE穩(wěn)定性的改進(jìn)做了很多嘗試，在一定程度上改善了MOT算法的時(shí)間后期穩(wěn)定性。這里基于時(shí)域電場(chǎng)積分方程的隱式算法[1-3]，提出了一種新的時(shí)間平均方法，并與時(shí)間基函數(shù)(BLIFs)相結(jié)合[4]，進(jìn)一步的推遲了振蕩產(chǎn)生時(shí)刻;其后將本文方法與B.P.Rynne的三步平均法和吳振森的先五步后三步的平均方法[5]做了比較，證明此方法的優(yōu)越性。

1 TDEFIE的隱式MOT迭代過(guò)程

由理想導(dǎo)體表面切向電場(chǎng)為零的邊界條件，得到時(shí)域電場(chǎng)積分方程(TDEFIE):

[A(chǔ)t+φ]tan=[Ei]tan(1)

其中矢量位為:

A=μ4π∫c∫∞z′=-∞J(ρ ′，t-R/c)Rdz′dc′(2)

標(biāo)量位為:

φ=14πε∫c∫∞z′=-∞qs(ρ ′，t-R/c)Rdz′dc′(3)

由電流連續(xù)性定理s#8226;J=-qst，可將電流J與電荷q聯(lián)系起來(lái)。為了利用時(shí)間步進(jìn)(MOT)算法求解(1)式，將散射體橫截面的周界上感應(yīng)電流J(r，t)用空間基函數(shù)fn(r)和時(shí)間基函數(shù)BLIFs展開(kāi)為:

J(r，t)=∑Nk=1∑Tl=1Ik(t)fk(r，t)Tl(t)(4)

式中:N為邊界的離散段數(shù);T為總的時(shí)間迭代步數(shù)。將式(4)代入式(2)和式(3)可得離散化后的矢量位A和標(biāo)量位φ:

A(rm，tn)=μ4π∑Nk=1∑∞l=-∞Ik(tn-Rmkl/c)kk，la(5)

φ(rm，tn)=-14πε∑Nk=1∑∞l=-∞∫tn-Rmkl/ct′=0Ik(t′)dt′#8226;

∫k，l∫patchfk/τRds′a(6)

式中:a是第m段矢量位單位向量，

kk，l=ΔτΔz|rm-rk|2+z2l，fork≠morl≠0

4Δτkln(1+2)，fork=m l=0

空間基函數(shù)采用常用的脈沖基函數(shù)，這里要注意由于拐角處的電流連續(xù)性問(wèn)題脈沖基函數(shù)的分段區(qū)間在TM波和TE波激勵(lì)下是不同的:

TM:fm(ρ)=1，ρ∈ρm-1 to ρm

0，otherwise

TE:fm(ρ)=1，ρ∈ρm-1/2 to ρm+1/2

0，otherwise

具體剖分如圖1所示。

圖1 剖分示意圖

以脈沖函數(shù)為檢驗(yàn)函數(shù)，運(yùn)用矩量法(MOM)得到:

fmaτ，[A(chǔ)t+φ]=[fmaτ，Ei]

其中內(nèi)積被定義為:

[a，b]=∫Ca#8226;bdC′

這里采用隱式算法，矢量位A對(duì)時(shí)間的一階導(dǎo)數(shù)用的是后向差分近似。

由于隱式算法的時(shí)間步長(zhǎng)Δt>Rmin/c，則必然有落在待求時(shí)刻tn與前一時(shí)刻tn-1之間的點(diǎn)，此時(shí)刻的電流未知，所以要把矢量位A與標(biāo)量位φ分為兩個(gè)部分，時(shí)刻落在tn-1

2 新平均法和BLIFs時(shí)間基函數(shù)

當(dāng)前已經(jīng)提出的平均法有三步平均法和先五步后三步平均法等，這些平均法運(yùn)用在本文二維情況下顯得有些力不從心，其對(duì)穩(wěn)定性的改進(jìn)效果很差。具體表達(dá)式如下:

三步平均法:

Im(tj)=14[Im(tj+1)+2Im(tj)+Im(tj-1)]

先五步后三步平均法:

Im(tj)=18[Im(tj+2)+Im(tj+1)+4Im(tj)+

Im(tj-1)+Im(tj-2)]Im(tj+1)

=14[Im(tj+2)+2Im(tj+1)+Im(tj)]

本文提出的平均法是在三步平均法的基礎(chǔ)上加以改進(jìn)，其效果遠(yuǎn)勝于前兩種平均法，具體計(jì)算式如下:

Im(tj)=1M[0.5Im(tj+1)+(M-1)Im(tj)+

0.5Im(tj-1)]

式中:M為正整數(shù)，本文算例取M=5。

這里所用的時(shí)間插值基函數(shù)的表達(dá)式為:

T(t)=sin(πt/Δt)πt/Δt#8226;sin[a(t/NΔt)2-1]sinh(a)(t/NΔt)2-1，

|t|≤(N+0.5)Δt

0，|t|>(N+0.5)Δt

式中:a=πNδ;Δt=π(1-δ)/0;δ的取值范圍為[0，1];N為正整數(shù)，本文算例取N=7，Δt=1。

3 數(shù)值計(jì)算

本文采用Matlab軟件進(jìn)行仿真，并計(jì)算了4個(gè)算例，分別是TM波和TE波激勵(lì)情況下的無(wú)限長(zhǎng)方柱和圓柱。TM波圓柱半徑取1.0 m，方柱邊長(zhǎng)取1.0 m;TE波圓柱半徑取0.5 m，方柱取邊長(zhǎng)0.8 m。

文獻(xiàn)[1]的結(jié)果是指僅用時(shí)間基函數(shù)BILFs代替文獻(xiàn)[1]中的三角基函數(shù)，不經(jīng)過(guò)任何平均處理后得到的結(jié)果。

假定入射波為高斯平面波，其中TM波和TE波形式各如下:

3.1 TM波

電場(chǎng)只有z方向，Eiz  = E0 4Tπe-γ2，γ=4T(ct-cto-r#8226;k)。式中:k為波傳播方向，本節(jié)算例中入射角取π，即k=[-1，0]，E0= 120π，c為真空中光速c=2.9×108，ct0=6.0，T=4.0。由于z方向無(wú)限長(zhǎng)，所以Φ#8226;az=Φz(mì)=0，從而電場(chǎng)積分方程可簡(jiǎn)寫(xiě)為:[A(chǔ)t]tan=[Ei]tan，即Atan=∫t0Eitandt。

3.2 TE波

電場(chǎng)只有橫向分量，磁場(chǎng)只有z方向，Hiz = H04Tπe-γ2，γ=4T(ct-cto-r#8226;k)。其中H0=1.0。其他各參數(shù)同TM波。圖2為BLFS時(shí)間基準(zhǔn)函數(shù)圖，圖3為高斯平面波的時(shí)域圖，圖4為高斯平面波頻域圖。

圖2 BLIFs時(shí)間基函數(shù)N=7，Δt=1，δ=0.5

圖3 高斯平面波的時(shí)域圖

圖4 高斯平面波頻域圖

圖5為T(mén)M高斯平面波激勵(lì)下半徑α=1.0 m圓柱，截面圓心在z軸，離散段數(shù)N=24，時(shí)間步長(zhǎng)Δt=1.350 3×10-9電流密度隨時(shí)間變化關(guān)系，觀察點(diǎn)(1 m，0°，90°)。

圖5 電流密度隨時(shí)間的變化關(guān)系(一)

圖6為T(mén)M高斯平面波激勵(lì)下邊長(zhǎng)L=1.0 m方柱，截面中心在z軸，離散段數(shù)N=24，時(shí)間步長(zhǎng)Δt=6.095 7×10-10電流密度隨時(shí)間變化關(guān)系，觀察點(diǎn)(0.5 m，0°，90°)。

圖7為T(mén)E高斯平面波激勵(lì)下半徑α=0.5 m圓柱，截面圓心在z軸，離散段數(shù)N=18，時(shí)間步長(zhǎng)Δt=7.524 8×10-10電流密度隨時(shí)間變化關(guān)系，觀察點(diǎn)(0.5 m，0°，90°)。

圖6 電流密度隨時(shí)間關(guān)系(二)

圖7 電流密度隨時(shí)間關(guān)系(三)

圖8為T(mén)E高斯平面波激勵(lì)下邊長(zhǎng)L=0.8 m方柱，截面中心在z軸，離散段數(shù)N=24，時(shí)間步長(zhǎng)Δt=5.777 8×10-10電流密度隨時(shí)間變化關(guān)系，觀察點(diǎn)(0.4 m，0°，90°)(橫向坐標(biāo)為時(shí)間，縱向坐標(biāo)為切向電流密度)。

圖8 電流密度隨時(shí)間關(guān)系(四)

4 結(jié) 語(yǔ)

由Matlab數(shù)值仿真結(jié)果可知，在TM高斯平面波激勵(lì)下，三步平均法和先五步后三步平均都更提前了振蕩時(shí)刻，效果更差，而本文平均法在最后時(shí)刻仍沒(méi)有發(fā)現(xiàn)振蕩現(xiàn)象;對(duì)于TE高斯平面波激勵(lì)，雖然前兩種舊的平均方法都推遲了振蕩時(shí)刻但效果都不是很明顯，本文方法卻可以將振蕩時(shí)刻近一步延遲，從而證明了本文新平均方法的優(yōu)越性。

參 考 文 獻(xiàn)

[1]Rao S M.Time Domain Electromagnetics [M].New York:Academic Press，1999.

[2]Rao S M，Tapan K.Sarkar.Time-domain Modeling of Two-dimensional Conducting Cylinders utilizing an Implicit Scheme-TM INCIDENCE[Z].Microwave and Optical Technology Letters，1997.

[3]Rao S M，Douglas A Vechinski，Tapan K Sarkar.Transient Scattering by Conducting Cylinders-implicit Solution for the Transverse Electric Case[Z].Microwave and Optical Technology Letters，1999.

[4]Pisharody G，Weile D S.A Robust Solution to Time-domain Integral Equations for Perfect Electric Conductors using Loop-tree Decomposition and Bandlimited Extrapolation[A].IEEE Antennas and Propagation Society AP-S International Symposium[C].2004(4):4 204-4 207.

[5]吳振森，泰三團(tuán)，代少玉.一種改善TDIE時(shí)間步進(jìn)算法穩(wěn)定性的新平均方法[J].電波科學(xué)學(xué)報(bào)，2004，19(Z1):102-105.

[6]Hu J L，Chan C H.An Improved Temporal Basic Function for the Time Domain Electric Field Equation Method[Z].Electronics Letters，1999.

[7]A.Al-Jarro.Time Domain Integral Equation Techniques:1 D and 3 D Models[J].Electrical and Electronic Engineering，Nottingham:University of Nottingham，2004.

[8]趙延文，聶在平，徐建華，等.時(shí)域電場(chǎng)積分方程的穩(wěn)定求解[J].電波科學(xué)學(xué)報(bào)，2004，19(2):148-152.

[9]逯貴禎.一種新的迭代算法研究[J].電波科學(xué)學(xué)報(bào)，2001，16(Z1):40-42.

[10]文柯一.電磁理論的新進(jìn)展[M].北京:國(guó)防工業(yè)出版社，1999.