劉燕婻，潘小敏，盛新慶

(北京理工大學(xué) 信息與電子學(xué)院，北京 100081)

?

一種矩量法精度分析驗(yàn)證方法研究

劉燕婻，潘小敏，盛新慶

(北京理工大學(xué) 信息與電子學(xué)院，北京 100081)

提出一種針對任意形狀目標(biāo)散射矩量法計算的精度驗(yàn)證方法. 根據(jù)等效原理，建立介質(zhì)目標(biāo)積分方程. 令介質(zhì)分界面兩側(cè)的介質(zhì)參數(shù)相同，推導(dǎo)出入射等效源在介質(zhì)體外產(chǎn)生的場為0. 利用這種性質(zhì)，對任意形狀目標(biāo)散射問題矩量法計算精度進(jìn)行驗(yàn)證. 對球和正8面體等算例進(jìn)行的數(shù)值實(shí)驗(yàn)，驗(yàn)證了不同奇異點(diǎn)處理方法下介質(zhì)目標(biāo)積分方程的精度. 數(shù)值實(shí)驗(yàn)證明了轉(zhuǎn)移法和坐標(biāo)變換法的精度優(yōu)于提取法. 結(jié)果同時表明所提出驗(yàn)證方法的有效性和通用性.

矩量法；精度；奇異點(diǎn)

基于積分方程的矩量法廣泛應(yīng)用于電磁計算. 理論上講，積分方程的建立過程沒有任何近似，這保證了矩量法的精度. 從算法設(shè)計過程看，使用局域基函數(shù)來離散積分方程需要對目標(biāo)進(jìn)行剖分. 對面積分方程，剖分單元既可以是平面三角形、四邊形及其組合，也可以是曲面的三角形、四邊形及其組合. 無論采用哪種剖分單元，都會在矩量法的離散過程中引入一些誤差. 更重要的是，當(dāng)場源重合時，矩陣元素的計算會出現(xiàn)奇異積分，奇異點(diǎn)的處理對矩量法的精度有重要影響[1-3]. 如何進(jìn)行精度分析，是矩量法研究的一個非常重要的方面. 對于有解析解的目標(biāo)，例如球，可以將計算結(jié)果與解析解進(jìn)行比較，驗(yàn)證矩量法的精度. 而對于復(fù)雜目標(biāo)，則沒有一種統(tǒng)一方法對矩量法精度進(jìn)行驗(yàn)證. 目前普遍的方法是將數(shù)值解與實(shí)驗(yàn)測量值，或其他算法及商業(yè)軟件計算結(jié)果進(jìn)行比較. 顯然，這種方法本身存在很多局限. 近年來，Ergül等[4]提到使用入射等效源不輻射的性質(zhì)來驗(yàn)證算子的精度. 在這一思想的基礎(chǔ)上，本文提出一種驗(yàn)證任意目標(biāo)矩量法計算精度的方法，并將其用于驗(yàn)證不同奇異點(diǎn)處理方式的精度.

1 驗(yàn)證原理和方法

均勻媒質(zhì)中的矢量波動方程為

(1)

(2)

并矢格林函數(shù)滿足

××(r，r′)-k2(r，r′)=δ(r-r′).

(3)

假定某一均勻媒質(zhì)中，k為常數(shù). 區(qū)域V2被區(qū)域V1包圍，S為V2的邊界，V1和V2內(nèi)部各有一個電流源J1和J2，如圖1所示.

(4)

其中

(5)

(6)

對式(4)左邊使用第二類矢量格林定理

(7)

(8)

V1中的總場E為V1中的電流源J1和V2中的電流源J2分別所產(chǎn)生的場的疊加. 而式(8)表明，E也可以表示為V1中的電流源J1與邊界S上的一組等效源(Js，Ms)產(chǎn)生的場的疊加. 這樣(Js，Ms)產(chǎn)生的場恰好代表了J2在V1中的場，這被稱作是等效原理.

求解介質(zhì)目標(biāo)的散射，可以使用等效原理來建立方程. 設(shè)V1中介質(zhì)參數(shù)為(ε1，μ1) ，V2中為(ε2，μ2).V1中的總場由入射場(E1，H1)與散射場(Es，Hs)疊加而成，即

(9)

(10)

則等效源(Js，Ms)可以分成兩部分，分別稱為入射等效源和散射等效源.

(11)

(12)

由此，式(8)化為

(13)

式(13)表明入射等效源在介質(zhì)體外產(chǎn)生的場為0. 使用這一性質(zhì)，就可以考察矩量法離散所產(chǎn)生的誤差.

以均勻介質(zhì)體的面積分方程CTF為例：

(14)

(15)

其中L和K算子定義為

