本發(fā)明屬于計(jì)算電磁學(xué)技術(shù)領(lǐng)域，具體涉及一種一維高精度迭代的磁化等離子體中的實(shí)現(xiàn)方法。

背景技術(shù)：

時(shí)域有限差分(finite-differencetime-domain，fdtd)方法已經(jīng)成功地應(yīng)用于求解磁化等離子體中的電磁問(wèn)題。由于它是一種顯式的數(shù)值計(jì)算的方法，其時(shí)間步長(zhǎng)受柯西穩(wěn)定性條件的限制，不能取的太大，方法時(shí)效性不好，計(jì)算速度較慢，尤其在磁化等離子體中的電磁計(jì)算中，時(shí)間步長(zhǎng)不僅僅受柯西穩(wěn)定性條件的限制，而且還與等離子體參數(shù)和差分方案有關(guān)，如此更進(jìn)一步的限制了時(shí)間步長(zhǎng)的選取，其計(jì)算速度和計(jì)算效率較空氣中的計(jì)算更差。為了消除柯西穩(wěn)定性條件的限制，人們提出了無(wú)條件穩(wěn)定時(shí)域有限差分方法，比如：交替方向隱式(alternating-direction-implicit，adi)的時(shí)域有限差分(adi-fdtd)方法和基于加權(quán)拉蓋爾多項(xiàng)式的時(shí)域有限差分(finite-differencetime-domainwithweighted-laguerre-polynomials，wlp-fdtd)方法。在這些方法中，adi-fdtd方法在使用較大的時(shí)間步長(zhǎng)時(shí)會(huì)產(chǎn)生很大的色散誤差，而wlp-fdtd方法既能消除柯西穩(wěn)定性條件的限制，又能解決adi-fdtd方法在使用較大的時(shí)間步長(zhǎng)時(shí)會(huì)產(chǎn)生很大的色散誤差這個(gè)難題，因此wlp-fdtd方法可以高效的求解磁化等離子體中的電磁問(wèn)題。然而，這種wlp-fdtd方法在求解電磁場(chǎng)過(guò)程中，會(huì)產(chǎn)生一個(gè)大型的稀疏矩陣方程，直接求解此方程會(huì)使得計(jì)算較復(fù)雜，內(nèi)存消耗較大，于是提出了一種因式分裂的wlp-fdtd方法，該方法在計(jì)算速度和計(jì)算效率上得到了較大的提高，但是由于該方法是通過(guò)添加微擾項(xiàng)，進(jìn)行因式分裂得來(lái)的，計(jì)算中會(huì)產(chǎn)生分裂誤差，為了減小分裂誤差、提高計(jì)算精度，同時(shí)保證較快的計(jì)算速度，提出了一種高精度迭代的wlp-fdtd方法。

另外由于計(jì)算機(jī)容量的限制，電磁場(chǎng)的計(jì)算只能在有限區(qū)域進(jìn)行。為了能模擬開(kāi)域電磁波傳播過(guò)程，必須在計(jì)算區(qū)域的截?cái)噙吔缣幗o出吸收邊界條件。有人提出了完全匹配層(perfectlymatchedlayer，pml)吸收邊界，后來(lái)pml被廣泛應(yīng)用于計(jì)算區(qū)域的截?cái)?，而且被證明是非常有效的，但是研究發(fā)現(xiàn)這種傳統(tǒng)pml對(duì)低頻以及凋落波的吸收效果并不理想；使用帶有復(fù)頻率偏移(complexfrequencyshift，cfs)因子的擴(kuò)展坐標(biāo)的pml吸收邊界可以有效地改善傳統(tǒng)pml對(duì)低頻，凋落波的吸收效果。

技術(shù)實(shí)現(xiàn)要素：

本發(fā)明的目的是提供一種一維高精度迭代的磁化等離子體中的實(shí)現(xiàn)方法，計(jì)算精度高、計(jì)算速度快，且對(duì)于低頻和凋落波具有很好的吸收效果。

本發(fā)明所采用的技術(shù)方案是，一種一維高精度迭代的磁化等離子體中的實(shí)現(xiàn)方法，具體按照以下步驟實(shí)施：

