CN111090958A - 一种基于亚网格技术的电磁波时域高效数值混合算法 - Google Patents

Abstract

本发明提出了一种基于亚网格技术的电磁波时域高效数值混合算法，分别在二维TE模式和TM模式下，将基于亚网格技术的时域有限差分方法与时域精细积分方法相结合，其中电小尺寸结构使用局部细网格剖分，其它计算区域保持较大尺寸的粗网格。然后在粗网格区域采用对内存需求较小的时域有限差分方法进行计算，在细网格处使用与粗网格时间步长相匹配的具有较大时间步长的时域精细积分方法进行计算。最后采用二阶泰勒展开式插值来实现粗网格区域和细网格区域计算信息的交换，从而分别发挥时域有限差分方法内存需求小以及时域精细积分方法可以选取较大时间步长的优势，减少迭代步数，降低计算机内存需求和CPU执行时间，提高计算效率，避免计算资源的浪费。

Description

一种基于亚网格技术的电磁波时域高效数值混合算法

技术领域

本发明属于计算电磁学领域，具体涉及一种电磁波时域高效数值混合算法。

背景技术

对于大部分复杂电磁波数值问题，特别是包含电小尺寸的结构或单元的电大尺寸(即计算区域为电大尺寸的，内部包含电小尺寸的结构或单元)的电磁波数值问题，传统的处理方式是以电小尺寸为参考剖分空间网格，但这样剖分的网格尺寸较小，会使网格数目十分庞大，内存需求较高。同时，对于电小尺寸结构以外的区域，由于剖分得过于密集，实际上是对计算资源的一种浪费。因此，对于内存需求较大的时域精细积分方法(Precise-Integration Time-Domain Method，PITD)很难实现。为了解决这个问题，采用亚网格技术对空间网格的进行剖分，在电小尺寸结构使用局部细网格，而在其它计算区域则保持较大尺寸的粗网格，这样可以有效的减少网格数目。然而，对于传统基于亚网格技术的时域有限差分方法(Finite-Difference Time-Domain Method，FDTD)，由于受到Courant稳定性条件的限制，时间步长必须根据细网格来选取，因此时间步长比较小，导致CPU执行时间仍然很长。

发明内容

为了克服现有技术的不足，本发明提出了一种基于亚网格技术的电磁波时域高效数值混合算法，分别在二维TE模式和TM模式下，将基于亚网格技术的时域有限差分方法与时域精细积分方法相结合，首先采用亚网格技术将计算区域按照细网格与粗网格进行剖分，其中电小尺寸结构使用局部细网格剖分，其它计算区域保持较大尺寸的粗网格。然后在粗网格区域采用对内存需求较小的时域有限差分方法进行计算，在细网格处使用与粗网格时间步长相匹配的具有较大时间步长的时域精细积分方法进行计算。最后采用二阶泰勒展开式插值来实现粗网格区域和细网格区域计算信息的交换，从而分别发挥时域有限差分方法内存需求小以及时域精细积分方法可以选取较大时间步长的优势，使得整个计算过程能够采用统一的较大时间步长进行计算。

为达到上述目的，本发明提出了一种基于亚网格技术的电磁波时域高效数值混合算法，包括以下步骤：

步骤1：数值区域内的空间网格剖分计算

针对包含有电小尺寸结构或单元的电大尺寸电磁波数值问题，采用亚网格技术，对电小尺寸计算区域部分使用细网格进行剖分，对其它计算区域使用粗网格进行剖分；其中，粗网格的尺寸为所计算电磁波数值问题的电磁波波长的1/10，细网格尺寸小于等于粗网格尺寸的1/2；

在粗网格与细网格区域内，电磁场量都满足如下式的麦克斯韦旋度方程：

其中，E为电场强度矢量，H为磁场强度矢量，ε为计算区域介质的介电常数，μ为计算区域介质的磁导率，σ为计算区域介质的电导率；

步骤2：二维TE模式下混合数值算法求解

