李磊，張昕

(哈爾濱工程大學(xué)信息與通信工程學(xué)院，黑龍江哈爾濱150001)

隨著電磁脈沖武器性能的逐步提高，對裝有靈敏電子器件的現(xiàn)代電子設(shè)備造成毀傷的可能性也越來越大，因而電子設(shè)備對電磁脈沖的響應(yīng)成為國內(nèi)外研究的熱點［1-3］。前門耦合和后門耦合是電磁脈沖能量進入電子設(shè)備的2種途徑［4-5］。目前國內(nèi)外對于后門耦合的研究較多，而研究前門耦合的相對較少。雖然天線是主要的前門耦合途徑，但對其的研究主要集中于較為簡單的線天線。微帶天線由于具有重量輕、體積小、剖面薄、易于共形、電性能多樣化等優(yōu)點，被廣泛應(yīng)用于軍事和民用通信設(shè)備中［6-8］。電磁脈沖的研究屬于瞬態(tài)電磁學(xué)范疇，因而時域有限差分法在處理該問題時具有很大的優(yōu)勢。但時域有限差分法在處理微帶天線結(jié)構(gòu)時，由于天線輻射片邊緣的電磁場具有奇異性，因而在建模時需要特別處理，否則會導(dǎo)致局部的截斷誤差［9］。Taflove提出了亞網(wǎng)格技術(shù)來細化局部網(wǎng)格劃分［10］，但這樣減小時間間隔，降低計算效率。Railton等提出了修改材料電參數(shù)的方法［11］，提高了計算效率，但對導(dǎo)體嵌入面積大于網(wǎng)格1/2時的情況沒有處理方法。

1 基本原理

1.1 漸變非均勻網(wǎng)格原理

由于FDTD采用了電場和磁場在空間中離散排布的Yee元胞，把麥克斯韋旋度方程轉(zhuǎn)化為顯式差分方程，使得電場和磁場可以在時間上迭代求解，而不需要進行矩陣的求逆［12］，因而在利用FDTD分析電磁問題時，一個主要影響計算精度的因素就是空間離散間隔，即空間網(wǎng)格尺寸。盡管采用更小的網(wǎng)格大致可以保證需要的計算精度，但同時也使得元胞數(shù)量成倍地增長，極大地增加了計算量，延長了計算時間［13］。目前主要有2種情況必須要采用很小的空間網(wǎng)格才能保證計算精度，一種情況是在進行建模時，有不能忽略的細小結(jié)構(gòu)需要精確建模才能保證計算結(jié)果的可靠性;另一種情況是計算空間中，某一個區(qū)域存在劇烈變化的場，需要很小的空間網(wǎng)格才能準確描述。微帶天線同時符合以上2種情況，首先微帶天線的基底厚度與其長寬相比較小，其次微帶天線的不連續(xù)結(jié)構(gòu)周圍存在劇烈變化的電磁場。因此在計算微帶天線時，基底厚度通常至少需要3個空間網(wǎng)格，網(wǎng)格尺寸一般小于1 mm，如果采用均勻網(wǎng)格建模微帶天線，就需要很大的計算空間，進而導(dǎo)致較大的內(nèi)存消耗和較長的計算時間，計算效率極低。

漸變非均勻網(wǎng)格(variable mesh scheme，VMS)FDTD由choi等［14］提出，其在電磁場變化劇烈處采用小尺寸網(wǎng)格，在電磁場變化平緩處采用大尺寸網(wǎng)格，這種靈活的網(wǎng)格劃分方法較好地解決了計算效率與計算精度的矛盾。此外，除空間離散間隔Δx、Δy和Δz的大小在計算空間是變化的之外，VMSFDTD的電磁場更新方程與常規(guī)的均勻網(wǎng)格FDTD的相同。

例如，圖1所示的Ex和Hx在自由空間的更新方程為:

式中:Δx(i)、Δy(j)和Δz(k)為電場或磁場所處元胞的空間離散間隔。其他方向的電場和磁場的更新方程與式(1)、(2)類似。由圖1可見，電場Ez不處于磁場網(wǎng)格的中心，且2個相鄰的電場網(wǎng)格大小也不相等，因而在計算Ez時，磁場網(wǎng)格的大小被視為兩個相鄰元胞尺寸和的一半。為保證算法的穩(wěn)定性，Δt需要滿足:

