FDTD 算法的FPGA 实现

徐奇澎，郭裕顺

(杭州电子科技大学电子信息学院，杭州 310018)

FDTD 法即时域有限差分方法(Finite Difference in Time Domain)，是把Maxwell 方程式在时间和空间域进行差分，通过电场和磁场的交替计算，来模拟电磁波传播的一种数值计算方法［1］。自上世纪六十年代出现这种算法以来，人们已对它作了大量研究，是目前电磁波数值仿真的主要方法，被广泛地应用于科学与工程中各种电磁场与电磁波问题。

FDTD 算法应用中存在的一个主要问题是:当计算区域较大，即空间差分产生较多网格时，需要耗费大量的计算资源与时间［2］。为此降低计算代价，减少计算时间，一直是重要的研究课题。以往采取的措施不外乎两个方面:一是改进算法，设法提高计算效率;二是在通常的计算平台上采用并行计算。

FPGA 的发展为解决这类大规模计算问题提供了一种新途径。FDTD 计算的特点是在固定的一个空间网格中进行重复的算术运算，数据流规则，控制流较为简单，因此可利用FPGA 可编程的特点，设计专门的算术运算电路，实现FDTD 的计算。一块FPGA 芯片可同时实现多个乘法器、加法器，轻易实现并行运算。这样针对FDTD 算法专门设计的运算电路，可大大提高计算速度。最早用FPGA 实现FDTD 算法的是Schneider，他用定点整数格式实现了FDTD 算法的一维运算，但网格单元仅10个［3］。Gandhi、Durbano 分别用浮点数实现了FDTD 算法的二维、三维运算，但受限于当时FPGA 的性能，这些设计工作频率较低，性能不高，未能体现硬件计算的优势［4］。后来Wu Shuguang 在Gandhi 的研究基础上把数据存储位置由主机系统内存改为FPGA 片上SRAM，大大减少了数据传输耗费的时间［5］。最近Hideki Kawaguchi 等采用分配式并行存储使数据读写速度进一步加快，当FPGA 工作频率在50 MHz时，计算速度比软件工作在2.8 GHz 的个人电脑上快了50倍［6－7］。因此可以预计随着FPGA 的进一步发展，用FPGA 实现FDTD 算法的优势将更加明显。

1 二维FDTD 算法

本文考虑二维FDTD 算法的实现。二维情形下所有电磁量与Z 坐标无关［8］，即∂/∂z=0。这时电磁场的直角分量可以分为独立的两组，即以Hx、Hy、Ez为一组的TM 波，及以Ex、Ey、Hz为一组的TE波［9］。以TM 波为例，二维TM 波的Yee 元胞如图1所示，电场值根据周围4个磁场值完成更新，磁场值根据周围2个电场值完成更新，更新关系满足右手螺旋法则。

图1 二维TM 波的Yee 元胞

为了便于硬件实现，对上述关系中的1/2 按照一定规则做出处理［10］。处理后TM 波的计算公式如下:

下面介绍用FPGA 来完成上述计算的方法。

2 算法的硬件实现

用硬件实现式(1)～式(3)的计算首先要确定数字的表示方法。本文采用符合IEEE-754 32 位标准的单精度浮点数格式，即将一个数表示为V=(－1)S×1.M×2E-BIAS，其中最高位S 是符号位，随后的8 位E是指数位，后面的23 位M 是尾数位，BIAS=28－1=127。

图2 电路结构

图3 计算电路

实现一次FDTD 计算的工作流程如图4所示，分为以下几步:

(1)初始化:确定要计算的网格大小，所要计算的时间步，计算式系数。

(2)加入激励后，整个计算电路开始工作。本文中因为激励信号是添加在RAM Ez中，所以先计算Ez值。根据计数器计算出相应的地址值，从RAM 中取得相邻两个单元的磁场值，Ez的计算电路完成一个点的Ez值计算，E值计数器会加1，然后检查其是否数到最大值。如果计数器没有达到最大值则再次计算Ez值，进行RAM Ez中下一个点的更新。如果是计数器达到最大值则表明RAM Ez中所有的点都完成一个时间步的更新，转换更新算式开始进行磁场值的更新。

(3)Hx和Hy的更新同时进行，当一个点的Hx和Hy值更新完成之后，H值计数器会加1，然后检查其是否数到最大值。如果计数器没有达到最大值则再次计算Hx和Hy值，进行RAM Hx和RAM Hy中下一个点的更新。如果计数器达到最大值则表明RAM Hx和RAM Hy中所有的点都完成一个时间步的更新，同时也意味着网格中所有的点都完成了一个时间步的电场值和磁场值的更新。

