CN109657288A - 一种三维显隐时域电磁学数值方法 - Google Patents

Abstract

本发明属于三维时域电磁学数值求解技术领域，涉及一种三维显隐时域电磁学数值方法。本发明通过一种新的隐式时域杂交间断伽辽金来替换现有方法中采用的显隐通用的隐式求解算法，以获得较少的全局未知量，大大减少隐式待求全局矩阵的维数。最终利用这种显隐算法以高效地求解大规模多尺度的三维时域Maxwell方程组，在保证高阶精度的前提下，提升计算性能。

Description

一种三维显隐时域电磁学数值方法

技术领域

本发明属于三维时域电磁学数值求解技术领域，涉及两种数值算法：三维时域电磁学间断伽辽金数值方法和时域杂交间断伽辽金数值方法，通过这两种方法的有效结合，形成一种三维显隐时域电磁学数值方法，特别地，非常适用于求解大规模多尺度问题。

背景技术

在复杂电磁环境中，采用时域方法能够更加快速、直观地反映出多尺度装备的电磁响应特性。在时域电磁仿真中，最消耗时间和内存的步骤是电磁场的求解，而电真空器件显著的大规模多尺度特征使得要实现电磁场的快速求解极其困难。目前在时域中出现一种新型的具有很大优势的算法，即时域间断伽辽金方法，该方法是一种基于有限体积方法和有限元方法而发展出来的一种方法，支持非结构、非共形网格，支持区域分解技术，对基函数的选取非常灵活。

尽管时域间断伽辽金方法具有灵活的空间和时间离散方式，然而由于多尺度网格、时间稳定性条件的影响，要实现电真空器件的快速电磁场求解仍然效率很低。由于电真空器件的显著多尺度特征，四面体网格区域的网格尺寸会有很大差异。此时若采用常规的显式时域间断伽辽金方法，受稳定性条件限制，大量存在的精细网格将会导致极小的时间步长，以致求解效率极低。若采用常规的隐式时域间断伽辽金方法，虽然隐式时间格式具有无条件稳定性，能得到更大的时间步长，缩减计算时间，但是隐式时间格式需要求解一个全局线性矩阵，对矩阵处理的要求会进一步增大。因此，对于这种大规模多尺度问题，目前性能较高的是一种显隐混合算法，即将整体求解区域按照网格尺寸划分为粗细网格，在粗网格上采用显式方法，在细网格上采用隐式算法。目前在电磁中应用较多的是基于时域间断伽辽金的显隐混合算法。

然而，时域间断伽辽金方法有一个致命的缺点，即单元交界面上的未知量是重复的，特别地，对于多尺度大规模问题，随着模型网格数量和基函数阶数的提高，隐式时间格式所生成的全局矩阵的维数会更大，并且随着模型的复杂化，该矩阵极有可能是高度病态的，而这进一步加剧了矩阵求解的难度。因此，常规的基于时域间断伽辽金的显隐混合算法，随着大规模多尺度问题的复杂化，低阶基函数的精度已经不能满足设计者的要求，而高阶基函数在仿真中导致内存消耗和迭代时间更多。因此，迫切需要一种高效的时域电磁数值方法来求解这种大规模多尺度的问题，在能保证高阶精度的前提下,提升计算性能。

发明内容

针对上述存在问题或不足，为实现高效时域电磁数值法求解大规模多尺度问题时，在保证高阶精度的前提下，提升计算性能。本发明提供了一种三维显隐时域电磁学数值方法，通过结合显式的时域间断伽辽金和隐式的时域杂交间断伽辽金这两种数值方法,以高效地求解大规模多尺度的三维时域Maxwell方程组，且具有较少的全局未知量、较高的精度以及显著的计算性能。

一种三维显隐时域电磁学数值方法，包括以下步骤：

步骤A、根据目标电子器件的物理结构，结合工作环境与边界条件对其仿真建模；

