(空軍工程大學 導彈學院，陜西 三原 713800)

1 引 言

在采用時域有限差分(Finite Difference Time Domain,FDTD)方法對模型進行電磁場的數(shù)值計算時，先必須對所考慮的模型進行剖分。網(wǎng)格剖分的好壞，不僅影響到計算結果的精度，而且對計算過程的收斂程度起著決定性的作用，因而一種好的網(wǎng)格剖分算法對整個FDTD計算來說有著舉足輕重的意義。網(wǎng)格的描述應符合電磁場的變化，如果在每個計算域中都用細網(wǎng)格來描述，勢必要占用計算機過多的資源和較長的運行時間；如果在每個計算域中采用粗網(wǎng)格來描述，則將會給FDTD的計算帶來較大的誤差。為了解決這一矛盾，通常采用非均勻FDTD法[1-6]。文獻[4]中介紹了一種使網(wǎng)格尺寸在空間按一定的規(guī)律擴展的非均勻網(wǎng)格法，其基本思想是使網(wǎng)格密度沿軸向按等比級數(shù)增大。這種網(wǎng)格生成的方法，由于在整個網(wǎng)格空間都是用等比級數(shù)的方式增長得到的，使得6個場分量的迭代式的一致性較好。李響等人[2]在擴展網(wǎng)格算法的基礎上改進了分區(qū)擴展網(wǎng)格的算法，并在其編寫的非均勻剖分軟件中用到了分區(qū)遞變的非均勻網(wǎng)格算法并實現(xiàn)了設置和剖分結果的可視化。周國祥等人[3]同樣也以可視化的形式實現(xiàn)了非均勻網(wǎng)格剖分，提出了一種很具通用性的較合理的非均勻網(wǎng)格剖分算法，即分區(qū)遞變結合粗網(wǎng)格的非均勻網(wǎng)格算法。