式中:Δx(i)min、Δy(j)min和 Δz(k)min分別為 x、y和z方向的最小離散間。

圖1 漸變非均勻網(wǎng)格Fig.1 VMS

1.2 導(dǎo)體邊緣奇異場處理技術(shù)

FDTD在分析不連續(xù)結(jié)構(gòu)周圍存在劇烈變化的電磁場一般采取更為精細的網(wǎng)格劃分方式，但隨之會導(dǎo)致計算效率的降低。針對上述問題，為了在大網(wǎng)格劃分的基礎(chǔ)上保證計算的精度，本文將文獻［15］得到的導(dǎo)體邊緣場分布函數(shù)引入到FDTD中，得到了導(dǎo)體邊緣場的FDTD處理方法。

如圖2所示的導(dǎo)體邊緣場可表示為

國家環(huán)保模范城市考核指標(biāo)包括社會經(jīng)濟、環(huán)境質(zhì)量、環(huán)境建設(shè)、環(huán)境管理四部分?？紤]到該地區(qū)正在創(chuàng)建模范生態(tài)城市，礦區(qū)規(guī)劃環(huán)評在指標(biāo)選取上主要選擇了生態(tài)環(huán)境質(zhì)量及環(huán)境建設(shè)的部分指標(biāo)。

式中:r為導(dǎo)體邊緣到場點的距離，φ為二者之間的夾角，A和B為未知系數(shù)。由式(4)～(7)可見，由于1/項的存在使得橫向場具有奇異性。

圖2 導(dǎo)體邊緣橫向電磁場Fig.2 Transverse electric and magnetic field components around am et al edge

圖3 導(dǎo)體嵌入網(wǎng)格平面圖Fig.3 Conductors protrude intomesh

現(xiàn)以圖3為例說明導(dǎo)體邊緣場的奇異性處理方法。由圖3(a)可見，導(dǎo)體嵌入網(wǎng)格長度為 Δ，且Δ＞Δx/2，如果按照傳統(tǒng) Yee網(wǎng)格劃分方式，則Ex(i，j，k)、Ex(i，j，k+1)和 Hy(i，j，k)均為零，相當(dāng)于增大了導(dǎo)體嵌入網(wǎng)格的長度，這將導(dǎo)致計算誤差，因而本文將 Ex(i，j，k)、Ex(i，j，k+1)和 Hy(i，j，k)3個節(jié)點的位置沿x軸移動Δ/2，此時r=(Δx- Δ)/2 ，φ =180o，因而 Hy(i，j，k)=-Hφ，Ex(i，j，k)=Er，由式(4)和(7)可得

將式(8)、(9)分別代入式(4)、(7)可得導(dǎo)體邊緣橫向場的分布為

由法拉第定律可得

式中:S1為單元網(wǎng)格中未嵌入導(dǎo)體的區(qū)域。由式(12)可得

由安培定律計算圖3(b)中的 Ez(i，j，k)時，將磁場的線積分分為 ab、bc、cd、de、ef和 fa 共6 段進行計算，則

假設(shè)Ez在計算網(wǎng)格內(nèi)均勻分布，則結(jié)合式(7)可得

其中，A1、A2和A3可由式(9)類似方法求得

其他場量均可按類似方法處理。

2 天線建模與算法驗證

計算模型采用圖4(a)所示的微帶天線，基板材料采用高頻介質(zhì)板，相對介電常數(shù)εr1=2.2，通過同軸線進行饋電，同軸線填充介質(zhì)為聚四氟乙烯，其相對介電常數(shù)εr2=2.1?？臻g離散間隔是影響FDTD計算精度的主要因素之一，本文采用漸變非均勻網(wǎng)格技術(shù)，在天線周圍結(jié)構(gòu)不連續(xù)處采用較小的離散間隔以模擬劇烈變化的電磁場;另一方面，為了提高計算效率，降低內(nèi)存消耗，可以在電磁場平緩變化處采用較大的離散間隔。結(jié)合微帶天線結(jié)構(gòu)，建模時在xoy平面采用 4種離散間隔，即 Δl、1.59Δl、2.52Δl和4Δl，其中最小離散間隔 Δl=0.4 mm(l為x或y)。時間離散間隔Δt=4.514 2×10-13s。微帶天線的網(wǎng)格劃分如圖4(b)所示。