步驟1：輸入模型文件；

步驟2：初始化參數(shù)和設(shè)置參數(shù)；

步驟3：添加場(chǎng)源到x方向上的電場(chǎng)分量系數(shù)中，設(shè)置磁場(chǎng)分量系數(shù)和分別記為初始場(chǎng)值和其中，添加的場(chǎng)源的表達(dá)式為：

jx(t)＝exp(-(t-t0)²/τ²)

其中，t表示時(shí)間，單位為秒，t0,τ為場(chǎng)源參數(shù)；

步驟4：更新計(jì)算整個(gè)計(jì)算區(qū)域的x方向上電場(chǎng)分量系數(shù)

步驟5：更新計(jì)算整個(gè)計(jì)算區(qū)域的y方向上電場(chǎng)分量系數(shù)

步驟6：更新計(jì)算整個(gè)計(jì)算區(qū)域的磁場(chǎng)分量系數(shù)

步驟7：更新計(jì)算整個(gè)計(jì)算區(qū)域的磁場(chǎng)分量系數(shù)

步驟8：將k+1賦值給k，并判斷迭代次數(shù)k是否達(dá)到預(yù)設(shè)值，若未達(dá)到預(yù)設(shè)值，則返回步驟4，若達(dá)到預(yù)設(shè)值，則執(zhí)行步驟9；

步驟9：更新計(jì)算整個(gè)計(jì)算區(qū)域的極化電流密度系數(shù)

步驟10：更新計(jì)算整個(gè)計(jì)算區(qū)域的電磁場(chǎng)分量系數(shù)的輔助變量；

步驟11：更新計(jì)算觀測(cè)點(diǎn)處的電磁場(chǎng)分量；

步驟12：將q+1賦值給q，并判斷拉蓋爾多項(xiàng)式的階數(shù)q是否達(dá)到預(yù)設(shè)值，若未達(dá)到預(yù)設(shè)值，則返回步驟3，若達(dá)到預(yù)設(shè)值，則結(jié)束。

本發(fā)明的特點(diǎn)還在于：

步驟1具體為：

計(jì)算區(qū)域大小nz，nz為z軸方向上網(wǎng)格數(shù)；空間步長(zhǎng)δz；時(shí)間步長(zhǎng)δt；真空中的電導(dǎo)率σ，磁導(dǎo)率μ0，介電常數(shù)ε0；等離子體碰撞頻率υ，等離子體頻率ωp，電子回旋頻率ωb；等離子體在計(jì)算區(qū)域中的位置；吸收邊界層數(shù)npml與相關(guān)參數(shù)κzmax，αzmax，σzmax；κzmax取整數(shù)，κzmax取值范圍為[1，60]；αzmax取值范圍為[01)；σzmax/σopt取值范圍為(0，12]，σopt＝(m+1)/150πδz，m取值范圍為[1，20]，δz取值范圍為λ為源的波長(zhǎng)；仿真計(jì)算時(shí)長(zhǎng)tf；迭代次數(shù)k，k≥0且為整數(shù)；加權(quán)拉蓋爾多項(xiàng)式的階數(shù)q，q≥0且為整數(shù)；時(shí)間尺度因子s，s取值范圍為[10⁹，10¹³]；觀測(cè)點(diǎn)；場(chǎng)源參數(shù)。

步驟2初始化具體為：

將整個(gè)計(jì)算區(qū)域的電磁場(chǎng)分量系數(shù)整個(gè)計(jì)算區(qū)域的極化電流密度系數(shù)整個(gè)計(jì)算區(qū)域的電磁場(chǎng)分量系數(shù)的和整個(gè)計(jì)算區(qū)域的極化電流密度系數(shù)的和整個(gè)計(jì)算區(qū)域的輔助變量和其中表示輔助變量表示輔助變量拉蓋爾多項(xiàng)式其中全部初始化為零；磁化等離子體參數(shù)(p1,p2,p3)初始化為p1＝0，p2＝0，pml系數(shù)(c1z，c2z，c3)初始化為c1z＝2/ε0s，c2z＝c1z，c3＝ε0/μ0，其中，ε0是空氣中的介電常數(shù)，μ0是空氣中的磁導(dǎo)率；