分區(qū)遞變的非均勻網(wǎng)格剖分方法與一般的非均勻網(wǎng)格相比有了較大的改進，但也存在著一定的局限性：網(wǎng)格的最小尺寸和比例常數(shù)要人為設定；在算法處理的過程中，基本上是通過調整最小網(wǎng)格尺寸Δxmin及各網(wǎng)格的對應比例Rx使其能滿足要求，此過程不但結果不好控制，而且可能會產(chǎn)生一些意想不到的結果。對于Δxmin及Rx如何進行適當?shù)恼{整才能得到比較理想的結果，文獻中并沒有給出合理的方法?；诖怂惴ǖ幕舅枷耄疚奶岢隽艘环N更加合理的網(wǎng)格自動生成的算法。

2 改進的分區(qū)遞變非均勻網(wǎng)格剖分算法

首先，通過AutoCAD、3DMax等第三方繪圖軟件畫出模型的立體圖，以特定的格式將其導出；然后，通過程序將該文件讀入，就可以得到此模型的結構，同時也可以得到模型在3個軸向上的分界面及不連續(xù)點的坐標。

以x軸方向為例，假設通過模型的讀入，得到在x軸方向的一系列不連續(xù)點的坐標從小到大分別為x0,x1,x2,…，xn。依據(jù)媒質特性及所計算的頻率范圍，我們可以得到空間最大網(wǎng)格的尺寸Δxmax。依據(jù)Δxmax，取Δxmin=Δxmax/5，此精度一般可以滿足計算中不連續(xù)點處的需要，過小則由于同時還要考慮到時間步的關系，使得時間步過小，要使得內部的場穩(wěn)定必然增加所需要的時間步數(shù)，此值也依據(jù)實際情況手動設定。要保證計算空間場的穩(wěn)定性，相鄰的網(wǎng)格間的尺寸應該平緩地變化而不是突變，防止數(shù)值色散，為此，設定比例系數(shù)Rx≤1.3。

由于x0與xn處于計算空間的邊緣，其網(wǎng)格將延展到PML層中，此處的網(wǎng)格作為第三種形式予以考慮。而中間由x1到xn-1之間的網(wǎng)格均屬于由細網(wǎng)格到粗網(wǎng)格再到細網(wǎng)格的變化過程。

定義靠近分界點的網(wǎng)格均為最小網(wǎng)格Δxmin，以最大的比例系數(shù)1.3增大，直到達到最大網(wǎng)格長度Δxmax，由此可以求得以此規(guī)律增大可以得到的最大網(wǎng)格數(shù)ΔNmax及空間長度ΔLmax。而ΔNmax滿足如下關系式：

(1)

即其增長所達到的最大尺寸小于保證在計算空間內收斂的最大值Δxmax，若繼續(xù)增大將大于這一收斂條件。由此可以得到：

(2)

式中，LowI()表示對所求得的數(shù)向下取整，如LowI(9.2)=9。在得到ΔNmax后，可以求得按此規(guī)律得到的最大空間長度：

(3)

令ΔL=xi-xi-1(1

2.1 中間網(wǎng)格值小于Δxmax

設當半?yún)^(qū)間網(wǎng)格數(shù)為n時，網(wǎng)格剛好達到或超過半空間長度ΔL/2，即：

(4)

(5)

(6)

n=LowI(ΔL/Δxmin)

(7)

在此區(qū)間內可取Δx=ΔL/n。

2.2 中間網(wǎng)格值等于Δxmax

定義除去兩邊的遞增網(wǎng)格長度，中間剩余部分的總長度為ΔL′=ΔL-2ΔLmax。中間的網(wǎng)格長度均為Δxmax，網(wǎng)格個數(shù)為n，則有：

n=UpI(ΔL′/Δxmax)

(8)

2.3 邊緣網(wǎng)格

邊緣的網(wǎng)格由于此處不存在不連續(xù)性，故[x0,x1]區(qū)間為從x1向x0網(wǎng)格尺寸遞增，而在[xn-1,xn]區(qū)間為由xn-1向xn網(wǎng)格尺寸遞增。此區(qū)間一般大于3個網(wǎng)格，由于是單調遞增，不存在雙邊網(wǎng)格數(shù)的問題，故此處只需考慮大于3個網(wǎng)格的情況。此處如果ΔL<ΔLmax，則按照中間網(wǎng)格小于Δxmax的方式進行單邊帶處理；如果ΔL≥ΔLmax，則按照中間網(wǎng)格大于Δxmax的方式進行單邊帶處理。

由以上3種情況可以得到整個x軸方向分區(qū)間在以上3種形式下的剖分算法，在y軸方向和z軸方向可以依據(jù)此算法進行同樣的剖分。

3 剖分實例

根據(jù)以上算法，應用VC++編程對一個微帶貼片天線和一個低通濾波器模型進行剖分，結果如圖1～4所示。由圖可以看出，在不同媒質的交界面和金屬邊界條件處網(wǎng)格都采用最小尺寸，這樣在得到最小網(wǎng)格數(shù)的同時，最大限度地避免了網(wǎng)格尺寸所帶來的數(shù)值色散。

圖1 微帶貼片天線Fig.1 Microstrip patched antenna

圖2 剖分后的微帶貼片天線網(wǎng)格Fig.2 The divided grids for microstrip patched antenna

圖3 低通濾波器Fig.3 Lower-passed filter

圖4 剖分后xy平面的網(wǎng)格Fig.4 The divided grids in xy plane

4 結 論

采用改進算法所得空間網(wǎng)格的尺寸，在整個空間區(qū)域內并不是簡單的遞增關系，這一特性使得空間網(wǎng)格不但能夠避免數(shù)值色散，而且由于與一般算法相比其最大網(wǎng)格尺寸可取得更大一些，故其總的網(wǎng)格數(shù)相對其它算法來說也可以取得更少一些。同時，由于算法采用分區(qū)處理及遞增網(wǎng)格與最大網(wǎng)格相結合的方法，既使小網(wǎng)格出現(xiàn)在場變化劇烈的區(qū)域，又保證了場變化緩慢的區(qū)域均為大網(wǎng)格。與普通的剖分算法相比，這種剖分算法可產(chǎn)生更少的網(wǎng)格數(shù)，大大縮短了計算時間。由于這種剖分過程完全可由程序自動完成，因而該算法可以很好地應用于各種非均勻FDTD網(wǎng)格剖分程序。

參考文獻：

[1] 高本慶. 時域有限差分法[M]. 北京: 國防工業(yè)出版社, 1995.

GAO Ben-qing.Finite-Difference Time-Domain Method[M].Beijing:National Defense Industry Press,1995.(in Chinese)

[2] 李響, 楊明武, 王偉吉,等. 時域有限差分電磁建模系統(tǒng)的研究與設計[J].信息與電子工程, 2005, 3(4): 257-260.

LI Xiang,YANG Ming-wu,WANG Wei-ji,et al.Study of FDTD Electromagnetism Modeling System[J].Information and Electronic Eegineering,2005,3(4):257-260.(in Chinese)

[3] 周國祥,程萍,蔣經(jīng)國,等.直角坐標系下非均勻FDTD網(wǎng)格生成系統(tǒng)[J].微波學報,2005,21(2):56-59.

ZHOU Guo-xiang,CHENG Ping,JIANG Jing-guo,et al. A Non-uniform FDTD Mesh Generation System in Cartesian Coordinate[J]. Journal of Microwaves, 2005,21(2): 56-59. (in Chinese)

[4] 姚廣鋒,王積勤, 劉剛. 采用PML的FDTD方法對矩形微帶天線的研究[J].現(xiàn)代雷達,2003,25(11):36-38.

YAO Guang-feng,WANG Ji-qin,LIU Gang. A Study on the PML Absorbing Boundary Condition of the FDTD Method for Rectangle Microstrip Antenna[J].Modern Radar, 2003, 25(11):36-38. (in Chinese)

[5] Reineix A,Jecko B.Analysis of Microstrip Patch Antennas Using Finite Difference Time Domain Method[J].IEEE Transactions on Antennas and Propagation,1989,37(11):1361-1369.

[6] David M Sheen, Sami M Ali, Mohamed D Abouzahra,et al.Application of the Three-Dimensional Finite-Difference Time-Domain Method to the Analysis of Planar Microstrip Circuits[J].IEEE Transactions on Microwave Theory and Techniques,1990,38(7):849-857.