時(shí)域有限差分法連接邊界電磁泄漏

姬金祖

(北京航空航天大學(xué) 航空科學(xué)與工程學(xué)院，北京100191)

黃大慶

(中國(guó)航空工業(yè)集團(tuán)公司北京航空材料研究院，北京100095)

黃沛霖 魯振毅

(北京航空航天大學(xué) 航空科學(xué)與工程學(xué)院，北京100191)

飛行器隱身技術(shù)是現(xiàn)代戰(zhàn)爭(zhēng)中的一項(xiàng)重要軍事技術(shù)，對(duì)于提高飛行器的戰(zhàn)場(chǎng)生存力意義重大［1－2］.時(shí)域有限差分法(FDTD，F(xiàn)inite－Difference Time－Domain)是計(jì)算金屬和復(fù)雜介質(zhì)目標(biāo)電磁散射、傳播特性的一種重要方法，在隱身飛機(jī)概念設(shè)計(jì)、吸波材料設(shè)計(jì)等方面具有重要作用［3－5］.

用FDTD計(jì)算目標(biāo)電磁散射的過程中，一般要把計(jì)算區(qū)域分成總場(chǎng)區(qū)和散射場(chǎng)區(qū)兩個(gè)區(qū)域［6］.總場(chǎng)區(qū)的場(chǎng)量包括入射場(chǎng)和散射場(chǎng)，散射場(chǎng)區(qū)的場(chǎng)量只有散射場(chǎng)，沒有入射場(chǎng)，因此由散射場(chǎng)區(qū)進(jìn)行遠(yuǎn)場(chǎng)外推即可得到目標(biāo)的雷達(dá)散射截面RCS(Radar Cross Section)［7］.總場(chǎng)區(qū)和散射場(chǎng)區(qū)的邊界為連接邊界，入射場(chǎng)通過連接邊界引入到總場(chǎng)區(qū).

散射體位于總場(chǎng)區(qū)內(nèi)，當(dāng)總場(chǎng)區(qū)不存在散射目標(biāo)時(shí)，理想情況下散射場(chǎng)區(qū)的場(chǎng)量為0，總場(chǎng)區(qū)只有入射場(chǎng).但在實(shí)際計(jì)算過程中，即使總場(chǎng)區(qū)沒有散射目標(biāo)，散射場(chǎng)區(qū)的場(chǎng)一般也不會(huì)為0，這是因?yàn)樵谶B接邊界引入入射場(chǎng)時(shí)的數(shù)值計(jì)算產(chǎn)生了誤差，導(dǎo)致電磁泄漏.

文獻(xiàn)［8］采用分裂平面波時(shí)域有限差分法(SP－FDTD，Splitting Plane Wave FDTD)引入入射場(chǎng)，在一維網(wǎng)格中將場(chǎng)量分裂，建立一維輔助網(wǎng)格和三維網(wǎng)格的色散關(guān)系.文獻(xiàn)［9］在電磁散射計(jì)算中使用交錯(cuò)布置的方法引入入射波，該方法基于分裂場(chǎng)麥克斯韋方程組，擴(kuò)展了多點(diǎn)輔助入射場(chǎng)的作用.文獻(xiàn)［10］采用一維輔助多點(diǎn)模型，通過CN(Crank－Nicholson)差分方法引入入射波.該方法將二維色散方程投影到一維，使電磁泄漏在時(shí)間步長(zhǎng)取CFL(Courant Friedrichs Levy)時(shí)達(dá)到很小.

實(shí)際計(jì)算表明，時(shí)間步長(zhǎng)以及電磁波入射方向都對(duì)電磁泄漏有一定影響，本文研究這些參數(shù)對(duì)電磁泄漏的影響的大小.

1 FDTD計(jì)算區(qū)域劃分

計(jì)算目標(biāo)電磁散射過程中，F(xiàn)DTD計(jì)算區(qū)域一般劃分為總場(chǎng)區(qū)和散射場(chǎng)區(qū)，在散射場(chǎng)區(qū)還設(shè)置輻射邊界，用于將近場(chǎng)外推至遠(yuǎn)場(chǎng)，進(jìn)一步計(jì)算RCS.FDTD計(jì)算區(qū)域劃分如圖1所示.

計(jì)算過程中，在連接邊界處引入入射場(chǎng)，在此邊界上的時(shí)間迭代公式要做一些修改，以保證該點(diǎn)周圍的環(huán)量均為散射場(chǎng)或均為總場(chǎng).時(shí)間迭代結(jié)束后，散射場(chǎng)區(qū)只有散射場(chǎng).在輻射邊界根據(jù)等效原理，將電磁場(chǎng)外推至遠(yuǎn)場(chǎng)，得到目標(biāo)的RCS.

圖1 FDTD計(jì)算區(qū)域劃分