設(shè)置的參數(shù)具體為：

設(shè)置帶有cfs因子的sc-pml吸收邊界的參數(shù)σz,κz,αz；具體為：

σz＝σzmax|z-z0|^m/d^m；

κz＝1+(κzmax-1)|z-z0|^m/d^m；

αz＝αzmax；

式中z0為pml層與非pml截面位置，d是pml吸收邊界的厚度，κzmax取整數(shù)，κzmax取值范圍為[1，60]；αzmax取值范圍為[0，1)；σzmax根據(jù)σopt來(lái)設(shè)置，σzmax/σopt取值范圍為(0，12]；σopt＝(m+1)/150πδz，m取值范圍為[1，20]，δz取值范圍為λ為源的波長(zhǎng)；

設(shè)置pml系數(shù)c1z，c2z；具體為：

c1z＝1/(κzαz+σz+0.5κzε0s)，c2z＝(2αz/ε0s+1)c1z；

設(shè)置等離子體參數(shù)為：

步驟4具體為：

步驟4.1：電場(chǎng)分量系數(shù)在計(jì)算區(qū)域的方程為：

其中，k表示迭代次數(shù)，n表示z軸上第n個(gè)計(jì)算網(wǎng)格處的位置，n+1/2表示z軸上第n個(gè)半網(wǎng)格處的位置，c2z|n表示系數(shù)c2z在z軸方向上第n個(gè)網(wǎng)格處的值，c2z|n-1/2表示系數(shù)c2z在z軸方向上第n-1個(gè)半網(wǎng)格處的值，δz|n表示z軸方向上第n個(gè)網(wǎng)格的大小，第表示第k+1次迭代第n-1個(gè)網(wǎng)格處的電場(chǎng)分量系數(shù)的值，

步驟4.2：使用追趕法求解步驟4.1的方程，得到整個(gè)計(jì)算區(qū)域的電場(chǎng)分量系數(shù)

步驟5具體為：

步驟5.1：電場(chǎng)分量系數(shù)在計(jì)算區(qū)域的方程為：

其中，表示第k+1次迭代第n-1個(gè)網(wǎng)格處的電場(chǎng)分量系數(shù)的值，

步驟5.2：使用追趕法求解步驟5.1的方程，得到整個(gè)計(jì)算區(qū)域的y方向上電場(chǎng)分量系數(shù)

步驟6的更新方程具體為：

其中

步驟7的更新方程具體為：

其中

步驟9更新公式具體為：

步驟10更新公式具體為：

其中fζz代表exz,eyz,hxz,hyz，因此表示輔助變量表示電磁場(chǎng)分量系數(shù)的和對(duì)z的求導(dǎo)，表示輔助變量的和，表示輔助變量

步驟11的更新公式具體為：

上式中u表示電磁場(chǎng)分量ex,ey,hx,hy，r表示電磁場(chǎng)分量的位置，u^q表示q階電磁場(chǎng)分量系數(shù)，是q階加權(quán)拉蓋爾多項(xiàng)式，是帶有時(shí)間尺度因子s>0的擴(kuò)展時(shí)間，是q階拉蓋爾多項(xiàng)式。

本發(fā)明的有益效果是：

1).本發(fā)明一種一維高精度迭代的磁化等離子體中的實(shí)現(xiàn)方法，在直角坐標(biāo)系下，通過(guò)用加權(quán)拉蓋爾多項(xiàng)式表示電磁場(chǎng)分量，來(lái)解時(shí)域麥克斯韋方程，使得在更新計(jì)算整個(gè)計(jì)算區(qū)域的電磁場(chǎng)分量系數(shù)時(shí)不涉及到時(shí)間步長(zhǎng)，只是在最后計(jì)算觀測(cè)點(diǎn)處的電磁場(chǎng)分量時(shí)用到時(shí)間步長(zhǎng)，因此計(jì)算過(guò)程中時(shí)間步長(zhǎng)可以取得比柯西穩(wěn)定性條件限制的時(shí)間步長(zhǎng)更大；

