計(jì)算電磁學(xué)-第5章-時(shí)域有限差分法1

第5章時(shí)域有限差分法(1)

計(jì)算電磁學(xué)－時(shí)域有限差分法（FDTD）張洪欣電子工程學(xué)院

1主要內(nèi)容Maxwell偏微分方程的離散問(wèn)題Yee元胞FDTD差分格式蛙跳格式FDTD穩(wěn)定性和數(shù)值色散FDTD截?cái)噙吔鐥l件FDTD激勵(lì)源一維仿真實(shí)例二維仿真實(shí)例23FDTD介紹時(shí)域有限差分法（Finite-DifferenceTime-Domain,FDTD）是直接離散時(shí)域Maxwell偏微分方程組的方法.計(jì)算電磁學(xué)領(lǐng)域的時(shí)域有限差分法起源于20世紀(jì)60年代美籍華人K.S.Yee提出Yee離散格式。K.S.YeeNumericalsolutionofinitialboundaryvalueproblemsinvolvingMaxwell'sequationsinisotropicmedia[K.S.Yee,IEEETrans.AntennasPropagat.Page(s):302-307,1966,Volume:AP-14]FDTD方法提出之后，隨著計(jì)算技術(shù)，特別是電子計(jì)算機(jī)技術(shù)的發(fā)展，F(xiàn)DTD方法得到了長(zhǎng)足的發(fā)展，在電磁學(xué)，電子學(xué)，光學(xué)等領(lǐng)域都得到了廣泛的應(yīng)用4為求解由偏微分方程定解問(wèn)題所構(gòu)造的數(shù)學(xué)模型，有限差分法是將定解區(qū)域（場(chǎng)區(qū)）離散化為網(wǎng)格離散節(jié)點(diǎn)的集合。并以各離散點(diǎn)上函數(shù)的差商來(lái)近似該點(diǎn)的偏導(dǎo)數(shù)，使待求的偏微分方程定解問(wèn)題轉(zhuǎn)化為一組相應(yīng)的差分方程。根據(jù)差分方程組解出各離散點(diǎn)處的待求函數(shù)值—離散解。5

能直接給出非常豐富的電磁場(chǎng)問(wèn)題的時(shí)域信息。如果需要頻域信息，則只需對(duì)時(shí)域信息進(jìn)行Fourier變換。也是一種頻域分析方法。

為獲得寬頻帶的信息，只需在寬頻譜的脈沖激勵(lì)下進(jìn)行一次計(jì)算。67

任何問(wèn)題只要能正確地對(duì)源和結(jié)構(gòu)進(jìn)行模擬，時(shí)域有限差分法就能夠給出正確的解答；不管是散射、輻射、傳輸、透入或吸收中的哪一種，也不論是瞬態(tài)問(wèn)題還是穩(wěn)態(tài)問(wèn)題。8計(jì)算電磁學(xué)方法的比較低頻方法：MoM、FEM，諧振點(diǎn)以下的問(wèn)題；中頻方法：FDTD，第一諧振點(diǎn)附近4個(gè)量級(jí)；高頻：GTD幾何繞射，遠(yuǎn)高于諧振點(diǎn)的問(wèn)題；混合方法：例如，GTD＋FDTDFDTD最適于分析瞬態(tài)響應(yīng)問(wèn)題。FDTD用于分析低頻響應(yīng)問(wèn)題時(shí)計(jì)算時(shí)間很長(zhǎng)，例如電力線傳輸；而矩量法在分析高頻響應(yīng)時(shí)往往誤差過(guò)大，例如封閉金屬體內(nèi)接近諧振點(diǎn)問(wèn)題，對(duì)于低頻響應(yīng)問(wèn)題則存在優(yōu)勢(shì)。91.

Maxwell偏微分方程的離散問(wèn)題時(shí)域Maxwell偏微分方程組（理想介質(zhì)）Faraday’sEquationAmpère’sEquation5.1自由空間麥克斯韋方程的差分格式10Maxwell方程組各分量組成的6個(gè)方程

