矢量波動(dòng)方程的新瀑布型多重網(wǎng)格方法

陸康梅，李郴良，曹艷斌

(桂林電子科技大學(xué)數(shù)學(xué)與計(jì)算機(jī)科學(xué)學(xué)院，廣西桂林541004)

在微波理論和技術(shù)中，諧振腔本征值問題是最基本的問題之一．很多微波部件和系統(tǒng)的分析與最優(yōu)化設(shè)計(jì)又往往以該問題的求解為基礎(chǔ)．矢量有限元方法是近10年來在電磁場(chǎng)計(jì)算中應(yīng)用比較廣泛的一種方法．只要選擇了合適的矢量基，所考慮結(jié)構(gòu)的內(nèi)部和外部邊界都能夠從數(shù)學(xué)上自然滿足，就能夠很好地解決偽解問題;又因?yàn)橛邢拊椒ǖ木W(wǎng)格劃分能很好地模擬實(shí)際結(jié)構(gòu)，因此我們選擇了矢量有限元對(duì)諧振腔進(jìn)行離散計(jì)算．

多重網(wǎng)格方法，對(duì)于求解由微分方程離散化得到的方程組來說，是目前最快速高效的方法之一，它的求解的工作量可降為O(n)或O(nlnn)．因此它在計(jì)算電磁場(chǎng)問題上得到了廣泛的關(guān)注，最近一些學(xué)者發(fā)展了基于棱邊元的多重網(wǎng)格算法［1－6］．瀑布型多重網(wǎng)格無(wú)需粗網(wǎng)格校正，它比一般的多重網(wǎng)格方法計(jì)算量減少了，從而提高了計(jì)算效率，然而基于棱邊元上的瀑布型多重網(wǎng)格方法的研究還是比較少的，特別是矢量場(chǎng)的本征問題，目前研究的也較少．

外推算法是由林群，陳傳淼等［7－8］引入到有限元求解偏微分方程，楊一都在文獻(xiàn)［9］中引進(jìn)了本征有限元外推的一個(gè)新技術(shù)，李永明等［10］將有限元外推應(yīng)用到了波導(dǎo)本征問題．本文在其基礎(chǔ)上，將其推廣應(yīng)用到矢量波動(dòng)本征問題，具體的算例表明其可行性和高精度性．

本文結(jié)合外推技術(shù)，提出了一種基于矢量場(chǎng)本征問題的瀑布型多重網(wǎng)格方法．?dāng)?shù)值算例結(jié)果說明該方法的有效性和實(shí)用性．

1 數(shù)學(xué)基礎(chǔ)

對(duì)于一個(gè)填充相對(duì)介電常數(shù)為εr和相對(duì)磁導(dǎo)率為μr介質(zhì)的封閉諧振腔體，對(duì)應(yīng)的矢量波動(dòng)方程為

其合適的變分公式為

其中

這里S為諧振腔的面積．為了離散F，我們把S細(xì)分為小的矩形，每個(gè)小塊占有面積Se(e=1，2，3，…，M)，其中M是單元總數(shù)．在每個(gè)單元中，電場(chǎng)近似為

其中

Se表示單元面積．在進(jìn)行求和并使用全局編碼后，(4)式可寫

應(yīng)用里茲方法，可得到本征方程組

在強(qiáng)加狄利克雷邊界條件以使得腔體壁上的棱邊場(chǎng)為0后，由上式可解出本征值k20和對(duì)應(yīng)的本征向量．

2 求解矢量場(chǎng)本征問題的瀑布型多重網(wǎng)格法

對(duì)于不同的剖分步長(zhǎng)，矢量有限元離散方程(1)產(chǎn)生的相應(yīng)的一系離散方程組

需要求解最細(xì)網(wǎng)格層上的本征方程，即Alul=λlBlul．本文將提出一種求解矢量場(chǎng)本征問題的瀑布型多重網(wǎng)格方法．該方法是一種基于棱邊上的線性插值多重網(wǎng)格方法，并通過共軛梯度法來改善多重網(wǎng)格迭代的收斂性，并結(jié)合了外推技術(shù)提高解的精度．本文還將文獻(xiàn)［10］中的基于節(jié)點(diǎn)元上的外推方法推廣至棱邊元上，具體的算例表明算法具有高精度性．

2．1 棱邊上的插值算子

在求得粗網(wǎng)格層上的本征向量ul－1后，我們通過線性插值可以得到細(xì)網(wǎng)格層上的本征向量ul，

其中P為線性插值算子，矩陣P的構(gòu)造可有下式得到［1］

再根據(jù)Rayleigh商公式可以得到細(xì)網(wǎng)格上的近似本征值．

2．2 光滑過程

由上述的(5)式和(6)式，可以得到細(xì)網(wǎng)格上本征問題

的近似本征值?λl和本征向量?ul，我們可以以它們作為初始的迭代特征對(duì)，通過迭代過程可以得到細(xì)網(wǎng)格上的更好的近似解λl，ul．共軛梯度法是求解稀疏線性方程組，特別是對(duì)稱正定方程組的最有效的方法之一，并且矩陣和向量相乘的次數(shù)也比較少，所以本文的迭代過程采用了文獻(xiàn)［11］提出的求解本征問題的共軛梯度迭代．

令方程(7)的右端項(xiàng)為

則求解方程(7)的共軛梯度法的迭代過程如下：

算法 1［11］(NCGM)

其中ε為精度控制參數(shù)．以上迭代過程簡(jiǎn)記為(λj，uj)=

2．3 外推技術(shù)

