林智参，班涛
（华南师范大学，广东广州510006）
摘要：提出一种新型的分裂步长时域有限差分（NSS?FDTD）法，并对其数值色散进行分析。该方法基于Split?Step方案和Crank?Nicolson方案，采用新的矩阵分解形式，与传统的FDTD算法、SS?FDTD算法相比，减少了计算复杂度。新型算法的推导程序简单，且具有良好的数值色散特性，还加入了一阶Mur吸收边界条件，给出一阶Mur吸收边界差分方程。将数值实验的结果和传统FDTD方法及理论值进行比较，数值结果一致性较好。
教育期刊网 http://www.jyqkw.com
关键词 ：时域有限差分法；分裂步长；Split?Step方案；数值色散
中图分类号：TN802?34；O441 文献标识码：A 文章编号：1004?373X（2015）15?0117?03
收稿日期：2015?02?10
0 引言
时域有限差分法（Finite Difference Time Domain，FDTD）是一种简单直观的全波分析时域算法[1?3]，该方法以Yee氏立体网格作为电磁场离散单元，将麦克斯韦方程转化为差分方程，能够方便有效地结合计算机技术处理复杂的电磁场问题，目前已经在电磁学的各个领域中得到了广泛应用。然而，应用传统的FDTD方法也明显体现出其不足之处，为减小差分近似带来的数值色散[4]，空间网格尺寸必须远小于波长，这样反而增加运算负担，因而，出现了多钟FDTD的改进方法[5?7]。本文以TM波为例，提出一种基于Split?Step方案[8]和Crank?Nicolson方案[9]的时域有限差分法，所提出的算法采用新的矩阵分解方式，简化计算复杂度、减少差分近似所带来的数值色散。
1 NSS?FDTD 算法理论推导
考虑空间一无源区域，在均匀无耗、各向同性介质中，介电常数为ε，磁导率为μ，二维TM波Maxwell微分方程组的矩阵形式如下：
将矩阵M 分解成四个子矩阵，分别记为A 2，B 2，C 2，D 2，矩阵形式分别如下：
公式（1）可以写为：
利用Crank?Nicolson方案，对矩阵方程式（3）~式（6）的右端进行近似，进一步化简得到以下形式：
式中：I 为3 × 3 的单位矩阵。四个分步内需要求解的方程式以第（1）步为例化简如下：
2 数值色散分析
利用Fourier方法[10]在第n 个时间步内，场分量在空间区域内的表达形式如下：
将式（15）代入到式（7）~（10）中，并将式（7）~（10）进一步整理为如下矩阵形式：
将式（16）~（18）代入到式（19）中，得到一个完整时间步长内的矩阵方程式：
可以求得Ω 的特征值，其结果如下：
其中：
利用von Neumann方法[11]，假定一角频率为ω 的电磁波产生的电磁场满足：
将式（22）代入到式（20）中，得：
其中，Un 与初始值U0 相关，其具体关系式如下：
为了使式（23）中Un 有非零解，Un 的系数行列式的值应为零，即：
由式（21）得Ω 的特征值，可以得到NSS?FDTD 算法的数值色散表达式，其形式如下：
为了方便分析NSS?FDTD的数值色散特性，作如下定义：S = cΔt Δx,Δx = Δy,N = λ Δx,λ代表波长，N 表示单位波长元胞数，传输角度为θ。图1为S = 0.5，N = 8 时归一化数值相位速度随传输角度θ的变化曲线。
从图1 可以看出，NSS?FDTD 算法的归一化相位速度大于传统的FDTD 算法和传统SS?FDTD 算法的归一化数值相位速度，并接近1，且曲线变化值在0.98~0.995范围内，比较平缓，归一化数值相位速度各向异性误差也比较小。
3 一阶Mur 吸收边界条件
在分步1 内，电场分量Ez 在i = 1，i = Imax，j = 1 和j = Jmax 上的一阶Mur吸收边界差分方程式为：
在分步2，3和4内，电场分量Ez 在i = 1，i = Imax，j = 1 和j = Jmax 上的一阶Mur吸收边界差分方程式与分步1内的差分方程式类似，此处不再一一进行展开说明。
4 计算结果
将NSS?FDTD算法运算于尺寸101 cm×101 cm的自由空间，以一阶Mur为边界条件，以二维TM 波为例，在中心区域Ez 场分量上加正弦波激励源sin(2πft), 其中，f = 1.5 GHz。网格尺寸为Δx = Δy = 1 cm，为激励源最高频率对应波长的1 20，计算网格数为101×101。观察点在中心区域和吸收边界之间的中心位置。NSS?FDTD数值计算结果仿真如图2~图4所示。其中图2，图3两图为NSS?FDTD 数值计算过程中Ez 场分量的空间分布图，图4为传统FDTD算法和NSS?FDTD算法比较图。
由图可见，NSS?FDTD 算法的计算结果与传统的FDTD 算法的计算结果吻合的很好，且符合电磁场理论，从而证实新型分裂步长时域有限差分法的可行性。
5 结论
本文基于Split?Step 和Crank?Nicolson 方案提出了一种新型的二维FDTD算法，改进算法采用新的分解形式，与传统的SS?FDTD算法相比较，减少了计算复杂度，优化了计算公式，使推导过程更简单。结合算例使用Matlab对NSS?FDTD算法进行编程分析，结果表明，该算法具有良好的预期效果。
教育期刊网 http://www.jyqkw.com
参考文献
[1] 葛德彪，闫玉波.电磁波时域有限差分方法[M].西安：西安电子科技大学出版社，2002.
[2] 王秉中.计算电磁学[M].北京：科学出版社，2002.
[3] YEE K S. Numerical solution of initial boundary value prob?lems involving Maxwell′ s equations in isotropic media [J].IEEE Transactions on Antennas and Propagation， 1966， 14（3）：302?307.
[4] TAFLOVE A，HAGNESS S C. Computational electrodynamics：the finite?difference time?domain method [M]. 2nd ed. Boston：Artech House，2000.
[5] 徐利军，衰乃昌.高阶ADI?FDTD算法的数值色散分析[J].电子与信息学报，2005，27（10）：1662?1665.
[6] 党涛，郑宏兴.关于二维ADI?FDTD方法的数值色散分析[J].中国民航学院学报，2004，22（2）：42?46.
[7] 夏冬，党涛，郑宏兴.一维ADI?FDTD 方法的数值色散分析[J].中国民航学院学报，2005，23（2）：38?41.
[8] LEE J，FORNBERG B. A split step approach for the 3?D Max?well′ s equations [J]. Journal of Computational and Applied Mathematics，2003，158（2）：485?505.
[9] SMITH G D. Numerical solution of partial differential equa?tions：finite difference methods [M]. 3rd ed. Oxford：Oxford University Press，1986.
[10] SMITH G D. Numerical solution of partial differential equa?tions [M]. Oxford：Oxford University Press，1978.
[11] PEREDA J A，VIELVA L A，VEGAS A，et al. Analyzing the stability of the FDTD technique by combining the Von Neumann method with the Routh?Hurwitz criterion [J]. IEEE Transactions on Microwave Theory and Techniques，2001，49（2）：377?381.
作者简介：林智参（1986—），男，广东吴川人，硕士。研究方向为电磁场与微波技术。