一種任意比例LTS技術(shù)在DG-FETD中的應(yīng)用

李 敏，沈微微，葉曉東，楊逸婷

(1.宿遷學(xué)院信息工程學(xué)院，江蘇 宿遷 223800；2.南京郵電大學(xué) 電子與光學(xué)工程學(xué)院，南京 210033；3.南京理工大學(xué) 電子工程與光電技術(shù)學(xué)院，南京210094)

一種任意比例LTS技術(shù)在DG-FETD中的應(yīng)用

李 敏1,2，沈微微2，葉曉東3，楊逸婷3

(1.宿遷學(xué)院信息工程學(xué)院，江蘇 宿遷 223800；2.南京郵電大學(xué) 電子與光學(xué)工程學(xué)院，南京 210033；3.南京理工大學(xué) 電子工程與光電技術(shù)學(xué)院，南京210094)

不連續(xù)伽遼金時(shí)域有限元法(discontinuous Galerkin-finite element time domain，DG-FETD)便于處理多尺度電磁問題，但是由于精細(xì)結(jié)構(gòu)或高介電參數(shù)媒質(zhì)的存在，考慮到穩(wěn)定性條件，整體時(shí)間步長的選取會受到最小剖分尺寸的影響，導(dǎo)致計(jì)算效率降低。針對該問題，在傳統(tǒng)局部時(shí)間步長(local time stepping，LTS)技術(shù)的基礎(chǔ)上提出一種基于蛙跳格式的任意比例LTS技術(shù)，該方法在求解多尺度問題時(shí)，減少了迭代所需時(shí)間，提高了不同時(shí)間步長選取的靈活性，同時(shí)隨著未知量增加，其優(yōu)勢更加明顯，結(jié)合數(shù)值算例，驗(yàn)證了該方法的正確性和有效性。

不連續(xù)伽遼金；多尺度；局部時(shí)間步長

0 引 言

近幾十年來，時(shí)域有限元法(finite element method)因其便于模擬復(fù)雜外形及復(fù)雜媒質(zhì)的特點(diǎn)，在計(jì)算電磁學(xué)中得到了越來越多的關(guān)注[1-2]。作為一種時(shí)域方法，時(shí)域有限元法具有可進(jìn)行瞬時(shí)分析和寬帶特性分析的優(yōu)勢，然而，其在每個(gè)時(shí)間步長都需要直接或迭代求解一個(gè)大型的矩陣方程組，尤其針對于多尺度電磁問題，是十分消耗計(jì)算資源和時(shí)間的，因此，時(shí)域不連續(xù)伽遼金法(discontinuous Galerkin time domain method)被提出，以便提高其計(jì)算效率和計(jì)算能力[3-4]。1973年不連續(xù)伽遼金有限元方法(discontinuous Galerkin finite element method, DG-FEM)最早由美國國家重點(diǎn)實(shí)驗(yàn)室的W.H.Reed和T.R.Hill[3]提出，它既可以像有限元法一樣對目標(biāo)模型進(jìn)行非結(jié)構(gòu)網(wǎng)格剖分，又可像有限體積分法一樣單元之間采用數(shù)值通量[4](numerical flux)來傳遞。而不連續(xù)伽遼金時(shí)域有限元方法(discontinuous Galerkin-finite element time domain，DG-FETD)可以獲得塊對角特性的質(zhì)量矩陣，對大型稀疏矩陣可直接求逆運(yùn)算，同時(shí)易于并行計(jì)算，可處理非共形網(wǎng)格和使用分段連續(xù)基函數(shù)的特點(diǎn)。1999年，M.Y.Hussaini討論分析了不連續(xù)伽遼金波傳播問題；2014年，QIANG Ren和QING Huoliu等[5-6]討論分析了基于EB的SETD-FETD結(jié)合的不連續(xù)伽遼金方法。

針對多尺度問題中DG-FETD穩(wěn)定性受最小剖分尺寸影響導(dǎo)致整體時(shí)間步長的選取受到限制的問題，文獻(xiàn)[7]提出基于蛙跳的局部時(shí)間步長(local time-stepping，LTS)迭代格式，即不同的區(qū)域選取不同的時(shí)間步長，此時(shí)，各區(qū)域時(shí)間步長只受各自區(qū)域網(wǎng)格的影響，在時(shí)間步長的選取上有一定的限制，本文給出任意m∶n型LTS迭代格式，其中，m，n為任意正整數(shù)，可根據(jù)實(shí)際問題更加方便地選取不同時(shí)間迭代步長，并結(jié)合實(shí)際算例驗(yàn)證了該格式的正確性和有效性。