步骤B、采用四面体单元剖分三维求解区域，且根据网格尺寸分为粗网格集合与细网格集合面离散和体积离散必须相容。

步骤C、在粗网格集合上，采用显式时域间断伽辽金方法得到半离散格式；在粗网格集合上，定义电场和磁场[E_e,H_e]^Τ的待求系数为

结合(1)式

采用时域间断伽辽金方法，得到半离散格式如下：

其中ε为介质相对介电常数，μ为介质的相对磁导率，和是质量矩阵，是刚度矩阵，是通量矩阵。

步骤D、而细网格的网格尺寸较小，在细网格集合上，采用隐式时域杂交间断伽辽金方法得到半离散格式；

在细网格集合上，定义电场和磁场W_i＝[E_i,H_i]^Τ的待求系数为时域杂交间断伽辽金法是通过在单元交界面上引入杂交量Λ_i来定义电磁场的数值通量，并需要增加第三个守恒方程来保证全局系统的稳定性。其具体的公式推导与频域杂交间断伽辽金法类似，是一种公知过程，这里不再阐述。下面只给出三维时域杂交间断伽辽金数值方法的半离散形式：

其中n是在细网格集合的计算区域边界上的外法向单位矢量，τ＞0是局部稳定系数。表示吸收边界条件的边界方程，和是测试函数，是切向磁场，其中间断有限元函数空间与有限元迹空间是一种公知的空间，因此其具体形式不再详细描述。

步骤E.整个求解区域的时间离散，总共分为三步，从而得到完全离散格式；

根据时域显隐算法，在整个计算区域中将时域电磁场值(E,H)进一步划分为 定义E的待求系数为H的待求系数为定义时间步长Δt，将总仿真时间[0,T]离散为等间隔的时间步其中n的最大值为N_t，即时间迭代总次数。令t_n时刻的电磁场(E(t_n),H(t_n))＝(Eⁿ,Hⁿ)，t_n+1时刻的电磁场(E(t_n+1),H(t_n+1))＝(Eⁿ⁺¹,Hⁿ⁺¹)。

由于粗网格集合的网格尺寸较大，因此采用显式时间格式的时域间断伽辽金方法。而细网格集合的网格尺寸较小，如果仍然采用显式时间格式，这时会受到稳定性条件的限制而大大增加迭代次数，最终导致计算时间很大。为了减少隐式时间格式产生的全局矩阵维数，本发明在上采用时域杂交间断伽辽金数值方法，通过单元交界面上的杂交量，导出只与杂交量有关的全局矩阵。一旦求解杂交量，则每个单元的电磁场值就能获得。下面从时间离散来概述显隐混合算法的三部曲：

(1)在粗网格集合上，电磁场值的更新时间是t_n到t_n+1/2，时间步长为Δt/2，采用二阶显式的蛙跳格式；

(2)在细网格集合上，电磁场值的更新时间是t_n到t_n+1，时间步长为Δt，采用二阶隐式的Crank-Nicolson时间格式；

(3)在粗网格集合上，电磁场值的更新时间是t_n+1/2到t_n+1，时间步长为Δt/2，采用二阶显式的蛙跳格式。

步骤F.根据步骤E的完全离散格式，进行迭代，分别求出每个时刻的电场和磁场。

根据步骤E的完全离散格式，在计算区域上，按照粗细网格的划分，采用不同的方法。即粗网格采用显式时域间断伽辽金方法，细网格采用隐式时域杂交间断伽辽金方法。特别的，在每一次时间迭代上，电磁场值的更新，需要按照步骤E分为三部曲，先是显式，然后是隐式，最后是显式。最终以求出每个时刻的电场和磁场。

相比现有的时域显隐算法，本发明提供了一种新的隐式时域杂交间断伽辽金来替换现有方法中采用的显隐通用的隐式求解算法，以获得较少的全局未知量，大大减少隐式待求全局矩阵的维数。利用这种显隐算法以高效地求解大规模多尺度的三维时域Maxwell方程组，在保证高阶精度的前提下，提升计算性能。