各場(chǎng)分量關(guān)于時(shí)間、空間的一階偏導(dǎo)的方程組，并且各場(chǎng)量相互嵌套。場(chǎng)分量是時(shí)間、空間的四個(gè)變量函數(shù)。11離散取樣根據(jù)有限差分法，對(duì)包含時(shí)間t、空間(x,y,z)的場(chǎng)分量偏導(dǎo)數(shù)采用中心差分法，進(jìn)行差分離散。首先對(duì)連續(xù)場(chǎng)量進(jìn)行離散，包含時(shí)間離散和空間離散。時(shí)間離散：假設(shè)在時(shí)間軸上均勻離散，離散的時(shí)間步長(zhǎng)為Δt

，用字符

n分別表示時(shí)間

nΔt的時(shí)刻標(biāo)示。用un上標(biāo)時(shí)間取樣點(diǎn)。

12離散取樣空間離散：假設(shè)在各方向上均勻離散，網(wǎng)格步長(zhǎng)

Δx,Δy,Δz，用字符

i,j,k分別表示x,y,z方向上的網(wǎng)格標(biāo)示。這樣連續(xù)的空間

(x,y,z)離散為用(i,j,k)表示的離散空間點(diǎn)——空間取樣點(diǎn)。

13場(chǎng)量時(shí)空離散場(chǎng)量u(x,y,z,t)的時(shí)空離散通常表示為在一定體積內(nèi)和一段時(shí)間上對(duì)連續(xù)電磁場(chǎng)的數(shù)據(jù)離散取樣。

14場(chǎng)量偏導(dǎo)的中心差分場(chǎng)量u(x,y,z,t)對(duì)空間一階偏導(dǎo)的差分格式場(chǎng)量u(x,y,z,t)對(duì)時(shí)間一階偏導(dǎo)的差分格式

15場(chǎng)量方程的離散問(wèn)題？是否可以把上面對(duì)空間、時(shí)間偏導(dǎo)的的差分格式代入6個(gè)Maxwell偏微分方程組？？？？6個(gè)場(chǎng)分量如何離散？？？？如果6個(gè)分量都在一個(gè)空間點(diǎn)、一個(gè)時(shí)間點(diǎn)上，有什么問(wèn)題？？？？

162.Yee元胞在一個(gè)空間點(diǎn)上的離散比較復(fù)雜。有沒(méi)有更為簡(jiǎn)便的離散方式？Yee元胞17Yee元胞場(chǎng)量分布6個(gè)場(chǎng)分量在Yee元胞的表面上進(jìn)行離散。在空間上，各電場(chǎng)分量

Ex,Ey,Ez在Yee元胞的棱邊中間取樣，方向與棱邊一致；各磁場(chǎng)分量

Hx,Hy,Hz在Yee元胞表面的中間取樣，方向垂直元胞面。在時(shí)間上，電場(chǎng)分量在棱邊上，在整數(shù)網(wǎng)格線上，電場(chǎng)分量在整時(shí)刻離散；磁場(chǎng)分量在元胞面中間，在半網(wǎng)格位置，磁場(chǎng)分量在半時(shí)刻離散。（或者磁場(chǎng)整數(shù)取樣，電場(chǎng)半時(shí)刻取樣）18YEE元胞特點(diǎn)每一磁場(chǎng)分量由四個(gè)電場(chǎng)分量環(huán)繞；同樣，每個(gè)電場(chǎng)分量由四個(gè)磁場(chǎng)分量環(huán)繞.19YEE元胞特點(diǎn)這種空間取樣方式符合Faraday感應(yīng)定律和Ampère環(huán)路定律的自然結(jié)構(gòu)。能描述電生磁、磁生電的電磁場(chǎng)傳播特性202122電磁場(chǎng)量在YEE氏網(wǎng)格上的定義法：電場(chǎng)在棱上，磁場(chǎng)在面上：電場(chǎng)分量在相對(duì)應(yīng)的坐標(biāo)方向上少一個(gè)（例如Ex在x方向少一個(gè)），在其它兩個(gè)方向上多一個(gè)（即Ey，Ez在x方向上比Ex多一個(gè)）。磁場(chǎng)分量在各個(gè)方向上對(duì)應(yīng)的個(gè)數(shù)正好相反。對(duì)于理想金屬表面場(chǎng)的模擬，在面上能直接滿足切向電場(chǎng)為零，法向磁場(chǎng)為零的條件。這是比較方便的。23磁場(chǎng)在棱上，電場(chǎng)在面上。磁場(chǎng)分量在相對(duì)應(yīng)方向上少一個(gè)，在其它兩個(gè)方向上多一個(gè)。電場(chǎng)分量在各個(gè)方向上對(duì)應(yīng)的個(gè)數(shù)正好相反。對(duì)于理想金屬表面的模擬，切向電場(chǎng)為零，法向磁場(chǎng)為零，要考慮半個(gè)空間的步長(zhǎng)。這種情況對(duì)于計(jì)算表面電流比較方便。24YEE元胞網(wǎng)格編號(hào)Ex(i+1/2,j,k)Ey(i,j+1/2,k)Ez(i,j,k+1/2)Hx(i,j+1/2,k+1/2)Hy(i+1/2,j,k+1/2)Hz(i+1/2,j+1/2,k)半個(gè)網(wǎng)格編號(hào)問(wèn)題?253.Yee元胞的Maxwell差分格式分析Maxwell的一偏微分方程26偏微分方程離散以場(chǎng)分量Hx的Yee元胞面為例，偏微分方程在YEE元胞上的離散方程離散問(wèn)題，空間上以Hx在點(diǎn)(i,j+1/2,k+1/2)，時(shí)間整時(shí)刻n取樣27偏微分方程離散寫(xiě)成場(chǎng)量的離散格式由中心差分法28偏微分方程離散整理后得同理，可以推導(dǎo)其他偏微分方程的差分格式2930315.2有耗媒質(zhì)麥克斯韋方程的差分格式假定媒質(zhì)為線性，各向同性，有耗媒質(zhì)，于是無(wú)源區(qū)內(nèi)Maxwell方程組中的兩個(gè)旋度方程為：式中

