- 1 -
中国科技论文在线
使用 Jury 准则对时域有限差分算法进行稳
定性分析#
肖飞,唐小宏,王玲,吴涛*
基金项目:高等学校博士学科点专项科研基金(新教师基金)(20070614034)
作者简介:肖飞,(1975-),男,副教授,主要研究方向为计算电磁学和无源微波器件
(电子科技大学电子工程学院极高频复杂系统国防重点实验室,成都 611731)
摘要:稳定性分析是时域有限差分算法应用中重要的内容。然而,在分析过程中要找出对应
特征多项式的所有的根是相当困难的,尤其是对于高阶算法而言。为了解决这个问题,控制
领域中的 Jury 准则将被引入,它仅仅需要使用特征多项式的系数而不需要直接求解整个方
程获得所有的根。因此,使用它可以简化算法的稳定性分析。为了验证它的有效性,文中给
出了两个例子。
关键词:时域有限差分算法;稳定性分析;Jury准则
中图分类号:O44
Simplified stability analysis of the FDTD methods by using
the Jury test
Xiao Fei, Tang Xiaohong, Wang Ling, Wu Tao
(The EHF Key Laboratory of Fundamental Science, School of Electronic Engineering,University
of Electronic Science and Technology of China, ChengDu 611731)
Abstract: Stability analysis is very important for the application of finite-difference time-domain
(FDTD) methods. However, to find out all roots of a characteristic polynomial in stability analysis of
the FDTD methods is a very difficult task, especially for high-order accurate FDTD methods. In this
paper, the Jury test in automatic control theory is introduced in detail, which only uses those
coefficients of the characteristic polynomial for stability analysis instead of directly solving it. Such
treatment can simplify the procedure of stability analysis. As demonstrated in this paper, the Jury test is
successfully applied in stability analysis of two typical FDTD methods.
Keywords:finite-difference time-domain method; stability analysis; Jury test
0 引言
在计算电磁学中,时域有限差分算法已经被成功用于求解麦克斯韦方程之中[1,2]。到目
前为止,已经涌现出众多的相关算法,它们使用不同的时间步进格式或者空间有限差分格式。
这些时间步进格式都有自己的优缺点。例如,对于以 Yee-FDTD 算法为代表所采用的蛙跳时
间步进格式[1],它仅是二阶精度。这既是说,这些算法内在的数值色散将导致数值相速变成
频率和传播角的函数。另外一个重要的时间步进格式是交替方向隐式格式(ADI),基于它的
算法能够具有无条件稳定的特点 [3,4]. 从本质上来说,ADI-FDTD 算法仍然是二阶精度的。
在另一方面,针对高阶精度的时域有限差分算法正引起研究者的兴趣。到目前为止,一
些高阶精度的时间步进格式包括 Runge-Kutta 格式[5,7],Symplectic Integrator Propagator 格式
[8-9],Stagger Backward Differentiation 格式[10-11]和 Taylor Expansion 格式[12],已经被成功用于
构造高阶精度时域有限差分算法。这些高阶算法可以在一定程度上改善数值色散。然而,这
些高阶算法的稳定性分析却是比较困难的事情。按照稳定性要求,算法对应的特征多项式的
所有根都必须等于或者小于 1[13]。由于这些特征多项式的阶数一般都很高,因此很难找出其
所有的根来一一判别。举一个例子来说,在文献[12]中提出了一个显式四阶交错时域有限差
- 2 -
中国科技论文在线
分算法,其稳定性分析中将得到一个包含放大因子的高阶多项式,相对而言就很难找出其所
有对应的根。
在控制理论中,Jury 准则可以通过分析特征多项式的系数来判断线性离散时间系统的稳
定性。它是 Routh-Hurwitz 准则针对离散时间的应用 [14-16]。Jury 准则要求系统的极点都位于
以原点为中心的单位圆内。在本文中,Jury 准则将被详细得介绍。然后,它将被用于时域有
限差分算法的稳定性分析之中,从而极大得简化了分析过程。最后,Jury 准则将被用到两个
典型时域有限差分算法的稳定性分析中以验证它的有效性。
1 Jury 准则
在自动控制理论中,Jury 准则被用于离散控制系统的稳定性分析之中[14-15]。假设离散控
制系统的特征多项式是
( ) 0012211 =+++++= −−−− azazazazazR nnnnnn L , 0>na (1)
Jury 准则通过判断一系列的不等式来完成。如果,任何一个不等式不成立,则系统就是
不稳定的。如果所有的不等式都成立,则系统是稳定的。整个过程如下所示:
第一步:
J1: ( ) 01 >R
J2: ( ) ( ) 011 >−− Rn
第二步:
如果在第一步中的不等式都成立,则接下来需要构建 Jury 表,如表 1 所示。注意,
kn
kn
k aa
aa
a
b −= 0
0
1
(其中 0=k 到 1−n )和
jn
jn
j bb
bb
b
c
1
10
0
1
−
−−= ( 0=j 和 2−n ),其它类似。
表 1 Jury表
Table 1 Jury table
M
L
L
L
L
L
L
M
032
210
0121
1210
0121
1210
00
0
0
Entries
6
5
4
3
2
1
Row
ccc
ccc
bbbb
bbbb
aaaaa
aaaaa
nn
n
nn
nn
nn
nnn
−−
−
−−
−−
−
−−
如果奇数行的第一个元素都为正,即 00 >a 、 00 >b 和 00 >c 等等,则多项式 ( )zR 就
是稳定的多项式。
很显然,Jury 准则的优点就在于它仅仅使用多项式的系数来进行稳定性判断,而不需要
求解整个方程。另外,与 Routh-Hurwitz 准则不同的是,Jury 准则不需要双线性变换。
- 3 -
中国科技论文在线
2 使用 Jury 准则对两个典型时域有限差分算法进行稳定性分析
使用 Jury 准则对显式四阶精度交错时域有限差分算法进行稳定性分析
在文献[12]中的算法,以 TE 极化情况为研究对象。首先,在离散化过程中,它采用泰
勒展开方法引入了一个三阶时间导数,如下所示:
( ) ( )( )4213322
1
1
24
tO
t
Et
t
E
t
EE
n
z
n
z
n
z
n
z Δ+∂
∂Δ+∂
∂=Δ
− +++
(2)
其中, tΔ 式时间步长。在使用麦克斯韦方程将三阶时间导数降为二阶导数之后,可以
得到下面的方程。
( )
⎟⎟
⎟
⎠
⎞
⎜⎜
⎜
⎝
⎛
−∂
∂−∂
∂⎥⎦
⎤⎢⎣
⎡
∂
∂Δ+=Δ
− ++++
2
1
,
2
1
,
2
1
,
2
22
,
,
1
,
,24
11 n jiz
E
n
ji
x
n
ji
y
ji
n
jiz
n
jiz E
y
H
x
H
t
t
t
EE
jiσε (3)
接着,可以使用一个四阶差分格式去近似空间导数。例如,以 x 方向为例,得到
⎟⎠
⎞⎜⎝
⎛ −+−Δ=∂
∂
++−− jijijiji
ji
uuuu
xx
u
,23,21,21,23
, 24
1
8
9
8
9
24
11
(4)
其中, xΔ 是 x 方向上的空间步长。
最终,可以得到显示四阶交错时域有限差分算法[12]。这个算法需要保存前四个时间点
上的场量值。与 Yee 算法相比,这个四阶算法需要额外的计算存储,但是它并不导致太多
的计算量。
现在,Von-Neumann 方法将被用来进行稳定性分析。设
( )yjkxikJn
z
y
x
n
jiz
n
jiy
n
jix
yxe
E
H
H
E
H
H
Δ+Δ−
⎟⎟
⎟
⎠
⎞
⎜⎜
⎜
⎝
⎛
=
⎟⎟
⎟⎟
⎠
⎞
⎜⎜
⎜⎜
⎝
⎛
ξ
0
0
0
,
,
,
(5)
其中, 1−=J , θcoskkx = , θsinkkx = , k 是波数,θ 是传播角。而ξ 是放大
因子,它的幅度决定了算法的稳定性。把(5)带入到场分量的更新方程之中,由于系数矩阵
的行列式必须等于零,于是可以得到下面的一个高阶特征多项式。
( ) ( ) ( )[ ]( ) 0145261576 2232225 =−+−++− ξξξξξ yyxx QPQP (6)
其中, ⎟⎠
⎞⎜⎝
⎛ ΔΔ
Δ= vk
v
tcP vv 2
1sin , ⎥⎦
⎤⎢⎣
⎡ ⎟⎠
⎞⎜⎝
⎛ Δ+= vkQ vv 2
1sin6
3
1 2 ,且 yxv ,= 。由于稳定性条
件要求(6)中的所有的根都必须等于或者小于 1[13],即 1≤iξ ,而(6)的阶数很高,如果直接
分析它将会很困难。
现在使用 Jury 准则来分析(6)。首先,把(6)展开为
( ) ( )
082692
2335762601152676576
23
4567
=+⋅−⋅+⋅−
⋅++−+−+
WWWW
WWW
ξξξ
ξξξξ
(7)
其中, ( ) ( )22 yyxx QPQPW += 。接着,将 Jury 准则直接应用到(7)。分析表明,下列的
- 4 -
中国科技论文在线
不等式必须被满足,即
( ) ( ) 01296230411 7 >−=−− WR (8)
另外,(8)的上边界也满足稳定性条件,所以也必须包括进来。最后,可以得到
9
16≤W (9)
因为
( ) ( )
⎥⎥⎦
⎤
⎢⎢⎣
⎡ ⎟⎠
⎞⎜⎝
⎛
Δ
Δ+⎟⎠
⎞⎜⎝
⎛
Δ
Δ≤
⎭⎬
⎫
⎩⎨
⎧
⎥⎦
⎤⎢⎣
⎡ ⎟⎠
⎞⎜⎝
⎛ Δ+⋅⎟⎠
⎞⎜⎝
⎛ ΔΔ
Δ+
⎭⎬
⎫
⎩⎨
⎧
⎥⎦
⎤⎢⎣
⎡ ⎟⎠
⎞⎜⎝
⎛ Δ+⋅⎟⎠
⎞⎜⎝
⎛ ΔΔ
Δ=
+=
22
2
2
2
2
22
3
7
2
1sin6
3
1
2
1sin
2
1sin6
3
1
2
1sin
y
tc
x
tc
ykyk
y
tcxkxk
x
tc
QPQPW
yyxx
yyxx
(10)
所以,同时考虑(9)和(10),可以得到下面的稳定性条件。
7
4
22
≤⎟⎟⎠
⎞⎜⎜⎝
⎛
Δ
Δ+⎟⎠
⎞⎜⎝
⎛
Δ
Δ
y
tc
x
tc (11)
这个结果与文献[12]中的结果一致。然而,采用 Jury 准则的分析过程却要容易得多。
使用 Jury 准则分析 HIE-FDTD 算法的稳定性
最近,一个混合隐式-显式时域有限差分算法(HIE-FDTD)在文献[17-18]中提出。该算法
是弱条件稳定的。仿照上节的处理方法,首先仍然使用 Von-Neumann 方法获得特征多项式
如下
( ) ( ) ( ) ( )[ ] 01122211 22222222 =++−++++− yzyxy PPPPP ξξξ (12)
其中
⎟⎠
⎞⎜⎝
⎛ ΔΔ
Δ= vk
v
tcP vv 2
1sin , zyxv ,,= (13)
很明显得,(12)中的两个根等于 1。为了保证其它四个根都能够满足稳定性的要求,Jury
准则被用来到处相应的稳定性条件。于是,通过 Jury 准则可以发现下面的不等式必须被满
足,即
( ) ( ) 044411 222 ≥−−=−− zx PPR (14)
最后,可以得到 HIE-FDTD 的稳定性条件如下
1
22
≤⎟⎠
⎞⎜⎝
⎛
Δ
Δ+⎟⎠
⎞⎜⎝
⎛
Δ
Δ
z
tc
x
tc (15)
3 结论
在本文中,介绍了如何将 Jury 准则用于时域有限差分算法的稳定性分析。这种方法仅
仅使用特征多项式的系数就能完成稳定性判断,从而避免了直接求解特征多项式来获得其所
有的根进行判断,使整个分析更加简化。与 Routh-Hurwitz 准则相比,由于 Jury 准则不需要
进行双线性变换,所以在一定程度上 Jury 准则要简单一些。
- 5 -
中国科技论文在线
[参考文献] (References)
[1] K. S. Yee. Numerical solution of initial boundary value problems involving Maxwell’s equations in isotropic
media[J]. IEEE Trans. Antennas Propagat., 1966, 14 (3):302-307.
[2] A. in Computational Electrodynamics: The Finite-Difference Time-Domain Method[M].
Boston: MA: Artech, 1998.
[3] T. Namiki. A new FDTD algorithm based on alternating-direction implicit method[J]. IEEE Trans. Microw.
Theory Tech., 1999, 47 (10):2003-2007.
[4] S. W. Staker, C. L. Holloway, A. U. Bhobe, M. Piket-May. Alternating-direction implicit (ADI) formulation
of the finite-difference time-domain (FDTD) method: algorithm and material dispersion implementation[J].
IEEE Trans. Electromagn. Compat., 2003, 45 (2):156-166.
[5] J. L. Young, D. Gaitonde, J. S. Shang. Toward the construction of a fourth-order difference scheme for
transient EM wave simulation: Staggered Grid Approach[J]. IEEE Trans. Antennas Propagat., 1997, 45
(11):1573-1580.
[6] V. Vegh, . Turner, H. Zhao. Effective cell-centred time-domain Maxwell’s equations numerical solvers[J].
Appl. Math. Model., 2005, 29 (5):411–438.
[7] Q. Cao, R. Kanapady, F. Reitich. High-order Runge–Kutta multiresolution time-domain methods for
computational electromagnetics[J]. IEEE Trans. Microw. Theory Tech., 2006, 54 (8):3316-3326.
[8] T. Hirono, W. Lui, S. Seki, Y. Yoshikuni. A three-dimensional fourth-order finite-difference time-domain
scheme using a symplectic integrator propagator[J]. IEEE Trans. Microw. Theory Tech., 2001,
49(9):1640-1648.
[9] Z. X. Huang, W. Sha, X. L. Wu, M. S. Chen. A Novel High-order time-domain scheme for three-dimensional
Maxwell’s equations[J]. Microw. Opt. Techn. Let., 2006, 48 (6):1123-1125.
[10] M. Ghrist, B. Fornberg, T. Driscoll. Staggered time integrators for wave equations[J]. SIAM J. Numer. Anal.,
2000, 38 (3):718–741.
[11] K. P. Hwang. A fourth-order accurate FDTD scheme with long-time stability[J]. IEEE Microw. Wireless
Compon. Lett., 2005, 15 (4):271–273.
[12] Z. Xie, C. H. Chan, B. Zhang. An explicit fourth-order staggered finite-difference time-domain method for
Maxwell’s equations[J]. J. Comput. Appl. Math., 2002, 147 (1):75–98, 2002.
[13] J. W. Thomas. Numerical partial differential equations: finite difference methods[M] . Springer-Verlag, 1997.
[14] C. D. Johnson, H. A. Malki. Control systems technology[M]. New Jersey : Upper Saddle River, 2002.
[15] Jamshidi, Mohammad. Linear control systems[M]. Great Britain : A. Wheaton Co. Ltd., 1986.
[16] Fei Xiao, Guohong Dua, and Xiaohong Tang. Application of the Jury test in stability analysis of an explicit
fourth-order staggered FDTD method[Z]. Asia-Pacific Microwave Conference, 2008.
[17] J. Chen and J. Wang. A 3D hybrid implicit-explicit FDTD scheme with weakly conditional stability[J].
Microw. Opt. Technol. Lett., 2006, 48 (11):2291–2294.
[18] J. Chen and J. Wang. A three-dimensional semi-implicit FDTD scheme for calculation of shielding
effectiveness of enclosure with thin slots[J]. IEEE Trans. Electromagn. Compat., 2007, 49 (2):354-360.