附图说明

图1是本发明的流程图。

具体实施方式

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

参照图1，一种三维显隐时域电磁学数值方法，包括以下步骤：

步骤A、根据目标电子器件的物理结构，结合工作环境与边界条件对其仿真建模；

本发明以无源三维时域Maxwell方程组为例，进行说明。首先给出无源三维时域Maxwell方程组如下所示：

其中，T表示计算电磁学时域中模型仿真计算的最终时间，E(x,y,z,t)和H(x,y,z,t)分别是电场强度和磁场强度矢量。计算区域Ω为三维，即

步骤B、采用四面体单元剖分三维求解区域，且根据网格尺寸分为粗网格集合与细网格集合

采用四面体单元剖分三维求解区域是一种公知过程，因此本步骤不再详细描述。需要注意的是，面离散和体积离散必须相容。本实施将计算区域Ω划分成N_h个四面体网格的集合其中每个体单元用τ_i(i＝1,2,3,…,N_h)表示，即

由于本发明采用显隐时域数值方法，这里需要根据实际模型的网格尺寸将整个求解区域进一步划分为粗网格集合与细网格集合即且由于粗网格的网格尺寸较大，因此在四面体集合上，采用显式时域间断伽辽金方法；而细网格的网格尺寸较小，因此在四面体集合上，采用隐式杂交时域间断伽辽金方法。隐式杂交时域间断伽辽金方法是通过在单元交界面上引入新的杂交量来定义电磁场的数值通量，使得全局矩阵只由该杂交量导出，由于杂交量只存在于单元交界面上，并保持单值，因此相比常规的隐式时域间断伽辽金方法，本发明所采用的隐式杂交时域间断伽辽金方法能够极大地减少需要求解的全局矩阵维数，自由度大大降低，求解效率更高。下面我们定义中的面集合是由N_f个三角形面单元D^f组成，即

步骤C、在粗网格集合上，采用时域间断伽辽金方法得到半离散格式；

在粗网格集合上，定义电场和磁场[E_e,H_e]^Τ的待求系数为结合(1)式，采用时域间断伽辽金方法，得到半离散格式如下：

其中ε为介质相对介电常数，μ为介质的相对磁导率，和是质量矩阵，是刚度矩阵，是通量矩阵。由于半离散格式的推导是一种公知过程，因此本步骤不再详细描述。

步骤D、在细网格集合上，采用时域杂交间断伽辽金方法得到半离散格式；

在细网格集合上，定义电场和磁场W_i＝[E_i,H_i]^Τ的待求系数为时域杂交间断伽辽金法是通过在单元交界面上引入杂交量Λ_i来定义电磁场的数值通量，并需要增加第三个守恒方程来保证全局系统的稳定性。其具体的公式推导与频域杂交间断伽辽金法类似，是一种公知过程，这里不再阐述。下面只给出三维时域杂交间断伽辽金数值方法的半离散形式：

其中n是在细网格集合的计算区域边界上的外法向单位矢量，τ＞0是局部稳定系数。表示吸收边界条件的边界方程，和是测试函数，是切向磁场，其中间断有限元函数空间与有限元迹空间是一种公知的空间，因此其具体形式不再详细描述。

步骤E.整个求解区域的时间离散，总共分为三步，从而得到完全离散格式；

根据时域显隐算法，在整个计算区域中将时域电磁场值(E,H)进一步划分为 定义E的待求系数为H的待求系数为定义时间步长Δt，将总仿真时间[0,T]离散为等间隔的时间步t_n＝nΔt,其中n的最大值为N_t，即时间迭代总次数。令t_n时刻的电磁场(E(t_n),H(t_n))＝(Eⁿ,Hⁿ)，t_n+1时刻的电磁场(E(t_n+1),H(t_n+1))＝(Eⁿ⁺¹,Hⁿ⁺¹)。由于粗网格集合的网格尺寸较大，采用简单的显式时间格式会使得计算更快，并且结合时域间断伽辽金方法更容易实现并行性。而细网格集合的网格尺寸较小，如果仍然采用显式时间格式，这时会受到稳定性条件的限制而大大增加迭代次数，最终导致计算时间很大。因此，对于这种细网格，采用无条件稳定的隐式时间格式是最合适的。为了减少隐式时间格式产生的全局矩阵维数，本发明在上采用时域杂交间断伽辽金数值方法，通过单元交界面上的杂交量，导出只与杂交量有关的全局矩阵。一旦求解杂交量，则每个单元的电磁场值就能获得。下面从时间离散来概述显隐混合算法的三部曲：