步骤2-1：二维TE模式下，将麦克斯韦旋度方程式(1)在二维直角坐标系展开，如下式：

其中，E_x和E_y分别为二维TE模式下电场强度矢量在x方向和y方向的分量，H_z为二维TE模式下磁场强度矢量在z方向的分量；

步骤2-2：在粗网格区域内，应用时域有限差分方法进行计算：对式(2)中的空间偏微分算子和时间偏微分算子进行二阶精度中心差分离散，得到时域有限差分方法的离散迭代递推公式，如下式：

其中，

为第n+1个时间迭代步时电场强度在x方向分量的值，

为第n个时间迭代步时电场强度在x方向分量的值，

为第n+1个时间迭代步时电场强度在y方向分量的值，

为第n个时间迭代步时电场强度在y方向分量的值，i和j分别表示空间网格节点的行和列的序号，

为第n+1/2个时间迭代步时磁场强度在z方向分量的值，

为第n-1/2个时间迭代步时磁场强度在z方向分量的值，△x、△y为粗网格在x方向和y方向的空间步长，△t为离散迭代时间步长，n为迭代时间步数，ca、cb为求解系数：

步骤2-3：在细网格区域内，应用时域精细积分方法进行计算：对式(2)中的空间偏微分算子进行二阶精度中心差分离散，时间偏微分算子不做处理，得到时域精细积分方法的常微分方程组，如下式：

其中，E_xf和E_yf分别为二维TE模式下细网格区域内电场强度矢量在x方向和y方向的分量，H_zf为二维TE模式下细网格区域内磁场强度矢量在z方向的分量；

将常微分方程组式(6)-(8)统一写成矩阵形式，如下式：

其中，Y为包含计算区域内所有电磁场量的一维列向量，M为由空间步长和介质参数决定且不随时间变化的系数矩阵，f(t)为由激励源引入的一维列向量；

根据常微分方程理论，得到式(9)解的时域离散递推公式，如下式：

其中，k为第k个时间迭代步，t_k＝k△t(k＝0,1,2,...)，t_k+1＝(k+1)△t(k＝0,1,2,...)，Y_k+1为t_k+1时刻电磁场量的值Y((k+1)△t)，Y_k为t_k时刻电磁场量的值Y(k△t)，T为系数矩阵M的指数矩阵，T^k+1为指数矩阵T的k+1次幂，s为被积变量；

利用高斯积分公式对式(10)的右端积分项

进行近似，应用两点插值,构造高斯积分公式，得到离散迭代递推公式，如下式：

步骤2-4：粗网格区域与细网格区域计算信息交换；

粗网格与细网格边界处包含粗网格的电场强度矢量在y方向的分量E_y1和细网格的电场强度矢量在y方向的分量E_yf1、E_yf2；

E_y1的迭代公式为：

H_zf1、H_zf2、H_zf3和H_zf4分别为边界一侧的4个细网格处的磁场强度矢量在z方向的分量；

粗网格的电场强度矢量在x方向的分量E_x1的迭代公式为：

E_yf1的迭代公式为：

其中，E_y2、E_y3分别为E_y1上侧和下侧网格的电场强度分量，n_f为粗网格步长与细网格步长的比值；

E_yf2的迭代公式为：

步骤3：二维TM模式下混合数值算法求解

步骤3-1：在二维TM模式下，将麦克斯韦旋度方程式(1)在二维直角坐标系展开，如下式：

其中，H_x和H_y分别为二维TM模式下磁场强度矢量在x和y方向的分量，E_z为二维TM模式下电场强度矢量在z方向的分量。

步骤3-2：在粗网格区域内，应用时域有限差分方法进行计算：对式(16)所示偏微分方程组的空间偏微分算子和时间片微分算子进行二阶精度中心差分离散，得到时域有限差分方法的离散迭代递推公式，如下式：

其中，

为第n+1/2个时间迭代步时磁场强度在x方向分量的值，