1 不連續(xù)伽遼金時(shí)域有限元法

1.1 DG-FETD基本原理

本文中采用四面體的網(wǎng)格離散目標(biāo)模型。不連續(xù)伽遼金時(shí)域有限元法用于求解麥克斯韋方程時(shí)，必須保證磁場和電場在分界面上的切向連續(xù)性，因此，我們引入了連續(xù)性條件，已知的連續(xù)性條件有2種，中心通量(central flux)和迎風(fēng)通量(upwind flux)[8]，本文需要通過中心通量來加強(qiáng)相鄰元胞在公共面上的場量的連續(xù)性，其具體的公式推導(dǎo)過程為

(1)

(2)

(3)—(4)式中：E+，H+表示相鄰體在該面上電場和磁場；n表示元胞面的法向量。由(1)式和(2)式，采用伽遼金測試可得以下矩陣形式的半離散格式

(5)

(6)

(5)—(6)式中：

(7)

1.2 LTS迭代格式

對1.1節(jié)的半離散格式，在時(shí)間上采用蛙跳迭代格式為

(8)—(9)式中：Δt為選取的迭代時(shí)間步長；hn+1/2，hn-1/2分別表示n+1/2，n-1/2時(shí)刻的磁場強(qiáng)度；en,en-1表示n及n-1時(shí)刻的電場強(qiáng)度。

本文使用的是蛙跳的時(shí)間差分格式，與隱式差分格式相比，它的效果更好。因?yàn)槠鋾r(shí)間步長的選取可以通過CFL(Courant-Friedrichs-Lewy)條件計(jì)算出來，為了保證穩(wěn)定性，必須使用最小剖分網(wǎng)格計(jì)算。針對多尺度問題中，DG-FETD穩(wěn)定性受最小剖分尺寸影響導(dǎo)致整體時(shí)間步長的選取受到限制的問題，文獻(xiàn)[9]提出基于蛙跳的局部時(shí)間步進(jìn)(local time-stepping，LTS)迭代格式，即不同的區(qū)域選取不同的時(shí)間步長，此時(shí)，各區(qū)域時(shí)間步長只受各自區(qū)域的網(wǎng)格影響，提高了計(jì)算效率。它包含以下4個(gè)步驟。

步驟1計(jì)算最小時(shí)間步長，根據(jù)CFL條件，利用已知最小剖分網(wǎng)格的尺寸，計(jì)算出能滿足整個(gè)區(qū)域穩(wěn)定性條件所需的時(shí)間步長。

步驟3將所有網(wǎng)格分為2部分，根據(jù)步驟2的分組找出最佳的劃分方法，達(dá)到最優(yōu)的加速。在上述例子中，最優(yōu)的劃分方法是在第2組和第3組中間劃分。

步驟4按照DG-FETD的方法對2個(gè)區(qū)域分別求解。

在LTS技術(shù)中最重要的一個(gè)部分就是對控制方程差分近似方法的改變。當(dāng)整個(gè)求解區(qū)域選取Δt1，Δt22種大小的時(shí)間步長時(shí)，其中Δt1∶Δt2=m∶n，各自區(qū)域分別采用蛙跳迭代格式，此時(shí)，需要考慮2個(gè)區(qū)域分界面上的場值傳遞，下面以Δt=3Δt1=2Δt2為例介紹其迭代過程。

(10)

大部分實(shí)際應(yīng)用問題都包含著更為復(fù)雜的結(jié)構(gòu)，這時(shí)，固定比例的LTS技術(shù)效果就顯得不夠突出。我們通過改進(jìn)，實(shí)現(xiàn)了任意比例的LTS技術(shù)，并將其應(yīng)用在了DG-FETD中，取得了顯著的效果。通過上面的描述，我們知道在LTS技術(shù)中，最重要的是時(shí)間上的近似，將需要用到的信息用已經(jīng)存在的時(shí)間步內(nèi)的信息近似，在每個(gè)時(shí)間上的銜接中加入嵌套循環(huán)，設(shè)1區(qū)域的未知量為e1和h1，2區(qū)域的未知量為e2和h2，若已知每個(gè)未知量n時(shí)刻的值，需求解n+1時(shí)刻的值，若Δt1∶Δt2=m2∶m1，即Δt=m1Δt1=m2Δt2，其中，Δt為實(shí)際的時(shí)間步長。1區(qū)域每個(gè)時(shí)間步應(yīng)迭代m次，因電磁的未知量和磁場未知量分開計(jì)算，則1區(qū)域應(yīng)計(jì)算2m1次，2區(qū)域應(yīng)計(jì)算2m2次，具體實(shí)現(xiàn)步驟如圖1所示。