(1)在粗网格集合上，将电磁场值的更新时间是t_n到t_n+1/2，时间步长为Δt/2，采用二阶显式的蛙跳格式；

(2)在细网格集合上，将电磁场值的更新时间是t_n到t_n+1，时间步长为Δt，采用二阶隐式的Crank-Nicolson时间格式；

(3)在粗网格集合上，将电磁场值的更新时间是t_n+1/2到t_n+1，时间步长为Δt/2，采用二阶显式的蛙跳格式。

对于(1)与(3)的完全离散格式，可以从(4)式出发，分别得到如下的结果：

其中是t_n+1/2时刻的电磁场的待求系数，是t_n+1时刻的电磁场的待求系数。

由于②式是基于时域杂交间断伽辽金方法，从半离散格式(5)式出发，结合二阶隐式的Crank-Nicolson时间格式，即 从而可以推导出全离散格式。需要注意的是：时域杂交间断伽辽金的全离散格式包括两个系统：局部线性系统和全局线性系统。具体过程如下：从(5)式的前两式出发，我们可以将细网格集合上的每个四面体单元的电场和磁场用杂交量线性表出，从而形成局部线性系统。考虑到杂交量只存在于单元交界面上，并且保持单值，因此依据(5)式的第三个守恒方程，将局部线性系统中的杂交量，按照面单元进行集成，从而得到最终的只含杂交量的全局线性系统，即这里是全局线性矩阵，Λ是杂交量Λ_i的待求未知系数，y是右端项。一旦求解出杂交量，通过全局面编号与体单元内的局部面编号的一一对应信息，可以得出每一个体单元内的第n+1时刻的四个面上的杂交量系数值，进一步由局部线性系统求出每个体单元的电场与磁场

以上分析是对于时域显隐算法中该四面体与之相邻的四面体同属于粗网格集合或者细网格集合特别地，对于显隐交界的四面体，例如某个四面体属于粗网格集合但是其相邻的四面体有一个属于细网格集合那么在数值通量矩阵与电磁场相乘时，需要使用相邻四面体的上一时刻的电磁场值，即(6)与(7)式等式右端不一定全是 也可能是含有

步骤F.根据步骤E的完全离散格式，在计算区域上，粗网格采用显式时域间断伽辽金方法，细网格采用隐式时域杂交间断伽辽金方法，进行迭代；在每一次时间迭代上，电磁场值的更新，按照步骤E分为三部曲，先是显式，然后是隐式，最后是显式。最终以求出每个时刻的电场和磁场。

以边长为1m的正立方体为模型进行分析其平面波时域传播特性。模型的网格剖分，是从立方体的中心处，局部加密，由2968个四面体单元构成整个模型。四面体网格最大与最小尺寸之比为110.9，符合多尺度的特性。通过仿真10个周期的一阶时域平面波，得到如下的计算结果：

表1四种算法仿真10个周期平面波在一阶的计算性能

	CFL	CPU(time)	加速比	峰值内存	L2范数误差
						exDGTD	0.3	6419s	—	98Mb	6.14e-02
imHDGTD	38.4	205s	31.3	163Mb	7.08e-02
						eximDGTD	38.4	226s	28.4	211Mb	7.12e-02
eximHDGTD	38.4	96s	66.9	94Mb	7.02e-02