为第n-1/2个时间迭代步时磁场强度在x方向分量的值，

为第n+1/2个时间迭代步时磁场强度在y方向分量的值，

为第n-1/2个时间迭代步时磁场强度在y方向分量的值，

为第n+1个时间迭代步时磁场强度在z方向分量的值，

为第n个时间迭代步时电场强度在z方向分量的值；

步骤3-3：在细网格区域内，应用时域精细积分方法进行计算：对式(16)所示偏微分方程组的空间偏微分算子进行二阶精度中心差分离散，时间偏微分算子不做处理，得到时域精细积分方法的常微分方程组，如下式：

其中，E_zf为二维TM模式下细网格区域内电场强度矢量在z方向的分量，H_xf和H_yf为二维TM模式下细网格区域内磁场强度矢量在x和y方向的分量；

利用高斯积分公式，最终得到离散迭代递推公式，如下式：

步骤3-4：粗网格区域与细网格区域计算信息的交换；

粗网格与细网格边界处包含粗网格的电场强度矢量在z方向的分量E_z1、E_z4和细网格的电场强度矢量在z方向的分量E_zf1；

E_z1的迭代公式为：

其中，H_x2、H_x3分别为E_z1上侧和下侧网格的磁场强度分量，H_y1为E_z1左侧网格的磁场强度分量；

E_z4的迭代公式为：

其中，H_y2、H_y3分别为E_z4左侧和右侧网格的磁场强度分量，H_x1为E_z4上侧网格的磁场强度分量，H_xf1、H_xf2分别为E_z4下侧的2个细网格上的磁场强度矢量在x方向的分量；

E_zf1的迭代公式为：

其中，E_z2、E_z3分别为E_z1上侧和下侧网格的电场强度分量。

本发明的有益效果是：由于采用了本发明的一种基于亚网格技术的电磁波时域高效数值混合算法，能够减少网格剖分数目，最大程度发挥时域有限差分方法与时域精细积分方法的优势，保证粗细网格区域内可使用统一的较大的时间步长进行计算，减少迭代步数，降低计算机内存需求和CPU执行时间，提高计算效率，避免计算资源的浪费；本发明中的时域精细积分方法采用两点高斯积分技术获得离散迭代递推公式，能够避免系数矩阵不可逆的问题，降低迭代过程的复杂度；本发明的基本思想与求解步骤具有普适性，适用于各类复杂结构电磁波问题的求解，实用性强。

附图说明

图1是本发明算法流程图。

图2是本发明二维TE模式下的空间网格剖分示意图。

图3是本发明二维TM模式下的空间网格剖分示意图。

图4是本发明在细网格时混合算法的计算流程图。

图5是本发明在粗网格时混合算法的计算流程图。

图6是本发明实施例1中二维TE模式下空间某点处电场强度E_y随时间的变化曲线图。

图7是本发明实施例1中二维TM模式下空间某点处电场强度E_y随时间的变化曲线图。

图8是本发明实施例2空间某点处电场强度E_y随时间的变化曲线图。

具体实施方式

下面结合附图和实施例对本发明进一步说明。

如图1所示，本发明提出了一种基于亚网格技术的电磁波时域高效数值混合算法，在二维TE模式和TM模式下，将基于亚网格技术的时域有限差分方法与时域精细积分方法相结合，首先采用亚网格技术将计算区域按照细网格与粗网格进行剖分，其中电小尺寸结构使用局部细网格剖分，其它计算区域保持较大尺寸的粗网格。然后在粗网格区域采用对内存需求较小的时域有限差分方法进行计算，在细网格处使用与粗网格时间步长相匹配的具有较大时间步长的时域精细积分方法进行计算。最后采用二阶泰勒展开式插值来实现粗网格区域和细网格区域计算信息的交换，从而分别发挥时域有限差分方法内存需求小以及时域精细积分方法可以选取较大时间步长的优势，使得整个计算过程能够采用统一的较大时间步长进行计算。

