FDTD原理及例子

DbHy
i,
j,
k

Ex
i,
j,
k

x
1, n
z
Ex
i,
j,
k,

n

H z i, j, k, n 1 DaHz i, j, kH z i, j, k, n
Ex i, j 1, k, n Ex i, j, k, n
求所有空间离散点上n+1时间步的磁场 求所有空间离散点上n+1时间步的电场
n=n+1
No
n>nmax
Yes 结束
三维Maxwell方程的Yee算法
介绍了求解矢量Maxwell方程的FDTD Yee算法，归纳起 来，Yee算法的主要特点有：
1）Yee算法采用耦合的Maxwell旋度方程，同时在时间和 空间求解电场和磁场，而不是采用波动方程只求解电场或 磁场。
2 2 2
E y n
E E y i, j1,k1
y i, j 1,k

2
2
O
z
2
可以看到有个2阶
z i, j 1 ,k 1
z
22
Ez n
E E n z i, j 1,k 1
n z i, j,k1

2
2
O
y
2
小量，这个要忽略 掉，因此是近似的
Ex 1 H y
t
0 r z
H y 1 Ex
t
0r z
利用一阶导数的二阶中心差分近似，上面的方程变为
E
n x
1
(k
)

E
n x
(k
)


1
H
n 1 / y
2
(k

1 2
)

H
n 1 / y
2
(k

1 2
)
t
0 r (k)
z
H
n1 y
/
2
(k

3k
H~
n1/ 2 y
(k

1) 2

H~
n 1 / y
2
(k

1) 2

ct zr (k

1 2
)
E
n x
(k

1)

E
n x
(k
)
E
n x
1
(k
)

E
n x
(k)

ct
z r (k )

H~
n 1 y
/
2
(k

1) 2

H~
n 1 y
/
2
(k

1 2
)
三维Maxwell方程的Yee算法
基础知识
麦克斯韦方程微分形式：
FDTD方式将时间进 行差分，并且磁场与 电场交替迭代更新
对于有耗媒质：
H D J t
E B t
B 0
D
时谐场形式：
H

D t

Je,
Je

E

E


B t

Jm,
Jm

sH
B m
ca[k] ct ;cb[k] ct
z r
(k

1 2
)
z r (k)
只要给定了所有空间点上电/磁场的初值，就可以一步一步地求出任 意时刻所有空间点上的电/磁场值。
一维Maxwell方程的Yee算法
Ex
Hy Ex
Hy Ex
0
1
n2 n 3/2
n 1 n 1/ 2
n0
2


H
z
i,
j

1,
k,
n
z
1
H
z
i,
j,
k, n

1

y

2 t CaEw 2 t Ew所在空间位置
2t
CbEw 2 t Ew所在空间位置
三维Maxwell方程的Yee算法
媒质参数赋值 在所有空间点给电磁场分量赋初值
z
k, n
Ez
i,
j, k,

n

y

三维Maxwell方程的Yee算法
同理，可以得到其他2个磁场分量的FDTD方程
H y i, j, k, n 1 DaHy i, j, kH y i, j, k, n
Ez i 1, j, k, n Ez i, j, k, n
差分近似所带来的数值色散
一维标量波动方程为例 2u c 2 2u
t 2
x 2
设在离散空间点xi ,tn ,离散行波解为
u
n i
u
xi , tn
e j ntk~ix
将上式代入差分方程
u in 1
ct x
H~
n1 y
/
2
(k

1) 2
ct
zr (k
1 2
)
E
n x
(k

1)

E
n x
(k
)
E
n x
1
(k
)

E
n x
(k
)

ct
z r (k )

H~
n y
1
/
2
(k

1) 2

H~
n y
1
/
2
(k

1 2
)
式中， c 1
为自由空间中的光速。
0 0
一维Maxwell方程的Yee算法
D e
H j( j )E
E jH
E(x, y, z, t) E0 (x, y, z,t)e jt H (x, y, z, t) H0 (x, y, z,t)e jt
一维Maxwell方程的Yee算法
一维Maxwell方程，介质参数和场量均与x，y无关(无损，电导率和 磁导率为0)
0
E x t
yEx


H zx H zy y
0
E y t
xEy

H zx H zy x
0
H zx t


* x
H
zx
E y x
换句话说，Berenger构造了一般新的非物 理媒质（称为PML媒质），在该媒质中 场满足左端的方程（并不一定是Maxwell 方程）
吸收边界条件
理论上说，求解空间是无限大的，但是由于计算的数据容量 问题，需要在有限空间的周围做特殊处理，使得向边界面行进 的波在边界处保持“外向行进”，无明显的反射现象，并且不 会使内部空间的场产生畸变。
1）Mur吸收边界条件 2）Berenger完美匹配层（PML） 3）各向异性匹配层吸收边界条件
1.每一个磁场分量由 四个电场分量环绕； 每一个电场分量由 四个磁场分量环绕
2.这种空间取样方式 符合法拉第感应定 律和安培环路定律
3.电场和磁场在时间 上交替抽样，抽样 时间彼此相差半个 时间步
Hz Ey
Ex Ez
Ex
Ey
Ez
Ez
Hy
Hx
z
Ey
Ex
( i, j, k )
y
x
FDTD离散中的Yee元胞
数值稳定性问题
二维的Yee算法数值Courant稳定性条件
t
1
1
1
c

x2 y2
三维的Yee算法数值Courant稳定性 条件
t
1
c
111
x2 y2 z2
便于理解，当Yee元胞 为立方体时，此两式表 明时间间隔必须等于或 小于波以光速通过yee 元胞对角线长度1/3或 1/2所需时间

1 x
cos
1
1


x ct
2
cost


1

我们想要k与没加元胞之前的k一样，但是差分近似后k变
的不一样
差分近似所带来的数值色散
将平面波带入差分方程所出现的稳定性、色 散、以及各向异性，这些特性并非由介质的 物理特性所引起，而是数值计算中的差分近 似所致，在FDTD数值计算中，稳定性、色 散、各向异性将影响计算精度。
Berenger完美匹配层（PML）
只考虑二维TE情况，对于二维TM和三维情况可采用类似方法 进行分析。 分析范围在PML层内
0
E x t
Ex

H z y
0
E y t
Ey

H z x
Hz*
0
H z t
*Hz

E x y

E y x
分分别别表表示示自自由由空空间间中中的的电电导导率率和和磁磁阻阻率率
Da Ca , Db Cb
，并注意到E与H在时间上差半个步长，
可以直接从磁场FDTD公式得到电场的FDTD公式。如：
Ex i, j, k, n 1 CaEx i, j, kEx i, j, k, n CbEx i, j, k
H y i, j, k 1, n 1 H y i, j, k, n 1
2
u
n i 1

2u
n i
uin1
2uin uin1
t2 c 2O x2 O t2
得： e jt
ct 2 e jk~x 2 e jk~x

2 e jt
x
最后得色散关 系
~ k
1 2
)

