《《计算电磁学》第四讲ppt课件.ppt》由会员分享，可在线阅读，更多相关《《计算电磁学》第四讲ppt课件.ppt（40页珍藏版）》请在taowenge.com淘文阁网|工程机械CAD图纸|机械工程制图|CAD装配图下载|SolidWorks_CaTia_CAD_UG_PROE_设计图分享下载上搜索。
1、8/7/2022 6:01 PM第四讲 FDTD:差分格式及解的稳定性Dr. Ping DU (杜平)E-mail: School of Electronic Science and Applied Physics， Hefei University of Technology (HFUT)8/7/2022 6:01 PMp FDTD基本原理u Yee差分算法 考虑一无源区域，其媒质的参数不随时间变化且各向同性，则Maxwell旋度方程可以写成 1tEHE(4.1)(4.4)其中，E为电场强度，H是磁场强度，是介电常数， 是电导率， 是磁导率。 在直角坐标系中，式(4.1)、（4.2）变为 1
2、yxzHEEtxz(4.3)(4.2)1yxzEHEtzy1t HE8/7/2022 6:01 PM1yxzEEHtyx1yxzxHEHEtyz1yxzyEHHEtzx1yxzzHHEEtxy(4.5)(4.6)(4.8)(4.7)式(4.3)-(4.8)是FDTD算法的基础。 1966年，美籍华人K. S. Yee对上述6个方程引入了一种差分格式Yee网格。其原理是，首先在空间建立矩形差分网格，网格节点与一组相应的整数标号相对应：( , , )(,)i j ki x j y k z(4.9)8/7/2022 6:01 PM该函数在时刻 的值可以表示为 n t( , , )(,)nFi j k
3、F i x j y k z n t(4.10)其中， , , 分别为矩形网格沿x, y, z方向的空间步长，为时间步长。Yee采用中心差分来代替对空间和时间坐标的微分，因而具有二阶精度xyz211(, , )(, , )( , , )22() )nnnFij kFij kFi j kOxxx11222( , , )( , , )( , , )() )nnnFi j kFi j kFi j kOttt(4.11)(4.12) 为获得(4.11)中的精度，并满足(4.3)-(4.8),Yee将空间任一矩形网格上的E和H的6个分量按如图4-1所示放置。8/7/2022 6:01 PM图4-1 Yee
4、差分网格(Yees cell) 8/7/2022 6:01 PM112211( ,1)( , )111122,( ,)112222,2211( , ,)( ,1,)22nnyynnxxnnzzEi jkEi jktHi jkHi jkzi jkEi j kEi jky(4.13)112211,1, ,111122,112222,2211,(1, )22nnxxnnzznnyyEijkEij ktHijkHijkyijkEi jkEijkx (4.14)(4.15)1122111, , ,111122, , ,112222, ,2211, ,(, ,1)22nnzznnyynnxxEij kEi
5、 j ktHij kHij kxij kEij kEij kz 8/7/2022 6:01 PM(4.16)111221, ,2112, ,1112, , ,112211, , ,1, ,2, ,2222112, ,211,22nnxxnnzzij ktij ktEij kEij kij ktij kij ktij kij kHijkH112211,221111, , ,2222nnyyijkyHij kHij kz8/7/2022 6:01 PM(4.17)111221,2112,1112,112211,1,2,2222112,211,22nnyynnxxi jkti jktEi jkEi
6、jki jkti jki jkti jki jkHi jkH112211,221111,2222nnzzi jkzHijkHijkx8/7/2022 6:01 PM(4.18)111221, ,2112, ,1112, , ,112211, , ,1, ,2, ,2222112, ,211, ,22nnzznnyyi j kti j ktEi j kEi j ki j kti j ki j kti j ki j kHij kH112211, ,221111,2222nnxxij kxHi jkHi jky8/7/2022 6:01 PM 以上通过对各网格点上电磁场交替进行计算，在执行到恰当的时
7、间步后，就能得到所需要的时域结果。这种差分格式成为蛙跳格式。式(4.13)-(4.18)的差分格式中，从第n层更新到第n+1层时，提供了逐点计算场分量的直接表达式，这被称为显式格式。于此相对的是隐式格式，它通常需要求解代数方程组，才能得到解。以后介绍的ADI-FDTD就是隐式格式的。 由式(4.13)-(4.18)可以看出，每个网格点上的各个场分量的新值由该点在前一个时间步长时刻的值和该点周围邻近点上另一场量的场分量早半个时间步长时刻的值所决定。我们可以一次算出网格点上的场分量。也可以采用多个并行处理器一次算出多个点上的场分量，所谓的并行FDTD就是这样。8/7/2022 6:01 PM111
8、SCddt 1DSHlu 环路积分解释 上面从Maxwell旋度方程出发，利用中点差分公式，导出Yee差分方程。其实，也能够从积分形式的Maxwell方程、Ampere定律和Farady定律推导出。为简化，仅考虑自由空间情形。 如图4-2所示，将Ampere定律用于环路 ，有 1C(4.19)图4-2 环路 1C8/7/2022 6:01 PM11111111, , ,22222222xyxyRHSHi jkxHij kyHi jkxHij ky 1011, , ,22nnzzx yLHSEi j kEi j kt 假设场分量在环路每边中点的值等于场分量在该边的平均值，于是 (4.20)再假设
9、 等于 在小面元S1的平均值。于是用中点差分代替对时间的偏导，1, ,2zEi j kzE中点取在12tnt ，得 (4.21)8/7/2022 6:01 PM111112222011, , ,2211111111, , ,22222222nnzznnnnyyxxEi j ki j kHij kHij kHi jkHi jktxyExE由式(4.20)、(4.21),得(4.22)相似地，可以对 、 ,利用Ampere定律导出相应的差分方程。 yE8/7/2022 6:01 PM同样地，可以将Farady定律用于图4-3所示的环路 2C22SCddt 22BSEl(4.23)8/7/2022
10、6:01 PM图4-3 环路 2C8/7/2022 6:01 PM1111, ,1,1,2222xyxyRHSE ij kx E ijky E ijkx E i jky 110221111,2222nnzzx yLHSHijkHijkt 假设场分量在环路每边中点的值等于场分量在该边的平均值，于是 (4.24)再假设 等于 在小面元S2的平均值。于是用中点差分代替对时间的偏导，11,22zHijkzHzH中点取在 ，得 tn t (4.25)由式(4.24)、(4.25),得 8/7/2022 6:01 PM1122011,1,1,111122,222211,1,22nnxxnnzznnyyE
11、ijkE ijktHijkHijkyE i jkE ijkx(4.26)xHyH相似地，可以对 、 ，利用Farady定律导出相应的差分方程。 8/7/2022 6:01 PM 从积分形式的Maxwell方程出发推导出的差分方程，对于处理细线、槽缝、弯曲表面等结构处理特别方便，可将环路路径选取来与弯曲表面、槽缝等结构共形匹配。 8/7/2022 6:01 PMp 解的稳定性及数值色散 在FDTD法中，时间步长 和空间步长 ， 和 并不是无关的。它们的取值必须满足一定关系，才能避免数值结果的不稳定。不稳定表现为，在求解显示差分方程时，随着步长步数的增加，计算结果将无限制地增加。 txyz 为了确
12、定数值解稳定的条件，必须考虑在FDTD算法中出现的数字波模，其基本方法是把有限差分算式分解为时间和空间的本征值问题。任何波都可以展开为平面波的叠加。因而，如果一种算法对平面波是不稳定的，则它对任何波都是不稳定的。因而，我们只需要考虑平面波本征模在数字空间中传播，这些模的本征值谱由数字空间微分方程来确定，并与由数字时间微分方程确定的稳定本征值谱比较。空间本征值谱必须包含在稳定区间，以确保这种算法中所有可能的数字波模是稳定的。 8/7/2022 6:01 PM为简单起见，近考虑无耗媒质空间( ，和 均为实数/张量)。 0考虑 、 的FDTD方程，分别为 zEzH11122112211, , ,22
13、1111, , ,22221111,2222nnzznnyynnxxEi j kEi j kHij kHij ktxHi jkHi jky (4.27)8/7/2022 6:01 PM11221111,2222111,2211,1, ,22nnzznnyynnxxHijkHijkEijkEi jktxEijkEij ky (4.28)对其他4个分量也可以求得类似的方程。 8/7/2022 6:01 PM 这两个方程左边可构成个相应场分量的时间本征值方程。用V表示各分量，则他们可统一写成 1122nnnVVVt(4.29)定义增长因子 12/nnqVV根据冯诺依曼稳定性条件，要求 | 1q 将其
14、代入到(4.29)，两边同除以 12nV可得 必须满足 q210qtq (4.30)8/7/2022 6:01 PM其解为 2122ttq(4.31)为满足 | 1q 只需要 Re( )022Im( )tt这就是说，为保证算法的稳定性，时间本征值必须落在虚轴的稳定区间。 (4.32)8/7/2022 6:01 PM另外，我们可以将V表示为如下的平面波本征模 0( , , )exp()xyzV i j kVj ikxjkykkz (4.33)代入到(4.27)和(4.28)，其中 ,可得到各分量之间的关系.01j 其中 和 的表达式为 zEzH02sinsin22yyxxzHkyjkxHExy0
15、2sinsin22yyxxzEkyjkxEHxy(4.34)(4.35)8/7/2022 6:01 PM其他4个分量的表达式为 02sinsin22yyzzxkyHjHkzEyz(4.36)02sinsin22xxzzyjHkxkzHEzz02sinsin22yyzzxkyEjEkzHyx02sinsin22xxzzyjEkxkzEHzx(4.37)(4.38)(4.39)8/7/2022 6:01 PM将方程组(5.34)-(5.39),用矩阵形式写为 A0 xyzxyzEEEHHH(4.40)该方程组有非零解的条件是 det(A)0(4.41)8/7/2022 6:01 PM解得 2222
16、222sinsinsin2422()()()yxzkykxkzxyz (4.42)对所有可能的 ， 满足 ,xyzkkk和Re( )0222111Im( )2()()()vxyz(4.43)其中， 1v1v8/7/2022 6:01 PM为保证数值稳定性，式(4.43)所示的区域必须落入到时间本征值的稳定区域内。 于是，由式(4.32)和(4.43)，可得 2221111()()()tvxyz (4.44) 这就是FDTD算法的数值稳定条件(CFL稳定性条件)。对非均匀区域，应选最大的v。 对二维问题(如场不随z变化)，可令z 对一维问题(如场不随y和z变化)，可分别令y z 8/7/2022
17、 6:01 PMu 数值色散0( , , )exp()nxyzVi j kVj ikxjkyzkznt 差分方法会在计算网格中引起所模拟波模的色散，即在FDTD网格中，数值波模的传播速度将随频率变化。这会引起脉冲波形畸变、人为的各向异性及虚假折射现象。(4.45)考虑一个单色平面波，其各分量可表示为 将其代入到(4.27)、(4.28)，得 因而,FDTD中必须考虑数值色散问题。8/7/2022 6:01 PMsinsinsin222yyxxzEkykxEttHxy(4.47)sinsinsin222yyxxzkyHHkxttEyx(4.46)8/7/2022 6:01 PM对其余4个分量可得
18、到类似的关系。它们构成齐次线性方程组 0 xyzxyzEEEBHHH (4.48)该齐次方程组有非零解的条件是det(B)0(4.49)8/7/2022 6:01 PM可解得 22222222sinsinsin2122sin2()()()yxzkykxkztv txyz(4.50)式(4.50)即为三维情况的数字色散关系式， 其中 分别为波矢量沿x，y，z方向的分量，为角频率，v是被模拟的均匀媒质中的光速。 ,xyzkkk和与数字色散相对应，在无耗介质中的平面波，其解析色散关系式为 2222xyzkkkv(4.51)8/7/2022 6:01 PM由式(4.50)可以看出，当t、x、y、z均趋
19、于零时，它变为式(4.51)。这说明数字色散可以降到任意小，只要空间步长和时间步长足够小。在实际中，步长只能取有限大。因而，这会造成数字色散。 为定量说明数字色散与FDTD网格关系，以二维TM波为例进行数值计算假定xy 波的传播方向与x轴夹角为。 cos,sinxykkkkk为波矢量的模。 数字色散关系为 2222cossinsinsinsin222tkkv t(4.52)8/7/2022 6:01 PM给定 t 、用牛顿迭代法求式(4.52)，可得到 0k进而求得相应的相速度 pkv例如取2222002sin2sin2sin2sin222 2kftc t进一步，取 0=/5 2202sin2
20、sin1010 2ko=45 ,/ 2v tc t 由式(4.52)可得8/7/2022 6:01 PM用牛顿迭代法求此方程，可得到 0k由下式求得 000220.9823pppvvvcfk /pvc数值实验表明，FDTD算法：(1) 存在各向异性；8/7/2022 6:01 PM(3)可能会导致伪折射现象。 在实际应用中，空间步长取信号最高频率对应波长的1/10-1/20。 (2) Yee差分格式有数字低通滤波特性;仿真时，可能出现高频拖尾；解决办法：适当选取脉冲的频谱分量和网格步长，使得主要的频谱分量的波长至少为10倍网格步长。 8/7/2022 6:01 PMn良导体中的FDTD差分格式
21、 在Yee差分格式中(4.16)-(4.18)中，已经给出含电导率的FDFD差分格式，但其会导致时间步进趋于不稳定。 为使解稳定，电流项 -E的差分格式应介于中心平均和前向近似之间 对导体媒质中的Maxwell方程 001yzyEHEtx (4.53)Luebbers不以中心平均近似 112( )( )( )2nnnyyyEiEiEi(4.54)u 差分格式8/7/2022 6:01 PM而是以介于中心平均近似和前向近似之间的半步长前向近似格式 112( )( )nnyyEiEi(4.55)代入式(4.53)。于是得差分方程 112210001122( )( )2nnzznnyyHiHitEiEittx (4.56)该方程对良导体是稳定的。 另外，我们还可以采用指数差分格式，它对良导体也是稳定的，在此不讲授。 8/7/2022 6:01 PMu 良导体中FDTD格式的步长选择在良导体内部，电磁波传播的衰减常数和相位常数相等。其中电磁波的波长为22其中 为导体的趋肤深度(skin depth)。 2 为准确模拟良导体内部的电磁场，FDTD的空间步长h应小于最高信号频率的趋肤深度 。min比如取 min/5h或更小。为保证解的稳定性，时间步长应满足CFL条件。 8/7/2022 6:01 PMThank you!