为达到上述目的，发明提出了一种基于亚网格技术的电磁波时域高效数值混合算法，包括以下步骤：

步骤1：数值区域内的空间网格剖分计算

针对包含有电小尺寸结构或单元的电大尺寸电磁波数值问题，采用亚网格技术，对电小尺寸计算区域部分使用细网格进行剖分，对其它计算区域使用粗网格进行剖分；其中，粗网格的尺寸为所计算电磁波数值问题的电磁波波长的1/10，细网格尺寸小于等于粗网格尺寸的1/2；

在粗网格与细网格区域内，电磁场量都满足如下式的麦克斯韦旋度方程：

其中，E为电场强度矢量，H为磁场强度矢量，ε为计算区域介质的介电常数，μ为计算区域介质的磁导率，σ为计算区域介质的电导率；

步骤2：二维TE模式下混合数值算法求解

步骤2-1：二维TE模式下，将麦克斯韦旋度方程式(1)在二维直角坐标系展开，如下式：

其中，E_x和E_y分别为二维TE模式下电场强度矢量在x方向和y方向的分量，H_z为二维TE模式下磁场强度矢量在z方向的分量；

步骤2-2：在粗网格区域内，应用时域有限差分方法进行计算：对式(2)中的空间偏微分算子和时间偏微分算子进行二阶精度中心差分离散，得到时域有限差分方法的离散迭代递推公式，如下式：

其中，

为第n+1个时间迭代步时电场强度在x方向分量的值，

为第n个时间迭代步时电场强度在x方向分量的值，

为第n+1个时间迭代步时电场强度在y方向分量的值，

为第n个时间迭代步时电场强度在y方向分量的值，i和j分别表示空间网格节点的行和列的序号，

为第n+1/2个时间迭代步时磁场强度在z方向分量的值，

为第n-1/2个时间迭代步时磁场强度在z方向分量的值，△x、△y为粗网格在x方向和y方向的空间步长，△t为离散迭代时间步长，n为迭代时间步数，ca、cb为求解系数：

步骤2-3：在细网格区域内，应用时域精细积分方法进行计算：对式(2)中的空间偏微分算子进行二阶精度中心差分离散，时间偏微分算子不做处理，得到时域精细积分方法的常微分方程组，如下式：

其中，E_xf和E_yf分别为二维TE模式下细网格区域内电场强度矢量在x方向和y方向的分量，H_zf为二维TE模式下细网格区域内磁场强度矢量在z方向的分量；

将常微分方程组式(6)-(8)统一写成矩阵形式，如下式：

其中，Y为包含计算区域内所有电磁场量的一维列向量，M为由空间步长和介质参数决定且不随时间变化的系数矩阵，f(t)为由激励源引入的一维列向量；

根据常微分方程理论，得到式(9)的解析解，如下式：

其中，Y(t)为Y在任意时刻t的值，Y(0)为Y在初始时刻的值，s为被积变量。

求解离散时间点t_k＝k△t(k＝0,1,2,...)上的电磁场量对Y_k进行计算，得到式(9)的时域离散递推公式，如下式：

其中，k为第k个时间迭代步，t_k＝k△t(k＝0,1,2,...)，t_k+1＝(k+1)△t(k＝0,1,2,...)，Y_k+1为t_k+1时刻电磁场量的值Y((k+1)△t)，Y_k为t_k时刻电磁场量的值Y(k△t)，T为系数矩阵M的指数矩阵，T^k+1为指数矩阵T的k+1次幂，s为被积变量；

指数矩阵T可以被重新写为下式：

T＝exp(M△t)＝[exp(Mτ)]^l (10-1)

其中，τ＝△t/l为子时间步长，l＝2^N，N为预定义的正整数；

预取正整数N的取值通常不低于5，子时间步长τ相比于△t是一个非常小的值，因此，exp(Mτ)可以使用4阶泰勒展开式进行高精度近似，如下式：

exp(Mτ)≈I+T_a (10-2)