當(dāng)總場(chǎng)區(qū)沒有目標(biāo)時(shí)，散射場(chǎng)區(qū)仍然會(huì)有電磁場(chǎng)，這是數(shù)值計(jì)算誤差引入的電磁泄漏.本文對(duì)一維和二維情形采取不同的研究方法.一維情形下，用散射場(chǎng)區(qū)電場(chǎng)的平方和來衡量電磁泄漏的大小.在二維情形下，將散射場(chǎng)區(qū)的電磁場(chǎng)通過遠(yuǎn)場(chǎng)外推得到“無目標(biāo)”下的遠(yuǎn)場(chǎng)RCS，用這種方法得到的RCS類似于微波暗室的背景電平.研究不同時(shí)間步長(zhǎng)、入射角度下的遠(yuǎn)場(chǎng)RCS來研究電磁泄漏的程度.

2 一維情形電磁泄漏

2.1 差分方程

一維情形下電磁場(chǎng)微分方程具有最簡(jiǎn)單的形式，設(shè)電磁波傳播方向?yàn)閦，電場(chǎng)方向?yàn)閤，磁場(chǎng)方向?yàn)閥，則電磁波在真空中傳播的微分方程如式(1)所示.

式中，Ex為電場(chǎng)強(qiáng)度;Hy為磁場(chǎng)強(qiáng)度;μ0為真空磁導(dǎo)率;ε0為真空介電常數(shù).設(shè)c為真空中光速，令，結(jié)合真空中光速公式，將上式改寫成較簡(jiǎn)單的形式，如式(2)所示.

式(2)顯得更加簡(jiǎn)單、對(duì)稱，本文涉及的電場(chǎng)、磁場(chǎng)也均采用上式簡(jiǎn)單、對(duì)稱的形式.

用一階Mur吸收邊界條件在兩個(gè)端點(diǎn)進(jìn)行截?cái)?采用一階差分格式，對(duì)微分方程進(jìn)行離散.將計(jì)算區(qū)域劃分成100個(gè)網(wǎng)格，在距邊界20個(gè)網(wǎng)格處設(shè)置連接邊界條件.將ex在網(wǎng)格離散點(diǎn)處采樣，hy在網(wǎng)格中點(diǎn)處采樣，則得到101個(gè)ex離散點(diǎn)和100個(gè) hy離散點(diǎn)，設(shè)離散點(diǎn)為 ex(i，n)和hy(i，n)，其中，i表示網(wǎng)格點(diǎn)編號(hào)，n表示時(shí)間步.設(shè)入射波為正諧波，波長(zhǎng)為λ.

設(shè)網(wǎng)格長(zhǎng)度為Δz，時(shí)間步長(zhǎng)為Δt，定義空間網(wǎng)格因子為 fs=z/λ，時(shí)間網(wǎng)格因子為 ft=Δz/(cΔt).根據(jù)色散條件和穩(wěn)定性條件的要求，fs應(yīng)大于15，一維情形下ft應(yīng)大于1.fs主要與色散誤差相關(guān)，而本文主要研究時(shí)間網(wǎng)格因子ft的影響，故fs一律取為20.

用fs、ft表示真空中電磁波傳播的差分格式，如式(3)所示.

該差分格式將介電常數(shù)、磁導(dǎo)率、網(wǎng)格尺寸、時(shí)間步長(zhǎng)等參數(shù)用空間網(wǎng)格因子和時(shí)間網(wǎng)格因子進(jìn)行歸一化，形式更加簡(jiǎn)潔，便于編寫程序和進(jìn)行相應(yīng)的計(jì)算.

2.2 時(shí)間網(wǎng)格因子影響

設(shè)電磁波傳播方向?yàn)閺淖笙蛴遥O(shè)置不同時(shí)間網(wǎng)格因子，對(duì)比計(jì)算結(jié)果.運(yùn)行1000個(gè)時(shí)間步后，兩種情形下電場(chǎng)分布對(duì)比如圖2所示.

圖2 不同時(shí)間網(wǎng)格因子電場(chǎng)分布

由圖2可見，ft的取值對(duì)電磁泄漏影響嚴(yán)重.ft=1時(shí)，電磁泄漏最小，散射場(chǎng)區(qū)電場(chǎng)強(qiáng)度基本為0.而ft=10時(shí)，電磁泄漏導(dǎo)致散射場(chǎng)區(qū)電場(chǎng)強(qiáng)度有較大起伏.

為定量研究ft對(duì)電磁泄漏的影響，以散射場(chǎng)區(qū)的電場(chǎng)平方和來表征泄漏的能量.取ft的值從1到10，間隔為1，分別計(jì)算其電磁泄漏能量.電磁泄漏能量隨時(shí)間網(wǎng)格因子的變化如圖3所示.