為介電常數(shù)（F/m）;μ為磁導(dǎo)率（H/m）;σ為電導(dǎo)率

；

σm為等效磁導(dǎo)率

，其引入是為使方程具有對(duì)稱性。32六個(gè)分量所滿足的方程33343536373839

任一網(wǎng)格點(diǎn)上的電場(chǎng)

分量只與其上一個(gè)時(shí)間步的值及周?chē)鷩@它的磁場(chǎng)分量有關(guān)；同樣，任一網(wǎng)格點(diǎn)上的磁場(chǎng)分量也只與其上一個(gè)時(shí)間步的值及周?chē)鷩@它的電場(chǎng)分量有關(guān)。40414243444546472.FDTD的蛙跳計(jì)算Maxwell兩個(gè)矢量方程48蛙跳格式FDTD的電場(chǎng)、磁場(chǎng)矢量是相差半個(gè)時(shí)間步長(zhǎng)，兩個(gè)方程式顯示格式。FDTD步進(jìn)計(jì)算過(guò)程如同青蛙上下跳躍向前步進(jìn)的過(guò)程491.影響FDTD的因素單元尺寸、時(shí)間步長(zhǎng)、入射場(chǎng)、散射體結(jié)構(gòu)、場(chǎng)計(jì)算、吸收邊界條件、資源等；精度0.1dB－120dB由空間步長(zhǎng)確定單元尺寸；粗網(wǎng)格、中網(wǎng)格、細(xì)網(wǎng)格。由穩(wěn)定性條件決定時(shí)間步長(zhǎng)；邊界條件：Mur,PML時(shí)間步數(shù)足夠多，得到穩(wěn)定的周期解。估算存儲(chǔ)空間和總花費(fèi)?！?.3FDTD的穩(wěn)定性50FDTD方程組是時(shí)間離散上的顯式格式，必須要滿足一定的穩(wěn)定條件才能計(jì)算，否則會(huì)發(fā)散如何分析時(shí)間偏導(dǎo)差分離散的穩(wěn)定性？任意場(chǎng)量u(x,y,z,t)

在時(shí)域上可以分解為各時(shí)諧場(chǎng)的疊加，因此這里考慮考慮時(shí)諧場(chǎng)情形一階偏導(dǎo)512.Courant穩(wěn)定條件中心差分離散定義數(shù)值時(shí)間增長(zhǎng)因子q得到時(shí)間增長(zhǎng)因子的方程52Courant穩(wěn)定條件增長(zhǎng)因子q滿足|q|<1,那么時(shí)間步進(jìn)不會(huì)發(fā)散一個(gè)場(chǎng)量對(duì)時(shí)間離散要滿足的穩(wěn)定性要求對(duì)于Maxwell方程組，包含6個(gè)場(chǎng)分量，如何分析穩(wěn)定性？通過(guò)等價(jià)的頻域波動(dòng)方程來(lái)穩(wěn)定性分析平面波解53Courant穩(wěn)定條件中心差分離散平面波代入上式54Courant穩(wěn)定條件得到波動(dòng)方程差分離散后的傳播常數(shù)k(kx,ky,kz)與頻率w