(4)完成一个时间步更新后，时间步计数器T加1，然后检查T 是否小于nstep，如果是则表示还未完成设定的时间步，回到算法第3 步，开始网格内所有节点下一个时间步的更新。如果不是则表示计算完成设定的时间步，整个算法完成。

图4 算法流程

在设计过程中我们可以根据并行原理［11－12］和FPGA 架构特点，设计多个计算电路，让这些计算电路同时运算，可以大大减少运算所需要的时间，提高算法的效率。

3 实验结果

根据上述方案，用Verilog HDL 编写了代码，以ALTERA 公司的EP3C55F484C6为目标器件进行了编译，用ModelSim 进行了仿真。设置计算网格大小为60 ×60，在网格中间加一个高斯脉冲g(t)=exp(－0.5(t－t0)2/δ2)充当激励源，观察3 600个网格点经过60个步长时间的运算结果和运算所需要的时间。图5 就是仿真所得部分结果，端口data_out 输出的是EZ RAM 中的值，图中可以看到输出的值是符合IEEE 754 标准的32 位单精度浮点数，由于采用了流水线设计，每个时钟输出一个运算结果。图6 是3 600个点在经过60 步后仿真结果。

为了比较FPGA 计算的结果，我们对同样计算用Visual C++6.0 编写了软件程序，结果如图7，两者误差不超过0.11%。FPGA 设计的最高工作频率可达118.39 MHz，完成整个运算过程需670 000个时钟周期。假设FPGA 工作频率为10 MHz，则完成全部运算需要0.067 s，而软件程序在主频2.2 GHz的个人电脑上运行时间是1.48 s，因此FPGA 计算要快几乎22倍。当采用2倍并行计算时，完成全部运算需368 310个时钟周期，加速比可达40倍。因此FPGA 计算对提高FDTD 算法速度的效果是十分明显的。

图5 ModelSim 仿真结果

图6 FPGA 运算结果

图7 C 语言运算结果

4 结论

本文对FDTD 算法的FPGA 实现进行了初步研究。以两维情形为例，编写了Verilog HDL 实现代码，在Altera FPGA 上进行了编译仿真，获得了初步试验结果。在设计中采用了流水线技术、并行运算等手段，同时采用了FPGA 片内RAM 作为存储单元，取得了较好的加速效果。随着制造工艺的发展，FPGA 片内资源会越来越丰富，基于FPGA 的计算可能成为提高FDTD 算法性能的一条有效途径。

［1］Yee K S.Numerical Solution of Initial Boundary Value Problems Involving Maxwell's Equation in Isotropic Media［J］.IEEE Trans Antennas Propagate，1966，(14):302－307.

［2］丁伟.时域有限差分法关键技术及其应用研究［D］.西安电子科技大学，2007.

［3］Schneider R N，Turner L E，Okoniewski M M.Application of FPGA Technology to Accelerate the Finite-Difference Time-Domain(FDTD)Method［C］//The 2002 ACM International Symposium on Field-Programmable Gate Arrays(FPGA'02)，2002:24－26.

［4］Durbano J P，Ortiz F E，Humphrey J R.Hardware Implementation of a Three-Dimensional Finite-Difference Time-Domain Algorithm［C］//IEEE Antennas and Wireless Propagation Letters，2003:2，54－57.

［5］Wu S.An FPGA Implementation of FDTD Codes for Reconfigurable High Performance Computing［D］.M.S.Thesis，Dept.of Electrical and Computer Engineering and Computer Science，University of Cincinnati，2005.

［6］Kawaguchi H，Fujita Y，Fujishima Y.Improved Architecture of FDTD/FIT Dedicated Computer for Higher Performance Computation［J］.IEEE Transactions on Magnetics，2008，44(6):1226－1229.

［7］Fujita Y，Kawaguchi H.Full Custom PCB Implementation of FDTD/FIT Dedicated Computer［J］.IEEE Transactions on Magnetics，2009，45(3):1100－1103.

［8］张煊，张明德.二维光子晶体波导及耦合器的特性模拟［J］.电子器件，2005，28(6):63－67.

［9］葛德彪，闰玉波.电磁波时域有限差分方法［M］.西安电子科技大学出版社，2002.

［10］Dennis M Sullivan.Electromagnetic Simulation Using the FDTD Method［M］.IEEE Microwave Theory and Techniques Society，2000.

［11］雷继兆.PC和服务器集群下的并行FDTD 算法及其应用研究［D］.西安电子科技大学，2009.

［12］刘瑜.FDTD 算法的网络并行研究及其电磁应用［D］.电子科技大学，2008.