圖4 微帶天線結(jié)構(gòu)和網(wǎng)格劃分Fig.4 The geometry and FDTD grid of microstrip antenna

按上述方法進行建模后，利用本文的算法計算微帶天線的回波損耗，并與有限元法(finite element method，F(xiàn)EM)、沒有用奇異性處理技術(shù)的常規(guī)的FDTD算法及文獻［11］算法的計算結(jié)果進行對比，如圖5所示。

圖5 微帶天線的回波損耗Fig.5 Return loss of microstrip antenna

對比可見，常規(guī)的FDTD算法在整個頻段的計算誤差均較大;文獻［11］的算法在5 GHz以下頻段的精度較高，而在5 GHz以上頻段隨著頻率增加，計算誤差逐漸增大;而本文算法和FEM的計算結(jié)果在整個頻段均保持較好的一致性，從而驗證本文算法的可行性和計算精度。之后利用近遠場變換得到了天線yoz面的遠場方向圖，如圖6所示。由計算結(jié)果可見，天線的中心頻率約為4.7 GHz，rEθ在θ=180o方向最大(θ為輻射方向與正z軸夾角，φ為輻射方向在xoy平面的投影與正x軸夾角)，可知天線極化方向與y軸一致。

圖6 微帶天線的方向圖Fig.6 Radiation pattern of microstrip antenna

3 電磁脈沖入射時的計算結(jié)果與分析

3.1 入射電磁脈沖描述

高功率電磁脈沖的電場可以用微分高斯脈沖近似，其數(shù)學(xué)表達式為

式中:E0=4.664×104V/m，t0=1.264×10-10s，τ=4.750×10-11s。其時域波形和頻譜如圖7。

采用上述電磁脈沖從不同角度入射，輻照目標(biāo)天線，如圖8所示，計算饋電同軸線中距離口徑100Δl處的瞬態(tài)響應(yīng)電壓Ut，并通過傅里葉變換得到其頻譜。其中θ為入射脈沖與正z軸夾角，φ為入射脈沖在xoy平面的投影與正x軸夾角，入射電磁脈沖采用平行極化波。

圖7 入射脈沖的波形和頻譜Fig.7 Waveform and spectrum of incident pulse

圖8 電磁脈沖入射示意圖Fig.8 Incidence of EMP

3.2 電磁脈沖沿不同φ入射時天線的響應(yīng)

如圖9(a)所示，高功率電磁脈沖分別沿著θ=π/4，不同的φ角入射，得到饋電同軸線中距離口徑100Δl處的瞬態(tài)響應(yīng)電壓Ut的時域波形，由于微帶天線的導(dǎo)體結(jié)構(gòu)對電磁脈沖的反射作用，使得電磁脈沖在微帶天線輻射片與接地板之間來回反彈產(chǎn)生了振蕩，振蕩的電磁場通過口徑的耦合進入到饋電同軸線內(nèi)，且隨著能量的耗散，振蕩幅度逐漸衰減。對比來看，φ=0時的響應(yīng)電壓峰值最高，約為8.7 V，φ=π/4時的響應(yīng)電壓峰值最低，約為6.3 V，這主要是由脈沖沿不同路徑傳播，導(dǎo)致微帶天線中場的分布模式也不同引起的。