之間，即數(shù)值色散關(guān)系。當(dāng)Δx,Δy,Δz為零時(shí)，由前面推導(dǎo)的時(shí)諧場(chǎng)量時(shí)間離散要求色散關(guān)系不等式（上式兩側(cè)同乘）55Courant穩(wěn)定條件式中得即時(shí)間、空間步長(zhǎng)的穩(wěn)定條件，Courant穩(wěn)定條件56Courant穩(wěn)定條件假設(shè)均勻網(wǎng)格，即三維FDTD二維FDTD一維FDTD57FDTD空間步長(zhǎng)要求Courant條件，是對(duì)FDTD離散的時(shí)間步長(zhǎng)與空間步長(zhǎng)之間的關(guān)系要求，那么空間步長(zhǎng)如何確定？空間步長(zhǎng)反映場(chǎng)量在空間上離散取樣的情況，由取樣定理可知，一個(gè)波長(zhǎng)必須取樣兩個(gè)點(diǎn)以上采用還原信號(hào)，空間離散越緊密，越能還原信號(hào)。對(duì)于電磁場(chǎng)量離散滿足什么樣的取樣？還是用波動(dòng)方程推導(dǎo)的色散關(guān)系進(jìn)行分析討論一維情況58FDTD空間步長(zhǎng)要求當(dāng)寫(xiě)為相速形式當(dāng)也就是把代入上式，得對(duì)于寬帶時(shí)域信號(hào)，滿足信號(hào)最小波長(zhǎng)(fmax)，那么其它小于最大頻率的波也都滿足，并且離散精度更高。593.數(shù)值色散物理色散：非線性媒質(zhì)中速度隨波長(zhǎng)的變化。數(shù)值色散：Yee網(wǎng)格抽樣的間斷性，使得不同頻率的電磁波在網(wǎng)格傳播中表現(xiàn)不同的速度，并且與傳播方向有關(guān)。量化時(shí)各波列采用固定的步長(zhǎng)，有的波長(zhǎng)為整數(shù)倍關(guān)系，有的不是，則造成跳躍式傳播，在速度和方向上都造成誤差；Yee網(wǎng)格在網(wǎng)格節(jié)點(diǎn)處才有意義，與物理空間相比存在缺陷，只能近似表達(dá)物理規(guī)律和過(guò)程；60數(shù)值色散規(guī)律?s=λ/5時(shí)，比要求的標(biāo)準(zhǔn)網(wǎng)格大；?s=λ/10時(shí)，基本滿足要求的標(biāo)準(zhǔn)網(wǎng)格；?s=λ/20時(shí)，比要求的標(biāo)準(zhǔn)網(wǎng)格細(xì)；61相速度最大值出現(xiàn)在45度，最小值出現(xiàn)在0和90度；當(dāng)網(wǎng)格空間步長(zhǎng)?s減小時(shí)，網(wǎng)格中的速度與實(shí)際物理空間的速度偏差減小；相速度隨空間步長(zhǎng)的增加而減小，直到為零；不同的入射角，相速度下降對(duì)應(yīng)的?s不同；在45度，急劇下降的臨界步長(zhǎng)值最大。一定頻率、入射角的平面波在Yee網(wǎng)格中有一空間步長(zhǎng)極限，超過(guò)該極限電磁波將不能在此空間傳播。一個(gè)給定的Yee網(wǎng)格相當(dāng)于一個(gè)低通濾波器，給寬頻脈沖的計(jì)算帶來(lái)困難。高頻分量被截至。

