《计算电磁学》第五讲(二)ppt课件.ppt
2022/11/7,1,Oct. 31, 2011,计算电磁学第五讲（二）,Dr. Ping DU (杜平),E-mail: ,School of Electronic Science and Applied Physics， Hefei University of Technology (HFUT),2022/11/7,2,Berenger完全匹配层,Gedney完全匹配层,Outline,2022/11/7,3,Berenger完全匹配层,PML媒质的定义,首先，以二维TE( )情形为例建立PML媒质的方程。,如图4所示，在直角坐标系中，电磁场不随z坐标变化，电场位于(x, y)平面。,图4 TE问题,2022/11/7,4,电磁场只有3个分量,Maxwell方程退为3个方程。媒质的介电常数、磁导率分别为,磁阻率为 .,Maxwell方程为,(5-42),(5-43),(5-44),2022/11/7,5,如果,成立，则该媒质的波阻抗与自由空间波阻抗相等。当波垂直入射到媒质-自由空间分界面时，不存在反射。,定义TE情形下的PML媒质。需要将,分裂为,这样，PML媒质有4个场分量,，其满足方程,两个分量。,(5-45),(5-46),(5-47),2022/11/7,6,其中，,为电导率，,为磁阻率。,(5-49),(5-48),2022/11/7,7,PML媒质中平面波的传播,考察正弦平面波在PML媒质中的传播。,令电场分量的振幅为,，与y 轴夹角为,，如图4所示。,分裂磁场分量,的振幅分别为,这些场分量分别可表示为,(5-50),(5-51),2022/11/7,8,是角频率，t是时间，,为复常数。,其中,已知，方程(5-50)-(5-53)中有4个待定量：,假设,将他们代入到方程(5-46)- (5-49)中，可得4个待定量的关系式为,(5-53),(5-52),(5-54),2022/11/7,9,消去,，可得,(5-57),(5-56),(5-55),(5-58),2022/11/7,10,我们可以解出,首先，,的比值为,(5-59),(5-60),2022/11/7,11,然后由式(5-60)和(5-59)可求得 。,由(5-60)、(5-58)可求得,由此，存在两组符号相反的,，代表着两个相反的传播方向,选择正的一组解，有,(5-62),(5-61),2022/11/7,12,其中，,于是求得了PML媒质中平面波的传播常数。,(5-65),(5-64),(5-63),2022/11/7,13,用,代表该平面波的任一场分量，其振幅为,。,由式(5-50)-(5-53)及式(5-61)-(5-62)，有,最后，将式(5-61)、(5-62)的值代入式(5-56)(5-57)，可解得,(5-66),(5-67),(5-68),2022/11/7,14,由式(5-63)-(5-65)，可得,电场与磁场的振幅比为,现在假设两组参数 和 均满足(5-45)。,此时，对任何频率 均等于1.,(5-69),(5-70),2022/11/7,15,场分量和波阻抗的结果为, (71)第一指数项表明，PML媒质中的波在垂直于电场的方向 以自由空间中的波速传播； 后两个指数项表明，此时波的振幅沿x和y方向按指数规律减小, 式(72)表明：此时PML媒质的波阻抗与自由空间波阻抗相同,(5-71),(5-72),2022/11/7,16,对式(5-66)，如果波沿y方向传播，并且 ，,反之，如果波沿x方向传播，并且 ，,则波沿x方向无衰减地传播，PML不能吸收波。,则波沿y方向无衰减地传播，PML不能吸收波。,观察波函数表达式(5-71)，如果 ，,则沿y方向的指数项等于1，波只在x传播方向被吸收；,反之，如果 ，则沿x方向的指数项等于1，波,只在y传播方向被吸收,2022/11/7,17,PML-PML媒质分界面处波的传播,考察电磁波从一种PML到另一种PML媒质的传播问题。,先考虑分界面垂直于x轴的情形，如图5所示。,图5 垂直于x轴的PML-PML媒质分界面,2022/11/7,18,设,分别表示入射、透射、反射电场,相对分界面的角度,和,媒质参数分别为,如果PML媒质满足匹配条件(5-45)，根据式(5-71)，,分别为相对分解面法向定义的入射角、透射角、反射角(见图5),假设分界面为无限大，入射、透射、反射波均为平面波。,记,振幅分别为,由式(5-66)可写出入射、反射、透射电场在x=0分界面的表达式为,2022/11/7,19,其中，,可由(5-64)和(5-65)算出,(5-75),(5-74),(5-73),(5-76),2022/11/7,20,由(5-69)和(5-70)式可知入射、反射、透射磁场在x=0分解面分别为,(5-77),(5-78),(5-79),在分界面处，切向电磁场分量,必须连续，因而有,(5-80),(5-81),2022/11/7,21,将式(73)-(75)代入到(80)，该式要对所有的y都成立，必须使式(73)-(75)中的指数因子相等.,由此可得,(5-82),这两式为PML-PML媒质分界面处的Snell定律。,(5-83),最后，式(5-80)和(5-81)变为,2022/11/7,22,(5-84),(5-85),定义反射系数为分界面处反射与入射电场切向分量之比。,由上面两式，可得TE情形的反射系数为,由式(5-70)，反射系数还可写成,其中,是,的函数，可由(5-76)计算,(5-86),(5-87),2022/11/7,23,考察特殊情况：两种媒质的,相同。,此时Snell定律(5-83)变为,将(5-88)代入到式(5-87)，有,(5-88),(5-89),对(5-89)两边平方，将,用式(5-76)代入，注意此时,得,(5-90),2022/11/7,24,由(5-89)和(5-90)，可得,(5-91),这表明：即使媒质不满足匹配条件式(5-45)，只要两种媒质具有相同的,反射系数将不随入射角,变化；,只通过式(64)、(65)随频率而改变,如果媒质是匹配的，有,式(5-88)、(5-89)简化为,(,),(5-92),(5-93),2022/11/7,25,这表明：如果媒质满足匹配条件式(5-45)，只要两种媒质具有相同的,任意入射角、任意频率的平面波将无反射地通过PML-PML媒质分界面,当PML-PML分界面垂直于y轴时，有类似结论,综上所述，可得出结论：, 在垂直于x轴的PML-PML媒质分界面，如果两种媒质的,则反射系数始终为零,相同, 在垂直于y轴的PML-PML媒质分界面，如果两种媒质的,相同,则反射系数始终为零,2022/11/7,26,用于FDTD的PML,图6为PML与FDTD网格结合的方案。FDTD仿真区域假设为自由空间，它被PML媒质包围，PML又被PEC包围。,图6 FDTD网格与PML的结合,2022/11/7,27,(5-94),(5-95),(5-96),在PML媒质中距离分界面为r的地方，外向平面波的幅度可写成,其中 是相对媒质分界面定义的入射角，,是,或,如果PML媒质厚度为 ，则其表面的反射系数为,电导率从自由空间-PML分界面的0到PML最外面逐渐增加。,设离开自由空间-PML分解面距离为r处的导电率为,则分界面反射系数为,，,2022/11/7,28,如果取,(5-97),(5-98),将其代入到(5-96)，可得,当,(垂直入射)时，,(5-99),2022/11/7,29,PML中差分格式,在PML中可以采用指数差分，也可以采用普通的中心差分公式,二者的计算复杂度相当，准确程度也差不多,对于PML媒质中的偏微分方程(5-46)-(5-49)，仍旧采用Yee差分网格,的取样位置不变；,都在标准Yee差分网格中,的取样位置取样(见图7),PML媒质内部，式(5-46)-(5-49)的指数差分格式分别为,2022/11/7,30,图7 右上角区域的FDTD网格,2022/11/7,31,(5-100),(5-101),2022/11/7,32,(5-102),(5-103),其中，,在左、右侧边及角处是 x(i)的函数，在上、下侧边为零,在上、下侧边及角处是 y(j)的函数，在左、右侧边为零,2022/11/7,33,对于位于分界面上的电场分量，其差分格式中要用到自由空间的磁场分量,在,的上分界面,(5-104),在,的右分界面,(5-105),2022/11/7,34,在公式中，,均取以取样点为中心的单元上的平均值,如,(5-106),(5-107),(5-108),(5-109),(5-110),2022/11/7,35,结论：, 在垂直于x轴的PML-PML媒质分界面，如果两种媒质的,相同,则反射系数始终为零, 在垂直于y轴的PML-PML媒质分界面，如果两种媒质的,相同,则反射系数始终为零,对TM模，可以像TE情形那样构造无反射PML吸收层。,三维情况的PML,此时，直角坐标系下的6个场分量均被分裂：,2022/11/7,36,PML媒质中的场分量满足下列方程,(5-111),(5-112),(5-113),(5-114),(5-115),2022/11/7,37,(5-116),(5-119),(5-118),(5-117),(5-120),(5-121),(5-122),2022/11/7,38,匹配条件为,当该条件满足时，波在i方向很快衰减。,(5-123),(5-124),(5-125),2022/11/7,39,对方程(5-111)-(5-122)进行差分化，得PML媒质中的差分格式,(5-126),(5-127),(5-128),2022/11/7,40,(5-129),(5-130),(5-131),2022/11/7,41,(5-134),(5-133),(5-132),2022/11/7,42,(5-135),(5-136),(5-137),2022/11/7,43,其中，系数A、B分别按不同差分格式分别给出：,(1)指数差分格式,(5-138),(5-139),(5-141),(5-140),2022/11/7,44,(2) 标准中心差分格式,(5-142),(5-143),(5-144),(5-145),2022/11/7,45,PML参数选择,PML层数N、电导率分布阶数n、PML表面反射系数R(0)。数值实验表明:,减小反射误差措施,误差来自PML外侧PEC引起的理论反射系数；数值网格引起的误差。当电导率分布阶数较小时，前者为主；当其较高时，后者为主。,减小误差的措施有：,分界面出参数合理取值；,1) 当PML厚度固定，减小R(0)可以使局部及总体误差都单调减小。但当其小于10-5时，这种现象不再存在。一般选取R(0)=1%.,2) 增加PML厚度可以使局部及总体误差都单调减小。但层数太多增加计算量。一般N为48层。,3) 电导率分布阶数的改变不会影响计算量、但会影响吸收边界效果。常数分布会带来较强的数值反射；线性分布可以改善结果；选择n=2可进一步改进吸收效果。,2022/11/7,46, 用吸收算子代替PML外侧的理想导体, 衰减模的吸收,1997年，南京理工大学 D. G. Fang (方大纲)教授等在亚太微波会议上提出了MPML(修正PML)来吸收衰减模。,2022/11/7,47,Gedney完全匹配层(UPML),与Berenger提出的PML在数学上是等效的,理论体系是Maxwell方程的,便于理解和高效数值实现, 既可吸收传播模，也能同时吸收凋落模。 初始Berenger的PML难以完成的,完全匹配单轴媒质,设一任意极化的时谐平面波,在各向同性媒质中传播，并照射到半无限大空间的单轴各向异性媒质,假设分界面为z=0平面,(146),2022/11/7,48,在单轴各向异性媒质中激起的是平面波，其满足Maxwell方程,传播常数为,这样能保证分界面处相位匹配,将平面波代入到Maxwell方程，得,where,are the relative permeability and the relative permittivity,of the isotropic media, respectively, which are,and,(148),(147),2022/11/7,49,(149)和(150)式中，假设了该各向异性媒质关于z轴旋转对称。,注：当时谐场的时间因子(time dependence)为,时，,这样，,k为波的传播矢量,(149),(150),2022/11/7,50,由(149)和(150)，可推导得波动方程,where,(151) can be expressed in matrix form as,where,and,(151),(152),2022/11/7,51,只有系数行列式为零时，才可以得到非零解，也就得到了该单轴媒质的色散关系。这样，我们可得到4个本征模解。,对其进行去耦分解，得到前传和后传的,模和,模,它们满足,计算分界面的反射系数。先考察,模入射。,上半部各向同性媒质中的场表示为,(153),(154),(155),2022/11/7,52,透射进入各向异性媒质半空间的波也是, 满足色散关系(153), 可表示为,where,and,are the reflection and transmission coefficients, respectively,利用边界切向分量连续条件，得,(156),(157),(158),(159),2022/11/7,53,我们需要使反射系数对任意入射角都为0.,也就是说，,由式(152)，可进一步表示为,with,显然，当,且,时，,成立.,对,模, 可相似地推导出无反射条件为,且,。,由此，可以得出：,当平面波照射到由式(149)和(150)描述的单轴各向异性媒质时，,如, 则平面波无反射.,且不随入射角、极化、频率变化而变化,另外，由式(153)可知，此时,模和,模的传播特性相同 .,(160),2022/11/7,54,为使单轴媒质有耗，可取,此时,已对,归一化,于是，相对介电常数和相对磁导率张量为,(161),2022/11/7,55,色散关系现在可表示为,(162),由此可解出,(163),单轴媒质中的,模场表达式为,(164),2022/11/7,56,(165),where,(166),是相对于z轴的入射角。透射波与入射波的相速相同。,透射波沿z轴是衰减的，衰减常数不随频率变化，但与入射角及单轴媒质的导电率有关。透射波的横向特征波阻抗与入射波相同。,因此，媒质是完全匹配的。,2022/11/7,57,如果上半部各向同性媒质是有耗的，相应的入射波是有衰减的平面波,Maxwell旋度方程为,where,.,is the complex relative dielectric constant.,is the conductivity of isotropic media.,The complex wavenumber becomes,(167),(169),(168),2022/11/7,58,仍然假设单轴各向异性媒质是关于z轴旋转对称的，,和,满足(149)和(150),此时，单轴各向异性媒质中的,模和,模满足色散关系,(170),模的反射系数为,(171),由匹配条件,要求,(172),2022/11/7,59,According to the dispersion relation, we have,(173),It can be observed that (172) holds when c=a and b=a-1.,For the mode, we can derive the condition with zero reflection, which is,可得出结论：当有耗平面波照射到由式(149)、(150)描述的单轴各向异性媒质时，如果 ,则平面波将全部透射到该单媒质中。,它不会随入射波的入射角、极化、频率的变化而变化。,2022/11/7,60,如果仍取,由,和色散关系可以推出,(174),当入射波主要特征是凋落波时，,观察(173)式发现单轴媒质几乎不提供附加衰减，电磁波在其中以与入射波相同的速度衰减。,To accelerate the attenuation, we can choose,(175),2022/11/7,61,where,Thus, we can obtain,(176),UPML媒质中FDTD差分格式(FDTD formulation in the UPML),如果电磁波入射区是无耗的各向同性媒质，以相对于z轴旋转对称的完全匹配单轴各向异性媒质为例，Maxwell旋度方程可写成,(177),2022/11/7,62,(178),式(176)中的第三式不是标准形式，必须单独处理.,定义,(179),采用两步方式进行FDTD仿真,2022/11/7,63,于是，(177)中的第三式可写成,(180),采用标准FDTD方式刷新,(181),2022/11/7,64,By (179), we have,(182),通过替换, 可得其时域表达式,(183),By using the central difference, we obtain,(184),于是，利用式(181)和(184)可以通过两步方式刷新法向电场分量,2022/11/7,65,对式(178)，记等效磁阻率为, 满足关系,(185),于是，式(178)可重写为,(186),式(186)中第一、二两式与横向场分量,相关，表达式同各向同性,有耗媒质中的公式完全一样，因而可使用标准的FDTD差分格式。,2022/11/7,66,式(186)中第三式不是标准形式，需要单独处理。,采用两步方式进行FDTD仿真。,定义,(187),于是，式(186)中第三式可写成,(188),采用标准FDTD方式刷新,2022/11/7,67,(189),由式(187)，有,(190),通过替换,，得到其时域表达式,2022/11/7,68,(191),用中心差分，得到其差分格式如下,(192),这样，利用式(189)和(192)可以通过两步方式刷新法向磁场分量,如果电磁波入射区是有耗导电的各向同性媒质，相对介电常数为复数,，完全匹配PML媒质中的Maxwell旋度方程可写为,2022/11/7,69,(193),(194),2022/11/7,70,对于(193)中的第一、第二两式，定义,(195),(196),于是，第一、第二式可以写为,(197),(198),2022/11/7,71,Transforming them into time-domain counterparts, we have,(199),(200),他们与各向同性有耗媒质中的公式完全相同，因此可使用标准的FDTD差分格式刷新得到 和 .,转到第二步，由辅助关系式(195)、(196)刷新,和,.,对各向同性的有耗媒质，,(201),2022/11/7,72,辅助关系(195)、(196)的时域表达式为,(202),(203),由中心差分可得到刷新,的公式为,(204),2022/11/7,73,刷新,的公式与此相似。,对于式(192)中的第三式，定义,(205),这样，式(192)的第三式，可写成,(206),对于各向同性的有耗媒质，,上式可按标准的FDTD公式刷新,2022/11/7,74,(207),从解的稳定性考虑，对良导体也可采用指数差分格式或Luebbers差分格式。,然后，转入第二步，由辅助关系(205)可得到时域方程,2022/11/7,75,(208),用中心差分可得到其差分格式如下,(209),于是，利用式(207)和(209)可以通过两步方式刷新法向电场分量,对于式(192)，采用等效磁阻率,的定义(185)，可将式(194) 重写为,2022/11/7,76,(210),进一步，记,(211),(212),2022/11/7,77,则式(210)可进一步写成,(213),该式与式(186)形式相同，差分格式也应该相同。这里不再叙述。,对于以,或,平面为分界面的单轴各向异性完全匹配层，,其公式也类似。,2022/11/7,78,角区的差分格式,在角区，电磁场满足Maxwell旋度方程,(214),(215),其中，,(216),分别只与x、y、z方向的分界面有关，且分别只沿x、y、z方向有变化,2022/11/7,79,(217),(218),(219),在角区之外，相应的,为零，式(216)退化为单一平面边界的结果。,角区的差分格式可按相似的过程推导出。以式(215)中的第三个方程为例，定义,(220),2022/11/7,80,该式中的第三个方程可写成,(221),将其转换到时域，可用标准的FDTD差分格式刷新,然后由式(220)，可得,(222),将其转换到时域，采用中心差分格式，可得到下列二阶精度的差分方程,(223),2022/11/7,81,由式(220)的标准FDTD公式和式(222)，按两步方式可以刷新角区的,其余分量的FDTD刷新公式可相似地得到。,如仿真区是有耗导电媒质，此时设在角区的电磁场满足Maxwell旋度方程,(224),(225),其中，,依然按照式(215)定义,(226),(227),(228),差分格式的推导与前面相似。,2022/11/7,82,PML的参数选择 Parameter choice for PML,为避免在PML表面引入反射，以,的PML表面为例，,、,应从表面起渐变增加，可采用如下形式,(229),(230),其中，,为PML媒质厚度。根据Gedney的经验，一般取,m=4,.,注意：PML媒质中的时间步长仍要满足CFL稳定条件。,(231),2022/11/7,83,Thank you very much!,