注：表1中的exDGTD表示，显式的时域间断伽辽金方法，imHDGTD表示隐式的时域杂交间断伽辽金方法，eximDGTD表示常规的显隐的时域间断伽辽金方法，eximHDGTD表示本发明中的显隐算法(即显式采用时域间断伽辽金方法，隐式采用时域杂交间断伽辽金方法)，CFL表示时间步长的稳定性条件的取值。

从表1中，可以看出eximHDGTD能够在相比显式的exDGTD，时间步长能够扩大128倍，同时具有较好的L2范数误差，同时其时间加速比能够达到66.9，计算性能最好。该实施例，进一步说明，本发明的方法能够在处理多尺度复杂问题中减少计算时间，大大提高计算性能。

Claims (2)

1.一种三维显隐时域电磁学数值方法，包括以下步骤：

步骤A、根据目标电子器件的物理结构，结合工作环境与边界条件对其仿真建模；

步骤B、采用四面体单元剖分三维求解区域，且根据网格尺寸分为粗网格集合与细网格集合面离散和体积离散必须相容；

步骤C、在粗网格集合上，采用显式时域间断伽辽金方法得到半离散格式；在粗网格集合上，定义电场和磁场[E_e,H_e]^Τ的待求系数为

结合(1)式

采用时域间断伽辽金方法，得到半离散格式如下：

其中ε为介质相对介电常数，μ为介质的相对磁导率，和是质量矩阵，是刚度矩阵，是通量矩阵；

步骤D、在细网格集合上，采用隐式时域杂交间断伽辽金方法得到半离散格式；

定义电场和磁场W_i＝[E_i,H_i]^Τ的待求系数为时域杂交间断伽辽金法是通过在单元交界面上引入杂交量Λ_i来定义电磁场的数值通量，并需要增加第三个守恒方程来保证全局系统的稳定性；三维时域杂交间断伽辽金数值方法的半离散形式：

其中n是在细网格集合的计算区域边界上的外法向单位矢量，τ＞0是局部稳定系数；表示吸收边界条件的边界方程，和是测试函数，是切向磁场，为间断有限元函数空间，为有限元迹空间；

步骤E.整个求解区域的时间离散，总共分为三步，从而得到完全离散格式；

根据时域显隐算法，在整个计算区域中将时域电磁场值(E,H)进一步划分为 定义E的待求系数为H的待求系数为定义时间步长Δt，将总仿真时间[0,T]离散为等间隔的时间步t_n＝nΔt,其中n的最大值为N_t，即时间迭代总次数；令t_n时刻的电磁场(E(t_n),H(t_n))＝(Eⁿ,Hⁿ)，t_n+1时刻的电磁场(E(t_n+1),H(t_n+1))＝(Eⁿ⁺¹,Hⁿ⁺¹)；

在上采用时域杂交间断伽辽金数值方法，通过单元交界面上的杂交量，导出只与杂交量有关的全局矩阵；

从时间离散来概述显隐混合算法的步骤：

(1)在粗网格集合上，电磁场值的更新时间是t_n到t_n+1/2，时间步长为Δt/2，采用二阶显式的蛙跳格式；

(2)在细网格集合上，电磁场值的更新时间是t_n到t_n+1，时间步长为Δt，采用二阶隐式的Crank-Nicolson时间格式；

(3)在粗网格集合上，电磁场值的更新时间是t_n+1/2到t_n+1，时间步长为Δt/2，采用二阶显式的蛙跳格式；

步骤F.根据步骤E的完全离散格式，在计算区域上，按照粗细网格的划分：粗网格采用显式时域间断伽辽金方法，细网格采用隐式时域杂交间断伽辽金方法，进行迭代，分别求出每个时刻的电场和磁场。

2.如权利要求1所述三维显隐时域电磁学数值方法，其特征在于：所述步骤F在每一次时间迭代上，电磁场值的更新，需要按照步骤E分为三步：先是显式，然后是隐式，最后是显式；最终求出每个时刻的电场和磁场。