如圖9(b)所示，高功率電磁脈沖分別沿著θ=π/4，不同的φ角入射，得到饋電同軸線中距離口徑100Δl處的瞬態(tài)響應(yīng)電壓Ut的頻譜，由圖可見，φ=0時的響應(yīng)電壓主要集中在4、8、9.6 GHz等頻段，由腔模理論可知［16］，其分別為相對于z軸的TM10、TM20和TM21(TMmn，m為x方向變化的半周期數(shù)，n為y方向變化的半周期數(shù))等模式的諧振頻率，由此可以得出，由于脈沖的激勵使得微帶天線中的電磁場存在上述分布模式。當(dāng)脈沖沿φ=π/4入射時，與φ=0相比，響應(yīng)電壓頻譜分布向低頻轉(zhuǎn)移，在4.7 GHz頻段響應(yīng)最大，其對應(yīng)于該天線工作的主模TM01的諧振頻率，此外在6.5、8、9.4 GHz頻段的響應(yīng)也較大，說明天線中存在較強的TM11、TM20和TM02等高次模。當(dāng)脈沖沿φ=π/2入射時，頻譜分布更集中于該天線主模TM01的諧振頻率4.7 GHz頻段，高次模頻段的響應(yīng)已相對較小。由此可見隨著入射φ角的變化，天線中的電磁場模式分布也不同，由腔模理論分析可知，這主要是由于脈沖沿不同φ角入射時，脈沖沿不同的傳播路徑傳播，等效于脈沖在不同尺寸的漏波空腔中傳播造成的。

圖9 脈沖沿θ=π/4，不同φ時的Ut波形和頻譜Fig.9 Waveform and spectrum of Ut for differentφ ，but withθ=π/4

3.3 電磁脈沖沿不同θ入射時天線的響應(yīng)

如圖10(a)所示，高功率電磁脈沖分別沿著φ=π/2，不同的θ角入射，得到饋電同軸線中距離口徑100Δl處的瞬態(tài)響應(yīng)電壓Ut的時域波形。由圖可見，入射脈沖在微帶天線內(nèi)產(chǎn)生振蕩，且振蕩幅度逐漸衰減，其原因如前所述。對比來看，θ=0時的響應(yīng)電壓峰值最高，約為10.4 V，θ=π/2時的響應(yīng)電壓峰值最低，約為3.1 V。對照該天線的輻射方向圖可知，當(dāng)入射脈沖沿φ=π/2，θ=0方向入射時，天線接收效能最好，因而此時的天線響應(yīng)電壓最大。

如圖10(b)所示，高功率電磁脈沖分別沿著φ=π/2，不同的θ角入射，得到饋電同軸線中距離口徑100Δl處的瞬態(tài)響應(yīng)電壓Ut的頻譜，由圖可見，脈沖沿θ=0入射時，響應(yīng)電壓主要集中在4.7 GHz頻段，這與天線主模TM01的諧振頻率一致。當(dāng)脈沖沿θ=π/4入射時，與 θ=0相比，在6.5、8、9.4 GHz頻段的響應(yīng)電壓增大，說明由于脈沖的激勵使得天線中出現(xiàn)了較強的TM11、TM20和 TM02等高次模。當(dāng)脈沖沿θ=π/2入射時，與θ=0和π/4相比，6.5、8、9.4 GHz頻段的響應(yīng)電壓不斷增大，說明隨著θ角增大，脈沖激勵起的高次模也越強。造成上述現(xiàn)象的原因與上節(jié)類似。

圖10 脈沖沿φ=π/2，不同θ時的Ut波形和頻譜Fig.10 W aveform and spectrum of Ut for differentθ，but withφ=π/2

4 結(jié)論

本文提出一種新的導(dǎo)體邊緣奇異性處理技術(shù)，并利用該算法計算了同軸饋電微帶天線在高功率電磁脈沖輻照下產(chǎn)生的瞬態(tài)響應(yīng)。通過計算結(jié)果可以看出:

1)電場強度為20 kV/m的高功率電磁脈沖可在微帶天線的同軸饋電線內(nèi)產(chǎn)生最大幅值超過10 V的響應(yīng)電壓。

2)入射高功率脈沖在微帶天線內(nèi)不但激勵起天線的主模 TM10，而且激勵起 TM11、TM20、TM21和 TM02等高次模，這些高次模電磁場通過口徑耦合使得同軸饋電線在一些頻段也產(chǎn)生了的較大響應(yīng)電壓。