(16)

(17)

(18)

令ε1=ε2，μ1=μ2，求解方程(14)和(15)得到入射等效源(Ji，Mi)的展開系數(shù)，并由此求得(Ji，Mi)在遠(yuǎn)處所產(chǎn)生的場值. 如取181個觀察點(diǎn)，遠(yuǎn)場的均方根(RMS)可以按照下式進(jìn)行計算為

(19)

從上面的分析中看到，在介質(zhì)體外部Ji和Mi不應(yīng)輻射，所以從式(19)中計算得到的RMS的值就直接反映了計算的均方根誤差.

2 奇異點(diǎn)處理方法及數(shù)值算例

2.1 奇異點(diǎn)處理方法

在矩量法中，當(dāng)場點(diǎn)和源點(diǎn)相距較近時，矩陣元素計算涉及到的被積函數(shù)變化劇烈，尤其是在場點(diǎn)和源點(diǎn)重合在一起的時候，被積函數(shù)是奇異函數(shù)，不能直接使用數(shù)值積分方法進(jìn)行計算.必須要采用一些特殊方法來處理奇異點(diǎn).

介質(zhì)目標(biāo)積分方程中包含L算子和K算子. 在本文中，只考慮L算子中的奇異點(diǎn). 使用RWG基函數(shù)[6]進(jìn)行離散，并使用伽略金法進(jìn)行檢驗(yàn). 對于積分

(20)

有以下幾種處理的方法：

方法1：提取法，從式(20)出發(fā)，對其中主值積分部分的G和G進(jìn)行泰勒展開，對含有1/Rn項的積分通過極坐標(biāo)變換進(jìn)行處理[2,7-8]. 由于積分核中含有的G具有較高的奇異性，其計算精度可能會相對較差.

方法2：轉(zhuǎn)移法，對式(20)使用矢量恒等式，將作用于格林函數(shù)G上的梯度算子轉(zhuǎn)移到試函數(shù)上[2,6]，得到

降低了積分的奇異性. 理論上其計算精度應(yīng)該比提取法高.

方法3：坐標(biāo)變換法，使用等邊參數(shù)空間變換，將奇異性分為共點(diǎn)、共線和共面3種情況分別處理，處理方法詳見文獻(xiàn)[3,9-10].

2.2 數(shù)值算例

以CTF為例，通過數(shù)值實(shí)驗(yàn)來研究和驗(yàn)證上述3種不同的L算子奇異點(diǎn)處理方法對矩量法計算精度的影響.下面實(shí)驗(yàn)中數(shù)據(jù)類型為雙精度，迭代器為GMRES，迭代精度為10-8.

首先選擇半徑為0.5波長(λ)的球和邊長為0.5λ的正8面體來驗(yàn)證算法. 入射波為沿-z方向入射，垂直極化的平面波，幅度為1.

分別計算了球和正8面體散射場的RMS. 圖2中，分別給出了誤差隨不同網(wǎng)格剖分密度的變化情況. 其中網(wǎng)格剖分大小分別為λ/10，λ/20，λ/30和λ/40，球?qū)?yīng)的未知數(shù)個數(shù)為972，4 332，10 092，18 252，正8面體對應(yīng)的未知數(shù)個數(shù)為789，2 946，7 029，12 537. 從圖中可以看到隨著剖分逐漸變密，RMS值也逐漸減小，也即矩量法的誤差逐漸減小. 同時，可以清楚地看到不同的奇異點(diǎn)處理方法下計算結(jié)果精度的差異. 方法1的誤差明顯大于另外兩種方法，這與前面的分析相吻合. 而方法2和方法3則達(dá)到了類似的計算精度. 但是在計算稍復(fù)雜的正8面體時，方法3的精度略優(yōu)于方法2.

為進(jìn)一步驗(yàn)證本文提出方法的有效性和不同奇異點(diǎn)處理方式的精度，計算了幾何構(gòu)圖如圖3所示的圓柱腔. 入射波為沿-z方向入射的垂直極化平面波，幅度為1，頻率為1.0 GHz. 圖4給出了誤差數(shù)據(jù)，未知數(shù)個數(shù)分別為3 039，6 705，16 998. 結(jié)果再次表明，方法1的精度差于方法2和方法3.

前面的算例中，使用入射等效源在外部不輻射的性質(zhì)，驗(yàn)證了不同奇異點(diǎn)處理方式對于矩量法散射問題計算精度的影響，并得到了一些結(jié)論. 為進(jìn)一步驗(yàn)證這些結(jié)論，計算了半徑為0.5λ，相對介電常數(shù)ε為2.0的介質(zhì)球散射場的雷達(dá)散射截面(RCS).和前面的算例相同，入射波為沿-z方向入射，垂直極化的平面波. 圖5給出了以Mie級數(shù)解析解為參照的RCS均方根誤差. 圖中數(shù)據(jù)表明，方法2和方法3的精度優(yōu)于方法1，即轉(zhuǎn)移法和坐標(biāo)變換法的精度優(yōu)于提取法，與前面的結(jié)論相吻合.

