倪中非，黃斌科，師振盛

（西安交通大學(xué)電子與信息工程學(xué)院，710049，西安）

波導(dǎo)的傳輸特性一直是微波工程領(lǐng)域的一項(xiàng)重要基礎(chǔ)研究課題。隨著微波傳輸系統(tǒng)的不斷發(fā)展，大量性能各異的非規(guī)則截面波導(dǎo)應(yīng)運(yùn)而生，由于結(jié)構(gòu)的復(fù)雜性，傳統(tǒng)方法在研究這類新型波導(dǎo)的傳輸特性問(wèn)題時(shí)面臨重重困難。目前，非規(guī)則波導(dǎo)問(wèn)題的求解大多是利用有限差分法和有限元法等數(shù)值方法直接進(jìn)行求解［1-2］。由于非規(guī)則波導(dǎo)的截面具有不規(guī)則性，為保證計(jì)算精度，往往需要大量的剖分網(wǎng)格，導(dǎo)致計(jì)算資源消耗大、處理效率低、編程實(shí)現(xiàn)復(fù)雜等問(wèn)題，因此如何高效、精確、方便地求解非規(guī)則波導(dǎo)的傳輸特性問(wèn)題仍是目前微波工程領(lǐng)域的研究熱點(diǎn)之一。

早期的傳統(tǒng)坐標(biāo)變換法以解析保角變換為基礎(chǔ)，其關(guān)鍵問(wèn)題在于如何確定具體的變換函數(shù)，然而這一過(guò)程通常比較困難，導(dǎo)致其在解決實(shí)際工程問(wèn)題上的局限性很大［3］。2006年，Pendry和Leonhardt提出的坐標(biāo)變換理論為新型微波器件的分析和設(shè)計(jì)提供了一條全新的思路［4］。該理論依據(jù)Maxwell方程組基于坐標(biāo)變換的形式不變性，將空間坐標(biāo)變換和媒質(zhì)參數(shù)變換有機(jī)聯(lián)系起來(lái)，使得以往復(fù)雜的分析和設(shè)計(jì)過(guò)程變得直觀高效。目前，該方法已被廣泛應(yīng)用于隨機(jī)粗糙面散射、理想匹配層（PML）吸收邊界、新型微波器件設(shè)計(jì)和多尺度復(fù)雜邊值問(wèn)題數(shù)值仿真等領(lǐng)域［5-7］。

本文將坐標(biāo)變換法引入到非規(guī)則截面波導(dǎo)的分析中，通過(guò)坐標(biāo)變換，將不規(guī)則邊界區(qū)域變換為規(guī)則邊界區(qū)域，而區(qū)域變換的影響則通過(guò)填充變換媒質(zhì)來(lái)等效，從而將原非規(guī)則波導(dǎo)轉(zhuǎn)變?yōu)樘畛渥儞Q媒質(zhì)的等效規(guī)則波導(dǎo)來(lái)處理，降低了邊界處理算法的復(fù)雜度和剖分區(qū)域的網(wǎng)格密度，有效提高了數(shù)值計(jì)算效率。本文利用有限元法結(jié)合坐標(biāo)變換法計(jì)算了橢圓波導(dǎo)和矩形脊波導(dǎo)的截止波數(shù)和場(chǎng)分布，并且與利用有限元法直接對(duì)原問(wèn)題求解的結(jié)果進(jìn)行比較，充分驗(yàn)證了本文方法的高效性。

1 坐標(biāo)變換法基本理論

坐標(biāo)變換法的理論基礎(chǔ)為Maxwell方程組的形式不變性［8］。首先，定義一個(gè)坐標(biāo)變換F：Ω→~Ω，將原空間域Ω中的所有點(diǎn)P映射到變換域~Ω中，記為點(diǎn)~P，即

當(dāng)原空間Ω中的媒質(zhì)無(wú)耗時(shí)，變換媒質(zhì)同樣為無(wú)耗媒質(zhì)。

本文在分析非規(guī)則截面波導(dǎo)模型時(shí)，利用坐標(biāo)變換法將非規(guī)則截面波導(dǎo)變換為規(guī)則截面波導(dǎo)進(jìn)行數(shù)值仿真。在此取非規(guī)則波導(dǎo)為原空間域Ω，經(jīng)坐標(biāo)變換所得的規(guī)則波導(dǎo)為變換域~Ω。在原非規(guī)則波導(dǎo)中電場(chǎng)矢量E（r）和磁場(chǎng)矢量H（r）滿足