其中，I为单位矩阵，T_a为泰勒展开式的1阶项到4阶项的求和；

将式(10-2)代入到式(10-1)中，得到：

T＝(I+T_a)^l (10-4)

对式(10-4)做如下分解得到：

同时，又有：

(I+T_a)²＝I+2T_a+T_a×T_a (10-5)

在实际计算过程中，将式(10-5)所示分解进行N次，只需对2T_a+T_a×T_a进行N次计算；计算过程由式(10-3)计算T_a开始，进行循环得到指数矩阵T；

利用高斯积分公式对式(10)的右端积分项

进行近似，应用两点插值,构造高斯积分公式，得到离散迭代递推公式，如下式：

步骤2-4：粗网格区域与细网格区域计算信息交换；

如图2所示，以n_f＝2为例，粗网格与细网格边界处包含粗网格的电场强度矢量在y方向的分量E_y1和细网格的电场强度矢量在y方向的分量E_yf1、E_yf2；

E_y1的迭代公式为：

H_zf1、H_zf2、H_zf3和H_zf4分别为边界一侧的4个细网格处的磁场强度矢量在z方向的分量；

粗网格的电场强度矢量在x方向的分量E_x1的迭代公式为：

采用二阶泰勒展开式插值求解E_yf1、E_yf2，得到E_yf1的迭代公式为：

其中，E_y2、E_y3分别为E_y1上侧和下侧网格的电场强度分量，n_f为粗网格步长与细网格步长的比值；

E_yf2的迭代公式为：

步骤3：二维TM模式下混合数值算法求解

步骤3-1：在二维TM模式下，将麦克斯韦旋度方程式(1)在二维直角坐标系展开，如下式：

其中，H_x和H_y分别为二维TM模式下磁场强度矢量在x和y方向的分量，E_z为二维TM模式下电场强度矢量在z方向的分量。

步骤3-2：在粗网格区域内，应用时域有限差分方法进行计算：对式(16)所示偏微分方程组的空间偏微分算子和时间片微分算子进行二阶精度中心差分离散，得到时域有限差分方法的离散迭代递推公式，如下式：

其中，

为第n+1/2个时间迭代步时磁场强度在x方向分量的值，

为第n-1/2个时间迭代步时磁场强度在x方向分量的值，

为第n+1/2个时间迭代步时磁场强度在y方向分量的值，

为第n-1/2个时间迭代步时磁场强度在y方向分量的值，

为第n+1个时间迭代步时磁场强度在z方向分量的值，

为第n个时间迭代步时电场强度在z方向分量的值；

步骤3-3：在细网格区域内，应用时域精细积分方法进行计算：对式(16)所示偏微分方程组的空间偏微分算子进行二阶精度中心差分离散，时间偏微分算子不做处理，得到时域精细积分方法的常微分方程组，如下式：

其中，E_zf为二维TM模式下细网格区域内电场强度矢量在z方向的分量，H_xf和H_yf为二维TM模式下细网格区域内磁场强度矢量在x和y方向的分量；

利用高斯积分公式，最终得到离散迭代递推公式，如下式：

步骤3-4：粗网格区域与细网格区域计算信息的交换；

如图3所示，以n_f＝2为例，粗网格与细网格边界处包含粗网格的电场强度矢量在z方向的分量E_z1、E_z4和细网格的电场强度矢量在z方向的分量E_zf1；

E_z1的迭代公式为：

其中，H_x2、H_x3分别为E_z1上侧和下侧网格的磁场强度分量，H_y1为E_z1左侧网格的磁场强度分量；

E_z4的迭代公式为：

其中，H_y2、H_y3分别为E_z4左侧和右侧网格的磁场强度分量，H_x1为E_z4上侧网格的磁场强度分量，H_xf1、H_xf2分别为E_z4下侧的2个细网格上的磁场强度矢量在x方向的分量；

E_zf1的迭代公式为：

