（凝聚态物理专业论文）移动坐标系fdtd方法在电磁脉冲传播中的应用研究.pdf

摘要 摘要 时域有限差分( f d t d ) 方法是一种在计算电磁学中得到广泛应用的数值计算方 法。由于f d t d 方法是一种全波方法，计算时需要对整个空间进行网格划分和计算推 进，在模拟区域很长时需要大量的计算资源，因此f d t d 方法难以对电磁脉冲的长距 离传播进行数值模拟。通过在f d t d 方法中引入移动坐标系，使得坐标系与脉冲相互 “冻结 ，坐标系伴随脉冲一起传播，让脉冲始终位于计算区域内，可以用有限的计 算区域实现对脉冲长距离传播问题的模拟，从而将f d t d 方法应用于脉冲传播问题。 本篇论文研究了移动坐标系f d t d 方法在电磁脉冲传播模拟中的应用。 本篇论文首先介绍了电磁脉冲的种类以及电磁脉冲大气传播的研究情况。 介绍了f d t d 方法，传统f d t d 方法在模拟电磁脉冲传播时面临的困难，移动 坐标系f d t d 方法概念的提出及其研究现状。研究了引入移动坐标系后的f d t d 方法在一维情况下的应用，分析了一维情况下该方法的数值稳定性和数值色 散性质，给出了稳定性条件，证明了该方法在一维时没有数值色散。采用该 方法分别计算了均匀介质及非均匀介质的一维算例。通过将计算结果与理论 解以及传统f d t d 方法得到的结果相比较，验证了该方法的正确性以及优于传 统f d t d 方法的数值色散性质。将该方法推广应用于三维情况，推导了三维情 况下的移动坐标系f d t d 的差分格式，将其应用于脉冲在方波导中的传播，计 算了多个频率的脉冲传播算例，并将结果与传统f d t d 方法得到的结果相比较。 计算结果与传统f d t d 方法得到的结果吻合，从而验证了该方法的正确性。 本篇论文验证了移动坐标系f d t d 方法的可行性，为电磁脉冲传播问题的模拟提 供了一个有效工具。 关键词：移动坐标系时域有限差分方法电磁脉冲传播 a b s 仃a c t a b s t r a c t f i n i t e d i f f e r e n c et i m e - d o m a i n ( f d t d ) m e t h o di sw i d e l yu s e di ne l e c t r o m a g n e t i c s i m u l a t i o n b e c a u s ef d t dm e t h o di saf u l l - w a v em e t h o d ，w h e ns i m u l a t i n gp r o p a g a t i o n o fe l e c t r o m a g n e t i cp u l s e so v e rl o n gd i s t a n c e ，t h ec o m p u t a t i o n a lb u r d e nq u i c k l yb e c o m e s t o ol a r g e t h a ti st h em a j o rd i f f i c u l t yw h yf d t dm e t h o dc a nn o tb ea p p l i e dt os i m u l a t i n g l o n gd i s t a n c ep r o p a g a t i o no fp u l s e sd i r e c t l y w h e nu s i n gm o v i n gc o o r d i n a t ef l a m ef d t d m e t h o d ，t h ec o o r d i n a t es y s t e mi s f r o z e n w i t hp u l s e s ，a n dp u l s e sa r ei ns m a l lc a l c u l a t i o n z o n e sw h e ns i m u l a t i n g t h ef d t dm e s h e s ( i nt h ed i m e n s i o na l o n gt h ep r o p a g a t i o np a t h ) n e e do n l yt ob el a r g ee n o u g ht oc o n t a i nt h ep u l s e sw h i l eu s i n gt h i sm e t h o dt os i m u l a t e p r o p a g a t i o no fp u l s e s t h ea p p l i c a t i o no fm o v i n gc o o r d i n a t ef r a m ef d t d i np r o p a g a t i o n o fp u l s e si ss t u d i e di nt h i sd i s s e r t a t i o n f i r s t l y , t h et y p e so fe l e c t r o m a g n e t i cp u l s e sa n dt h ec u r r e n ts t u d i e s o fn u m e r i c a l s i m u l a t i o no fe m pp r o p a g a t i o na r ei n t r o d u c e d t h e n ，f d t dm e t h o da n dt h ed i f f i c u l t i e s i ns i m u l a t i n gt h ep r o p a g a t i o no fp u l s e so v e rl o n gd i s t a n c eu s i n gf d t dm e t h o dd i r e c t l y a r ei n t r o d u c e d i no r d e rt oo v e r c o m et h e s ed i f f i c u l t i e s ，m o v i n gc o o r d i n a t ef r a m ef d t d m e t h o da r ed e v e l o p e d a n d ，t h es t a b i l i t ya n dn u m e r i c a ld i s p e r s i o no ft h i sm e t h o di no n e d i m e n s i o na r es t u d i e d ，s t a b i l i t yc o n d i t i o ni sd e d u c e da n dt h er e s u l t ss h o wt h a tm o v i n g c o o r d i n a t ef l a m ef d t dm e t h o di no l l ed i m e n s i o ni sf l e eo fn u m e r i c a ld i s p e r s i o n s i m u l a t i o nr e s u l t so fh o m o g e n e o u sa n di n h o m o g e n e o u sm e d i au s i n gb o t hf d t da n d m o v i n gc o o r d i n a t ef l a m ef d t da r ec o m p a r e d ，t h ea d v a n t a g eo fm o v i n gc o o r d i n a t e f r a m ef d t di sv a l i d a t e d t h i sm e t h o di se x t e n d e dt ot h r e ed i m e n s i o nc a s e t h e d i s c r e t i z e df o r mo fm o v i n gc o o r d i n a t ef l a m ef d t di nt h r e ed i m e n s i o n si sg i v e n ，a n di s u s e di np r o p a g a t i o no fw a v e so fd i f f e r e n tf r e q u e n c i e si nar e c t a n g u l a rw a v e g u i d e a f t e r c o m p a r i n gt h er e s u l t so fb o t hf d t da n dm o v i n gc o o r d i n a t ef r a m ef d t d ，t h i sm e t h o di s n u m e r i c a l l y v a l i d a t e d t h ef e a s i b i l i t yo fm o v i n gc o o r d i n a t ef d t da r ep r o v e di n t h i sd i s s e r t a t i o n t h i s m e t h o dc a l lb ew i d e l yu s e df o rs i m u l a t i n gp r o p a g a t i o no fe l e c t r o m a g n e t i cp u l s e s k e yw o r d ：m o v i n gc o o r d i n a t e f l a m ef d t dm e t h o de l e c t r o m a g n e t i cp u l s e p r o p a g a t i o n 中国科学技术大学学位论文原创性和授权使用声明 本人声明所呈交的学位论文，是本人在导师指导下进行研究工作所 取得的成果。除已特别加以标注和致谢的地方外，论文中不包含任何他 人已经发表或撰写过的研究成果。与我一同工作的同志对本研究所做的 贡献均已在论文中作了明确的说明。 本人授权中国科学技术大学拥有学位论文的部分使用权，即：学校 有权按有关规定向国家有关部门或机构送交论文的复印件和电子版，允 许论文被查阅和借阅，可以将学位论文编入有关数据库进行检索，可以 采用影印、缩印或扫描等复制手段保存、汇编学位论文。 保密的学位论文在解密后也遵守此规定。 作者签名：奎噩 2 0 0 9 年0 5 月6 日 第一章电磁脉冲大气传播的研究进展 第一章电磁脉冲大气传播的研究进展 1 1 电磁脉冲简介 电磁脉冲( e l e c t r o m a g n e t i cp u l s e ，e m p ) 特别是高功率电磁脉冲( h i g hp o w e r e l e c t r o m a g n e t i cp u l s e ，h p e m p ) ，因为其在雷达、通信和军事上的广泛应用，发展 很快。虽然，高功率电磁脉冲没有明确的定义，但是通常认为电场大于1 0 0 v m 的 脉冲为高功率电磁脉冲。高功率电磁脉冲主要包括核电磁脉冲( n u c l e a r e l e c t r o m a g n e t i cp u l s e ，n e m p ) 、高功率微波( h i g hp o w e rm i c r o w a v e ，h p m ) 、超宽 带脉冲( u l t r a w i d e b a n d ，u w b ) 以及雷电产生的电磁脉冲。研究电磁脉冲在大气中 的传播，对于电磁脉冲在通信和军事上的应用，以及电磁脉冲对电子器件系统的防 护，都有很大的指导意义。 本文研究了利用时域有限差分方法( f i n i t e d i f f e r e n c et i m e d o m a i nm e t h o d f d t d ) 在电磁脉冲传播中的应用，首先简要介绍一下电磁脉冲传播的相关情况。 1 2 影响电磁脉冲传播的因素 由于大气中各种因素的影响，电磁波会发生折射现象，传播路径发生改变。同 时，由于电磁波对大气的电离和击穿等效应，大气中气体分子的吸收，大气沉积物 ( 雨、云雾、冰雪和沙尘等) 的吸收和散射，电磁波在大气中传播时会发生衰减， 电磁波强度减小。下面将介绍上述几个方面的国内外相关研究情况。 1 2 1 折射率 电磁脉冲在真空中传播时，可以很容易求得电磁脉冲传播的方向和传播路径。 在大气中传播时，由于大气分布的不均匀性，大气的参数在空间上不一致，电磁脉 冲在大气中的传播速度小于真空中的光速，且传播速度随位置变化，从而传播过程 中发生偏转，发生大气折射效应。折射效应的显著程度主要依赖于沿路径的折射指 第一章电磁脉冲大气传播的研究进展 数变化梯度和射线初始仰角 1 c 2 。 通常用无线电折射率刀描述折射指数，以- 定义为： 以= 0 一1 ) x 1 0 6 其中1 1 是大气的折射率。 射频段的折射率近似公式为 n ，n， 刀= 5 8 2 ( 形名) + 2 8 x1 0 5 ( 哆名2 ) ( 1 2 ) 式中，p 为大气压强，单位是t o r r ，t 是绝对温度，单位是k ，砌是水汽分压，单位 为毫巴。从上式可以看出，折射率主要依赖与温度和水汽浓度。高度升高时，大气 压强和水蒸气的压强随之减小，温度上升，折射指数随着高度升高而减小，地球表 面折射指数的典型数值为1 0 0 0 3 ，随高度变化的递减率为4 x 1 0 q m 。 g o r d o n 和d t h a y e r 七十年代提出一种2 0 g h z 频率的折射率改进方程 3 ： 如5 8 2 ( ) 酊1 + 4 8 6 ( ) g ：1 + 2 8 3 3 1 0 5 ( 2 ) g ：1 ( 1 3 ) 式中，见为干燥空气压力，单位是t o r r ，虻1 和瓦1 分别为干燥空气和水汽对非理想 气体特性的修正。 分析大气折射指数对电磁脉冲传播的影响时，常常采用一个介电常数沿径向变 化的球壳来模拟大气。将球壳分成许多薄的介质层，每层内的介质均匀分布，可以 用平面界面求反射系数和透射系数的方法确定电磁脉冲在这些介质层内的传播。确 定了大气中折射率分布之后，电磁脉冲在大气中的传播路径可以采用抛物方程方法 等不同方式求出【4 】- 6 】，进一步可以确定出电磁脉冲在大气中传播路径上相应的介 电常数和电导率等参数分布情况。 1 2 2 大气电离击穿 电磁脉冲在大气中传播时，会与大气中的各种成分发生许多复杂的物理作用， 使得脉冲在传播过程中逐渐衰减。其中起最主要作用，人们最感兴趣的是脉冲激发 出的强电场对大气的击穿电离作用【7 】- i 1 5 】。 一般情况下，在大气除了中性的气体原子和分子外，还存在有一定数量的自由 电子，被称为种子电子。电磁脉冲传播过程中激起的电场使大气中自由电子获得能 量，速度增加，被加速的电子在运动过程中会与中性气体原子和分子发生碰撞。当 电磁脉冲激发出的电场强度超过了大气击穿阈值时，电子获得的能量大于气体原子 2 第一章电磁脉冲大气传播的研究进展 和分子的电离能，高速电子与中性分子或原子碰撞过程中会使其电离，产生更多的 自由电子，这过程重复下去，在大气中产生越来越多的自由电子，当自由电子的 浓度超过临界浓度时，大气被击穿。这种击穿机制被称为雪崩式击穿。 图i - i 大气电离后自由电子增多 如果大气发生雪崩式击穿，在高功率微波脉冲传播途径上会形成等离子体区域。 脉冲能量被部分吸收，转化成等离子体区的场能和粒子的热能及动能。电离过程进 一步发展，自由电子密度随之增加，从而使等离子体频率增加。当等离子体频率增 至与高功率微波频率相当时，由于受到等离子体的屏蔽作用，脉冲无法透过等离子 体区域继续传播，电磁脉冲将会被截断。由于等离子体形成需要一定时间，在等离 子体区发展到完全屏蔽高功率微波前，高功率微波前沿已经穿过了等离子体区，有 部分微波脉冲能量透过大气，耦合进入目标，所以脉冲的前沿衰减较少；等离子区 域完全形成后，脉冲无法透过，所以尾部衰减较多。脉冲波形结构发生变化，脉冲 变窄，出现尾蚀现象【9 】。 发生击穿后，由于部分能量被吸收和反射，电磁脉冲出现衰减，到达目标的脉 冲能量减小，电磁脉冲的有效传播距离与效能的发挥受到限制。由于电离击穿效应 在电磁脉冲传播上的巨大影响，国内外较早就开展了相应的研究。美国在上世纪七 八十年代已经开展了这方面的研究【7 】- 【9 】。 第一章电磁脉冲大气传播的研究进展 周l o 尾蚀现象 jhy e e 等人通过对一维平面电子流体方程和m a x w e l l 方程进行数值模拟给出 了大气击穿阈值随压强变化的规律以及h p m 大气传输的物理图像【7 】。s p k u o 等人 采用电磁能量方程、电子密度方程以及半经验的空气离化率公式，对大气传输进行 了数值模拟，得到了一组可以模拟脉冲长距离传输的方程 9 ，并研究了击穿闽值电 压和大气压强的关系。国内也有不少单位从事脉冲的大气击穿和能量传输等方面的 研究【1 0 】1 15 】。 通过分子动力学等方法对电子与气体分子以及离子的碰撞频率进行分析可咀 第一牵电磁脉冲大气传播的研究进展 导出电离状态下的大气介电常数和电导率 1 4 】，从而为电磁脉冲大气传播的数值计 算提供相关的数据。 1 2 3 其他影响电磁脉冲大气传播的因素 电磁波在大气中传播时，除了因大气电离和击穿而产生能量损失外，某些特定 频率的电磁波也会被大气中的气体分子也会吸收，传播受到影响。同时，电磁波也 会因为大气沉积物( 雨、云雾、冰雪和沙尘等) 的吸收和散射而产生衰减。 ( a ) 气体分子的吸收作用对电磁脉冲传播的影响 在晴朗天气下，大气分子吸收是影响电磁波传播的最主要因素。大气中的气 体分子在电磁波的作用下，会吸收电磁波的能量，从能量较低的能级跃迁到能量高 的能级，吸收的电磁波的能量等于跃迁前后能级的能量差。所以，气体分子对电磁 脉冲能量的吸收主要集中在某些特定频率。在通常研究的电磁波范围内，对电磁波 起吸收作用的主要是氧分子和水蒸气分子，其他分子的吸收作用很小，一般可以忽 略。 氧分子在微波段主要有两个吸收区，一个是波长2 5 3 m m 的单谱线吸收线， 另一个以波长5 m m 为中心，由若干吸收线组成。即频率在1 1 9 g h z 或者6 0 g h z 附 近的电磁波会被氧分子强烈吸收。 水蒸气分子的微波波段也有两个吸收区，中心波长分别为1 3 4 8 c m 和0 1 6 4 c m ， 对应的频率分别为2 2 g h z 和1 8 3 g h z 。此外在亚毫米和远红外波段还有一些吸收线。 经验公式表i l f j 2 ，氧气和水蒸气对微波的吸收与压强和温度都有依赖关系， 通常，在海拔5 k m 以下，1 0 3 5 g h z 波段内，水蒸气的吸收系数比氧气高约一个数 量级，此时的微波吸收主要由水蒸气引起。 ( b ) 雨对电磁脉冲传播的影响 对于1 0 1 0 0 g h z 的电磁波，降雨也会带来较大的衰减。降雨衰减率与降雨强度 有关，有拟合公式 1 7 】： 口，= 口霹 其中，孵是降雨衰减率，单位是d b k m ， 由实验可以给出经验关系： ( 1 4 ) b 为降雨强度。参数口，b 与频率相关， 第一章电磁脉冲大气传播的研究进展 a = 3 x ( f - 2 ) 2 - 2 x ( f - 2 ) x 1 0 6 = 【1 1 4 0 0 7 ( y 一2 ) 7 3 】 1 0 0 8 5 ( f 一3 5 ) x e x p ( 0 0 6 f 2 ) 】 可见降雨带来的衰减与降水量和频率有关。 降水强度的大小会影响到脉冲的形状，可能会导致失真。实验表明，对于 3 1 a s 的脉冲，当降水强度为2 5 m m h r 时，会导致脉冲尾部出现带尾态，当降水 强度达到2 5 m m h r ，脉宽展宽为1 0 郎，脉冲严重变形。此外，电磁波发射角的 改变也会使得降雨对电磁波衰减率变化。 ( 1 5 ) ( 1 6 ) ( c ) 云雾和冰雪对电磁脉冲传播的影响 云和雾都是由空气中的细小水滴形成的天气现象，所以，它们对电磁衰减率的 计算公式相同，都依赖于云雾中的含水量，并且不同频率的电磁波衰减系数不相同。 在一1 0 2 0 0 c 温度范围内，雨雾的衰减率经验公式为 1 6 ： ：4 “坐唑坐嘿掣幽m ( 1 7 ) 其中t 是雨雾中的温度，单位为k ，m 为含水量，单位为g m 3 ，名为波长，单位为 c m o 冰雪由水的结晶形成。由于雪花的形状，介电常数和尺寸分布依气象条件变化， 并且，不同的雪引起的衰减相差也很大，干雪与湿雪引起的衰减可能会相差数倍。 因此降雪引起的微波衰减难以分析。实验数据表明，相同降水量下，不含水的降雪 引起的微波衰减比降雨引起的衰减要小一个数量级，而湿雪或含水雪引起的衰减则 与相同降水量的降雨相当，甚至超过降雨。 ( d ) 沙尘对电磁脉冲传播的影响 沙暴是由大量沙粒在风力作用下沿地表飞散形成的。沙粒的半径一般在8 0 1 5 0 微米。尘暴则是由更小的固体微粒( 一般半径小于3 0 微米) 形成，尘暴时，固体颗 粒在风力作用下，可上升到几百至上千米高空，遮蔽太阳。沙暴和尘暴会出现在许 多地区，影响范围可达数百平方公里，持续时间可长达数小时。 电磁波在沙尘暴中的衰减系数依赖于两个因素：( 1 ) 沙尘系统的介电常数，( 2 ) 沙尘粒子的形状分布( 粒子的半径概率密度) 。 s i g h o b r i a l 和s m s h a r i e f 等人研究了电磁波在尘暴中的传播情况 1 8 1 。他们认 为尘暴随机系统的介电常数可以写成 6 第一章电磁脉冲大气传播的研究迸展 岛= - + v 丽6 - 6 0 c z = 1 ，2 ，3 ， c 。8 ， 式中占是尘暴的介电常数，是背景介质的介电常数，1 ，是沙尘暴所占的相对体积， 1 ，= 警，是能见度，单位为k m ，4 是由实验确定的系数。通过对介电常数 的测量，得出两个定性规律： 1 干燥地面电介质的实数部分和虚数部分随频率变化较小。 2 湿度较大时，介电常数与频率相关，实数部分和虚数部分变化不一致。 一般认为尘暴中粒子呈现椭球形状，平均三轴比为1 ：0 7 1 ：0 5 3 ，最长轴小于 6 0 微米。相应的研究表n j j 1 9 ，如果假定粒子大小一致，半径均为1 0 0 微米，在湿 度较高、能见度为1 0 0 m 的沙暴中，3 7 g h z 的电磁波衰减率为1 7 d b k m 。如果假定 粒子的形状服从典型分布，湿度较大时，3 7 g h z 频率电磁波的衰减率变为l d b k m 。 粒子的形状分布对电磁波的衰减有一定影响，湿度的影响更为显著，在低湿度条件 下，衰减效应可以忽略。 1 3 电磁脉冲大气传播的相关研究 上一节介绍了影响电磁脉冲大气传播的因素，并相应介绍了部分电磁脉冲大气 传播的研究情况。下面简单介绍一下电磁脉冲大气传播其他方面的相关研究。 自由空间中的电磁脉冲大气传播受实验条件的影响很大。早期，受微波源功率 的限制，大气传播的相关研究( 例如大气击穿等实验) 都是在谐振腔中进行。由于 微波在谐振腔内是驻波状态，所以这种方法不能准确地反映真实空间微波的传输条 件。在波导内，由于电磁波是以行波形式传播，并且对微波源功率要求较低，容易 准确测量，所以近些年来利用波导管开展了电磁脉冲大气传播的相关实验，特别是 大气击穿方面的实验 2 0 【2 l 】。 波导由金属壁围成，管壁会对实验产生一定影响，产生壁效应，其影响程度依 赖与电子在气体中平均自由程与波导尺寸的关系，平均自由程远小于波导尺寸时， 壁效应的影响可以忽略。电子在气体中自由程可以近似写成： 4 = 1 ( p m o - ) 7 ( 1 9 ) 第一章电磁脉冲大气传播的研究进展 其中，成为气体分子密度，= 3 5 4 1 0 uxp ，单位为c m 一，仃为电子与 气体分子的碰撞界面，与电子能量有关，当电子能量在几个到几十个电子伏范围内 时，可以近似为仃1 0 1 5 c m 2 。电子的自由程可以近似写成：以= 0 0 2 8 p ， 单位为c m 。波导尺寸一般为厘米量级，因此，除了在高真空条件下之外，一般的 实验条件，管壁的影响都- 7d a 忽略。 上图分别列出波导和自由空间中的大气击穿实验数据，可以看出，波导中得到的 实验数据与自由空间中的数据以及理论解基本一致，说明波导击穿实验是一种正确 的方法，可以得到可信数据。 除了实验研究外，在数值模拟方面，电磁脉冲大气传播的研究方法主要有抛物方 程和射线追踪法，胡绘斌、毛钧杰等人利用抛物方程方法计算了电磁不同折射率条 件下的传输损耗【4 】。官莉、顾松山等人利用射线追踪法，研究了大气波导的形成及 传播路径 2 2 】。但是采用这些方法难以得到电磁脉冲具体的场分布情况。 时域有限差分法( f i n i t e - d i f f e r e n c et i m e d o m a i nm e t h o d ，f d t d ) 是- - 种在计算电 磁学中得到广泛应用的数值计算方法。f d t d 方法具有优良的稳定性，收敛性，快 速性，所以，本文选择了f d t d 方法研究电磁脉冲在大气中的传播。采用f d t d 方 法，可以计算脉冲产生的电磁场在空间中分布及其推进演化。与抛物方程和射线追 踪法等方法相比，采用f d t d 方法可以得到电磁脉冲中电磁场强度的具体分布情况。 第一章电磁脉冲大气传播的研究进展 在f d t d 方法的计算过程中，需要输入计算区域内的许多参数例如：折射率、 电导率和介电常数等。利用抛物方程等方法可以解析求出电磁脉冲传播路径经过的 区域，电磁脉冲传播路径上的各自电磁参数可以由计算得到，或者通过实验测得， 在本章中第二节中介绍了相应的参数研究情况。传播路径上电磁参数的确定，是能 够采用f d t d 方法研究电磁脉冲传播的重要基础。 1 4 本文的研究任务 电磁脉冲的大气传播一直是人们很关心的问题，本课题研究了移动坐标系 f d t d 方法在电磁脉冲大气传播中的应用。主要工作包括两方面的内容，首先验证 了一维移动坐标系f d t d 方法的可行性，其次将移动坐标系f d t d 方法推广到三维， 并利用该方法研究了波导中脉冲的传播。 本篇论文分为四部分，第一部分是课题的背景，介绍了国内外电磁脉冲大气传 播的研究情况，介绍了影响电磁脉冲传播的相关因素；第二部分介绍移动坐标系 f d t d 方法及其数值色散和稳定性分析；第三部分是移动坐标系f d t d 方法的计算 实例，包括一维和三维情况下的例子。第四部分是论文的总结。 9 第二章移动坐标系f d t d 方法 第二章移动坐标系f d t d 方法 本章主要介绍移动坐标系f d t d 方法的基本原理，推导在f d t d 方法中引入移 动坐标系后相应的差分公式。在引入移动坐标系前，首先简单介绍一下f d t d 方法。 2 1f d t d 方法 时域有限差分法( f i n i t e d i f f e r e n c et i m e d o m a i nm e t h o d ，f d t d ) 是一种在计算电 磁学中得到广泛应用的数值计算方法，它由k s y e e 子1 9 6 6 年提l k l 2 3 。f d t d 方 法通过对电磁场的e 和h 分量在空间和时间上交替抽样进行离散，将麦克斯韦方程 转化为一组差分方程，并在时间轴上逐步推进，对空间电磁场进行求解。f d t d 方 法具有很好的稳定性和收敛性。f d t d 方法可以给出计算区域内的电磁场随时间演 化的过程，能够得到计算区域内电磁场的具体分布。随着电子计算机容量和速度的 飞速发展，f d t d 方法也得到了不断的发展和改进【2 4 】- 【2 6 】。吸收边界的引入 2 7 【2 9 】，使得f d t d 方法可以处理开放区域问题，大大拓展了f d t d 方法的应用 范围。f d t d 方法逐渐成为一种在多种领域都得到了广泛应用的标准算法，几乎可 以覆盖所有涉及电磁场问题的领域【2 4 】- 【2 6 】，如微波器件设计分析、电磁兼容问题、 电磁波传播与散射以及天线设计等。 2 1 1f d t d 方法的差分格式 麦克斯韦方程是描述电磁场传播的基本方程，f d t d 方法是从微分形式的麦克斯 韦旋度方程出发，进行差分离散，从而得到一组时域上的推进公式。 对于任意介质，麦克斯韦方程的旋度方程可以写成 v 一e ：一o b 一石 一o t ( 2 1 ) v x 耳：丝+ 了 l o 第二章移动坐标系f d t d 方法 式中面是电场强度，单位为v m d 是电通量密度，单位为c m 2 再是磁场强度，单位为a m 面是磁通量密度，单位为w b m 2 ，是电流量密度，单位为a m 2 石是磁流量密度，单位为v m 2 在各向同性介质中，本构关系为 至= 占三1 b = , u hi 三础一f 厶= d 。hj c 2 ) 其中，s 为介电常数，是磁导系数，盯为电导率，o m 为磁导率。仃和描述介 质的电损耗和磁损耗，在真空中两者均为0 。 在直角坐标系中，麦克斯韦旋度方程可以写成 c 等一争= 一警一风( ) v0 z0 l f a e x l 一 c i 2 f a e y l i c 等一警鲁+ 盯晟o vo z o l 翠a z 0 咖- 2 ) = 占鲁蚋a xa f ( 孕d x 一等等川 d 1 ，a r 采用y e e 网格离散后，e 和风的f d t d 差分格式可写成 ( 2 3 ) ( 2 4 ) ， z q 也 m m 吒 y 一 三 研一钟甜一钟 上 工 一 一 i i = ) ) 堡舐堡砂 第二章移动坐标系f d t d 方法 e + 1 ( f + 圭，j ，七) = c a ( m ) 群( “壬，j ，k ) + ( 咖型塾盟型喾竺噬止坦 il ， 郦卅圮( “ ，j + 如一坨( “ ，工七一圭) l 一一 心 j 月：佗( f ，+ ，| j + ) = c p ( ，行) 月：_ 1 胆( f + 圭，j + i ，七 嗍吣l 盟业等盟幽 髟( f ，+ ，七+ 1 ) 一髟( f ，- ，+ 士，i ) j 一一。 & i ( 2 6 ) 其中的下标i 、j 和k 标记h 和e 的空间位置，i - - _ 标n 标记h 和e 取值时间。式中 的系数分别为： 叫啪= 焉o ( m ) a t ， 1 c p ( m ) = 1 一等等 1 + 犁l , l a 筹l m j l 叫砌= 痞 c q ( 砌= 碡 ( 2 7 ) ( 2 8 ) 其他分量的f d t d 差分格式可以同样得到，文献【2 4 】中有各个场分量的具体表达式。 在一维情况下，可以简化为 ( 七) ：酬( 肌) e ( j j ) 一c b ( m ) 卫t t n + _ 1 1 2 i r , y 产r r n _ + 1 1 2 t 型l _ l 、 矿) = c p ( 吣珊咿c q ( 朋) 趔等塑 其中的系数c 么，c b ，c p ，与前文相同。在真空中，一维的差分格式可写成 e + ( 七) = 露( 七) 一老 彬州2 ( 七+ 圭) 一点1 7 2 ( 七一专) 郦川2 ( 七+ ) = 彬( | | + 上) 一盖 譬( 七+ 1 ) 一霹( _ | ) 2 1 2f d t d 方法中的约束条件 ( 2 9 ) ( 2 1 0 ) ( 2 1 1 ) ( 2 。1 2 ) 由于f d t d 方法是一种差分方法，在计算中为了保证数值不发散，需要满足一 1 2 第二章移动坐标系f d t d 方法 定的约束条件，其中最重要的是c o u r a n 稳定性条件( c o u r a n t - f r i e d r i c h l e w y ，c f l 稳定性条件) 以悠孺焉面 c f l 条件给出了空间和时间离散间隔之间应当满足的关系。 一维情况下，c f l 条件化为 c a t a x ，2 1 3 ) ( 2 1 4 ) 上式指出，一维情况下，时间间隔必须等于或者小于波以光速通过一维y e e 网 格的时间。 除了稳定性条件外，在f d t d 方法中由于用有限差分式代替了波动方程中 的导数项，使得电磁波的波数k 与角速度国之间依赖关系发生改变，由此产生 色散，称之为数值色散。数值色散是由对波动方程进行差分近似引起的，即使 介质本身是无色散的，离散之后也会带来数值色散。这种现象会对时域数值计 算带来误差，为了减小这一误差，需要限制空间离散间隔的取值。研究表明 2 4 】， 对于一维平面波，当空间间隔血与波长允满足 g x 三 1 2 ( 2 1 5 ) 数值色散带来的误差可以忽略。在模拟时域脉冲信号时，由于信号不是单色波， 可以采用信号带宽中对应的上限频率的波长丸访来代替波长允。 时间离散也需要同样的约束条件，时间间隔出与空间脉冲的周期r 要求满 足 岔i t ( 2 1 6 ) 1 2 、7 同样，对应非单色波情形时的时域脉冲信号，采用瓦i 。代替。 2 2f d t d 方法研究电磁脉冲大气传播面临的困难 由于f d t d 方法优良的稳定性、收敛性以及快速性，所以本文选择了f d t d 方 第二章移动坐标系f d t d 方法 法来研究电磁脉冲在大气中的传播，计算脉冲产生的电磁场在空间中的分布及其推 进演化。 由于f d t d 方法是一种全波方法，计算时需要对整个空间进行网格划分和计算 推进，在计算脉冲波包的长距离传播时，如果采用一般的f d t d 方法，对整个传播 路径划分网格，需要占用大量的计算机内存资源，也会使得计算需要的计算机时变 得巨大，在数值色散性要求下，也很难通过增大网格的间距来减少网格数目。由于 电磁脉冲的特性，脉冲能量集中于脉冲波包所在位置附近，传播路径上的其他地方 没有电磁场或电磁场非常弱，只在脉冲波包附近有较强的电磁场。如果采用一般的 f d t d 方法，对整个空间划分网格进行计算推进，由于脉冲能量分布的局域性，大 量内存和计算被用于没有电磁场的区域，这是一个极大的浪费。随着传播路径的增 加，计算所需要的内存和计算机时也会增加，当路径足够长时，传统的f d t d 方法 已经不能处理脉冲的传播问题。 2 3 移动坐标系f d t d 方法 2 3 1 移动坐标系f d t d 方法原理介绍 由于脉冲的局域性质，电磁脉冲在传播过程中能量集中分布在传播路径上的一 个很小区域内( 即脉冲波包所在的位置附近) ，与波包处相比，传播路径上其他位置 上的电磁场强度可以忽略不计。计算过程中，可以只考虑波包附近区域的电磁场， 忽略其他位置的电磁场，这样的处理可以极大降低计算量和需要的计算机内存，同 时不会带来很大的误差。 为了解决f d t d 方法在处理脉冲长距离传播问题时面临的困难，移动坐标系 f d t d 方法被提出并被应用于研究脉冲的长距离传播【3 2 】。这种方法的提出使得脉 冲长距离传播的求解成为n f l 皂 3 9 】，因为这种方法不像传统的f d t d 方法那样需要 对全部感兴趣的空间划分网格，进行离散，引入移动坐标系之后只需要对传播路径 上包含波包的那部分感兴趣区域划分网格，在f d t d 方法中引入移动坐标系后，计 算网格将会伴随电磁脉冲一同传播，只需要一个能包括脉冲波包的计算网格就可以 1 4 第二章移动坐标系f d t d 方法 模拟脉冲在整个路径上的传播情况，减小了对计算机内存的需求。由于每个时间步 需要计算更新的网格数目大为减少，也有效减少了计算所需时问。并且移动坐标系 下的波速可以被解析地求出，所以可以有效减小数值色散。 移动坐标系f d t d 方法的计算步骡【3 9 】： 1 构建一个足够大可以容纳脉冲的f d t d 网格。 2 用时域有限差分方法计算麦克斯韦方程组。 3 在每个时间步，通过确定脉冲中心的位置，确定脉冲在传播路径所处的位置。 根据传播路径上介质参数的分布情况在计算网格的前边界上加上对应的介质参数， 与此同时，在后边界上舍弃电磁波已经传过区域的介质参数，实现计算网格伴随脉 冲传播一起前进。 圈2 - 1f d t d 目格伴随着脉h 在传播略径上咎动 2 3 2 移动坐标系f d t d 方法研究进展 移动坐标系方法在计算物理的很多领域都有广泛的应用。在f d t d 方法中引入 移动坐标系的研究开始于上世纪九十年代3 1 1 。b f i d e l 、e h e y m a n 和r k a s m er ，等 人在i e e e 的一次会议上提出将f d t d 方法与移动坐标系相结合进行了一维情况 下的计算，利用该方法求解了一维情况脉冲长距离传播的问题。文献1 3 2 中详细分 析了一维情况的移动坐标系f d t d 方法的色散关系和数值稳定性。他们利用该方法 计算了脉冲传播1 0 0 0 个脉冲宽度的长距离传播情形，并将计算结果与一股f d t d 方法的计算结果相比较，验证了该方法的可行性和优点。e a k l e m a n 和i s e v g i 等进 第二章移动坐标系f d t d 方法 步研究了二维情况下该方法的应用 3 5 3 6 1 。 除了在电磁波传播方面得到广泛应用外，移动坐标系f d t d 方法也被应用于其 他方面的研究，例如声学波的传播【3 7 】【3 8 】，p e m p e r 等人利用该方法研究了声学波 脉冲在声波导中传播，计算了传播脉冲长距离后的波形，与理论值的比较，该方法 可以得到可信的数据。他们还迸一步编写了采用移动坐标系f d t d 方法计算声学波 脉冲传播的三维程序。 2 4 移动坐标系f d t d 差分格式：一维情形 2 4 1 非均匀介质中f d t d 方法的稳定性和色散关系 前文中介绍了移动坐标系f d t d 的基本原理，下面将推导引入移动坐标系后的 差分格式。为了便于理解，下文中研究的是在均匀无损耗介质中传播的例子。 一维情况下，在无损耗介质中，电磁场可以用下面的方程描述 堡= 一檠竽，丝= 一c ( z ) z o ) o _ h y 一= 一一一一= 一f ：i z ，厶iz o t z ( z 、8 z j o t 、。、 移z 其中的c ( z ) 是波速，z ( z ) 是波阻抗。 c ( z ) = 1 再，z ( z ) = 万 ( 2 1 7 ) ( 2 1 8 ) 一维情况下，可以将h 和e 的标记取向的下标x ，y 去掉。按照y e e 网格进行离散， 得到的差分格式为 hn + 1 + l 2 2 = _ i n + 1 ，- 1 7 2 2 一五a t c l + l 2 ( e i n 1 一霹) 霹+ 1 = g 一五tc i z i 1 _ l n + l 2 1 h i n - 1 + l 2 2 ) 式中e 和h 的上标n 标记时间步，下标f 标记空间网格。假定h 和e 是空间谐变 量，上式可以写成如下形式 第二章移动坐标系f d t d 方法 硝警= 弼盖陀e f = 弘喝 ( 2 2 0 ) 其中的七是平面波的波矢，花体字母鬈一和孑行表示电场和磁场的空间缓变部分， 满足 2 s i n 尬2 = 2 s i n k 2 2 s i n k a z 2 = 2 s i n k 2 , 兹州+ n + ：l 2 = 髯麓坨一f 2 竽竽孑8s i n 足2 a z 么“1 2 孑斛l = 孑”一f 2 塞c f 互兹旌门s i n k 2 2 2 2 中的两式合并后，可以得到 1 = 纠1 _ 4 ( 各去s i n 2 肌h 2 石a t 嘱兹拶s i 皑2 ( 2 2 1 ) ( 2 2 2 ) ( 2 2 3 ) 由于波速q 和阻抗z i 都是空间缓变量，所以采用近似q + l ，：c ，和z ，+ l ，2 z ，乙后， 可以将2 2 2 式以矩阵的形式写出 举h 髻彪 g 柳 t 是传递矩阵 z ：l 1 2 乃虿s i ，n k 2l ( 2 2 5 ) 一i - i 2 y t z ts i n k 2l 一4 y s i n 2 k 2 p 叫 其中扎= c , a t l a z 。 将传递矩阵对角化，可以得到矩阵的本征值 g = 1 2 疗s i n 2 ( k 2 ) f 2 以s i n ( k 2 ) 、j i 万i 丽 ( 2 2 6 ) 为了让方程能够稳定，g 的模必须要小于等于1 ，该条件只有在形s i n 2k 2 l 时 才能满足，此时g 的模等于1 。该条件对所有k 都成立，即可得到不均匀介质下的 c f l 稳定性条件 y ( 力量c ( z 石a t 尝 l ( 2 2 7 ) 1 7 第二章移动坐标系f d t d 方法 上式中的是介质中的最大波速。 进一步分析数值色散性质，假设h 和e 在时间和空间上均为谐变量，即 麟1 1 2 2 - - 一形。m n + 1 ，1 2 2e - i ( i ，2 一锄+ l 2 f = 彦一场) 2 瑚 其中的形? 和矿? 是缓变部分，满足 孑k 一者： 2 s i n k 2 ， 2 s i n o 2 ( 2 2 9 ) 式中的k = k a z ，q = c o a t 。 甚叫讣 s i n 2f ) 2 = y 2 ( z ) s i n 2k 2 由群速度的定义式，就可以得到群速度的表达式： (