由Maxwell方程組的形式不變性可知，變換后的等效規(guī)則波導(dǎo)中電場(chǎng)、磁場(chǎng)滿足

由式（5）～式（8）可知，原非規(guī)則波導(dǎo)和變換后等效規(guī)則波導(dǎo)中的電場(chǎng)、磁場(chǎng)之間存在如下關(guān)系

取波導(dǎo)軸向沿z方向，由于波導(dǎo)軸向的均勻性，可將三維模型簡(jiǎn)化為二維處理。在二維情況下，變換域中的縱向場(chǎng)方程為

由式（15）～式（17）可見(jiàn)，通過(guò)求解坐標(biāo)變換所得等效規(guī)則波導(dǎo)的截止波數(shù)和場(chǎng)分量），即可獲得原非規(guī)則波導(dǎo)的傳輸特性。方程（11）和（12）中的截止波數(shù)和場(chǎng)分量）可利用有限元法轉(zhuǎn)化為廣義本征值問(wèn)題來(lái)計(jì)算［9］

式中：A和B為有限元離散得到的稀疏矩陣。

本文在使用有限元方法計(jì)算時(shí)，剖分單元采用三角線性元，這樣式（4）的Jacobian矩陣就可方便地由剖分單元的基函數(shù)和節(jié)點(diǎn)坐標(biāo)值計(jì)算得到［10］。

2 仿真驗(yàn)證及分析

本文利用坐標(biāo)變換法分析了橢圓波導(dǎo)和矩形脊波導(dǎo)的傳輸特性。通過(guò)坐標(biāo)變換分別將橢圓波導(dǎo)和矩形脊波導(dǎo)轉(zhuǎn)化為等效的圓波導(dǎo)和矩形波導(dǎo)，然后利用有限元法分別計(jì)算了其截止波數(shù)和縱向場(chǎng)分量。截止波數(shù)和縱向場(chǎng)分量在原波導(dǎo)中記為k′c、E′z（x，y）（TM 模）、H′z（x，y）（TE模），在變換后的規(guī)則波導(dǎo)中記為經(jīng)坐標(biāo)反變換后所得原非規(guī)則波導(dǎo)橫截面上的縱向場(chǎng)記為Ez（x，y）和Hz（x，y），由等效波導(dǎo)求得的截止波數(shù)相對(duì)于原非規(guī)則波導(dǎo)截止波數(shù)的誤差系數(shù)記為δkc）×100%，縱向場(chǎng)的誤差系數(shù)記為δf）×100%，其中f′代表E′z（x，y）（TM 模）或者 H′z（x，y）（TE模），f″代表由Ez（x，y）或者 Hz（x，y）在原非規(guī)則波導(dǎo)剖分網(wǎng)格點(diǎn)上進(jìn)行加權(quán)插值得到的場(chǎng)值。本文結(jié)果由Matlab 2012計(jì)算所得，計(jì)算機(jī)配置：內(nèi)存為1GB，CPU為Celeron（R）。

2.1 橢圓波導(dǎo)傳輸特性

利用坐標(biāo)變換法將空氣填充的橢圓波導(dǎo)用填充變換媒質(zhì)的圓波導(dǎo)來(lái)等效，Ω和~Ω分別為橢圓波導(dǎo)截面和等效圓波導(dǎo)截面，其坐標(biāo)變換關(guān)系為F（r）：），即

式中：rα、rβ分別為空間~Ω和Ω中邊界點(diǎn)Pα和Pβ的位置矢量；r為空間Ω中點(diǎn)的位置矢量。橢圓波導(dǎo)到等效圓波導(dǎo)的坐標(biāo)變換過(guò)程如圖1所示。

圖1 橢圓波導(dǎo)到等效圓波導(dǎo)的坐標(biāo)變換示意圖