其中，E_z2、E_z3分别为E_z1上侧和下侧网格的电场强度分量。

在实际计算过程中，工作流程分为两种情况：当激励源在细网格内部时，工作流程如图4所示；当激励源在粗网格内部时，其工作流程如图5所示。

实施例1：在一个40×40的自由空间计算区域，空间步长为△x＝△y＝0.001m。计算区域中心处的4×4个网格按照1:2的比例剖分为细网格，在计算区域四周设置完全匹配层吸收边界条件，厚度为10层。

图6所示为二维TE模式下，分别应用亚网格时域有限差分方法、本发明提出的混合方法以及时域有限差分方法计算实施例1时，空间某点处电场强度E_y随时间的变化曲线图。

图7所示为二维TM模式下，分别应用亚网格时域有限差分方法、本发明提出的混合方法以及时域有限差分方法计算实施例1时，空间某点处电场强度E_y随时间的变化曲线图。

实施例2：以TE模式为例，在一个201×201的自由空间计算区域，空间步长为△x＝△y＝0.001m。计算区域中心处的1个网格按照1:10的比例剖分为细网格，在计算区域四周设置完全匹配层吸收边界条件，厚度为10层。

图8所示为分别应用亚网格时域有限差分方法与本发明提出的混合方法计算实施例2时，空间某点处电场强度E_y随时间的变化曲线图。

Claims (1)

1.一种基于亚网格技术的电磁波时域高效数值混合算法，其特征在于，包括以下步骤：

步骤1：数值区域内的空间网格剖分计算

针对包含有电小尺寸结构或单元的电大尺寸电磁波数值问题，采用亚网格技术，对电小尺寸计算区域部分使用细网格进行剖分，对其它计算区域使用粗网格进行剖分；其中，粗网格的尺寸为所计算电磁波数值问题的电磁波波长的1/10，细网格尺寸小于等于粗网格尺寸的1/2；

在粗网格与细网格区域内，电磁场量都满足如下式的麦克斯韦旋度方程：

其中，E为电场强度矢量，H为磁场强度矢量，ε为计算区域介质的介电常数，μ为计算区域介质的磁导率，σ为计算区域介质的电导率；

步骤2：二维TE模式下混合数值算法求解

步骤2-1：二维TE模式下，将麦克斯韦旋度方程式(1)在二维直角坐标系展开，如下式：

其中，E_x和E_y分别为二维TE模式下电场强度矢量在x方向和y方向的分量，H_z为二维TE模式下磁场强度矢量在z方向的分量；

步骤2-2：在粗网格区域内，应用时域有限差分方法进行计算：对式(2)中的空间偏微分算子和时间偏微分算子进行二阶精度中心差分离散，得到时域有限差分方法的离散迭代递推公式，如下式：

其中，

为第n+1个时间迭代步时电场强度在x方向分量的值，

为第n个时间迭代步时电场强度在x方向分量的值，

为第n+1个时间迭代步时电场强度在y方向分量的值，

为第n个时间迭代步时电场强度在y方向分量的值，i和j分别表示空间网格节点的行和列的序号，

为第n+1/2个时间迭代步时磁场强度在z方向分量的值，

为第n-1/2个时间迭代步时磁场强度在z方向分量的值，△x、△y为粗网格在x方向和y方向的空间步长，△t为离散迭代时间步长，n为迭代时间步数，ca、cb为求解系数：

步骤2-3：在细网格区域内，应用时域精细积分方法进行计算：对式(2)中的空间偏微分算子进行二阶精度中心差分离散，时间偏微分算子不做处理，得到时域精细积分方法的常微分方程组，如下式：

其中，E_xf和E_yf分别为二维TE模式下细网格区域内电场强度矢量在x方向和y方向的分量，H_zf为二维TE模式下细网格区域内磁场强度矢量在z方向的分量；

将常微分方程组式(6)-(8)统一写成矩阵形式，如下式：