2).本發(fā)明一種一維高精度迭代的磁化等離子體中的實(shí)現(xiàn)方法，在求解求解磁化等離子體中的電磁場(chǎng)分量系數(shù)時(shí)，使用迭代的因式分裂的方案將大型稀疏矩陣方程分裂成兩個(gè)迭代的三對(duì)角矩陣方程，使得它在計(jì)算時(shí)比wlp-fdtd方法更簡(jiǎn)單、計(jì)算速度更快、內(nèi)存消耗更少而且比無(wú)迭代的因式分裂的wlp-fdtd方法精度更高；

3).本發(fā)明一種一維高精度迭代的磁化等離子體中的實(shí)現(xiàn)方法，在設(shè)置pml系數(shù)時(shí)，由于采用了cfs因子，并且通過(guò)調(diào)整cfs因子中的參數(shù)，可以使得該吸收邊界對(duì)低頻與凋落波的吸收更加有效；

4).本發(fā)明一種一維高精度迭代的磁化等離子體中的實(shí)現(xiàn)方法，由于采用了復(fù)擴(kuò)展坐標(biāo)系，使得pml在實(shí)現(xiàn)時(shí)避免了場(chǎng)的分裂且與媒質(zhì)無(wú)關(guān)。

附圖說(shuō)明

圖1是本發(fā)明實(shí)現(xiàn)方法的流程示意圖；

圖2是使用本發(fā)明的方法與解析解、wlp-fdtd方法和無(wú)迭代的因式分裂的wlp-fdtd方法計(jì)算的左旋極化波的反射系數(shù)振幅圖；

圖3是使用本發(fā)明的方法與解析解、wlp-fdtd方法和無(wú)迭代的因式分裂的wlp-fdtd方法計(jì)算的左旋極化波的反射系數(shù)相位圖。

具體實(shí)施方式

下面結(jié)合附圖和具體實(shí)施方式對(duì)本發(fā)明進(jìn)行詳細(xì)說(shuō)明。

本發(fā)明的一種一維高精度迭代的磁化等離子體中的實(shí)現(xiàn)方法，原理為：首先導(dǎo)出一維磁化等離子體中，電磁場(chǎng)所滿足的復(fù)擴(kuò)展坐標(biāo)系下的麥克斯韋方程，然后使用一維磁化等離子體中迭代的因式分裂的wlp-fdtd方法推導(dǎo)出整個(gè)計(jì)算區(qū)域的電磁場(chǎng)分量系數(shù)更新方程，最后采用公式(21)求解觀測(cè)點(diǎn)處的電磁場(chǎng)分量。

一維磁化等離子體中，電磁場(chǎng)所滿足的復(fù)擴(kuò)展坐標(biāo)系下的麥克斯韋方程推導(dǎo)過(guò)程如下：

在各向異性色散介質(zhì)碰撞磁化等離子體中，麥克斯韋方程組和相關(guān)的聯(lián)立方程為

式中，h是磁場(chǎng)強(qiáng)度；e是電場(chǎng)強(qiáng)度；j是極化電流密度；ε0、μ0分別為真空中的介電常數(shù)和磁導(dǎo)率；υ是等離子體碰撞頻率；是等離子體角頻率的平方；ωb＝eb0/m是電子回旋頻率；b0是外磁場(chǎng)；e、m分別是電子的電量和質(zhì)量。

設(shè)外磁場(chǎng)的方向?yàn)?z軸，方程(3)可寫(xiě)為

應(yīng)用擴(kuò)展坐標(biāo)的cfs-pml，僅考慮二維tez的情況，上述麥克斯韋方程組和相關(guān)的聯(lián)立方程可化為：

式中是帶有cfs因子的坐標(biāo)伸縮因子，取為：

下面我們引入幾個(gè)輔助變量如下：

將(12)式分別代入(13)-(16)，然后做如下變化于是可得下面幾式：

我們可以知道