3)隨著入射脈沖在xoy面的投影與正x軸夾角的增大，時域上，響應(yīng)電壓均小于10 V;頻域上，高次模逐漸減小，其頻譜分布更集中于天線主模TM01的諧振頻率。

4)隨著入射脈沖與正z軸夾角的增大，時域上，響應(yīng)電壓逐漸減小;頻域上，高次模逐漸增大，其頻譜分布逐漸分散于 TM11、TM20、TM21和 TM02等高次模對應(yīng)的諧振頻率頻段。

［1］ZHANG Ka，SU Donglin，LIAO Yi，et al.Fast EMP response calculation of dipole antennas based on equivalent circuit［C］//10th International Symposium on Antennas，Propagation and EM Theory.Xi＇an，China，2012:1147-1149.

［2］LIN Jinyong，DONG Wenjie，ZHOU Chunmei.Investigation of injury effects to strong electromagnetic pulse on typical flying control system ［C］//2nd International Conference on Instrumentation and Measurement， Computer，Communication and Control.Guangzhou，China，2012:1532-1535.

［3］JIANG Zheng，YIN Wenyan，LIU Qifeng，et al.Transient responses of some antennas under the impact of an intentionally incident high-power electromagnetic pulse［J］.Progress in Electromagnetics Research，2010，105:365-381.

［4］王韶光，魏光輝，代平.超寬帶電磁脈沖及其對電子裝備的影響研究綜述［J］.軍械工程學(xué)院學(xué)報，2006，18(4):16-20.WANG Shaoguang，WEI Guanghui，DAIPing.Research of the impact of UWB electromagnetic pulse on electronic equipment［J］.Journal of Ordnance Engineering College，2006，18(4):16-20.

［5］CAMP M.Influence of the technology on the destruction effects of semiconductors by impact of EMP and UWB pulse［C］//Electromagnetic Compatibility Conference.Minneapolis，USA，2002，1:87-92.

［6］王泰春，賀云漢，王玉芝.電磁脈沖導(dǎo)論［M］.北京:國防工業(yè)出版社，2011:98-113.WANG Taichun，HE Yunhan，WANG Yuzhi.Introduction of electromagnetic pulse［M］.Beijing:National Defense Industry Press，2011:98-113.

［7］楊雨川.導(dǎo)體板上單極天線對電磁脈沖響應(yīng)特性的矩量法分析［J］.強激光與粒子束，2008，20(4):649-652.YANG Yuchuan.Response characteristics and protection techniques of monopole on conductive plane exposed to electromagnetic pulse［J］.High Power Laser and Particle Beams，2008，20(4):649-652.

［8］STUTZMAN W L，THIELE G A.Antenna theory and design［M］.New York:John Wiley＆ Sons，1981:197-202.

［9］MUR G.The modeling of singularities in the finite-difference approximation of the time-domain electroma gneticfields equations［J］.Electromagnetic Compatibility，1981，23(1):377-382.

［10］TAFLOVE A.The finite-difference time-domain method［M］.Boston:Artech House，1995:63-89.

［11］RAILTON J，PAUL D L.Analysis of structures containing sharp oblique metal edges in FDTD using MAMPs［J］.IEEE Transactions on Antennas and Propagation，2010，58(9):2954-2959.

［12］TEIXEIRA F.Summary review on 25 years of progress and future challenges in FDTD and FETD techniques［J］.Applied Computational Electro magnetics Society Journal，2010，25(1):1-14.

［13］ABDULHUSSEIN N E.3-D FDTD design and analysis of a three-band composite microstrip patch antenna for wireless communication applications［C］//ICNIT.Gyeongju，Korea，2011:104-108.

［14］CHOID H，HOEFER.A graded mesh FDTD algorithm for eigenvalue problems［C］//IEEE MTT-S International Microwave Symposium Digest.Warsaw，Poland，1987，151(6):413-417.

［15］BLADEL JV.Singular electromagnetic fields and sources［M］.Oxford:IEEE Press，1991:125-158.

［16］LO Y T，SOLOMAN D.Theory and experiment onmicrostrip antenna［J］.IEEE Transactions on Antennas and Propagation，1979，27(2):137-145.