（地球探测与信息技术专业论文）地质雷达二维时域有限差分正演.pdf

! 堕叁堂婴主堡苎一一j 堕 摘要 地质雷达是在高频微电子技术以及计算技术迅速发展的2 0 世纪后期，此项 技术才获得根本性的进展。由于它是一种非破坏性的原位探测技术，抗电磁干扰 能力强，故可在各种噪声环境下工作，环境干扰小。同时地质雷达还具有工程上 较满意的探测深度和分辨率，现场能直接提供实时剖面记录，图像清晰、直观； 再上使用便携式微机控制数据采集、存储和处理，工作效率高，重复性好。地质 雷达具有如此多的优点，故应用十分广泛；但是目前地质雷达亦还有待改进之处， 如雷达图像的识别主要还停留在人工的肉眼和手工的基础上；所以对常见探测对 象的物理正演模拟，有助于我们了解电磁波的传播规律，挖掘雷达波中更多有用 信息。 作者在本文中首先以麦克斯韦方程组为出发点，着重于对麦克斯韦的两个旋 度方程的分析，运用k s y 雎氏的空间网格模型，并利用中心差商代替微商，把 空间连续变量离散化，把电场进行规约化，得出二维空间的时域有限差分方程， 即我们所要推导的地质雷达正演模拟方程：为了保证算法的精确性，我们对二维 空间的t m 电磁波( 雷达波) 进行数值稳定性分析及数值频散问题的讨论，推导 出了数值稳定性条件表达式和理想频散关系；同时为了防止雷达模拟模型截断边 界引起波的超强反射，我们推导出了超吸收边界条件，并把该条件运用到自制的 用c 语言编写的程序中，取得到了良好的效果。 在本论文中我们还讨论了地质雷达的 d z t 文件格式，并着重分析了$ d z t 的头文件格式，在清楚地了解 d z t 文件格式的基础上，我们把头文件和用上 面自制的程序模拟得到的数据用c 语言的库函数f w r i t e 写入到 d z t 雷达格式 文件中，再将d z t 雷达格式文件导入到随机配套的雷达处理软件中，得到了 常见的几种规则模型的波形图。模拟结果的波形与实际探测的波形非常吻合，证 明了用时域有限差分法在地质雷达正演模拟方面具有良好效果并且用时域有限 差分法模拟地质雷达模型的正演模拟结果对工程实践亦具有重要的指导意义，再 加之时域有限差分法计算程序的对电磁场的广泛通用性，对不同的问题或不同的 计算对象只需修改有关的少部分便可运用，故我们编制的程序对电磁场的模拟工 作者亦有一定的参考价值。 关键词：地质雷达，时域有限差分法，正演模拟，超吸收边界条件 ! 堕查兰堡主堡苎一竺翌墨! ! 三 a b s t r a c t i np e r i o d1 a t eo f2 0c e m u r y ，t h et e c 上1 1 1 0 i o g yo fg r o u n dp e n e t r a t i n g r a d a r a c q u i r i n gb a s i c a l l yd e v e l o p m e n t ，j u s t a r e rt h e h i g hf e q u e n c y m i c r o e l e c t r o n i c s t e c h n i q u e a n dc a l c u l a t i o n t e c h n i q u eg e tq u i c k d e v e l o p m e n t b e c a u s ei t i sae x p l o r a t i o nt e c h n i q u eo fn o tb r e a k i n gt l e o r 垮n a ip r o p e r y ，a n di t sa b i l i t yo f 锄t i e l e c t r o m a g n e t i s mi m e r f b r e n c ei s s t r o n g ，t l e r e f o r e i tc a nb ew o r kw e l lu n d e rt h ee n v i r o m e n to ft h ea l l k l n d so fn o i s ei n t e r f e r e n c e ，也ee n v i r o n m e n ti n t e r f e r e n c oi ss f n a l l ，t o o ，k p o s s e s s e s m es a t i s f a c t o r y e x p l o r a t i o nd e p t h a n d d i s t i n g u i s h r a t ei n e n g i n e e r i n g ，c a no f 亿r i n g t 1 1 er e a lt i m es e c t i o nr e c o r di nt h e1 0 c a l e d i r e c t l y t h e d i a g r a m c l e a r s ， i n m i t i o n i s t t h ec o n v e n i e n t s c h l 印p i n g m i c r o c o m p u t e rc o n 订o l st h ed a t ac 0 1 】e c t e d ，s t o r a g ea 1 1 dh a n d l e ，a n dt h e w o r ke f f i c i e n c yi s h i g h ，r e p e t i t i o n i s g o o d g r o u n dp e n e t r a t i n gr a d a r ( g p r ) h a v es om a n ya d v a n t a g e ，m e r e f o r e 印p l i e dv e r ye x t e n s i v e b u t n o wi ts t e x i s ts o m e p r o p o r t i o n t o w a h i n g f o r i m p r o v e m e n t ，f o r e x a m p l e ，d i s t i n g u i s h i n g t 1 1 eg p r d i a 黟锄p r i m a r i l ys t a y s i nt h e f o u n d a t i o no fm ea r t i f i c i a ln a k e de y ed i s t i n g u i s ha n dh a n d c r a f ti d e n t i 一 t h e r e f o r ei ti sb e n e f i c i a lt 0u st ow l d e r s t a n dt i ee l e c t r o m a g n e t i cw a v e s p r e a d i n gr e g u l a t i o n ，a j l ds c o o po u tm u c hm o r eu s e m li n f o 肿a t i o nf r o m t 1 1 eg p r 、v a v ei n s i d e ，t h r o u 曲f o n v a r ds i m u l a t em e f a m i l i a rp e n e t r a t i n g o b j e c ti np h y s i c a l a t 丘r s t ，i nt h i s 渤【ta u t h o rs e to u t 8 n 越y s i s 蠹o mm a x w e ；le q u a t i o n g r o u p ，a n dp u tg r e a te m p h a s i so nm et w ov o r t i c e se q u a t i o no fm a ) ( w e e q u a t i o ng m u p a p p l i c a t i o nk s y e e ss p a c eg r i dm o d e l ，a 1 1 dt h eu t i l i z e c e n t e rd j 腩r e n c ec o e 佑c i e n ts u b s t i t u t ef o rt h ed i f r r e n t i a ic o e m c i e n t w e c a l l d i s p e r s et h es p a t i a lc o n t i n u o u sv a r i a b l e ，r e g u l 撕z a t i o nt h ee l e c t r i c f i e l d ，a n do b t a i nt h et w od i m e n s i o ns p a t i a jt j i n e d o m a i nf i n i t e d i f 话r e n c e e q u a t i o n ，w h i c h i st h eg p rf o r w a r d s i m u l a t e e q u a t i o nw er e q u e s t d e d u c i n g f o f g u a r a n t e e i n g 也ea r 三l h l e t i c a c c 硅r a t e ，w ea n a l v z et h e n u m e r i c a lv a l u e s s t a b i l i 哆a 1 1 dd i s c u s s i o nq u e s t i o no fn u m e r i c a lv a l u e n e q u e n c yd j s p e r s l o nt o m et w o d i m e n s i o n1 m e l e c t r o m a g n e t i cw a v e ( g p rw a v e ) ，d e d u c eo u tn u m e r i c a lv a i u e ss t a b i i 时e x p r e s s i o na n di d e a i f r e q u e n c yd i s p e r s i o nr e l a t i o n a tt h es 柚et i m ei no r d e rt op r e v e n tt h e r a d a rs l m u l a t l o nm o d e lw a v e a tt h et r u n c a t e d b o u n d a us t r o n 置r e 日e c t i 咖 ! 堕奎堂堡主堡苎 塑塑呈! ! ! 、ed e d u c eo u tm es u p e ra b s o r bb o u n d a 哆c o n d i t i o n ，a n d m a k eu s eo tt h e c o n d i t i o ni n t ot h es e l f r e g u l a t i n gp r o g r 锄w r i t eb yc1 a n g u a g e ，a n d o b t a i n e dt h e 旦o o de 舒皓c t i nt h i sm e s i sw ed i s c u s sg p r + d z tf i l ef d r m a t ，a n dp m g r e a t e m p h a s i sp o i n to na n a l y z er a d a r 枣d z t h e a d6 1 ef o h n a t o nb a s i co f c l e a r l vu n d e r s t a l l dt h e d z tf i l ef o 衄a t 1 1 1 eh e a df i l ea i l dt h ed a t ao b t a i n b ym n n i n ga n ds i m u l a t i n gt 1 1 es e j f r e g u l a t i n gp r o g r 锄a r eb e e n 、v r i t t e n i n t o 也ef i l eo f 车d z tr a d a rf o m l a tu s i n gt h em n c t i o n - b a s e 舳r i t ei nc 1 a n g u a g e a f t e n v 盯dw et r a n s m i t 恤e 母d z tf i l e o fr a d a rf o n l l a ti n t ot h e r a d a rd i s p o s es o r w a r et h a tt h em a c h i n ee q u i p ，a 1 1 dg e tm ew a v e d i a g r a m 行o ms e v e r a lf a m i i i a rf o 衄u l am o d e l t h ew a v e so fi m i t a t er e s u l ta n dt h e w a v e so fp r a c t i c a ip e n e t r a t e 礅w 以，l t 弘c r v e dt h a tt h e 矗n h e d l f f e r e n c e t i m e d o m a i n ( f d t d ) m e t h o dc a nr e a c hg o o de f 诧c ti n ( 狰rf o n v a r d s i m u l a t e ，b e s i d e st h er e s u l to fu s i n g 矗n i t e d i f f e r e n c et i m e d o m a i nf f d t d ) m e t h o df o n v a r ds i m u l a t et h eg p rm o d e lh a v e i m p o r t a n tg u i d a n c e m e a n i n g f o rt h e e n g i n e e r i n gp r a c t i c e 。a d d i t i o n a l l y 6 n i t e d i 恐r e n c e t i m e - d o m a i n ( f d t d ) m e t h o dc a l c u l a t i n gp r o g r 锄h a v em ee x t e n s i v e c u r r e n c yc h a r a c t e r i ne l e c t m m a g n e t i s mf i e l d t h ep r o g r a mc a nb ee x e r t e d o n l yn e e d i n gm d d i f i c a t i o n a1 i t t l eo fr e l e v a n t p a r t s a i ma td i f 话r e m p r o b l e m o rd i f 代r e n tc a l c u l a t i o n o b j e c t t h e r e f o r e t h ep r o g r a mo fu s a u m o r i z e dh a v eaf e wo fr e f e r e n c ev a l u et ot h ee l e c 仃o m a 2 n e t i s m 丘e l d s i m u l a t eo d e r a t o r k e yw o r d s ：g r o u l l d p e n e 仃a t i n gr a d 弛f j n i t e d i f f e r e n c et i m e d o m a i n m e t h o d ，f b 九v a r ds i m u l a t i o n ，s u p e ra b s o r bb o u n d a uc o n d i t i o n i i 硕。1 ：学位论文 j 一 日l j苗 在应用时域有限差分法对地质雷达的正演模拟之前，我们应首先了解一下地 质雷达进行探测的基本原理及相关理论。 地质雷达的基本原理是基于高频电磁波理论，工作方式是以宽频带短脉冲 tr 一1 。 电磁波形式，在地面通过发射天线向地下发 射，经地层或目的体反射后形成反射回波返 回地面，为另接收天线所接收，其示意图如 o 一1 所示。地质雷达的脉冲波双程旅行时为： z = 嘲 图。一1 反射雷达探测原理 当地下介质的波速v 为己知时，可根据精确得 到的走时f ，由上式求得目标体深度h 。式中z 值即收发距。地质雷达现场测量时，我们通常使发射天线和接收天线以固定的间 距沿测线同步移动，这种测量可反映地下介质同一深处的反射信号。其现场探测 示意图如o 一2 所示。 图o 2 雷达现场探测示意图 地质雷达的检测过程是雷达子波以宽频带短脉冲形式通过发射天线不断地 向地下介质发射的过程，雷达子波即高频电磁波，它遇到电性差异的地下地层或 硕士学位论文 前言 目标体便反射回地面，由接收天线接收。高频电磁波在传播时，其路径、波场强 度与波形将随所通过介质的电性及几何形态而变化，故通过对时域波形的采集、 处理与分析，可确定地下界面或地质体的空间位置及结构。3 。这便给我们对地质 雷达的模拟提供了理论依据。地质雷达测量方式示意图如0 3 所示。 天线 a 射波l 兀毒至黑色) 地表面 7 7 研f ? 77 、趔 酯吟芒骂 移动方向 天钱移动时反射波的记录 豸色悸 黑色至白色i 开斤一 乞- - - 自至黑厂一l 畸 l 产 墓抽l 一| ，冬肿 二= 二j ! = 丽”ir n 湾 黑至白 tp t 图0 3 地质雷达模拟示意图 地质雷达是根据发射雷达波的特征识别地下介质的性质，因此，能否获得可 分辨的地质雷达反射波是影响地质雷达应用效果的关键因素，它决定于雷达波的 强度、地质雷达仪器的分辨能力、波在地质界面上的反射特性以及波在地下介质 中行进时的衰减情况1 3 1 。而雷达波的强度以及地质雷达仪器的分辨能力由地质雷 达仪器本身决定，地质雷达波是高频电磁波，其传播可近似为平面波，在地质界 面的反射系数，为： 一酬 式中r 为反射系数，z 为界面上层介质的波阻抗，z ，为界面下层介质的波阻抗。 其定义为 扛j 篇 其中= = 了，珊= 2 万为角频率，为导磁率，c 为介质的介电系数，。为介 质的电导率( 以下同) 。由于地质雷达频率高，一般有o ，故反射系数可简 ，：邀 _ s 、+ 0 2 堡主堂堡丝茎! ! 旦 。，为介质的相对介电常数，角标r 表示波所在的介质”3 。反射系数反映了反射 能量占入射能量的比率，取决于界面上下介质的介电系数差异。即探测目标体与 周围介质的介电系数差异是地质雷达应用的物性前提，其差异的大下影响其探测 的效果。下面是常见介质的地质雷达勘探应用效果按电导率的分类： a 类：低电导率( 小于l os m ) ，如空气、干的花岗岩、灰岩、混凝土、沥 青等应用效果好的介质。 b 类：中电导率( 1 0 1 0 s m ) ，如纯净水、纯净水结成的冰、雪、砂、干粘 土、干的玄武岩、海冰等应用效果一般的介质。 c 类：高电导率( 大于l o s 【1 1 ) ，如湿粘土、湿页岩、海水等应用效果较差的 介质。 对于地质雷达的探测方式，它的分辨率也是一个必须了解的内容，地质雷达 的分辨率包括垂真分辨率和水平分辨率。地质雷达的垂直分辨率主要由地质雷达 的波长决定，从波的传播规律可知，可识别目标体的尺度一般需大于l 2 波长， 简记垂直最小可分辨的层的厚度为d m 。则它的计算式为： d 。= o s 丸= c 娩f 如 其中，f 为电磁波在真空中的传播速度。可见频率越高，介质的介电常数越大， 则屁越小，即垂直可分辨层的厚度越薄，垂直分辨率越高。 地质雷达的水平分辨率是指地质雷达在水平方向上所能分辨的最小异常体 的尺寸。根据波的干涉原理，法线反射波与第一f r e s e n e l 带外缘的反射波的旅 行差为丑2 。当反射波之间发生相长性干涉时，振幅增强；而第二f r e s e n e l 带 内的反射波发生相消性干涉时，振幅减弱，因此雷达波向下传播的区域是一个圆 锥体，其反射能量主要来自第一f r e s e n e l 带。设目标体埋深为h ，雷达波的波 长为五时，则第一f r e s e n e l 带的半径尼为： 耻f 胴+ 牙 - 由于每一雷达记录是第一f r e s e n e l 带内反射波的综合反映，因此尼是水平分辨 率的最小尺度”1 。从b 的计算公式可以看出，当目标体埋深越大，雷达波频率越 低，波长越长，则届越大，水平分辨率越低，反之，水平分辨率越高。同时， 水平分辨率与地质雷达的观测方式有关，对于剖面观测方式，采样点距越大，水 平分辨率越低，一般的，采样点距设为最小目标体水平尺度的1 5 。 通过了解地质雷达的探测原理及机制，则可以设计好地质模型，然后在地面 拟输入高频电磁脉冲波，如高斯脉冲，模拟完第一个点，待模拟数据稳定后，便 移动天线到下一个模拟点，一般发射天线移动的距离为定值，以此类推模拟出各 个点的波形，这些波形组合在一起便可得到雷达剖面，这样便可模拟出地质雷达 雯主兰焦笙苎旦堕， 的检测过程，对雷达波的正演能够有效地验证我们反演算法的正确性。 目前正演模拟技术是她质雷达理论研究的主要内容之一，也是地质雷达的研 究重点之一，其意义主要表现在两个方面“1 ：首先，通过分析地电模型的正演结 果，可以加深对地质雷达反射剖面的认识，提高解析精度；其次，利用已知模型 的正演结果进行反演，可验证反演算法的正确性。同时，正演技术的发展，也能 有效地推动地质雷达技术向前发展。正演模拟是研究高频电磁波在地下介质中传 播规律的有效途径，目前所用的正演方法主要有两类：一类是基于几何光学的射 线追踪法，另一类是以物理光学原理为基础的绕射迭加法。前者适于起伏变化较 为缓慢的，无电性突变点的双层或多层介质，后者同时考虑了波的动力学和运动 学特征，可用于合成曲率较大的目的体的雷达剖面，其结果也较为精确。中国地 质大学邓世坤”1 着重应用射线追踪法对地质雷达的正演合成的理论研究；东南大 学的沈飚”“”运用有限元方法对地质雷达的正演模拟亦取得了一定的成果；还有 中国科学院地质与地球物理研究所的底青云“；石油大学的何兵寿“o “”“”等在地 质雷达的数值模拟工作中取得了一定的成就。 作者在结合几位学者成果的基础上，对比了如时域有限差分法( f d t d ) 、有 限元法、矩量法等地质雷达具体的正演模拟方法的特点，基于时域有限差分法节 约存储空间和计算时间、计算程序的通用性强、适合并行运算、简单直观易掌握 等原因选择了时域有限差分法进行地质雷达的正演模拟。 f d t d 法是近几年发展起来的以物理光学原理为基础的高频电磁场模拟方法， 是一种较为理想的地质雷达正演模拟工具。针对地质雷达这一特定正演问题，首 要任务是解决吸收边界条件和数值频散这两个基本理论问题。本文通过分析雷达 波在y e e 氏网格中的传播规律，结合f d t d 算法的收敛性和稳定性条件，导出了 适于地质雷达正演的理想频散关系：通过理论分析，引入超吸收技术，提出了超 吸收边界条件，同时用数值试验论证了把超吸收边界条件用在地质雷达正演模拟 中，其精度能得到有效的保证。 4 堡主堂垡笙奎一塑二兰翌堑复堡叁坌塑 第一章时域有限差分法的概况 1 1 时域有限差分法的发展 f d t d 法是近几年发展起来的以物理光学原理为基础的高频电磁场模拟方法， 是一种较为理想的地质雷达正演模拟工具。下面我们便谈一谈时域有限差分法的 发展。早在1 9 6 6 年k a n es y e e 在他发表的著名论文“n u m e r i c a ls o l u t i o no f i n i t i a lb o u n d a r y v a l u ep r o b l e m si n v 0 1 v i n gm a x w e l l se q u a t i o n i n i s o t r o 口i c 矾e d i a ”中“，用后来被称为y e e 氏网格的空间离散方式，把带时间 变量的m a x w e l l 旋度方程转化为差分格式，并成功地模拟了电磁脉冲与理想导体 作用的时域响应。这就诞生了后来被称作时域有限差分法( f i n i t e d i f f e r e n c e t i 【n e d 0 1 i l a i nm e t h o d 或f d t d ) 的一种新的电磁场的时域计算方法。当然，从现 在的观点看y e e 氏在初创时所用方法还只是时域有限差分法的雏形，后来经过一 批科学家们的不断改进，经历了近2 0 年的发展才逐渐走向成熟的。对这一方法 的发展贡献最大的科学家除y e e 氏外，还有r h 0 1 1 a n d ，k s k u n z 和a t a f l o v e 等。在这近2 0 年的发展中主要解决了以下的一些问题“。 ( 1 ) 吸收边界条件的应用和不断改善。在y e e 氏的最初方法中使用的是硬 边界，把边界设置为理想导体( 令边界的切向电场和法向磁场分量为零) 。在这 种情况下只能模拟电磁脉冲在到达边界以前一段时间的电磁过程，否则由于边界 对电磁脉冲的反射，将改变电磁脉冲的传播路径，而不再是电磁脉冲在自由空间 中的实际散射过程。为了用有限的计算网格空间模拟无限大的物理空间，就要设 法消除电磁波在网格空间边界上的反射，也就是要能吸收掉到达边界上的电磁 波。为吸收电磁波采用了硬边界和软边界两种方法。所谓硬边界就是在边界处设 置一定的吸收材料，使传播到边界的电磁波受到衰减，从而减小反射。这方法 由于存在种种缺点而没有得到发展。所谓软边界是在边界处让电磁波满足定的 吸收边界条件，以消除电磁波从边界上的反射。曾提出过多种吸收边界条件的方 案，但得到广泛应用的是b e n g q u i s t 和a m a j d a 提出的单向波方程，后来g m u r 给出了一阶和二阶近似的差分形式，更促进了这一吸收边界条件的推广。 ( 2 ) 总场区和散射场区的划分。在y e e 氏的初始算法中没有场区的划分， 入射波和散射波只能靠时一空上的分离。这只有对窄脉冲才能实现，而且入射波 只能沿某个坐标轴传播。发展后的时域有限差分法利用连接条件把计算网格空间 分为内部的总场区和外部的散射场区，吸收边界条件用在散射场区。这样做的结 果带来了以下的好处： 堡主兰垡垒塞 篁二童堕蔓查堡董坌i ! 塑翌 ( i ) 可实现任意入射波的设计。在这种方法中入射平面波是在单独的一维 空间中传播，通过连接条件接入总场区，因此可以实现对入射波的入射方向、极 化角度及其随时间变化规律的任意设置，这大大扩展了这方法的适用范围。 ( i i ) 可以保证宽的计算动态范围。在某些结构的阴影区或腔内散射场值很 小处，若直接计算散射场可能由于相减噪声的存在而产生很大误差，只计算总场， 则可避免这种误差的影响，从而提高了计算的动态范围。 ( i i i ) 使散射体的设置变得简单。由于在时域有限差分法中在总场区内介 质交界面的边界条件是自动得到满足的，使得任意复杂的散射体结构的设置变得 简单，因为散射体总是设置在总场区中。如果用散射场编程，就需要在所有介质 交界面处计算入射场，以便使散射场和入射场之和的总场满足边界条件，而今只 需在连接边界上计算入射场。 ( i v ) 方便远场的计算。由于近区教射场已经获得，由此计算远区散射场已 经不成问题。而且在散射场区可规定一个规则的等效面，只要该等效面上的切向 电场和磁场分量已知，便可计算出远区散射场。这样可使远区散射场的计算程序 通用化，而与具体的散射体形式无关。 ( 3 ) 实现了稳态场的计算。由于有了以上两种技术，网格空间边界的反射 可降到很低的水平，而且入射波可以单独设置和总场区与散射场区的分离，都为 稳态电磁场问题的计算创造了条件。在解决了场的幅度和相位的检测后，计算远 区稳态散射场也成为可能。如此一来，直接时域方法和直接频域方法实现了直接 转化，当需要单频率或窄频带的信息时，时域有限差分法就可以用于直接频域计 算。 8 0 年代后期以来时域有限差分法进入一个新的发展阶段，即由成熟转入被广 泛接受和应用，在应用中又不断有新的发展。在8 0 年代中期以前，时域有限差 分法的研究和应用始终限于不大的一个圈子里，而在这之后的几年里发生了明显 的变化，大批科学家参加了进来，使得它的应用范围迅速扩大。随着应用范围的 扩大，不断提出新的要求，这就促使对时域有限差分法进行更深入的研究，使其 和到了进一步的发展。在这一阶段主要解决了以下几个问题。 ( 1 ) 回路积分法和变形网格。根据y e e 氏网格中电磁场的空间关系和积分 形式的f a r a d a y 定律和a 瑚p e r e 定律推导出了回路积分表示的差分格式，使之适 用于任意形态的网格结构，使得原来用阶梯法模拟的曲面结构得到更精确的模 拟。 ( 2 ) 亚网格技术。在传统的时域有限差分法中对物体模拟的最小尺度为一 个网格，为了模拟细微的结构只能把网格划细，这样会大大增加对存储空间和计 算时间的要求。在一些情况下这不仅意味着不经济，甚至可能超出计算机所具有 6 堡主堂垡笙茎 苎二童堕堕查里茎坌i ! ! 垦翌 的能力。所谓亚网格技术是一种用于模拟小于一个网格尺度的结构的方法，如小 于一个网格的导体窄缝或介质薄层等，这样就扩大了时域有限差分法的应用范 围。 ( 3 ) 广义正交曲线坐标系中的差分格式和非正交变形网格。这些方法把网 格的形式多样化，以适应于更精确地模拟各种形状的物体。正交曲线坐标系中的 网格体系适合于模拟各种外形与坐标面相重合的物体，由于网格体系的规则性， 使得编程比较容易。非正交变形网格则比较灵活，适于模拟形状复杂的结构。 ( 4 ) 适于色散介质的差分格式。在用时域有限差分法计算脉冲电磁场与色 散介质的相互作用时，由于脉冲含有宽的频率成分，差分格式中必须考虑介质的 色散性质，传统的差分格式已不再适用。近来用适当的迭代技术解决了卷积计算 问题。从而获得了适用于色散介质的差分格式，但色散关系必须有合适的数学模 型。 ( 5 ) 超吸收边界条件和色散吸收边界条件。传统的吸收边界条件只考虑边 界面上的切向电场或切向磁场，超吸收边界条件则既考虑电场需满足的条件，又 通过磁场的计算进行修正，从而进一步提高了相应吸收边界条件的效果。色散吸 收边界条件是为了适应不同频率时传播速度不同的情况。 时域有限差分法近期发展的另一个特点是迅速扩大了它的应用范围。在8 0 年代中期以前它还主要用于电磁散射问题，到8 0 年代中期首先成功地用到了生 物电磁剂量学问题的计算和电磁热疗的计算机模拟。到8 0 年代后期证明了时域 有限差分法用于微波电路的时域分析非常成功。进入9 0 年代以来又被用于天线 辐射特性的计算问题。随着新技术的不断提出，应用的范围和质量正在不断地扩 大和提高。 1 2 时域有限差分法的特点 作为一种电磁场的数值计算方法，时域有限差分法具有一些非常突出的特 点，也是它的优点。正是由于这些，使得越来越多的人对它产生了浓厚的兴趣， 并得到越来越广泛的应用。这些特点中最重要的是以下的几方面。 ( 1 ) 直接时域计算。时域有限差分法直接把含时间变量的】i i a x w e l l 旋度方 程在y e e 氏网格空间中转换为差分方程。在这种差分格式中每个网格点上的电场 ( 或磁场) 分量仅与它相邻的磁场( 或电场) 分量及上一时间步该点的场值有关。 在每时间步计算网格空间各点的电场和磁场分量，随着时间步的推进，即能直 接模拟电磁波的传播及其与物体的相互作用过程。时域有限差分法把各类问题都 硕士学位论文 第一章时域有限差分法概况 作为初值问题来处理，使电磁波的时域特性被直接反映出来。这一特点使它能直 接给出非常丰富的电磁场问题的时域信息，给复杂的物理过程描绘出清晰的物理 图像。 ( 2 ) 广泛的适用性。由于时域有限差分法的直接出发点是概括电磁场普遍规 律的m a x w e l l 方程，这就预示着这一方法应具有最广泛的适用性。从具体的算法 看，在时域有限差分法的格式中被模拟空间电磁性质的参数是按空间网格给出 的，因此，只需要设定相应空间以适当的参数，就可模拟各种复杂的电磁结构。 媒质的非均匀性、各向异性、色散特性和非线性等均能很容易地进行精确模拟。 ( 3 ) 节约存储空间和计算时间。在时域有限差分法中每个网格电场和磁场 的六个分量及其上时间步的值是必须存储的，此外还有描述各网格电磁性质的 参数以及吸收边界条件和连接条件的有关参量，它们一般是空间网格总数数 倍。所以时域有限差分法所需要的存储空间直接由所需的网格空间决定，与网格 总数成正比，在计算时，每个网格的电磁场都按同样的差分格式计算，故它所 需的主要计算时也是与网格总数成正比的。相比之下，若离散单元也是m 则 矩量法所需的存储空间与( 3 ) 2 成正比，而所需的c p u 时间则与( 3 ) 2 至( 3 ) 3 成 正比。当比较大时两者的差别是很明显的。 ( 4 ) 适合并行计算。当代电子计算机的发展方向是运用并行处理技术，以 进一步提高计算速度。时域有限差分法的计算特点是，每一网格点上的电场值( 或 磁场) 只与其周围相邻网格点处的磁场( 或电场) 及其上一时间步的场值有关， 这使得它特别适合并行运算。 ( 5 ) 计算程序的通用性。由于m a x w e l l 方程是时域有限差分法计算任何问 题的数学模型，因而它的基本差分方程对广泛的问题是不变的。因而一个基础的 时域有限差分计算程序，对广泛的电磁场问题具有通用性，对不同的问题或不同 的计算对象只需修改有关部分，而大部分是共同的。 ( 6 ) 简单、直观、容易掌握。由于时域有限差分法直接从方程出发，不需 要任何导出方程，这样避免了使用更多的数学工具，使得它成为所有电磁场的计 算方法中最简单的一种。其次，它能直接在时域中模拟电磁波的传播及其与物体 作用的物理过程，所以它又是非常直观的一种方法。这样，时域有限差法既简单 又直观很容易得到推广，并在很广泛的领域发挥作用。 本章小结 本章系统地介绍了时域有限差分法发展的几个阶段及它独特的优势。 硕士学位论文 第二章时域有限差分法的基本原理 第二章时域有限差分法的基本原理 2 1 y e e 氏算法 2 1 1 微商的差商近似 有限差分法是用变量离散的、含有限个未知数的差分方程近似地替代连续 变量的微分方程。因此，首要任务是构造合理的差分格式，使得它的解能保持原 问题的主要性质，并且有相当高的精确度。建立差分方程的基本步骤是把变量按 某种方式离散化，然后用差商近似地替代微分方程中的微商。为了表明用差商替 代微商的精确程度，我们用一元函数为例来加以说明。假设，( x ) 为x 的连续函 数，若在z 轴上每隔厅长取一个点，其中任一点用一表示，则在z 。点上的函数值 厂( x 。) 可通过t a y l o r 级数表示为 m 。m “，掣。等字。等挚l h ， 由此可得 掣= 掣l ，宝警l 等挚h 融l ，；21 知2 i z = - 3 1 舐3 = 警o ) ( 2 1 2 ) 而 m “h “， 掣卜等警卜蓦擎b 弦， 故 学叫做，。) 在一点的向前差商，而笪掣则叫做，( x ) 在一 点的向后差商。由式( 2 1 - 2 ) 和( 2 1 4 ) 可知，向前差商和向后差商与微商的差 均为变量离散步长 的一阶近视。 若把式( 2 1 1 ) 和式( 2 1 3 ) 相减，则可解得 9 挚 + 蛳 警 鬻 掣 堡主堂垡垒塞 笙三童堕堕查里茎坌婆堕苎兰! ! 塞里 丝! ! l 也= ，2 ：型i + 堡笪掣i + 笙掣+ 2 a xb - 3 1 叙3 l q 51 叙 ( 2 1 5 ) = 掣l 。+ 0 ( 确 丛墨盐掣曰q 做厂( 。) 在t 点的中心差商。由式( 2 1 5 ) 可知，中心差商与微 z ，z z 。： 、e ve 夕 陀 广 e ： ，- ， i ； e z e ： 卜b 乍 e 。 ，臣e ，f 7 r卢 ( i ，l k )e v 一 商之差为变量离散步长的二 阶近似。显然，就差商对微商 逼近的精度而言，在上面三种 差商形式中，中心差商的精度 最高。在时域有限差分法中正 是用中心差商代替微商由方 程来建立m “w e l l 差分方程式 的。 2 1 2y e e 氏网格 为了建立差分方程，首先要在变量空间把连续变量离散化。通常是用一定形 式的网格来划分变量空问，且只取网格结点上的未知量作为计算对象。这样，自 变量为离散的，以只在有限个点上计算未知量。当在每个离散点上用差商来替代 微商对，就把在一定空间解微分方程式的问题化为解有限个差分方程的问题。由 微分方程导出的差分方程式的一般通式，常常称为该方程的差分格式。 个逼近程度高的差分格式，不一定能给出好的近似解。这是因为一个合理 的差分格式还必须保持原总是的基本物理性质，所以在构造差分格式时常常从物 理定律出发，以便建立能给出高精度近似解的差分方程。 电磁场的最基本规律是m a x w e n 方程组，它们的一般形式是依赖时间变量的 旋度方程。从含有时间变量的m a x w e l l 旋度方程式出发，建立计算时域电磁场的 数值方法是很自然的。k s y e e 正是由此出发于1 9 6 6 年创立了计算电磁场的时 域有限差分法。 一般情况下，在时域计算电磁场要在包括时间在内的四维空间进行。如果 采用在限差分法，首先就要把问题的变量空间进行离散化，也就是要建立合适的 网格剖分体系。从m a x w e l l 方程出发建立差分方程式的复杂性在于，不仅要在四 维空间中进行，还要能同时计算电场和磁场的六个分量。在四维空间中合理地离 o 硕士学位论文 第二二章时域有限差分法的基本原理 散六个未知场量成为建立具有高精度的差分格式的关键问题ay e e 氏网格正是因 为提出了一个合理的网格体系，才成功地创立了时域有限差分法。我们把他所使 用的网格体系称为y e e 氏网格。在直角坐标系中的y e e 氏网格示于图2 1 。这个 网格体系的特点是，电场和磁场各分量在空间的取值点被交叉地放置，使得在每 个坐标平面上每个电场分量的四周由磁场分量环绕。这样的电磁场空间配置符合 电磁场的基本规律一一f a r a d a y 电磁感应定律和a m p e r e 环流定律，亦即 m a x w e l l 方程的基本要求，因而也符合电磁波在空间传播的规律。正是由于电磁 场分量在空间网格中的这种配置，使得用时域有限差分法在计算机的存储空间可 以模拟电磁波的传播及其散射体的相互作用过程。在这种电磁场的配置下，当空 问出现介质突变面时，可以使突变面上场分量的连续性条件自然得到满足，因而 为一些复杂结构的电磁场计算问题带来很大方便。这一点保证了时域有限差分法 应用的广泛性。 电磁场的计算与计算空间媒质的电磁性质有重要关系，在网格空间中除了规 定电磁场的离散取值点以外，还必须同时给出各离散点相应媒质的电磁参量，即 电场离散点处的介电常数和电导率以及磁场离散点处的磁导率和等效磁阻率。这 也说明，通过赋予空间点电磁参数的方法可在网格空间中模拟各种媒质及各种电 磁结构，这使得用时域有限差分法模拟电磁波与各种复杂的电磁结构的相互作用 变得比较容易。 在y e e 氏网格中，每个坐标轴方向上场分量间相距半个网格空间步长，因而 同一种场分量之间相隔正好为一个空间步长。在图2 一l 的网格单元中没有给出时 间的离散规则。为了保证计算的稳定性，时间离散的步长与空间步长必须满足一 定的关系，不能任意给定。由以后的分析可知，时间步长可选为电磁波传播一个 空间步长所需的时间的一半。这样，在实际运用时域有限差分法时，网格的空问 步长按选定后，时间变量的离散规则也就完全确定了。也就是说在选定了空间网 格结构后，就可根据差分近似的基本原则来建立所需的差分方程。 2 1 3 m a w e f l 方程 m a x w e l l 方程组概括了宏观电磁场的基本规律，它由两个旋度方程和两个散 度方程构成。两个旋度方程是f a r a d a y 电磁感应定律和a i i l p e r e 环流定律的微分 形式a 从本质上讲，m a x w e l l 方程组的四个方程中两个旋度方程是基本的，因为 两个散度方程可以由它们导出。所以，研究电磁场问题可以以两个旋度方程作为 出发点。时域有限差分法是在时域计算电磁场的一种数值方法，自然应该从含时 间变量的两个m a x w e l l 旋度方程出发。 在叙述基本原理阶段，我们把问题尽量简化，以便突出关键问题。因此，假 定我们暂时限定所研究的电磁场问题只涉及各向同性、线性且与时间无关的媒 堡主堂垡笙壅 笙三童堕垫查堡差坌鎏塑茎! ! 壁 质，但可以存在电的和磁的损耗。于是在无源区域，我们可把m a x w e l l 方程的 两个旋度方程表示为如下的形式： v e ：一掣一仃。h ( 2 1 6 a ) v 日= s 半+ d 。e ( 2 1 6 b ) 其中e 为电场强度，单位为伏米( v m ) ：日为磁场强度，单位为安米( a m ) 占为介电常数，单位为法米( f m ) ：为磁导率，单位为亨米( h m ) ：盯。为 电导率，单位为西门子米( s m ) ：仃。为等效磁阻率，单位为欧米( m ) 。这 里引进等效磁阻率的目的主要在于使方程具有对称性。 在