通過引入LTS技術(shù)，在求解多尺度電磁問題時(shí)，DG-FETD能減少迭代所需的計(jì)算時(shí)間，更大地發(fā)揮優(yōu)勢，隨著未知量的增加，這一優(yōu)勢會變得更加顯著。

2 數(shù)值算例

為了驗(yàn)證本文提出的方法的正確性和有效性，算例1對半填充介質(zhì)矩形金屬諧振腔進(jìn)行了分析，如圖2所示。諧振腔尺寸是1 m×0.1 m×1 m，介質(zhì)部分的相對介電常數(shù)是2，激勵(lì)源是中心頻率為300 MHz，帶寬為400 MHz的調(diào)制高斯脈沖。計(jì)算區(qū)域中時(shí)間步長分別選取為Δt1=5.57 ps, Δt2=8.35 ps。

圖3給出Δt=3Δt1=2Δt2與選取傳統(tǒng)全局時(shí)間步長Δt1時(shí)的時(shí)域電場，結(jié)果吻合較好，表明所提出方法的正確性。

圖1 LTS時(shí)間迭代流程圖Fig.1 Flowchart for the time of iteration

圖2 半填充介質(zhì)諧振腔Fig.2 Half filled dielectric resonator

圖3 時(shí)域電場比較Fig.3 Comparison of electric field in time domain

表1給出2種方法的計(jì)算時(shí)間比較，表明了所提出方法的有效性。

表1 計(jì)算時(shí)間比較Tab.1 Comparison of computation time

算例2是里面含有一個(gè)介質(zhì)球的金屬球，金屬球半徑為0.5 m，介質(zhì)球半徑為0.05 m，相對介電常數(shù)εr=9.0，在z=0 m處加入調(diào)制高斯脈沖，其中心頻率為10 GHz，脈寬相關(guān)量為9，觀察點(diǎn)為(0,0,0)。金屬剖分尺寸為0.06 m，介質(zhì)部分剖分尺寸為0.015 m，離散后得到22 989個(gè)四面體單元，未知量個(gè)數(shù)為267 260，分別用傳統(tǒng)的DG-FETD和加入了LTS技術(shù)的DG-FETD計(jì)算該模型，其中，后者采用Δt1∶Δt2=1∶2.5的時(shí)間步長。

2種方法的時(shí)域電場波形如圖4所示，吻合較好，可以證明此方法的正確性。本算例在四核8 GByte內(nèi)存,主頻為2 GHz，64位的臺式機(jī)上運(yùn)行，2種方法的計(jì)算資源對比如表2所示。在未知量個(gè)數(shù)相等的情況下，第2種方法在迭代時(shí)間上明顯優(yōu)于第1種方法，這也進(jìn)一步證明了此方法的有效性。

圖4 LTS迭代與蛙跳迭代時(shí)域波形對比Fig.4 Comparison of time domain waveform between LTS and Leap-frog method表2 LTS迭代的DG-FETD與蛙跳迭代 DG-FETD計(jì)算資源對比Tab.2 Comparison of computation resource between iterative LTS method and iterative Leap-frog method for DG-FETD

迭代方法解析解/MHz諧振頻率/MHz迭代步數(shù)迭代時(shí)間/s蛙跳迭代263.017262.460120001144LTS迭代(1∶2)263.017262.3176000645LTS迭代(2∶5)263.017262.2314800515

3 結(jié) 論

基于蛙跳的LTS技術(shù)的重點(diǎn)在于不同時(shí)間迭代步長區(qū)域交界面上的電磁場迭代處理方式，以便實(shí)現(xiàn)各區(qū)域時(shí)間步長只受各自區(qū)域的網(wǎng)格影響，并且保持各區(qū)域的顯式迭代特性。本文在傳統(tǒng)LTS技術(shù)的基礎(chǔ)上通過引入任意比例LTS技術(shù)，提高了不同時(shí)間步長選取的靈活性，在求解多尺度電磁問題時(shí)，DG-FETD能減少迭代所需的計(jì)算時(shí)間，更大地發(fā)揮優(yōu)勢，隨著未知量的增加，這一優(yōu)勢會變得更加顯著，數(shù)值算例驗(yàn)證了該方法的正確性和有效性。

