另一方面,为了更好的拟合研究对象的形状,克服台阶逼近带来的误差,D.E.Merewether[19]提出了柱坐标系下的网格剖分方法,R.Holland[20]提出了球坐标系下的网格剖分方法,P.Monk和E.Suli[12]提出了变网格步长方法,S.S.Zivanovic等[21]和P.Thoma等[22]提出了亚网格技术(即在一般区域采用粗网格,在电磁场快变区域采用精细网格)。利用这些技术,可以更精确地模拟各种复杂的结构,适应各种复杂的介质,提高了复杂介质中数值计算的精度。
时域模拟一般获得的是近场电磁信息,为了得到诸如天线方向图或散射体雷达散射截面之类的远场信息,必须获得计算区域以外的频域场或瞬态场。多位学者在这方面做了许多工作,发展了一种高效的时域近远场变换方法[23-26]。借助这种方法,可以实现由计算区域内近场数据到计算区域外远场数据的外推。目前,粗糙面散射的FDTD,传递函数在FDTD中的应用,周期介质、各向异性介质、色散介质和含有集中元件的FDTD,以及网络并行FDTD技术等方面也取得了很大进展。
FDTD在迅速发展的同时,也获得了非常广泛的应用。目前,它几乎被应用到了电磁场工程中的各个方面,例如:电磁散射、生物电磁计量学、辐射天线的分析、微波器件和导行波结构的研究、散射和雷达截面的计算、周期结构的分析、电子封装和电磁兼容的分析、核电磁脉冲传播和散射的分析、以及微光学元器件中光的传播和衍射特性的分析等。随着新技术的不断提出,其应用范围和成效正在迅速地扩大和提高。
第二章 时域有限差分法的基本原理
Maxwell方程是描述宏观电磁现象的一组基本方程。这组方程即可以写成微分形式,又可以写成积分形式。FDTD方法由Maxwell旋度方程的微分形式出发,利用二阶精度的中心差分近似,直接将微分运算转换为差分运算,这样达到了在一定体积内和一段时间上对连续电磁场数据的抽样压缩。
2.1 Maxwell方程和Yee氏算法
根据
[27]
中电磁场基本方程组的微分形式,若在无源空间,其空间中的媒质是
各向同性、线性和均匀的,即媒质的参数不随时间变化且各向同性,则Maxwell旋度方程可写成:
??H???E??E (2-1a) ?t?H??mH (2-1b) ?t??E???式中,E是电场强度,单位为伏/米(V/m);H是磁场强度,单位为安/米(A/m);?表示介质介电系数,单位为法拉/米(F/m); ?表示磁导系数,单位为亨利/米(H/m);?表示介质电导率,单位为西门子/米(S/m);?m表示导磁率,单位为欧姆/米(?/m)。
在直角坐标系中,(2-1)式可化为如下六个标量方程:
??Ex?Hz?Hy?????Ex??y?z?t???Ey?Hx?Hz??????Ey? (2-2) ?z?x?t??Hy?Hx??Ez?????Ez??x?y?t????Hx?Ez?Ey??????mHx??y?z?t???Hy?Ex?Ez???????mHy? (2-3) ?z?x?t??Ey?Ex??Hz??????mHz??x?y?t??这六个偏微分方程是FDTD算法的基础。
K.S.Yee[3]在1966年建立了如图2-1所示的空间网格,这就是著名的Yee氏元胞网格。
zExHzExEyEzHxEyoExyEzHyEzEyx
图2-1 Yee氏网格及其电磁场分量分布
并引入如下的差分近似方法对(2-2)、(2-3)式中的六个偏微分方程进行了差分离散。令F(x,y,z,t)代表E或H在直角坐标系中某一分量,在时间和空间域中的离散可记为
F(x,y,z,t)?F(i?x,j?y,k?z,n?t)?Fn(i,j,k) (2-4)
式中,?x、?y和?z分别是长方体网格沿x、y、z方向的空间步长,?t是时间步长,i、j、
k分别是沿x、y、z方向的网格编号,n是时间步数。对F(x,y,z,t)关于时间和空间的一阶偏导数取中心差分近似,具有二阶精度,即
?F(x,y,z,t)?xFn(i?12,j,k)?Fn(i?12,j,k)2????O?x (2-5a) x?i?x?x???F(x,y,z,t)?y?F(x,y,z,t)?z?F(x,y,z,t)?tFn(i,j?12,k)?Fn(i,j?12,k)2????O?y (2-5b) y?j?y?y??Fn(i,j,k?12)?Fn(i,j,k?12)2????O?z (2-5c) z?k?z?z??Fn?12(i,j,k)?Fn?12(i,j,k)2????O?t (2-5d) t?n?t?t?? 在FDTD中,空间上连续分布的电磁场物理量离散的空间排布如图2-1所示。由图可见,电场和磁场分量在空间交叉放置,使得在每个坐标平面上每个电场分量被磁场环绕,每个磁场分量也被电场环绕。这种电磁场的空间结构与电磁感应和电磁波传播的规律相符,在每一个网格单元都能满足法拉第感应定律和安培环流定律。各分量的空间相对位置也适合于Maxwell方程的差分计算,能够恰当地描述电磁场的传播特性。同时,电场和磁场在时间上交替抽样,抽样时间间隔相差半个时间步,使Maxwell旋度方程离散以后构成显式差分方程,从而可以在时间上迭代求解,而不需要进行矩阵求逆运算。因此,由给定相应电磁问题的初始条件,FDTD就可以逐步推进地求得以后各个时刻空间电磁场的分布。
2.2 FDTD的基本差分方程
根据上述原则,可将(2-2)、(2-3)式离散为如下的差分方程形式:
Exn?1(i?1,j,k)?CA(i?1,j,k)?Exn(i?1,j,k)?CB(i?1,j,k)
2222(2-6a)
n?1n?12211?H(i?,j?,k)?H(i?1,j?1,k)?y??? zz2222 ???n?1n?1221111?H(i?,j,k?)?H(i?,j,k?)?z??yy??2222????Eyn?1(i,j?1,k)?CA(i,j?1,k)?Eyn(i,j?1,k)?CB(i,j?1,k) 2222(2-6b)
n?1n?12211?H(i,j?,k?)?H(i,j?1,k?1)?z??? xx2222 ???n?1n?1221111?H(i?,j?,k)?H(i?,j?,k)?x??zz??2222????Ezn?1(i,j,k?1)?CA(i,j,k?1)?Ezn(i,j,k?1)?CB(i,j,k?1) 2222(2-6c)
n?1n?12211?H(i,?j,k?)?H(i?1,j,k?1)?x??? yy ??222211?n?n?221111?H(i,j?,k?)?H(i,j?,k?)?y??xx??2222????Hx2(i,j?1,k?1)?CP(i,j?1,k?1)?Hx2(i,j?1,k?1)?CQ(i,j?1,k?1) (2-6d)
22222222n?1n?1?Ezn(i,j?1,k?1)?Ezn(i,j,k?1)?y??? 22 ???nn11?E(i,j?,k?1)?E(i,j?,k)?z??yy??22????Hn?12yn?11111(i?,j,k?)?CP(i?,j,k?)?Hy2(i?1,j,k?1)?CQ(i?1,j,k?1) (2-6e) 22222222nn??Ex(i?12,j,k?1)?Ex(i?12,j,k)?z?? ??? nn11?Ez(i?1,j,k?)?Ez(i,j,k?)?x????22????Hzn?12n?11111(i?,j?,k)?CP(i?,j?,k)?Hz2(i?1,j?1,k)?CQ(i?1,j?1,k) 22222222nn??Ey(i?1,j?12,k)?Ey(i,j?12,k)?x?? ??? nn11?Ez(i?,j?1,k)?Ez(i?,j,k)?y????22?(2-6f)
???式中
CA(i,j,k)?2?t2?(i,j,k)??(i,j,k)?t, (2-7a) CB(i,j,k)?2?(i,j,k)??(i,j,k)?t2?(i,j,k)??(i,j,k)?t2?t2?(i,j,k)??m(i,j,k)?t, (2-7b) CQ(i,j,k)?CP(i,j,k)?2?(i,j,k)??m(i,j,k)?t2?(i,j,k)??m(i,j,k)?t(2-6)式就是FDTD的基本差分方程组。从式中可以看出,方程组中含有半个空间步和半个时间步,为了便于编程,可将(2-6)式改写成如下形式[28]:
Ex(i,j,k)?CA(i,j,k)Ex(i,j,k)?CB(i,j,k)?H(i,j,k)?HZ(i,j?1,k)Hy(i,j,k)?Hy(i,j,k?1)?(2-8a)
??Z???y?z??Ey(i,j,k)?CA(i,j,k)Ex(i,j,k)?CB(i,j,k)?H(i,j,k)?Hx(i,j,k?1)HZ(i,j,k)?HZ(i?1,j,k)?(2-8b)
??x???z?x??EZ(i,j,k)?CA(i,j,k)EZ(i,j,k)?CB(i,j,k)?Hy(i,j,k)?Hy(i?1,j,k)Hx(i,j,k)?Hx(i,j?1,k)?(2-8c)
?????x?y??