取橢圓波導(dǎo)截面的長(zhǎng)、短半軸分別為a=10mm，b=6mm，等效圓波導(dǎo)截面的半徑為r=6mm。為較好地逼近曲線邊界，采用三角網(wǎng)格剖分。表1列出了橢圓波導(dǎo)和其等效圓波導(dǎo)前6個(gè)模式的截止波數(shù)及誤差系數(shù)δkc和δf。表2對(duì)本文方法和傳統(tǒng)方法的計(jì)算量和計(jì)算時(shí)間進(jìn)行了對(duì)比。圖2所示為橢圓波導(dǎo)和其等效圓波導(dǎo)橫截面上主模TEc11模的縱向場(chǎng)分布。

表1 橢圓波導(dǎo)和其等效圓波導(dǎo)的截止波數(shù)及相對(duì)誤差

表2 本文方法和傳統(tǒng)方法的計(jì)算量和計(jì)算時(shí)間對(duì)比

由表1可見(jiàn)，經(jīng)坐標(biāo)變換所得的圓波導(dǎo)和原橢圓波導(dǎo)的傳輸特性完全等效，截止波數(shù)的最大相對(duì)誤差不超過(guò)3／10 000，場(chǎng)分布的最大相對(duì)誤差不超過(guò)8／10 000。由表2可見(jiàn)，本文方法在計(jì)算內(nèi)存和計(jì)算時(shí)間上分別比傳統(tǒng)方法減少了近60%和50%，相比于利用有限元法直接求解橢圓波導(dǎo)，經(jīng)坐標(biāo)變換所得的等效圓波導(dǎo)不僅邊界條件處理更為簡(jiǎn)單方便，而且能有效節(jié)省計(jì)算資源，這也是本文方法的主要優(yōu)勢(shì)所在。

圖2 橢圓波導(dǎo)和其等效圓波導(dǎo)橫截面上TEc11模的縱向場(chǎng)分布

2.2 矩形脊波導(dǎo)的傳輸特性

利用坐標(biāo)變換法，將矩形脊波導(dǎo)變換為局部填充變換媒質(zhì)的等效矩形波導(dǎo)，如圖3所示。圖3中坐標(biāo)變換將左圖中虛線框與脊所圍的區(qū)域Ω映射到右圖中的矩形陰影區(qū)域即變換域~Ω，坐標(biāo)變換關(guān)系為

式中：rα、rβ、rγ分別為邊界點(diǎn)Pα、Pβ、Pγ的位置矢量；r為區(qū)域Ω內(nèi)點(diǎn)的位置矢量；變換域~Ω中填充變換媒質(zhì)的參數(shù)由式（2）、式（3）計(jì)算得到。

圖3 矩形脊波導(dǎo)到等效矩形波導(dǎo)的坐標(biāo)變換示意圖

取矩形脊波導(dǎo)的長(zhǎng)l、寬w為100mm、50mm，脊高度h分別為5mm、10mm和15mm；等效矩形波導(dǎo)的長(zhǎng)l、寬w為100mm、50mm，填充變換媒質(zhì)區(qū)域的長(zhǎng)、寬取30mm、15mm。

表3給出了矩形脊波導(dǎo)和其等效矩形波導(dǎo)在不同脊高度下幾種模式的截止波數(shù)及誤差系數(shù)δkc和δf。同樣，表4對(duì)本文方法和傳統(tǒng)方法的計(jì)算量和計(jì)算時(shí)間進(jìn)行了對(duì)比。

由表3可以明顯看出，加脊后主模TE10模的截止波數(shù)降低，且隨脊高度的增加，主模截止波數(shù)進(jìn)一步減小，而TM模及TE01、TE20模的截止波數(shù)增加，從而使得矩形脊波導(dǎo)的工作頻帶變寬。由表4可知，相比于矩形脊波導(dǎo)復(fù)雜的導(dǎo)體邊界，等效矩形波導(dǎo)的邊界簡(jiǎn)單規(guī)則，可以降低空間網(wǎng)格的剖分密度，有效節(jié)省計(jì)算資源。

表3 3種脊高度下矩形脊波導(dǎo)和其等效矩形波導(dǎo)的截止波數(shù)及相對(duì)誤差

圖4所示為脊高度取10mm時(shí)，矩形脊波導(dǎo)和其等效矩形波導(dǎo)橫截面上主模TE10模的縱向場(chǎng)分布。由圖4可見(jiàn)，等效矩形波導(dǎo)中通過(guò)填充的變換媒質(zhì)控制電磁場(chǎng)分布，使得其內(nèi)部場(chǎng)分布和原矩形脊波導(dǎo)的近乎相同。