H
n 1 / y
2
(k

1 2
)


1
E
n x
(k

1)

E xn
(k)
t
0r
(k

1 2
)
z
一维Maxwell方程的Yee算法
采用归一化磁场
H~ 0 H 0
使得电场与归一化磁场有相同的数量级，于是可以得到FDTD
迭代公式为
H~
n 1 / y
2
(k

1) 2

2） Yee网格在三维空间这样安排E和H分量，使得每一个 E或H分量由四个H或E循环的分量所环绕。
3） Yee算法以蛙跳算法在时间上安排E和H分量。在某一 时刻，使用前一时刻的E数据计算所有H分量。然后，再 使用刚计算的H数据计算所有的E分量。如此循环，直至 完成时间步进过程。
数值稳定性问题
（1）FDTD计算中每一步都是有误差的，随着时间步进， 误差会不断积累。如果误差的积累不会造成总误差的增 加，就成FDTD法是稳定的，否则成为不稳定的。数值 不稳定性会造成计算结果随时间步进无限增加。 （2）FDTD法是有条件稳定的，即：时间步必须必须小 于一定值以避免数值不稳定性。 （3）数值稳定性分析方法是建立在Courant等人几十年 前提出的经典方法基础上。这种方法首先把有限差分算 法分解为相互分离的时间和空间本征值问题。
三维Maxwell方程的Yee算法
Maxwell旋度方程为

E


H


H
t
H x t

1 Ey


z
Ez y
ຫໍສະໝຸດ Baidu
H x

H y t

1 Ez
x
Ex z

H
y

H z t

1


E x y
E y x
Berenger为了引入规定损耗和阻抗匹配的新自由度，将 HHzz分
裂为两个分量H zx 和 H zy，即 H z H zx H zy 同时引入了新的电导
率
x , y
和磁损耗

* x
,
* y
，并规定TE情形的四个场分量（而不
是通常的3个）由下列方程耦合在一起：
Berenger完美匹配层（PML）

H z

H

E

E
t
E x t

1


H z y
H y z

Ex

E y t

1

H x z
H z x

E
y

E z t

1

H y

x
H x y

Ez

三维Maxwell方程的Yee算法
FDTD数值分析法
目录
1.麦克斯韦方程的基础知识 2.一维和三维Maxwell方程的Yee算法 3.数值稳定性分析 4.吸收边界条件 5.波源的设置 6.编程思路
麦克斯韦方程微分形式：
BPM 方式
基础知识
H D J t
E B t
B 0
D
用计算机语言表示的FDTD公式
Hy[k] Hy[k] ca[k]* Ex[k 1] Ex[k] Ex[k] Ex[k] cb[k]* Hy[k] Hy[k 1]
式中，时间变量已隐含在迭代公式中，以及
Hy[k ]

H~
y
(k

1 ); Ex[k] 2

E x [k ];
y i, j 1,k 1
y
22
H x i, j, k, n 1 DaHx i, j, kH x i, j, k, n
最后得：
E y i, j, k 1, n E y i, j, k, n

DbHx i,
j,
k

Ez
i,
j
1,
采用时间平均近似
n 1
n1
H x
t
n i, j1,k1

H H 2 x i, j1,k1
2 x i, j1,k 1
22
22
t

O

t
2
22
n
n
n 1
n1
2
2
H H n
x i, j1,k1
x i, j1,k1
22
22
2
H O t x i, j 1,k 1

DbHz
i,
j, k

Ey
i
1,
y
j, k, n
x
Ey
i,
j,
k,

n

2 t DaHw 2 t Hw所在空间位置
2t
DbHw 2 t Hw所在空间位置
三维Maxwell方程的Yee算法


利用对偶原理： H E，E H， , ,