式中u代表ex、ey、hx、hy，是加權(quán)拉蓋爾多項(xiàng)式，q≥0；t≥0是q階拉蓋爾多項(xiàng)式。

將(21)、(22)代入(17)-(20)，然后應(yīng)用加勒金的測(cè)試過(guò)程，可得

式中

s>0是時(shí)間尺度因子，q是加權(quán)拉蓋爾多項(xiàng)式的階數(shù)。

將(21)、(22)代入(6)-(11)，再應(yīng)用加勒金的測(cè)試過(guò)程得到：

式中，

將(36)代入(37)得到

將(37)代入(36)得到

將(39)代入(32)得到

將(38)代入(33)得到

將方程(40)、(41)、(34)和(35)放置一起

整理后得到

其中

將(43)式寫(xiě)成矩陣的形式

(i-a-b)x^q＝v^q-1(45)

式中

b＝p1p2/(1+p2),a1＝c2zdz/(1+p2),a2＝c2zc3dz

添加微擾項(xiàng)到(45)式得到

為了強(qiáng)調(diào)迭代的特點(diǎn)，將(47)式重寫(xiě)為

引進(jìn)中間變量于是(48)式可以分裂為

將上式展開(kāi)得到

上兩式解得

將展開(kāi)得到

將(52)、(53)和(54)聯(lián)合解得

將b＝p1p2/(1+p2),a1＝c2zdz/(1+p2),a2＝c2zc3dz代入上式然后進(jìn)行中心差分得到

上面三式中，k表示第k次迭代，n表示第n個(gè)網(wǎng)格，n+1/2表示第n個(gè)半網(wǎng)格；在整個(gè)計(jì)算區(qū)域上，(56)式和(57)式可以寫(xiě)成三對(duì)角矩陣差分方程，與wlp-fdtd方法相比，這種迭代的因式分解的wlp-fdtd方法將大型稀疏矩陣方程的求解轉(zhuǎn)變成兩個(gè)三對(duì)角矩陣方程的求解，于是可以使用追趕法，非常簡(jiǎn)單的解得整個(gè)計(jì)算區(qū)域電磁場(chǎng)分量系數(shù)，最后通過(guò)(21)式解得觀測(cè)點(diǎn)的電磁場(chǎng)分量，而且這種方法添加了迭代的方案，使得其計(jì)算精度比因式分裂的wlp-fdtd要高。

本發(fā)明一種一維高精度迭代的磁化等離子體中的實(shí)現(xiàn)方法，流程如圖1所示，具體按照以下步驟實(shí)施：

步驟1：輸入模型文件，具體為：

計(jì)算區(qū)域大小nz，nz為z軸方向上網(wǎng)格數(shù)；空間步長(zhǎng)δz；時(shí)間步長(zhǎng)δt；真空中的電導(dǎo)率σ，磁導(dǎo)率μ0，介電常數(shù)ε0；等離子體碰撞頻率υ，等離子體頻率ωp，電子回旋頻率ωb；等離子體在計(jì)算區(qū)域中的位置；吸收邊界層數(shù)npml與相關(guān)參數(shù)κzmax，αzmax，σzmax；κzmax取整數(shù)，κzmax取值范圍為[1，60]；αzmax取值范圍為[01)；σzmax/σopt取值范圍為(0，12]，σopt＝(m+1)/150πδz，m取值范圍為[1，20]，δz取值范圍為λ為源的波長(zhǎng)；仿真計(jì)算時(shí)長(zhǎng)tf；迭代次數(shù)k，k≥0且為整數(shù)；加權(quán)拉蓋爾多項(xiàng)式的階數(shù)q，q≥0且為整數(shù)；時(shí)間尺度因子s，s取值范圍為[10⁹，10¹³]；觀測(cè)點(diǎn)；場(chǎng)源參數(shù)。

步驟2：初始化參數(shù)和設(shè)置參數(shù)，具體為：

初始化的參數(shù)包括：