[1] JIN J M.The Finite Element Method in Electromagnetics[M].2nd ed. New York: Wiley, 2002.

[2] JIN J M，RILEY D J. Finite Element Analysis of Antennas and Arrays[J].Hoboken, NJ: Wiley, 2008: 56(8): 222-240

[3] HESTHAVENS J S, WARBURTON T.High-order Nodal methods on unstructured grids. I. Time-domain solution of Maxwell’s equations[J]. Comput Phys, 2002, 181(1): 186-221.

[4] COHEN G,FERRIE ′RES X,PERNET S.A spatial high order hexahedral Discontinuous Galerkin method to solve Maxwell’s equations in time domain[J].Comput Phys, 2006, 217(2): 340-363.

[5] REED W H, HILL T R.Triangular mesh methods for the neutron transport equation[M].Los Alamos: Los Alamos Scientific Laboratory Report LA-UR, 1973.

[6] CHEN J, LIU Q H.Discontinuous Galerkin Time-Domain Methods for Multiscale Electromagnetic Simulations: A Review[J]. Proceedings of the IEEE, 2013, 101(2):242-254.

[7] MONTSENY E, PERNET S, FERRIERES X,et al. Dissipative terms and local time-stepping improvements in a spatial high order discontinuous Galerkin scheme for the time-domain method for the time-domain Maxwell’s equations[J].Comput Phys，2008,227(14)：6795-6820.

[8] LI P, SHI Y, JIANG L J, et al. A Hybrid Time-Domain Discontinuous Galerkin-Boundary Integral Method for Electromagnetic Scattering Analysis[J].IEEE Transactions on Antennas & Propagation, 2014, 62(5):2841-2846.

[9] MONTSENY E, PERNET S, FERRIERES X, et al. Dissipative terms and local time-stepping improvements in a spatial high order discontinuous Galerkin scheme for the time-domain method for the time-domain Maxwell’s equations[J].Comput Phys，2008，227(14): 6795-6820.

[10] SCHOMANN S, G?DEL N, WARBURTON T, et al.Local Timestepping Techniques Using Taylor Expansion for Modeling Electromagnetic Wave Propagation With Discontinuous Galerkin-FEM[J]. IEEE Transactions on Magnetics,2010, 46(8):3504-3507.

The Suqian Social Development Project Foundation(S201410)

ApplicationoflocaltimesteppingstrategywitharbitraryratioforthediscontinuousGalerkinFETDmethod

LI Min1,2, SHEN Weiwei2, YE Xiaodong3, YANG Yiting3

(1.School of Information Engineering, Suqian College, Suqian 223800, P.R.China； 2.College of Electronic and Optical Engineering, Nanjing University of Posts and Telecommunications, Nanjing 210033, P.R.China； 3.School of Electronic and Optical Engineering,Nanjing University of Science and Technology, Nanjing 210094, P.R.China)

The discontinuous Galerkin-finite element time domain(DG-FETD) method has the ability to deal with the multiscale problems. The size of time step is limited by the smallest spatial discretization of the simulation domain according to the stability condition when fine structures or high permittivity media occur. To handle this kind of problem, the local time-stepping technique with arbitrary ratio is proposed in the paper, which makes the selection of time step size more flexible and reduces the time of iteration. Numerical results demonstrate the correctness and efficiency of the proposed method.

discontinuous Galerkin; multi-scale; local time stepping

10.3979/j.issn.1673-825X.2017.06.005

2016-10-31

2017-06-01

沈微微 1105183311@qq.com

宿遷市社會發(fā)展(指令性)項(xiàng)目(S201410)

TN925

A

1673-825X(2017)06-0739-05

李 敏(1986 -)，女，山東濟(jì)寧人，講師，碩士，主要研究方向?yàn)闀r(shí)域電磁計(jì)算及其快速算法,超寬帶輻射、散射理論與技術(shù)研究。E-mail:18800608557@163.com。

沈微微(1985 -)，女，江蘇徐州人，碩士，研究方向?yàn)閳D形圖像技術(shù)。E-mail: 1105183311@qq.com。

葉曉東(1967 -)，男，江蘇南京人，副教授，主要研究方向?yàn)殡姶艌雠c微波技術(shù)，散射理論與技術(shù)研究。E-mail:yexiaodong@njust.edu.cn。

(編輯：王敏琦)