由計(jì)算結(jié)果可得，當(dāng) ft=1時(shí)，泄漏能量為 1.5179×10－29，非常接近 0，可以認(rèn)為此時(shí)泄漏能量很小，可以忽略不計(jì).由圖3可得，泄漏能量基本隨ft的增加而增加.當(dāng)ft＜4時(shí)，泄漏能量增加很快，但當(dāng)ft＞4時(shí)，能量泄漏增加很慢，基本保持在0.06 到0.07 之間.

圖3 電磁泄漏能量隨時(shí)間網(wǎng)格因子的變化

下面研究ft在1到2之間變化時(shí)，電磁泄漏的變化情況.ft步長(zhǎng)取為0.1，計(jì)算結(jié)果如圖4所示.

圖4 電磁泄漏隨時(shí)間網(wǎng)格因子的變化

由圖4可見，時(shí)間網(wǎng)格因子在1到2之間變化時(shí)，電磁泄漏能量基本上線性增加.

可見，在一維FDTD中，時(shí)間網(wǎng)格因子在滿足穩(wěn)定性條件的基礎(chǔ)上，取值應(yīng)盡量小，即接近于1.

3 二維情形電磁泄漏

3.1 二維情形電磁泄漏比較

下面研究二維情形下的電磁泄漏.為保證穩(wěn)定性，二維情形ft的下限是計(jì)算區(qū)域網(wǎng)格數(shù)為160×160，連接邊界到吸收邊界距離為20個(gè)網(wǎng)格，外推邊界到吸收邊界為10個(gè)網(wǎng)格.fs與一維情形相同，取為20.電磁波在真空中傳播，無散射體時(shí)，不涉及到邊界條件，故橫磁波 TM(Transverse Magnetic)和橫電波TE(Transverse Electric)沒有本質(zhì)區(qū)別.本文以TM波為例來進(jìn)行分析研究.

入射角度定義為入射波來波方向與水平線夾角，逆時(shí)針為正，如圖5所示.

ft取1.5 和10 兩個(gè)值，φ 取0°和45°，計(jì)算不存在散射目標(biāo)時(shí)的遠(yuǎn)場(chǎng)雙站RCS，得到4條曲線.雙站角為0°～359°，步長(zhǎng)1°.將4條曲線在一張圖上進(jìn)行對(duì)比，如圖6所示.二維RCS單位為m，在圖中換算為對(duì)數(shù)值，單位為dBm

圖5 入射角度示意圖

圖6 無目標(biāo)下遠(yuǎn)場(chǎng)雙站RCS

上面計(jì)算的RCS可以理解為類似微波暗室的“背景電平”.由上圖可見，斜向45°照射時(shí)雙站RCS明顯小于正入射時(shí)雙站RCS，一定程度上反映了斜向入射情況下電磁泄漏要小于正入射.將360個(gè)雙站RCS值進(jìn)行算術(shù)平均，RCS均值表如表1所示.

表1 雙站RCS全向均值統(tǒng)計(jì) dBsm

由表1可知，入射角45°下RCS遠(yuǎn)低于入射角0°，ft=1.5時(shí) RCS遠(yuǎn)低于 ft=10.可見，ft與入射角對(duì)電磁泄漏影響很大，且由表1可知其差別甚至接近30 dB.

3.2 時(shí)間網(wǎng)格因子的影響

下面研究ft對(duì)電磁泄漏的影響.計(jì)算過程中取入射角45°，ft取值從 1.5 到10，步長(zhǎng) 1.5，其他參數(shù)同前.RCS均值隨ft的變化曲線如圖7所示.

圖7 RCS均值隨ft變化

由圖7可知，ft從1.5 變化到2.5，RCS均值迅速增加，而ft超過2.5之后，RCS均值的變化區(qū)域平緩，逐漸趨于一個(gè)固定的值，最后維持在約－18 dB左右.ft取值接近穩(wěn)定性條件，可使遠(yuǎn)場(chǎng)RCS降低約30 dB.

3.3 入射角度影響

二維情形下電磁泄漏與入射角度也密切相關(guān).下面研究隨入射角度變化的影響，網(wǎng)格劃分同上.由于對(duì)稱性，入射角取為0°～45°就能夠代表所有的情形.入射角步長(zhǎng)取為5°.

用遠(yuǎn)場(chǎng)外推方法計(jì)算雙站RCS，雙站RCS均值隨入射角的變化如圖8所示.

圖8 雙站RCS均值隨入射角度的變化

由圖8可見，雙站RCS均值隨入射角度增大而減小.在入射角較小時(shí)，雙站RCS均值隨入射角改變比較緩慢，入射角大于15°后，雙站RCS均值改變逐漸加快.當(dāng)入射角為45°時(shí)，達(dá)到最小值，約－37 dB.不同的入射角雙站RCS均值差別可達(dá)30 dB.

4 結(jié)論