3 結(jié) 論

對于任意形狀目標(biāo)，矩量法計算精度的驗(yàn)證非常重要. 根據(jù)等效原理，給出了一種具有明確參考基準(zhǔn)，且適用于任意形狀目標(biāo)的精度驗(yàn)證方法. 使用這種方法對不同奇異點(diǎn)處理方法的精度進(jìn)行了驗(yàn)證. 對球和正8面體等算例進(jìn)行的數(shù)值實(shí)驗(yàn)，表明了所提出驗(yàn)證方法的有效性和通用性.

[1] Polimeridis A G，Vipiana F，Mosig J R，et al. DIRECTFN: fully numerical algorithms for high precision computation of singular integrals in galerkin sie methods[J]. IEEE Trans on Antennas Propagat，2013，61(6):3112-3122.

[2] 盛新慶.計算電磁學(xué)要論[M].2版.合肥:中國科學(xué)技術(shù)大學(xué)出版社，2008：23-24.

Sheng Xinqing. A brief treatise on computational electromagnetics[M]. 2nd ed. Hefei: Press of University of Science and Technology of China，2008:23-24. (in Chinese)

[3] Polimeridis A G，Yioltsus T V. On the direct evaluation of weakly singular integrals in galerkin mixed potential integral equation formulations[J]. IEEE Trans on Antennas Propagat，2008，56(9):3011-3019.

[4] Ergul ?，Gurel L. Discretization error due to the identity operator in surface integral equations[J]. Computer Physics Communications，2009，180:1746-1752.

[5] Chew W C，Tong M S ，Hu B. Integral equation methods for electromagnetic and elastic waves[M]. [S.l.]: Morgan & Claypool，2009.

[6] Rao S M，Wilton D R，Glission A W. Electromagnetic scattering by surfaces of arbitrary shape[J]. IEEE Trans on Antennas Propagat，1982，30(3):409-418.

[7] Wilton D R ，Rao S M ，Glisson A W，et al. Potential integrals for uniform and linear source distributions on polygonal and polyhedral domains[J]. IEEE Trans. on Antennas Propag，1984，32(3):276-281.

[8] Graglia R D. On the numerical integration of the linear shape functions times the 3-D green’s function or its gradient on a plane triangle[J]. IEEE Trans on Antennas Propagat，1993，41(10):1448-1455.

[9] Polimeridis A G，Mosig J R. On the direct evaluation of surface integral equation impedance matrix elements involving point singularities[J]. IEEE Antennas and Wireless Propagation Letters，2011，10:599-602.

[10] Polimeridis A G，Tamayo J M，Rius J M，et al. Fast and accurate computation of hypersingular integrals in garlerkin surface integral equation formulations via the direct evaluation method[J]. IEEE Trans on Antennas Propagat，2011，59(6):2329.

(責(zé)任編輯：劉芳)

A New Method for Accuracy Analysis in Method of Moments

LIU Yan-nan, PAN Xiao-min, SHENG Xin-qing

(School of Information and Electronics, Beijing Institute of Technology, Beijing 100081, China)

A method to validate the accuracy of method of moments (MoM) for arbitrarily shaped targets was proposed. Based on equivalence principle，integral equations for dielectric objects were established. The non-radiating property of incident sources was derived by setting the medium on both sides of the interface the same. According to this property，the accuracy of MoM for arbitrarily shaped targets can be investigated. Numerical experiments were conducted on sphere and octahedron targets. It is validated that among three singularity treatments，nabla operator transfer and coordinate transformation perform better than extraction in terms of accuracy. The results also show that the proposed method is possessed of applicability and generality.

method of moments (MoM); accuracy; singularity

2014-03-06

國家教育部新世紀(jì)優(yōu)秀人才支持計劃(NCET-12-0045)；國家重點(diǎn)基礎(chǔ)研究發(fā)展計劃項目(2012CB720702)；國家自然科學(xué)基金資助項目(60901005，61371002)；北京市優(yōu)秀人才培養(yǎng)資助項目(2012D0009011000002)

劉燕楠(1990—)，女，博士生，E-mail:ynliu1990@gmail.com.

潘小敏(1978—)，男，副研究員，E-mail:xman@bit.edu.cn.

O 441.4

A

1001-0645(2016)07-0723-04

10.15918/j.tbit1001-0645.2016.07.012