时域有限差分方法-林志立

Beihang University
编程举例1：一维FDTD问题（续）
% Initialize field and material arrays 3.初始化场量和介质参量阵列 Ceze = zeros(nx+1,1); Cezhy = zeros(nx+1,1); Cezj = zeros(nx+1,1); Ez = zeros(nx+1,1); 电场电流部分 Jz = zeros(nx+1,1); eps_r_z = ones (nx+1,1); % free space sigma_e_z = zeros(nx+1,1); % free space Chyh = zeros(nx,1); Chyez = zeros(nx,1); Chym = zeros(nx,1); Hy = zeros(nx,1); My = zeros(nx,1); mu_r_y = ones (nx,1); % free space sigma_m_y = zeros(nx,1); % free space

a0 a1 e j a2 e 2 j , j 2 j b0 b1 e b2 e
Beihang University
色散介质的FDTD模拟
模拟Lorentz色散介质的不同方法：
{a0 , a1 , a2 , b0 , b1 , b2 }
Error r ( ) r ()
E n [(b0 D n b1 D n 1 b0 D n 2 ) / 0 (a1 E n 1 a2 E n 2 )]/ a0
MSE approach :
Beihang University
FDTD编程流程
初始化 输出结果
主循环
Beihang University
J E
（电极化） （磁化） （欧姆定律）
一切电磁场和电磁波问题均可由以上方程，以及各类 具体的边界条件所决定！
Beihang University
麦克斯韦方程组中的运算符
散度（Divergence) 连续函数的偏微分运算
旋度（Curl)
Beihang University
FDTD的基本思想
－时域和空间域的离散化 －连续偏微分的有限阶近似
(1)电场在时间上取整数倍的 Δt; t=n *Δt;
(2)磁场在时间上取(整数＋ 1/2)倍的Δt; t=(n +1/2)*Δt;
Beihang University
麦克斯韦方程的离散化近似
以Hz为例：
H z 1 E E ( x y) t z y x
H zn1/ 2 ( i, j ) H zn1/ 2 ( i, j ) t
Beihang University
FDTD空间偏微分的近似
以Hz为例：
H z 1 Ex E y ( ) t z y x
E x E x ( i, j 1) E x ( i, j ) y y 2( ) 2 E y E y ( i 1, j ) E y ( i, j ) x x 2( ) 2
Beihang University
电磁学基本方程

麦克斯韦方程组( Maxwell Equations)
(安培环路定律) (高斯定律－电场)

(法拉第感应定律) (高斯定律－磁场) James Clerk Maxwell (1831–1879)
物质本构关系 (Constitutive Relations)
Beihang University
计算电磁学中的
时域有限差分方法
The Finite-Difference Time-Domain Method for Electromagnetics
林志立
zllin2008@ 北航仪器光电学院光电工程系
Beihang University
电磁学中几种重要的数值计算方法
Beihang University
FDTD空间域的离散化
(1)空间域的分割离散化
节点
Yee元胞（Δx, Δy, Δz)
Ex分量的空间离散分布图
Hx分量的空间离散分布图
Beihang University
FDTD空间域的离散化
YEE 元胞
例如：
H z 1 Ex E y ( ) t z y x
类似地，可实现各电磁场分量的 空间偏微分计算。
Beihang University
FDTD时间偏微分的近似
以Hz为例：
t=(n+1/2)Δt
H z 1 E E ( x y) t z y x
t=nΔt t=(n-1/2)Δt
n n E y H zn1/ 2 H zn1/ 2 1 E x ( ) t z y x 2( ) 2
以Lorentz 介质为例： r (t ) 2 (t )
时域：D(t ) 0 r (t ) * E 02
e0t sin( 02 02 t )U (t )
D() 0 r () E() 频域：
编程举例1：一维FDTD问题
x
基本旋度方程：
X向电导率 X向电流
X向磁电导率 X向磁流
Beihang University
编程举例1：一维FDTD问题
Matlab程序代码：
1.定义物理常量 % Define initial constants eps_0 = 8.854187817e-12; % permittivity of free space mu_0 = 4*pi*1e-7; % permeability of free space c = 1/sqrt(mu_0*eps_0); % speed of light % Define problem geometry and parameters 2.定义问题的参量和结构尺寸 domain_size = 1; % 1D problem space length in meters dx = 1e-3; % cell size in meters dt = 3e-12; % duration of time step in seconds number_of_time_steps = 2000; % number of iterations nx = round(domain_size/dx); % number of cells in 1D problem space source_position = 0.5; % position of the current source Jz
以减小数值色散。
sin 2 (
数值色散方程：
t
2
)
t 2 ( ) 2

c0 2
x , y , z

sin 2 (
k ) 2 2 ( ) 2
理想色散方程：
2
要求：
k 0 2

2 c0
2 2 (k x ky k z2 )
例如，取 xmax , ymax and zmax
1 min 10
t
2
0
1 例如，取 t 10 fmax
Beihang University
FDTD的离散参数的稳定性条件
★ 时间步长：Courant 稳定性条件
Ｚ域数值色散方程：（von Neumann method ）
4sin 2 ( k / 2) Z 1/2 Z 1/2 2 ( ) r (Z ) 0, 2 c0 t x, y,z
设有一连续函数 f ( x ) , 现欲求 f '( x ) 。
二阶中心差分近似表达式：
f '( x ) f ( x x ) f ( x x ) 2 x
当 x越小时，上式的近似程度 越高。
f ( x)
实际上：
f ( x x ) f ( x x ) x 2 f '( x ) f '''( x ) ... 2x 6

有限差分法 Finite Difference Method
– 静电场、静磁场的有限差分法；
– 时域行波的电磁场的时域有限差分法；

有限元法 (Finite Element Method)
– 数值求解各类独立的偏微分方程；（电磁学、材料力
学、工程热力学、声学等等）

矩量法 (Method of Moments) 适合于细线、平面形状结构的电磁场问题
n n n n 1 Ex ( i, j 1) E x ( i , j ) E y ( i 1, j ) E y ( i , j ) ( ) z y x
H zn1/ 2 ( i, j ) H zn1/ 2 ( i , j )
n n n n t E x ( i, j 1) E x ( i , j ) E y ( i 1, j ) E y ( i , j ) ( ) z y x
为了保持稳定性，该方程的所有解的模必须小于1。
对于非色散介质，时间步长不能大于以下表达式：
t c
1 1 1 1 x 2 y2 z 2
Beihang University
介质电磁参量的设定
H y H x Ez ( i, j ) 1 ( ) t z ( i, j ) x y
上式即为Hz的更新方程，由前一时刻的磁场和前半时刻 的临近空间格点的电场即可求出最新时刻的磁场。
Beihang University
麦克斯韦方程的离散化近似
采取类似的步骤，可以推导出其它场量的更新表达式：
例如，对于Ez：
Ez 1 H y H x ( ) t z x y