將整個(gè)計(jì)算區(qū)域的電磁場(chǎng)分量系數(shù)整個(gè)計(jì)算區(qū)域的極化電流密度系數(shù)整個(gè)計(jì)算區(qū)域的電磁場(chǎng)分量系數(shù)的和整個(gè)計(jì)算區(qū)域的極化電流密度系數(shù)的和整個(gè)計(jì)算區(qū)域的輔助變量和其中表示輔助變量表示輔助變量拉蓋爾多項(xiàng)式其中全部初始化為零；磁化等離子體參數(shù)(p1,p2,p3)初始化為p1＝0，p2＝0，pml系數(shù)(c1z，c2z，c3)初始化為c1z＝2/ε0s，c2z＝c1z，c3＝ε0/μ0，其中，ε0是空氣中的介電常數(shù)，μ0是空氣中的磁導(dǎo)率。

設(shè)置的參數(shù)具體為：

設(shè)置帶有cfs因子的sc-pml吸收邊界的參數(shù)σz,κz,αz；具體為：

σz＝σzmax|z-z0|^m/d^m；

κz＝1+(κzmax-1)|z-z0|^m/d^m；

αz＝αzmax；

式中z0為pml層與非pml截面位置，d是pml吸收邊界的厚度，κzmax取整數(shù)，κzmax取值范圍為[1，60]；αzmax取值范圍為[0，1)；σzmax根據(jù)σopt來(lái)設(shè)置，σzmax/σopt取值范圍為(0，12]；σopt＝(m+1)/150πδz，m取值范圍為[1，20]，δz取值范圍為λ為源的波長(zhǎng)；設(shè)置pml系數(shù)c1z，c2z；具體為：

c1z＝1/(κzαz+σz+0.5κzε0s)，c2z＝(2αz/ε0s+1)c1z；

設(shè)置等離子體參數(shù)為：

步驟3：添加場(chǎng)源到x方向上的電場(chǎng)分量系數(shù)中，設(shè)置磁場(chǎng)分量系數(shù)和分別記為初始場(chǎng)值和

其中，所添加的場(chǎng)源的表達(dá)式為：

jx(t)＝exp(-(t-t0)²/τ²)

其中，t表示時(shí)間，單位為秒；t0,τ為場(chǎng)源參數(shù)。

步驟4：更新計(jì)算整個(gè)計(jì)算區(qū)域的x方向上電場(chǎng)分量系數(shù)具體為：

步驟4.1：電場(chǎng)分量系數(shù)在計(jì)算區(qū)域的方程為：

其中，k表示迭代次數(shù)，n表示z軸上第n個(gè)計(jì)算網(wǎng)格處的位置，n+1/2表示z軸上第n個(gè)半網(wǎng)格處的位置，c2z|n表示系數(shù)c2z在z軸方向上第n個(gè)網(wǎng)格處的值，c2z|n-1/2表示系數(shù)c2z在z軸方向上第n-1個(gè)半網(wǎng)格處的值，δz|n表示z軸方向上第n個(gè)網(wǎng)格的大小，第表示第k+1次迭代第n-1個(gè)網(wǎng)格處的電場(chǎng)分量系數(shù)的值，

步驟4.2：使用追趕法求解步驟4.1的方程，得到整個(gè)計(jì)算區(qū)域的電場(chǎng)分量系數(shù)

步驟5：更新計(jì)算整個(gè)計(jì)算區(qū)域的y方向上電場(chǎng)分量系數(shù)具體為：

步驟5.1：電場(chǎng)分量系數(shù)在計(jì)算區(qū)域的方程為：

其中，表示第k+1次迭代第n-1個(gè)網(wǎng)格處的電場(chǎng)分量系數(shù)的值，

步驟5.2：使用追趕法求解步驟5.1的方程，得到整個(gè)計(jì)算區(qū)域的y方向上電場(chǎng)分量系數(shù)

步驟6：更新計(jì)算整個(gè)計(jì)算區(qū)域的磁場(chǎng)分量系數(shù)更新方程具體為：

其中

步驟7：更新計(jì)算整個(gè)計(jì)算區(qū)域的磁場(chǎng)分量系數(shù)更新方程具體為：

其中