本文用不同的衡量方法對(duì)一維、二維情形下連接邊界導(dǎo)致的電磁泄漏進(jìn)行了研究，主要結(jié)論如下:

1)為減小電磁泄漏，時(shí)間步長(zhǎng)在滿足穩(wěn)定性要求的前提下，應(yīng)越小越好.在二維情形下，時(shí)間步長(zhǎng)過大引起的電磁泄漏可比時(shí)間步長(zhǎng)接近穩(wěn)定性條件時(shí)大30 dB.電磁泄漏隨時(shí)間步長(zhǎng)變化的趨勢(shì)先慢后快，時(shí)間網(wǎng)格因子小于2時(shí)，基本呈線性關(guān)系.時(shí)間網(wǎng)格因子大于3后，電磁泄漏隨時(shí)間步長(zhǎng)的增加緩慢變化，趨于穩(wěn)定

2)入射角對(duì)電磁泄漏也有很大影響.當(dāng)入射波斜入射到連接邊界時(shí)，不同的入射角對(duì)電磁泄漏也有很大影響.當(dāng)入射角與計(jì)算邊界成斜角時(shí)，電磁泄漏較小.對(duì)于正方形區(qū)域，入射角度45°時(shí)電磁泄漏最小.垂直入射和斜入射時(shí)，電磁泄漏隨入射角度變化較為緩和，其他角度變化較劇烈.

References)

［1］鄭奎松，葛德彪，魏兵.導(dǎo)彈目標(biāo)的FDTD建模與RCS計(jì)算［J］.系統(tǒng)工程與電子技術(shù)，2004，26(7):896－899 Zheng Kuisong，Ge Debiao，Wei Bing.FDTD modeling of missile target and RCS computation［J］.Systems Engineering and Elec－tronics，2004，26(7):896－899(in Chinese)

［2］李軍，武振波，武哲.穩(wěn)定因子對(duì)FDTD數(shù)值計(jì)算的影響［J］.北京航空航天大學(xué)學(xué)報(bào)，2004，30(1):70－73 Li Jun，Wu Zhenbo，Wu Zhe.Effects of stability factor on FDTD numerical calculation［J］.Journal of Beijing University of Aeronautics and Astronautics，2004，30(1):70－73(in Chinese)

［3］黃志祥，沙威，吳先良，等.辛FDTD 算法［J］.系統(tǒng)工程與電子技術(shù)，2009，31(2):456－458 Huang Zhixiang，Sha Wei，Wu Xianliang，et al.Scheme of symplectic FDTD［J］.Systems Engineering and Electronics，2009，31(2):456－458(in Chinese)

［4］胡曉娟，盧兆林，葛德彪，等.基于三角面元的涂層目標(biāo)FDTD共形網(wǎng)格生成技術(shù)［J］.系統(tǒng)工程與電子技術(shù)，2010，32(9):1884－1888 Hu Xiaojuan，Lu Zhaolin，Ge Debiao，et al.Conformal FDTD mesh－generating scheme for coated targets based on trianglepatch［J］.Systems Engineering and Electronics，2010，32(9):1884－1888(in Chinese)

［5］葛德彪，閆玉波.電磁波時(shí)域有限差分方法［M］.2版.西安:西安電子科技大學(xué)出版社，2005 Ge Debiao，Yan Yubo.Finite－difference time－domain method for electromagnetic waves［M］.2nd ed.Xi′an:Xidian University Press，2005(in Chinese)

［6］王長(zhǎng)清.現(xiàn)代計(jì)算電磁學(xué)基礎(chǔ)［M］.北京:北京大學(xué)出版社，2004 Wang Changqing.Basic of modern computational electromagnetics［M］.Beijing:Peking University Press，2004(in Chinese)

［7］姜彥南，葛德彪，魏兵.時(shí)域有限差分并行算法中的吸收邊界研究［J］.系統(tǒng)工程與電子技術(shù)，2008，30(9):1636－1640 Jiang Yannan，Ge Debiao，Wei Bing.Study on absorbing boundary condition in parallel FDTD algorithm［J］.Systems Engineering and Electronics，2008，30(9):1636－1640(in Chinese)

［8］Wang Hui，Huang Zhixiang，Wu Xianliang，et al.Perfect plane wave injection into 3D FDTD(2，4)scheme［C］//Cross Strait Quad－Regional Radio Science and Wireless Technology Conference(CSQRWC)Volume 1.Harbin:IEEE，2011:40－43

［9］Hadi M F.A versatile split－field 1－D propagator for perfect FDTD plane wave injection［J］.IEEE Transactions on Antennas and Propagation，2009，57(9):2691－2697

［10］Huang Z，Pan G，Pan H K.Perfect plane wave injection for Crank－Nicholson time－domain method［J］.Microwaves，Antennas ＆ Propagation，IET，2010，4(11):1855－1862