李 超，何思遠(yuǎn)＊，朱國強(qiáng)，鄧方順，陶 華，肖柏勛

1 武漢大學(xué)電子信息學(xué)院，武漢 430079

2 中國船舶重工集團(tuán)第722研究所，武漢 430079

3 中國人民解放軍92474部隊(duì)

4 長江大學(xué)工程地球物理研究中心，湖北荊州 434023

1 引 言

目標(biāo)與粗糙面復(fù)合電磁散射研究在雷達(dá)目標(biāo)探測、目標(biāo)識(shí)別以及微波遙感等領(lǐng)域有著廣泛的應(yīng)用[1-20］.求解粗糙面上方或下方目標(biāo)的電磁散射問題必須充分考慮目標(biāo)與粗糙面之間的相互作用.Li等[1］提出用廣義前后向迭代算法與譜加速算法結(jié)合（GFBM／SAA）計(jì)算海洋粗糙面上方目標(biāo)的電磁散射場.Deng等[2］借助傳統(tǒng)UV算法，結(jié)合分層思想，發(fā)展得到適用于三維矢量波散射問題的MLUV算法，并成功將其應(yīng)用到目標(biāo)位于粗糙面上方和目標(biāo)位于粗糙面上時(shí)的復(fù)合散射問題.對于埋地問題，相關(guān)文獻(xiàn)[3-4］用分層媒質(zhì)并矢格林函數(shù)分析了平面分層大地中埋入目標(biāo)的散射特性.在此基礎(chǔ)上，Zhang等[5］引入離散復(fù)鏡像技術(shù)（DCIM），加快了分層媒質(zhì)并矢格林函數(shù)的計(jì)算，提高了計(jì)算效率.徐利明等[6］則利用積分方程以及半空間并矢格林函數(shù)的快速算法對埋地目標(biāo)矢量電磁散射進(jìn)行了計(jì)算.在一些埋入體散射問題中，需要計(jì)及大地表面起伏特性對散射場的影響.對于呈復(fù)雜曲面的地表面，埋入體的散射問題將變得更為復(fù)雜.為求解該類問題，一些文獻(xiàn)采用基于矩量法的快速算法計(jì)算粗糙地面埋入體的電磁散射問題，如El-Shenawee等[7-8］采用最陡下降的快速多級子方法（SDFMM），Zhang等[9］則采用了稀疏矩陣正則網(wǎng)格（SMCG）算法.

本文基于矩量法，在粗糙面和埋地目標(biāo)分別建立電場積分方程和PMCHW方程，結(jié)合MLUV算法，求解粗糙地面下方具有任意表面金屬目標(biāo)的電磁散射問題.三維多層UV方法[2］適用于三維矢量波散射問題的精確數(shù)值算法，其基于電磁相互作用的排序抽樣算法有效解決了因積分核震蕩導(dǎo)致傳統(tǒng)UV分解失效的問題，可實(shí)現(xiàn)任意形狀三維導(dǎo)電和介質(zhì)體目標(biāo)散射問題的求解.且通過矩陣快速填充技術(shù)[10］，在內(nèi)存和計(jì)算時(shí)間上的復(fù)雜度均可達(dá)到O（NlogN），其中N為未知量個(gè)數(shù).

2 耦合積分方程

假設(shè)三維金屬目標(biāo)埋于二維粗糙面下方，如圖1（截面圖）所示.上半空間為空氣，標(biāo)記為區(qū)域1，下半空間為大地，標(biāo)記為區(qū)域2.地表面為隨機(jī)粗糙面.E1，H1為區(qū)域1中的電磁場，E2，H2為區(qū)域2中的電磁場，Ei，Hi表示入射場.區(qū)域1、2的電磁參數(shù)分別為ε1，μ1，ε2，μ2，粗糙面和目標(biāo)表面法向向量分別為ns，no.下文中，時(shí)諧因子e-jωt省略.

圖1 三維目標(biāo)淺埋于二維粗糙面下方截面圖Fig.1 Cross section of 3-D object buried under the rough surface ground

區(qū)域1中任意場點(diǎn)r處的總電磁場為

區(qū)域2中任意場點(diǎn)r處的總電磁場為

算子L和K形式如下：

在PEC目標(biāo)表面利用電場邊界條件：

可得表面積分方程如下：

在粗糙面表面：

在PEC目標(biāo)表面：

使用 RWG 基函數(shù)f（r）[21］將表面等效電流、磁流離散，假定粗糙面和目標(biāo)表面離散后的未知量個(gè)數(shù)分別為N和P，其展開式如下：

將式（9）分別帶入積分方程（6～8），并使用伽略金（Galerkin）方法，可得到離散后的矩陣方程組：

其中I1，I2為分別存在于粗糙面和目標(biāo)上的未知量，兩者對應(yīng)的維數(shù)為2 N×1，P×1.

方程（10）將采用Bicgstable（ell）[22］迭代算法進(jìn)行求解.

3 MLUV算法

對未知量為N的矩陣方程，傳統(tǒng)矩量法對內(nèi)存的需求為O（N2），其中N為未知量個(gè)數(shù)，求解矩陣方程的計(jì)算復(fù)雜度為O（N3）.多層UV方法基于低秩矩陣分解技術(shù)，通過將原有矩陣分解為具有強(qiáng)相互作用的強(qiáng)子塊矩陣和弱相互作用的弱子塊矩陣，并根據(jù)矩陣元素的Green函數(shù)的性質(zhì)，將每個(gè)弱子塊矩陣分解為兩個(gè)小矩陣相乘的形式.這里以一個(gè)大小為m×n的遠(yuǎn)區(qū)弱子塊矩陣A為例說明UV方法的具體實(shí)施過程：

（1）使用文獻(xiàn)[10］提出的矩陣快速填充方法計(jì)算遠(yuǎn)區(qū)弱子塊矩陣A，并計(jì)算其行和列對應(yīng)的二階范數(shù)向量Rl2和Cl2；

（4）根據(jù)列索引SC，從中抽取相應(yīng)的s列構(gòu)成列矩陣使用揭示秩的矩陣正交三角分解（QR分解）對ACm×s作 QR 分 解，即其中為矩陣的秩；

（5）根據(jù)行索引Sr，從和中抽取相應(yīng)的s行構(gòu)成行矩陣，求解矩陣方程，

通過以上的步驟，矩陣A可以分解成A=Um×r·Vr×n.從這里可以看出，通過UV分解處理后，實(shí)際需要存儲(chǔ)的矩陣元素為r×（m+n），此外，在迭代過程，通過先將矩陣Vr×n與迭代中間值相乘后再與Um×r矩陣相乘實(shí)現(xiàn)矩陣向量積（Matrix Vector Multiplication），這樣即可減小計(jì)算復(fù)雜度.

4 數(shù)值計(jì)算與討論

取粗糙面為高斯譜粗糙面，相關(guān)長度lx，ly，均方根高度h，長度為Lx×Ly，金屬球位于粗糙面下方，半徑為R，埋地深度為dp，球體在xy面內(nèi)投影坐標(biāo)為 （xt，yt）.

為了消除數(shù)值計(jì)算中截取有限大小粗糙面產(chǎn)生的邊緣效應(yīng)，選用錐形平面波作為入射波.

其中，

算法的驗(yàn)證：取 Lx=Ly=8λ0，lx=ly=0.5λ0，h=0.02λ0，粗糙面相對介電常數(shù)εr=2-0.2j，其中λ0為真空中電磁波波長.金屬球半徑0.3λ0，埋地深度為0.6λ0，xt=yt=0.入射波采用水平極化，入射俯仰角θi=20°，錐型波因子g=Lx／4.方位角取為0°，36°，…，324°共10個(gè)，分為10次實(shí)現(xiàn)，并采用算數(shù)平均值統(tǒng)計(jì).圖2為3DMLUV算法和SDFMM算法及MOM計(jì)算結(jié)果比較，驗(yàn)證了3DMLUV的精度.

圖2 埋地金屬球雙站RCS結(jié)果Fig.2 Bistatic RCS of PEC sphere buried under rough surface

在表1中，記錄了在一次實(shí)現(xiàn)情況下，程序耗用的計(jì)算時(shí)間和內(nèi)存.由表可知，在3DMLUV框架下，不僅計(jì)算時(shí)間得到縮減，程序所需內(nèi)存也大幅度降低.計(jì)算平臺(tái)CPU為I7，主頻為2.8GHz，內(nèi)存2G.

表1 一次實(shí)現(xiàn)時(shí)內(nèi)存和計(jì)算時(shí)間的比較Table 1 Memory and time requirements for one realization

埋入球體目標(biāo)電尺寸參數(shù)變化對雙站RCS的影響.取粗糙面參數(shù)為Lx=Ly=10λ0，h=0.02λ0，lx=ly=0.5λ0，粗糙面相對介電常數(shù)εr=2-0.2j，其中λ0為真空中電磁波波長.入射波采用水平極化，入射俯仰角θi=20°.所有計(jì)算結(jié)果均是對50個(gè)樣本取算數(shù)平均值得到的.金屬球半徑分別為0.8λ0，1.2λ0，1.8λ0，埋地深度均為2λ0，xt=yt=0.相應(yīng)的E面內(nèi)雙站RCS如圖3所示.由圖可知，隨著圓柱半徑的增大，目標(biāo)與粗糙面的相互作用增強(qiáng)，雙站RCS增大，尤其在鏡面反射方向附近以外的很大范圍內(nèi)這一結(jié)論表現(xiàn)的尤為明顯.

圖3 粗糙面下不同半徑球體對應(yīng)的雙站RCSFig.3 Bistatic RCS of PEC Sphere with different radius buried under rough surface

粗糙面參數(shù)對雙站RCS的影響.其它計(jì)算參數(shù)保持不變，取粗糙面均方根高度分別為0.02λ0，0.03λ0，0.04λ0，此時(shí)，埋入球體半徑固定為1.2λ0.圖4為對應(yīng)的雙站RCS結(jié)果，其中圖4a為E面內(nèi)全角度雙站RCS，圖4b為鏡向附近角度內(nèi)雙站RCS.由圖可知，隨著均方根高度的增大，除鏡面反射方向附近一個(gè)小范圍內(nèi)RCS減小以外，其它散射角處RCS均增大.這是因?yàn)?，隨著粗糙度的增加，粗糙面漫反射現(xiàn)象越來越明顯.

埋地深度變化對雙站RCS的影響.其它參數(shù)保持不變，令球體埋藏深度取1.5λ0，2λ0，4λ0.對應(yīng)的雙站RCS如圖5所示，由圖可知，隨著目標(biāo)埋藏深度的增加，目標(biāo)與粗糙面的距離增大，相互作用減弱，鏡向反射方向外的RCS明顯變小.

最后，圖6為埋地目標(biāo)為立方體時(shí)的雙站RCS，立方體邊長為2λ0，圖7為埋地目標(biāo)為圓柱時(shí)的雙站RCS，立方體底面直徑為2λ0，高度為3λ0，其它參數(shù)設(shè)置均保持不變.計(jì)算結(jié)果表明了算法的通用性.

5 結(jié) 語

圖4 粗糙面下方埋有球體時(shí)不同粗糙面粗糙度的雙站RCS（a）E面內(nèi)雙站RCS；（b）鏡面附近角度內(nèi)雙站RCS.Fig.4 Bistatic RCS of rough surface with different surface height while buried object is a sphere（a）Bistatic RCS in E-plane；（b）Bistatic RCS in specular angle.

本文基于MLUV算法計(jì)算了二維粗糙面下方任意形狀金屬目標(biāo)的雙站RCS.并通過分析金屬球埋于粗糙面下方時(shí)，埋地目標(biāo)電尺寸，粗糙面粗糙度以及埋地深度對雙站RCS的影響.結(jié)果表明，隨著目標(biāo)尺寸增大和埋地深度的減小，目標(biāo)與粗糙面的相互作用增強(qiáng)，RCS明顯增大，反之則RCS減小.同時(shí)，粗糙面的粗糙度對RCS也有明顯影響.最后，通過立方體和圓柱的算例表明了算法的通用性.

圖6 埋地目標(biāo)為立方體時(shí)的雙站RCSFig.6 Bistatic RCS with cube buried

（References）

[1]Li Z X，Jin Y Q.Bistatic scattering from a fractal dynamic rough sea surface with a ship presence at low grazing angle incidence using the GFBM／SAA.Microw.Opt.Technol.Lett.，2006，31（2）：146-151.

[2]Deng F S，He S Y，Chen H T，et al.Numerical simulation of vector wave scattering from the target and rough surface composite model with 3-D multilevel UV method.IEEE Trans.Antennas Propagat，2010，58（5）：1625-1634.

[3]Cui T J， Wiesbeck W，Herschlein A.Electromagnetic scattering by multiple three-dimensional scatterers buried under multilayered media.II.Numerical implementations and results.IEEE Trans.Geosci.Remote Sensing，1998，36（2）：535-546.

[4]Cui T J，Chew W C.Fast evaluation of Sommerfeld integrals for EM scattering and radiation by three-dimensional buried objects.IEEE Trans.Geosci.Remote Sensing，1999，37（2）：887-900.

[5]Zhang Y H，Xiao B X，Zhu G Q.An improved weak-form BCGS-FFT combined with DCIM for analyzing electromagnetic scattering by 3-D objects in planarly layered media.IEEE Trans.Geosci.Remote Sensing，2006，44（12）：3540-3546.

[6]徐利明，聶在平.埋地目標(biāo)體矢量電磁散射的一種快速正演算法.地球物理學(xué)報(bào)，2005，48（1）：209-215.Xu L M，Nie Z P.A fast forward algorithm for modeling vector electromagnetic scattering from buried dielectric objects.Chinese J.Geophys.（in Chinese），2005，48（1）：209-215.

[7]El-Shenawee M，Rappaport C，Mille M L，et al.Threedimensional subsurface analysis of electromagnetic scattering from penetrable／PEC objects buried under rough surfaces：Use of the steepest descent fast multipole method.IEEE Trans.Geosci.Remote Sensing，2001，39（6）：1174-1182.

[8]El-Shenawee M.The multiple interaction model for nonshallow scatterers buried beneath 2-D random rough surfaces.IEEE Trans.Geosci.Remote Sensing，2002，40（4）：982-987.

[9]Zhang G F，Tsang L，Pak K.Angular correlation function and scattering coefficient of electromagnetic waves scattered by a buried object under a two-dimensional rough surface.J.Opt.Soc.Am.A.，1998，15（12）：2995-3002.

[10]Chen H T，Luo J X，Zhu G Q.Using UV technique to accelerate the MM-PO method for three-dimensional radiation and scattering problem.Microw.Opt.Technol.Lett.，2006，48（8）：1615-1618.

[11]雷銀照，馬信山.均勻半空間導(dǎo)體的并矢格林函數(shù)解析式.地球物理學(xué)報(bào)，1997，40（2）：265-271.Lei Y Z，Ma X S.An analytical formljl of dyadic Green′s function for homogeneous half-space conductor.Chinese J.Geophys.（in Chinese），1997，40（2）：265-271.

[12]魯來玉，張碧星，鮑光淑.電阻率隨位置線性變化時(shí)的三維大地電磁模擬.地球物理學(xué)報(bào)，2003，46（4）：568-575.Lu L Y，Zhang B X，Bao G S.Modeling of three-dimensional magnetotelluric response for a linear Earth.Chinese J.Geophys.（in Chinese），2003，46（4）：568-575.

[13]閆沛文，童創(chuàng)明.基于FGMRES-FBTG算法的介質(zhì)粗糙面散射特性的模擬.電波科學(xué)學(xué)報(bào)，2009，24（1）：115-119.Yan P W，Tong C M.Simulations of scattering characteristic of lossy dielectric surfaces based on FGMRES-PBTG method.Chinese Journal of Radio Science （in Chinese），2009，24（1）：115-119.

[14]Vitebskiy S，Carin L.Moment-method modeling of shortpulse scattering from and the resonances of a wire buried inside a lossy，dispersive half-space.IEEE Trans.Antennas Propagat.，1995，43（11）：1303-1312.

[15]Boix R R，Alexopóulos N G，Horno M.Efficient numerical computation of the spectral transverse dyadic Green′s function in stratified anisotropic media.J.Electromagn.Waves Appl.，1995，10（8）：1047-1083.

[16]Tran P，Maradudin A A.Scattering of a scalar beam from a two-dimensional randomly rough hard wall：Enhanced backscattering.Phy.Rev.B，1992，45（7）：3936-3939.

[17]Xiong Z H，Tripp A C.3-D electromagnetic modeling for near-surface targets using integral equations.Geophysics，1997，62（4）：1097-1106.

[18]Jandhyala V，Shanker B， Michielssen E，et al.Fast algorithm for the analysis of scattering by dielectric rough surfaces.J.Opt.Soc.Am.A，1998，15（7）：1877-1885.

[19]Chen Y H，Chew W C，Oristaglio M L.Application of perfectly matched layers to the transient modeling of subsurface EM problems.Geophysics，1997，62（6）：1730-1736.

[20]Thorsos E.The validity of the Kirchhoff approximation for rough surface scattering using a Gaussian roughness spectrum.J.Acoust.Soc.Am.，1988，83（1）：78-92.

[21]Rao S M，Wilton D R，Glisson A W.Electromagnetic scattering by surfaces of arbitrary shape.IEEE Trans.Antennas Propagat.，1982，30（3）：409-418.

[22]Topsakal E，Kindt R，Sertel K，et al.Evaluation of the BiCGSTAB（l）algorithm for the finite-element／boundaryintegral method.IEEE Trans.Antennas Propagation Magazine，2001，43（6）：124-131.