表4 本文方法和傳統(tǒng)方法的計(jì)算量和計(jì)算時(shí)間對(duì)比

圖4 矩形脊波導(dǎo)和其等效矩形波導(dǎo)橫截面上TE10模的縱向場(chǎng)分布

值得強(qiáng)調(diào)的是，有限元求解時(shí)，對(duì)于取不同脊高度的情形，如果直接對(duì)矩形脊波導(dǎo)進(jìn)行求解，每種情形都需重新建模剖分，而采用等效矩形波導(dǎo)時(shí)，僅建模剖分一次，每次只需重新計(jì)算填充變換媒質(zhì)的參數(shù)即可，有效降低了處理過(guò)程的復(fù)雜度。

3 結(jié) 論

本文將坐標(biāo)變換法引入到非規(guī)則波導(dǎo)傳輸特性問(wèn)題的分析中，通過(guò)坐標(biāo)變換，將非規(guī)則波導(dǎo)轉(zhuǎn)變?yōu)樘畛渥儞Q媒質(zhì)的等效規(guī)則波導(dǎo)。仿真結(jié)果表明，與傳統(tǒng)直接對(duì)非規(guī)則波導(dǎo)進(jìn)行數(shù)值求解相比，本文方法實(shí)現(xiàn)簡(jiǎn)單，在保證高精度的同時(shí)，可有效降低數(shù)值仿真的復(fù)雜度。本文方法的研究結(jié)果可應(yīng)用于曲面邊界、粗糙邊界等復(fù)雜邊值電磁問(wèn)題的高效數(shù)值仿真中。

［1］ 喬海亮，王玥，陳再高，等.全矢量有限差分法分析任意截面波導(dǎo)模式 ［J］.物理學(xué)報(bào)，2013，62（7）：1-8.

QIAO Hailiang，WANG Yue，CHEN Zaigao，et al.Full-vectorial finite-difference analysis of modes in waveguide with arbitrary shape ［J］.Acta Physica Sinica，2013，62（7）：1-8.

［2］ 周平.脊加載橢圓波導(dǎo)傳輸特性的矢量有限元分析［J］.電波科學(xué)學(xué)報(bào)，2009，24（6）：1164-1167.

ZHOU Ping.Propagation characteristics analysis for ridged elliptical waveguide by vector finite element method［J］.Chinese Journal of Radio Science，2009，24（6）：1164-1167.

［3］ 周海京，林為干，丁武，等.新型橢圓函數(shù)波導(dǎo)族的保角變換有限差分解法 ［J］.電波科學(xué)學(xué)報(bào)，1999，14（2）：121-128.

ZHOU Haijing，LIN Weigan，DING Wu，et al，Conformal mapping finite difference method for analyzing a novel family of elliptic function waveguides［J］.Chinese Journal of Radio Science，1999，14（2）：121-128.

［4］ PENDRY J B，SCHURIG D，SMITH D R.Controlling electromagnetic fields［J］.Science，2006，312（5781）：1780-1782.

［5］ OZGUN O，KUZUOGLU M.A transformation media based approach for efficient Monte Carlo analysis of scattering from rough surfaces with objects［J］.IEEE Transactions on Antennas and Propagation，2013，61（3）：1352-1362.

［6］ OZGUN O，KUZUOGLU M.Non-Maxwellian locallyconformal PML absorbers for finite element mesh truncation［J］.IEEE Transactions on Antennas and Propagation，2007，55（3）：931-937.

［7］ KWON D H， WERNER D H.Transformation electromagnetics：an overview of the theory and applications ［J］.IEEE Antennas and Propagation Magazine，2010，52（1）：24-45.

［8］ MILTON G W，BRIANE M，WILLIS J R.On cloaking for elasticity and physical equations with a transformation invariant form ［J］.New Journal of Physics，2006，8（10）：248-267.

［9］ JIN Jianming.The finite element method in electromagnetics［M］.New York，USA：Wiley，2002：233-236.

［10］OZGUN O，KUZUOGLU M.Form invariance of Maxwell’s equations：the pathway to novel metamaterial specifications for electromagnetic reshaping ［J］.IEEE Antennas and Propagation Magazine，2010，52（3）：51-65.