其中，Y为包含计算区域内所有电磁场量的一维列向量，M为由空间步长和介质参数决定且不随时间变化的系数矩阵，f(t)为由激励源引入的一维列向量；

根据常微分方程理论，得到式(9)解的时域离散递推公式，如下式：

其中，k为第k个时间迭代步，t_k＝k△t(k＝0,1,2,...)，t_k+1＝(k+1)△t(k＝0,1,2,...)，Y_k+1为t_k+1时刻电磁场量的值Y((k+1)△t)，Y_k为t_k时刻电磁场量的值Y(k△t)，T为系数矩阵M的指数矩阵，T^k+1为指数矩阵T的k+1次幂，s为被积变量；

利用高斯积分公式对式(10)的右端积分项

进行近似，应用两点插值,构造高斯积分公式，得到离散迭代递推公式，如下式：

步骤2-4：粗网格区域与细网格区域计算信息交换；

粗网格与细网格边界处包含粗网格的电场强度矢量在y方向的分量E_y1和细网格的电场强度矢量在y方向的分量E_yf1、E_yf2；

E_y1的迭代公式为：

H_zf1、H_zf2、H_zf3和H_zf4分别为边界一侧的4个细网格处的磁场强度矢量在z方向的分量；

粗网格的电场强度矢量在x方向的分量E_x1的迭代公式为：

E_yf1的迭代公式为：

其中，E_y2、E_y3分别为E_y1上侧和下侧网格的电场强度分量，n_f为粗网格步长与细网格步长的比值；

E_yf2的迭代公式为：

步骤3：二维TM模式下混合数值算法求解

步骤3-1：在二维TM模式下，将麦克斯韦旋度方程式(1)在二维直角坐标系展开，如下式：

其中，H_x和H_y分别为二维TM模式下磁场强度矢量在x和y方向的分量，E_z为二维TM模式下电场强度矢量在z方向的分量；

步骤3-2：在粗网格区域内，应用时域有限差分方法进行计算：对式(16)所示偏微分方程组的空间偏微分算子和时间片微分算子进行二阶精度中心差分离散，得到时域有限差分方法的离散迭代递推公式，如下式：

其中，

为第n+1/2个时间迭代步时磁场强度在x方向分量的值，

为第n-1/2个时间迭代步时磁场强度在x方向分量的值，

为第n+1/2个时间迭代步时磁场强度在y方向分量的值，

为第n-1/2个时间迭代步时磁场强度在y方向分量的值，

为第n+1个时间迭代步时磁场强度在z方向分量的值，

为第n个时间迭代步时电场强度在z方向分量的值；

步骤3-3：在细网格区域内，应用时域精细积分方法进行计算：对式(16)所示偏微分方程组的空间偏微分算子进行二阶精度中心差分离散，时间偏微分算子不做处理，得到时域精细积分方法的常微分方程组，如下式：

其中，E_zf为二维TM模式下细网格区域内电场强度矢量在z方向的分量，H_xf和H_yf为二维TM模式下细网格区域内磁场强度矢量在x和y方向的分量；

利用高斯积分公式，最终得到离散迭代递推公式，如下式：

步骤3-4：粗网格区域与细网格区域计算信息的交换；

粗网格与细网格边界处包含粗网格的电场强度矢量在z方向的分量E_z1、E_z4和细网格的电场强度矢量在z方向的分量E_zf1；

E_z1的迭代公式为：

其中，H_x2、H_x3分别为E_z1上侧和下侧网格的磁场强度分量，H_y1为E_z1左侧网格的磁场强度分量；

E_z4的迭代公式为：

其中，H_y2、H_y3分别为E_z4左侧和右侧网格的磁场强度分量，H_x1为E_z4上侧网格的磁场强度分量，H_xf1、H_xf2分别为E_z4下侧的2个细网格上的磁场强度矢量在x方向的分量；

E_zf1的迭代公式为：

其中，E_z2、E_z3分别为E_z1上侧和下侧网格的电场强度分量。

Images