取?s=λ/10非均勻網(wǎng)格還會(huì)導(dǎo)致折射效應(yīng)。62獲得理想色散關(guān)系的條件（1）：（2）：波沿網(wǎng)格的對(duì)角線方向傳播。63YEE元胞網(wǎng)格編號(hào)Ex(i+1/2,j,k)Ey(i,j+1/2,k)Ez(i,j,k+1/2)Hx(i,j+1/2,k+1/2)Hy(i+1/2,j,k+1/2)Hz(i+1/2,j+1/2,k)半個(gè)網(wǎng)格編號(hào)問(wèn)題?644.FDTD的截?cái)噙吔鐥l件Maxwell方程組在連續(xù)空間上離散成一個(gè)個(gè)Yee元胞，但實(shí)際上，仿真的電磁空間不可能把整個(gè)空間離散出來(lái)，這樣計(jì)算機(jī)無(wú)法計(jì)算。因此在仿真空間引入截?cái)噙吔鐏?lái)減少計(jì)算資源，這就是截?cái)噙吔鐥l件。65PEC邊界條件PEC邊界：邊界上的切向電場(chǎng)分量設(shè)置為零對(duì)于Xmin截?cái)噙吔缙渌较虻腜EC邊界相似66PMC邊界條件PMC邊界：邊界上的切向磁場(chǎng)分量設(shè)置為零。但截?cái)噙吔缟蠜](méi)有切向磁場(chǎng)分量，如何處理？虛擬磁邊界。對(duì)于Xmin截?cái)噙吔?75.FDTD的激勵(lì)源FDTD是對(duì)時(shí)域電磁波在空間中傳輸情況進(jìn)行仿真的數(shù)值算法，能直觀地仿真電磁波的時(shí)空信息。有源的激勵(lì)才能有電磁波的傳輸。1、時(shí)諧場(chǎng)源68激勵(lì)源2、脈沖源高斯脈沖升余弦脈沖69高斯脈沖源分析：在計(jì)算機(jī)中，高斯脈沖源受到截?cái)?。需要確定β、α、τ以保證帶寬。7071首先確定β。

β為從截?cái)帱c(diǎn)到高斯脈沖峰值的時(shí)間步數(shù)。即由τ=0持續(xù)到τ=2

β?t，峰值在τ=β?t處。截?cái)帱c(diǎn)：由τ=0持續(xù)到τ=2

β?t得選取α，避免引入無(wú)關(guān)的高頻浪費(fèi)時(shí)間。取α隨β變化。72截?cái)帱c(diǎn)處振幅衰減－140dB信號(hào)源問(wèn)題： Yee網(wǎng)格的抽樣和信號(hào)抽樣滿足頻譜的要求。例如：Yee網(wǎng)格?s=1cm。取β=32,脈沖如下：73當(dāng)波長(zhǎng)取10?s,λ=10cm，f=10G;當(dāng)波長(zhǎng)取4?s,λ=4cm，f=7.5G;均能滿足精度要求。當(dāng)存在介質(zhì)時(shí)，速度減小，例如εr=4，

f=3.0G,?s應(yīng)取0.5cm74激勵(lì)源微分高斯脈沖截?cái)嗳嘞颐}沖75激勵(lì)源雙指數(shù)脈沖Bell波形，應(yīng)用在核電脈沖和雷電脈沖的研究766.一維FDTD仿真仿真平行板傳輸線主模TEM7778方程簡(jiǎn)化79參數(shù)設(shè)置結(jié)構(gòu)參數(shù)激勵(lì)源參數(shù)（高斯源）

80仿真參數(shù)計(jì)算

81仿真參數(shù)計(jì)算仿真最大頻率－－＞仿真最小波長(zhǎng)數(shù)值色散的條件－－＞網(wǎng)格步長(zhǎng)－－＞仿真區(qū)域網(wǎng)格分布穩(wěn)定性條件－－＞時(shí)間步長(zhǎng)

82邊界設(shè)置1、終端短路截?cái)噙吔鐥l件2、終端匹配（吸收邊界條件）83場(chǎng)量觀察設(shè)置觀測(cè)面上點(diǎn)的電壓波形觀測(cè)仿真過(guò)程中的面上電場(chǎng)波動(dòng)情況84Hxi[0]=ez_lowm2; ez_lowm2=ez_lowm1; ez_lowm1=Hxi[1]; Ezi[JE-1]=ez_higm2; ez_higm2=ez_higm1; ez_higm1=Ezi[JE-2];Hxi[j]=Hxi[j]-dt/(m0*mr)*((Ezi[j+1]-Ezi[j])/dy);Ezi[j]=(2.0*e0*erdt*sigma)/(2.0*e0*er+dt*sigma)*Ezi[j]-2.0*