步驟8：將k+1賦值給k，并判斷迭代次數(shù)k是否達(dá)到預(yù)設(shè)值，若未達(dá)到預(yù)設(shè)值，則返回步驟4，若達(dá)到預(yù)設(shè)值，則執(zhí)行步驟9；

步驟9：更新計(jì)算整個(gè)計(jì)算區(qū)域的極化電流密度系數(shù)更新公式具體為：

步驟10：更新計(jì)算整個(gè)計(jì)算區(qū)域的電磁場(chǎng)分量系數(shù)的輔助變量，更新公式具體為：

其中fζz代表exz,eyz,hxz,hyz，因此表示輔助變量表示電磁場(chǎng)分量系數(shù)的和對(duì)z的求導(dǎo)，表示輔助變量的和，表示輔助變量

步驟11：更新計(jì)算觀測(cè)點(diǎn)處的電磁場(chǎng)分量，更新公式具體為：

上式中u表示電磁場(chǎng)分量ex,ey,hx,hy，r表示電磁場(chǎng)分量的位置，u^q表示q階電磁場(chǎng)分量系數(shù)，是q階加權(quán)拉蓋爾多項(xiàng)式，是帶有時(shí)間尺度因子s>0的擴(kuò)展時(shí)間，是q階拉蓋爾多項(xiàng)式。

步驟12：將q+1賦值給q，并判斷拉蓋爾多項(xiàng)式的階數(shù)q是否達(dá)到預(yù)設(shè)值，若未達(dá)到預(yù)設(shè)值，則返回步驟3，若達(dá)到預(yù)設(shè)值，則結(jié)束。

下面通過(guò)實(shí)驗(yàn)對(duì)本發(fā)明的效果進(jìn)行說(shuō)明：

為了檢驗(yàn)本發(fā)明方法的正確性、高效性以及高精度，我們計(jì)算了9mm厚碰撞磁化等離子體對(duì)垂直入射電磁波的反射和透射系數(shù)。入射電磁波為高斯脈沖，加到ex上，其表達(dá)式為exp(-(t-t0)²/τ²)，式中t0＝20ps，τ＝5ps。計(jì)算區(qū)域?yàn)?50個(gè)網(wǎng)格，每個(gè)網(wǎng)格大小為75微米，磁化等離子體占據(jù)第201到320的網(wǎng)格，其參數(shù)為ωp＝(2π)·50×10⁹rad/s，ωb＝3.0×10¹¹rad/s，υ＝2.0×10¹⁰hz，其余為真空。完全匹配層吸收邊界放置在計(jì)算區(qū)域的兩端，均為10個(gè)網(wǎng)格。時(shí)間尺度因子為s＝1.885×10¹²，時(shí)間步長(zhǎng)為0.25ps，仿真時(shí)間為tf＝1ns，階數(shù)為200，迭代次數(shù)k為2。pml吸收邊界參數(shù)κzmax＝1，σzmax＝σopt，αzmax＝0。采用本發(fā)明方法、wlp-fdfd方法和因式分裂的wlp-fdtd方法以及解析解計(jì)算反射系數(shù)振幅和反射系數(shù)相位，計(jì)算結(jié)果參見(jiàn)圖2和圖3。從圖中可見(jiàn)，本發(fā)明方法和wlp-fdtd方法、解析解計(jì)算結(jié)果一致，驗(yàn)證了本發(fā)明方法的正確性，且本發(fā)明方法的計(jì)算精度相比于因式分裂的wlp-fdtd方法的計(jì)算精度要高。在運(yùn)行計(jì)算時(shí)間上，wlp-fdtd方法需要0.34秒，因式分裂的wlp-fdtd方法需要0.15秒，本發(fā)明專利中的計(jì)算方法需要0.23秒，由此可知本發(fā)明專利中的計(jì)算方法相對(duì)于wlp-fdtd方法計(jì)算速度得到了提高，雖然計(jì)算速度比因式分裂的wlp-fdtd方法要稍慢點(diǎn)但精度比因式分裂的wlp-fdtd方法要高，尤其是在高頻段，精度提高更明顯。