定理1 設(shè)剖分是強(qiáng)正規(guī)的，又設(shè)λ為對(duì)應(yīng)于波導(dǎo)本征問題的簡(jiǎn)單本征值，相應(yīng)的規(guī)格化本征函數(shù)u∈C4(Ω)，則本征值外推估計(jì)為：

詳細(xì)證明見文獻(xiàn)［9］．

類似的，設(shè)某矩形諧振腔在粗矩形剖分單元下計(jì)算得到的本征值為λh，隨后在此粗剖分基礎(chǔ)上加密剖分1次，再計(jì)算得到的本征值為λh/2，則根據(jù)定理1，此矩形諧振腔本征值可由外推記為 λhw=E(λh/2，λh))，得到更精確的解．

2．4 矢量波動(dòng)本征問題的瀑布型多重網(wǎng)格算法

對(duì)于矢量有限元離散產(chǎn)生的一系列本征方程組Aiui=λiBiui(i=0，1，2，…，l)，若要求解最細(xì)網(wǎng)格層l層上本征方程Alul=λlBlul，我們結(jié)合外推技術(shù)給出如下求解該本征問題的新瀑布型多重網(wǎng)格算法．

算法2GWCMG

表1 最細(xì)網(wǎng)格(128×128)層上本文的GWCMG方法的數(shù)值結(jié)果

3 數(shù)值算例及分析

矢量波動(dòng)方程一般出現(xiàn)在電磁學(xué)中的腔體諧振問題以及波在封閉或開放結(jié)構(gòu)中的傳播問題．而在這些問題中，人們關(guān)注的是確定對(duì)應(yīng)于本征值的諧振頻率或轉(zhuǎn)播常數(shù)和對(duì)應(yīng)的諧振模式或轉(zhuǎn)播模式，它們一般都可以由特征值來確定，因此關(guān)鍵是特征值的求解．下面我們以求解矩形諧振腔本征值(諧振頻率)問題為例，其尺寸為(a×b)=2m×2m，計(jì)算矩形諧振腔的一些模式的本征值，并將以上算法得到的解λGWCMG、一般瀑布型多重網(wǎng)格方法解λCMG、一般有限元解λFEM和真解λ進(jìn)行比較，并分別列出其與真解的相對(duì)誤差，分別表示為eGW、eCM和eF．

表2 最細(xì)網(wǎng)格(128×128)層上一般瀑布型多重網(wǎng)格法的數(shù)值結(jié)果

表3 最細(xì)網(wǎng)格(128×128)層上一般有限元法的數(shù)值結(jié)果

由表1，表2和表3可以看到本文的算法使計(jì)算精度得到了很大的提高．

由圖2，我們也可以看出隨著網(wǎng)格的加密計(jì)算的精度越來越高，從而驗(yàn)證了算法的快速收斂性，且和一般的瀑布型多重網(wǎng)格方法比較計(jì)算精度也提高了很多．

4 結(jié)語(yǔ)

從上述的算法公式，以及矩形諧振腔的諧振頻率計(jì)算結(jié)果的分析可以看到，利用瀑布型多重網(wǎng)格方法大大減少了計(jì)算量，將外推推廣應(yīng)用到諧振腔本征問題求解精度得到了很大的提高，因此在比較少的單元剖分下便可以較大的提高計(jì)算精度．計(jì)算過程也比較簡(jiǎn)單，易于編程，因此該算法在現(xiàn)代數(shù)值計(jì)算中是一種十分實(shí)用、簡(jiǎn)單、高效的新方法．

［1］ WATANABE K ，IGARASHI H ．Robustness of nested multigrid method for edge－based finite element analysis［J］．IEEE Transactions on Magnetics，2009，45(3)：1 088 －1 091．

［2］ TSAI C L ，WANG W S ．An improved multigrid tcchnique for quasi－TEM analysis of a microstrip embedded in an inhomogeneous anisotropic medium［J］．IEEE Trans．Microwave Theory Tech，1997，45(5)：678 －686．

［3］ WEISS B ，BIRO O ．Edge element multigrid solution of nonlinear magnetostatic problems［J］．Comple，2001，20(2)：357 －365．

［4］ SCHINNERL M ，SHOBERL J ，KALTENBACHER M．Nested multigrid methods for the fast numerical computation of 3D magnetic fields［J］．IEEE Trans Magn，2000，36：1 539 －1 542．

［5］ WATANABE K ，IGARASHI H ，HINMA T．Comparison of geometric and algebraic multigrid methods in edge－based finite－element analysis［J］．IEEE Transactions on Magnetics，2005，41(5)：1 672 －1 675．

［6］金建銘，王建國(guó)，葛德彪．電磁場(chǎng)有限元方法［M］．西安：西安電子科技大學(xué)出版社，2001：165－186．

［7］陳傳淼，黃云清．有限元高精度理論［M］．長(zhǎng)沙：湖南科學(xué)技術(shù)出版社，1995：451 －492．

［8］李郴良，陳傳淼，許學(xué)軍．基于超收斂和外推方法的一類新的瀑布型多重網(wǎng)格方法［J］．計(jì)算數(shù)學(xué)，2007，37(9)：1 083 －1 098．

［9］楊一都．本征值有限元外推的一個(gè)新技術(shù)［J］．貴州大學(xué)學(xué)報(bào)：自然科學(xué)版，1989，6(3)：6 －11．

［10］李永明，俞集輝．有限元外推法在波導(dǎo)本征值問題中的應(yīng)用［J］．重慶大學(xué)學(xué)報(bào)：自然科學(xué)版，1999，22(4)：78－81．

［11］匡前義，李郴良，李明．一種求解廣義特征值的瀑布型多重網(wǎng)格方法［J］．云南民族大學(xué)學(xué)報(bào)：自然科學(xué)版，2009，18(3)：202 －205．