CN106991222B - 一种基于叠层矩阵分解的低频电磁特性仿真方法 - Google Patents

一种基于叠层矩阵分解的低频电磁特性仿真方法 Download PDF

Info

Publication number
CN106991222B
CN106991222B CN201710186894.1A CN201710186894A CN106991222B CN 106991222 B CN106991222 B CN 106991222B CN 201710186894 A CN201710186894 A CN 201710186894A CN 106991222 B CN106991222 B CN 106991222B
Authority
CN
China
Prior art keywords
matrix
decomposition
low
charge
laminated
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN201710186894.1A
Other languages
English (en)
Other versions
CN106991222A (zh
Inventor
宛汀
洪涛
时艳玲
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Nanjing University of Posts and Telecommunications
Original Assignee
Nanjing University of Posts and Telecommunications
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Nanjing University of Posts and Telecommunications filed Critical Nanjing University of Posts and Telecommunications
Priority to CN201710186894.1A priority Critical patent/CN106991222B/zh
Publication of CN106991222A publication Critical patent/CN106991222A/zh
Application granted granted Critical
Publication of CN106991222B publication Critical patent/CN106991222B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/16—Matrix or vector computation, e.g. matrix-matrix or matrix-vector multiplication, matrix factorization
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2113/00—Details relating to the application field
    • G06F2113/18—Chip packaging

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Mathematical Physics (AREA)
  • Computational Mathematics (AREA)
  • Mathematical Optimization (AREA)
  • Data Mining & Analysis (AREA)
  • Pure & Applied Mathematics (AREA)
  • General Engineering & Computer Science (AREA)
  • Mathematical Analysis (AREA)
  • Computing Systems (AREA)
  • Computer Hardware Design (AREA)
  • Evolutionary Computation (AREA)
  • Algebra (AREA)
  • Geometry (AREA)
  • Databases & Information Systems (AREA)
  • Software Systems (AREA)
  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)

Abstract

本发明公开了一种基于叠层矩阵分解的低频电磁特性仿真方法，该方法首先采用增量型电场积分方程方法来消除传统电场积分方程方法面临的低频崩溃问题；然后针对增量型电场积分方程的系数矩阵，采用低频多层快速多极子技术来构建其叠层矩阵表达式；进而采用一种重压缩方法对该叠层矩阵进一步压缩，去除冗余信息；最后采用叠层矩阵格式算法对压缩后的叠层矩阵进行上下三角分解，该叠层矩阵分解能够将计算复杂度和内存消耗分别降低到几乎线性，不仅为迭代解法构造了一种预条件技术，而且还构造了一种直接解法。本发明所述方法具有求解速度快、内存消耗低和精度可控的特点，可用于各种低频电磁特性分析。

Description

一种基于叠层矩阵分解的低频电磁特性仿真方法
技术领域
本发明涉及电磁仿真技术领域，特别是一种基于叠层矩阵分解的低频电磁特性仿真方法。
背景技术
随着微电子技术和制造工艺的发展，在计算机芯片封装设计、高速集成电路设计、天线设计、目标电磁散射特性分析、电磁兼容分析等领域中，出现了越来越多的具有复杂精细结构的目标。这些目标的尺寸远小于波长或几何建模剖分尺寸远小于波长，均属于低频问题的范畴。对于低频问题的分析，准静态方法随着工作频率的升高的误差将越来越大。因此，电磁数值仿真方法在低频问题的分析中发挥着越来越重要的作用。表面积分方程方法是一类具有代表性的电磁数值仿真方法，它只需离散目标表面，因此未知量少、求解效率高。然而，传统的表面积分方程方法在分析低频问题时会出现数值不稳定的情况，让低频电磁仿真无法获得精确解，这一现象被称为“低频崩溃”。许多方法致力于解决低频崩溃问题，例如：基函数Helmholtz分解方法和Calderón预条件方法等。
其中，增量型电场积分方程方法(Augmented Electric Field IntegralEquation,A-EFIE)是一种极其有效的消除低频崩溃的方法(Z.G.Qian and W.C.Chew.Anaugmented EFIE for high speed interconnect analysis.Microwave and OpticalTechnology Letter,vol.50,no.10,pp.2658-2662,2008)。A-EFIE通过在电场积分方程中引入电流连续性方程，对电场积分方程中不平衡的谱分支进行了归一化，从而有效解决了低频崩溃问题。尽管A-EFIE分离了电流矢量和电荷标量，A-EFIE系统矩阵仍然是高度稠密的。迭代解法是求解A-EFIE系统矩阵方程的常用方法，原因在于，迭代解法的主要操作是矩阵-矢量乘，其计算复杂度低且能够通过引入快速算法(如：快速多极子算法)进行加速。然而，迭代解法在面临病态的系统矩阵时，迭代收敛所需步数增多；预条件技术能够在一定程度上提高迭代解法的收敛速度，其中，约束(Constraint)预条件是低频电磁仿真中最为常用的一种预条件技术(Z.G.Qian and W.C.Chew.Fast full-wave surface integralequation solver for multiscale structure modeling.IEEE Transactions onAntennas and Propagation,vol.50,pp.3594-3601,2009)。此外，迭代解法在求解多激励问题时，对不同的右边向量均需重新计算，上述情况都会导致迭代解法效率低下。与迭代解法不同，直接解法不仅不存在收敛速度的问题，而且对于多右边向量问题无需重复计算。然而，由于直接解法巨大的计算消耗，采用直接解法来求解大规模低频问题通常被认为是不可行的。因此，发展一种针对低频电磁问题的高效求解方法具有重要的应用价值。
发明内容
本发明所要解决的技术问题是克服现有技术的不足而提供一种基于叠层矩阵分解的低频电磁特性仿真方法，该方法不仅消除了低频崩溃现象，而且具有求解速度快、内存消耗低和精度可调的特点。
本发明为解决上述技术问题采用以下技术方案：
根据本发明提出的一种基于叠层矩阵分解的低频电磁特性仿真方法，包括以下步骤：
第1步，针对目标表面建立电场积分方程，并引入电流连续性方程，采用三角形单元离散目标表面，用RWG基函数展开表面电流，用冲击基函数展开表面电荷，生成增量型电场积分方程的系统矩阵方程；
第2步，利用矩阵收缩方法或者缩减电荷自由度方法对该系统矩阵方程进行修正，消除电荷中性原理在极低频率时造成的矩阵奇异性；
第3步，建立增量型电场积分方程系统矩阵的叠层矩阵；
第4步，采用叠层矩阵格式的运算法则对该叠层矩阵进行上下三角分解；
第5步，将叠层矩阵的上下三角分解与前后向回代结合，采用预条件的迭代解法或者直接解法进行求解，获得目标表面的电流及电荷分布，通过计算提取所需的电磁特性参数。
作为本发明所述的一种基于叠层矩阵分解的低频电磁特性仿真方法进一步优化方案，第3步中，建立增量型电场积分方程系统矩阵的叠层矩阵的步骤如下：
第(1)步，分别针对定义在棱边上的RWG基函数簇和定义在面片上的冲击基函数簇，采用几何递归二分的方法，构造电流簇树和电荷簇树；
第(2)步，电流簇树和电荷簇树相互作用，构造电流-电流块树、电流-电荷块树、电荷-电流块树和电荷-电荷块树，并引入容许条件来划分容许块和非容许块，其中，非容许块以满阵形式表示，容许块以低秩分解矩阵形式表示；
第(3)步，将系统矩阵中各子矩阵的非零元素分别填入对应的块树中，其中，对于电流子矩阵和电荷子矩阵，采用矩量法构造非容许块，采用低频多层快速多极子技术构造容许块，并利用基于QR分解和奇异值分解的重压缩方法进一步对容许块进行压缩；对于其余子矩阵，将其非零元素对应填入非容许块或容许块，生成系统矩阵的叠层矩阵。
作为本发明所述的一种基于叠层矩阵分解的低频电磁特性仿真方法进一步优化方案，第(3)步中，采用低频多层快速多极子技术构造容许块具体包括，提取低频多层快速多极子的聚合、转移和配置矩阵，将容许块表示成低秩分解矩阵的形式，并通过调节多极子数目来控制低秩分解的精度。
作为本发明所述的一种基于叠层矩阵分解的低频电磁特性仿真方法进一步优化方案，第4步中，叠层矩阵格式的运算法则中涉及低秩分解矩阵运算，通过调节低秩分解矩阵运算中的截断精度来控制叠层矩阵上下三角分解的精度，高精度的分解用于构造直接解法，低精度的分解用于构造预条件来加速迭代解法的收敛，且不同精度的分解具有不同的加速效果。
作为本发明所述的一种基于叠层矩阵分解的低频电磁特性仿真方法进一步优化方案，叠层矩阵上下三角分解的精度通过分解前后矩阵之差的二范数和分解前矩阵的二范数之比来定义，高精度指达到实际仿真精度要求的精度，低精度指低于实际仿真精度要求的精度。
本发明采用以上技术方案与现有技术相比，具有以下技术效果：
(1)本发明提出的基于叠层矩阵分解的低频电磁特性仿真方法不仅为低频电磁问题的迭代解法提供一种新型预条件技术，而且为低频电磁问题的求解提供一种直接解法，显著提高现有低频电磁特性仿真方法的性能；
(2)与现有流行的约束预条件技术相比，本发明提出的基于叠层矩阵分解的预条件技术能够显著提高迭代收敛速度，有效解决在系统矩阵性态差时迭代收敛慢甚至不收敛的问题；
(3)本发明提出的基于叠层矩阵分解的直接解法能够将计算复杂度降低到O(Nlog2N)，将内存消耗降低到O(NlogN)(其中，N为未知量的数目)，为低频电磁问题的直接求解提供一种可行途径。
附图说明
图1是基于叠层矩阵分解的低频电磁特性仿真方法流程图。
图2是递归二分定界盒子方法示意图。
图3是电流簇树示意图。
图4是两个簇树相互作用示意图。
图5是电流-电流块树示意图。
图6a是低频单层快速多极子算法转移过程示意图。
图6b是低频多层快速多极子算法转移过程示意图。
图7是飞机目标模型及电流分布图。
图8是不同预条件技术的迭代收敛曲线图。
图9是叠层矩阵分解的计算复杂度曲线图。
图10是叠层矩阵分解的内存消耗曲线图。
具体实施方式
下面结合附图对本发明的技术方案做进一步的详细说明：
图1给出了本发明所述方法的主要流程，下面结合图1对本发明所述方法的步骤作进一步详细描述：
第1步，针对目标表面建立电场积分方程，
Figure BDA0001255066610000031
其中，i为根号负一，ω为角频率，μr为相对磁导率，εr为相对介电常数，
Figure BDA0001255066610000032
和
Figure BDA0001255066610000033
分别表示源点和观察点的位置矢量，
Figure BDA0001255066610000041
为格林函数，
Figure BDA0001255066610000042
表示观察点处的表面电流，S'表示目标表面，▽'为梯度算子，
Figure BDA0001255066610000043
表示入射电场矢量。然后，引入描述表面电流和表面电荷之间关系的电流连续性方程，
Figure BDA0001255066610000044
其中，
Figure BDA0001255066610000045
表示观察点处的表面电荷。
采用三角形单元离散目标表面，用RWG基函数来展开表面电流，用冲击基函数来展开表面电荷，采用伽辽金方法测试后，生成增量型电场积分方程的系统矩阵方程：
Figure BDA0001255066610000046
其中，
Figure BDA0001255066610000047
表示单位阵，上标T表示矩阵转置，k0表示自由空间波数，c0表示光速，η0表示自由空间波阻抗，jm和ρm分别表示未知的电流系数和电荷系数，矩阵
Figure BDA0001255066610000048
的元素
Figure BDA0001255066610000049
的元素
Figure BDA00012550666100000410
和右边向量
Figure BDA00012550666100000411
的元素bm具有如下形式：
Figure BDA00012550666100000412
Figure BDA00012550666100000413
Figure BDA00012550666100000414
其中，下标m和n分别表示矩阵的行编号和列编号，
Figure BDA00012550666100000415
和
Figure BDA00012550666100000416
分别表示第m个和第n个RWG基函数，hm(·)和hn(·)分别表示第m个和第n个冲击基函数，bm为入射波和第m个RWG基函数作用生成的右边向量元素，S表示目标的积分表面。矩阵
Figure BDA00012550666100000417
是一个稀疏的关联矩阵，描述了单元和棱边之间的相互关系，其矩阵元素
Figure BDA00012550666100000418
的表达式如下：
Figure BDA00012550666100000419
第2步，利用矩阵收缩方法或者缩减电荷自由度方法对该系统矩阵方程进行修正，消除电荷中性原理在极低频率时造成的矩阵奇异性。
对于矩阵收缩方法，首先将系统矩阵方程变换为：
Figure BDA00012550666100000420
令方程(8)的系数矩阵为
Figure BDA00012550666100000421
进而将
Figure BDA00012550666100000422
修正为
Figure BDA00012550666100000423
Figure BDA00012550666100000424
其中，
Figure BDA00012550666100000425
trace表示矩阵的迹，e表示单元离散生成的内边数目，p表示单元离散生成的三角形面片数目。；
Figure BDA0001255066610000051
是一个归一化的向量：
Figure BDA0001255066610000052
其中，
Figure BDA0001255066610000053
的前e个元素是0，后p个元素是
Figure BDA0001255066610000054
这样，基于矩阵收缩方法的增量型电场积分方程的系统可写为：
Figure BDA0001255066610000055
对于缩减电荷自由度方法，只需在目标的每个独立部分中去除一个电荷未知量，即在每个独立部分让表示电荷的面片数量从p缩减为p-1，从而消除矩阵的奇异性。
第3步，建立增量型电场积分方程系统矩阵的叠层矩阵表达式。进一步地，该过程分三步进行：
第(1)步，构造电流簇树和电荷簇树。簇树TI一般通过递归地细分一组有限指标簇I＝{1,2,…,N}来构造。对于电流簇树，I表示定义在棱边自由度上的RWG基函数簇Ij；对于电荷簇树，I表示定义在面片自由度上的冲击基函数簇Iρ。采用基于定界盒子的递归二分方法来构造具有二叉树结构的电流簇树
Figure BDA0001255066610000056
和电荷簇树
Figure BDA0001255066610000057
以构造电流簇树
Figure BDA0001255066610000058
为例，如图2所示，首先以定界盒子包围整个目标表面，然后不断沿坐标轴三个方向对定界盒子进行递归二分，这样包含在定界盒子中的RWG基函数簇Ij也被递归二分，这一过程直至达到终止条件时结束。终止条件定义为，定界盒子达到预设尺寸或者定界盒子中的指标个数小于预设值。假定目标在三角形离散后生成25条棱边，生成的电流簇树
Figure BDA0001255066610000059
如图3所示。类似地，电荷簇树
Figure BDA00012550666100000510
也通过这种方式来构造，不同的是递归二分的不是棱边簇，而是面片簇。
第(2)步，构造电流-电流块树、电流-电荷块树、电荷-电流块树和电荷-电荷块树。块树TI×J由两个簇树TI和TJ相互作用生成。在伽辽金方法中，TI×J由两个相同的簇树TI和TJ生成，其中，TI可视作原始基函数簇树，TJ可视作测试基函数簇树。TI×J通过递归细分块I×J来构造，细分过程直至块t×s∈TI×J满足如下容许条件：
min{diam(Bt),diam(Bs)}≤ηdist(Bt,Bs) (12)
其中，t表示簇树TI中的任意簇，s表示簇树TJ中的任意簇，Bt表示包围t的最小定界盒子，Bs表示包围s的最小定界盒子，diam和dist分别表示定界盒子的欧几里得直径和距离；常数η＞0控制容许条件的严格程度，η越小，容许条件越严格。
以构造电流-电流块树
Figure BDA00012550666100000511
为例，原始基函数簇的电流簇树
Figure BDA00012550666100000512
和表示测试基函数簇的电流簇树
Figure BDA00012550666100000513
层层相互作用，如图4所示，从而生成电流-电流块树
Figure BDA00012550666100000514
如图5所示。
容许条件将块树中的所有块划分成两类：容许块和非容许块。满足容许条件的块称为容许块，对于任一容许块
Figure BDA00012550666100000515
具有如下低秩分解的矩阵形式：
Figure BDA00012550666100000516
其中，
Figure BDA0001255066610000061
和
Figure BDA0001255066610000062
均为满阵，上标p,q和k表示矩阵维度，符号
Figure BDA0001255066610000063
表示复数域。容许块之外的其余块称为非容许块，以满阵形式表示。如图4所示，在
Figure BDA0001255066610000064
中，白色块为容许块，黑色块为非容许块。
类似地，电流-电荷块树
Figure BDA0001255066610000065
电荷-电流块树
Figure BDA0001255066610000066
和电荷-电荷块树
Figure BDA0001255066610000067
均可通过相应的簇树之间的相互作用并选用适当的容许条件来构造。需要指出的是，
Figure BDA0001255066610000068
和
Figure BDA0001255066610000069
是方树，即行簇树和列簇树相同；但
Figure BDA00012550666100000610
和
Figure BDA00012550666100000611
不是方树，原因在于行簇树和列簇树不同。
第(3)步，将系统矩阵中各子矩阵的非零元素分别填入对应的块树中，生成系统矩阵的叠层矩阵表达式。根据(9)式，基于矩阵收缩方法的增量型电场积分方程的系统矩阵
Figure BDA00012550666100000612
包含五个子矩阵：
Figure BDA00012550666100000613
和
Figure BDA00012550666100000614
因此，需要分别构造它们的叠层矩阵表达式
Figure BDA00012550666100000615
Figure BDA00012550666100000616
和
Figure BDA00012550666100000617
来构造出
Figure BDA00012550666100000618
的叠层矩阵表达式
Figure BDA00012550666100000619
对于
Figure BDA00012550666100000620
和
Figure BDA00012550666100000621
的构造，需要将
Figure BDA00012550666100000622
中的元素填入电流-电流块树
Figure BDA00012550666100000623
中，
Figure BDA00012550666100000624
中的元素填入电荷-电荷块树
Figure BDA00012550666100000625
中。其中，以满阵形式存储的非容许块采用矩量法直接构造，以低秩分解矩阵形式存储的容许块则采用低频多层快速多极子算法(Low-frequency MultilevelFast Multipole Algorithm,LF-MLFMA)来构造。在簇树层面，传统LF-MLFMA基于八叉树结构，而叠层矩阵算法通常基于二叉树结构，因此需要将八叉树结构转化为二叉树结构。由于八叉树的一层相当于二叉树的三层，因此，当一个n层的八叉树建立时，就相当于建立了一个3n层的二叉树。在块树层面，传统LF-MLFMA和叠层矩阵中，容许块和非容许块的划分方式必须一致，这要求容许条件的选取必须一致。根据传统LF-MLFMA的近远场划分方式可知，通过选取在(12)式所示的容许条件中选取η＝1能够实现LF-MLFMA和叠层矩阵块树的统一。这样，便可采用LF-MLFMA来构造
Figure BDA00012550666100000626
和
Figure BDA00012550666100000627
了。
三维LF-MLFMA的核心转移因子具有如下形式：
Figure BDA00012550666100000628
其中，
Figure BDA00012550666100000629
父层转移因子，
Figure BDA00012550666100000630
表示子层配置因子，
Figure BDA00012550666100000631
表示子层转移因子，
Figure BDA00012550666100000632
表示子层聚合因子。Li为多极子模式数。这里，
Figure BDA00012550666100000633
其中P表示多极子的数目，用来控制展开形式和转移因子的精度。(14)式可写成矩阵形式如下：
Figure BDA00012550666100000634
其中，
Figure BDA00012550666100000635
和
Figure BDA00012550666100000636
分别表示父层转移矩阵、子层配置矩阵、子层转移矩阵和子层聚合矩阵。考虑多层的情况，如图6a、图6b所示，每一层的矩阵形式如下：
第l层：
Figure BDA00012550666100000637
第l-1层：
Figure BDA00012550666100000638
对于满足容许条件(12)的两个簇t和s，根据(15)式集合其中所有指标j∈t和i∈s,即可构造出具有如下低秩分解形式的容许块
Figure BDA0001255066610000071
Figure BDA0001255066610000072
其中，
Figure BDA0001255066610000073
为配置矩阵，
Figure BDA0001255066610000074
为聚合矩阵，
Figure BDA0001255066610000075
为转移矩阵，下标#t、#s和r表示矩阵的维度。这里，符号“#”用来表示簇里的指标个数。对于
Figure BDA0001255066610000076
中的容许块，#t表示观察组(簇)t中RWG基函数的个数，#s则表示源组(簇)s中RWG基函数的个数；r＝3×(P+1)2，其中(P+1)2表示Li的总数，常数项“3”描述的是电流矢量的三个方向。对于
Figure BDA0001255066610000077
中的容许块，#t和#s分别表示簇t和s中冲击基函数的个数，由于电荷是标量，因此有r＝(P+1)2。通过调节多极子的数目P，可以控制低秩分解的精度。
尽管LF-MLFMA构造的容许块已经具有(16)式所示的低秩分解矩阵形式，该低秩分解矩阵仍然存在冗余信息，可被进一步压缩，生成更为紧凑的叠层矩阵表达式。下面采用一种基于QR分解和奇异值分解的重压缩方法对容许块进一步压缩，具体过程如下：
1.对
Figure BDA0001255066610000078
执行QR分解，生成矩阵Qt和Rt：
Figure BDA0001255066610000079
(
Figure BDA00012550666100000710
这里符号
Figure BDA00012550666100000711
表示复数域，上标#t和r表示矩阵维度)；
2.对
Figure BDA00012550666100000712
执行QR分解，生成矩阵Qs和Rs：
Figure BDA00012550666100000713
(
Figure BDA00012550666100000714
这里符号
Figure BDA00012550666100000715
表示复数域，上标#s和r表示矩阵维度)；
3.执行
Figure BDA00012550666100000716
的操作，并对
Figure BDA00012550666100000717
进行奇异值分解，生成矩阵U，∑和V：
Figure BDA00012550666100000718
其中∑为对角阵；
4.从矩阵∑的对角元素中提取部分元素构造出矩阵
Figure BDA00012550666100000719
这里
Figure BDA00012550666100000720
∑11,∑22…,∑kk分别表示矩阵∑的第1,2…,k个对角元素，且∑(k+1)(k+1)≤εrec∑11＜∑kk，其中，参数εrec表示相对截断误差，用来控制重压缩的精度；
5.分别从U和V中提取前k列元素获得矩阵
Figure BDA00012550666100000721
和
Figure BDA00012550666100000722
即：
Figure BDA00012550666100000723
其中，Uk和Vk分别表示U和V的第k列元素；
6.获得
Figure BDA00012550666100000724
其中，
Figure BDA00012550666100000725
和
Figure BDA00012550666100000726
对于
Figure BDA00012550666100000727
的构造，由于
Figure BDA00012550666100000728
是稀疏矩阵，因此只需将
Figure BDA00012550666100000729
中的非零元素填入电荷-电流块树
Figure BDA00012550666100000730
中。
Figure BDA00012550666100000731
中的非零元素均可被填入非容许块中，而不会被填入容许块中，原因在于，满足容许条件(12)的两个簇t和s必然离开一定的距离，而根据(7)式，离开一定距离的i∈t,j∈s生成的矩阵元素
Figure BDA00012550666100000732
因此，
Figure BDA00012550666100000733
可由
Figure BDA00012550666100000734
无损地表达。
Figure BDA00012550666100000735
只是
Figure BDA00012550666100000736
的转置。此外，
Figure BDA00012550666100000737
通过将电荷-电荷块树
Figure BDA00012550666100000738
的主对角块的主对角线元素赋值为1即可构造出来。
对于
Figure BDA00012550666100000739
的构造，只需将
Figure BDA00012550666100000740
中的非零元素填入电荷-电荷块树
Figure BDA00012550666100000741
中，这里，
Figure BDA00012550666100000742
表示提取向量
Figure BDA00012550666100000750
中与电荷相关的部分向量，即：
Figure BDA00012550666100000743
与
Figure BDA00012550666100000744
的构造过程不同,
Figure BDA00012550666100000745
中的容许块极易生成，这是因为
Figure BDA00012550666100000746
天然地具有低秩分解矩阵的形式，换句话说，
Figure BDA00012550666100000747
中的任一容许块
Figure BDA00012550666100000748
均可直接从
Figure BDA00012550666100000749
中提取，并具有秩k＝1的低秩分解矩阵形式：
Figure BDA0001255066610000081
Figure BDA0001255066610000082
表示向量
Figure BDA0001255066610000083
中与簇t相关的部分向量，其维度为#t×1；
Figure BDA0001255066610000084
表示向量
Figure BDA0001255066610000085
中与簇s相关的部分向量，其维度为#s×1；下标#t和#s分别表示簇t和s中基函数的个数。需要指出的是，当分析目标包含多个独立物体时，k的取值等于独立物体的数目。
Figure BDA0001255066610000086
中的非容许块通过直接计算
Figure BDA0001255066610000087
可方便地获得。
至此，系统矩阵
Figure BDA0001255066610000088
的五个子矩阵的叠层矩阵表达式
Figure BDA0001255066610000089
和
Figure BDA00012550666100000810
均构造完成，从而获得了系统矩阵
Figure BDA00012550666100000811
的叠层矩阵表达式
Figure BDA00012550666100000812
需要说明的是，上述第(1)步到第(3)步构造过程针对的是矩阵收缩方法。如果采用缩减电荷自由度的方法，只需在第(1)步中针对目标的每个独立物体去掉一个面片自由度即可，而无需加入修正矩阵
Figure BDA00012550666100000813
第4步，采用叠层矩阵格式的运算法则对该叠层矩阵进行上下三角分解。
根据(9)式，
Figure BDA00012550666100000814
可表示为一个2×2的分块矩阵形式：
Figure BDA00012550666100000815
其中，
Figure BDA00012550666100000816
及
Figure BDA00012550666100000817
这里，
Figure BDA00012550666100000818
表示
Figure BDA00012550666100000819
中与电荷相关的子矩阵。
首先，通过计算
Figure BDA00012550666100000820
来获得
Figure BDA00012550666100000821
这里，采用叠层矩阵格式的乘法
Figure BDA00012550666100000822
取代传统矩阵乘法来完成
Figure BDA00012550666100000823
的操作。此外，采用叠层矩阵格式的加法
Figure BDA00012550666100000824
取代传统矩阵加法来完成
Figure BDA00012550666100000825
的操作。然后，从(18)式所示分块矩阵的上下三角分解出发，
Figure BDA00012550666100000826
这里，
Figure BDA00012550666100000827
和
Figure BDA00012550666100000828
分别表示对
Figure BDA00012550666100000829
进行上下三角分解后获得的矩阵块，通过层层递归执行如下步骤来完成叠层矩阵的上下三角分解：
1.对
Figure BDA00012550666100000830
执行上下三角分解：
Figure BDA00012550666100000831
获得
Figure BDA00012550666100000832
和
Figure BDA00012550666100000833
2.求解下三角方程：
Figure BDA00012550666100000834
获得
Figure BDA00012550666100000835
3.求解上三角方程：
Figure BDA00012550666100000836
获得
Figure BDA00012550666100000837
4.计算
Figure BDA00012550666100000838
然后对其执行上下三角分解：
Figure BDA00012550666100000839
获得
Figure BDA00012550666100000840
和
Figure BDA00012550666100000841
上述步骤中涉及到的矩阵加法和乘法均采用叠层矩阵格式的加法和乘法
Figure BDA00012550666100000842
来完成。叠层矩阵格式的加法和乘法中涉及低秩分解矩阵的运算，采用基于QR分解和奇异值分解的截断算子
Figure BDA00012550666100000843
来定义，即：
Figure BDA00012550666100000844
和
Figure BDA00012550666100000845
引入自适应截断误差εt来控制截断算子的精度，εt的含义是用低秩分解矩阵与满阵之间的误差。叠层矩阵上下三角分解的精度也由εt来调控。叠层矩阵上下三角分解在低频电磁仿真时具有O(Nlog2N)计算复杂度和O(NlogN)内存消耗，其中，N为未知量的数目。完成叠层矩阵的上下三角分解之后，获得叠层矩阵格式的上下三角因子，结合叠层矩阵格式的前后向回代，即可用于完成对增量型电场积分方程的系统矩阵方程的求解。叠层矩阵格式的前后向回代在低频电磁仿真时的计算复杂度为O(NlogN)。
第5步，将叠层矩阵的上下三角分解与前后向回代结合，采用预条件的迭代解法或者直接解法进行求解。将系统矩阵方程简写作：
Figure BDA0001255066610000091
其中，
Figure BDA0001255066610000092
表示待求的未知系数向量，
Figure BDA0001255066610000093
表示右边向量。
采用预条件的迭代解法时，将方程(20)左右两端同时乘以一个预条件矩阵
Figure BDA0001255066610000094
让
Figure BDA0001255066610000095
具有比
Figure BDA0001255066610000096
更小的条件数。叠层矩阵的上下三角分解用作预条件时，
Figure BDA0001255066610000097
由于叠层矩阵分解之后，
Figure BDA0001255066610000098
和
Figure BDA0001255066610000099
均已获得，只需执行叠层矩阵格式的前后向回代即可完成一次
Figure BDA00012550666100000910
的矩阵-向量乘操作，参与到迭代求解方程(20)的每一次矩阵-向量乘操作中。
采用直接解法时，由于
Figure BDA00012550666100000911
方程(20)可写为：
Figure BDA00012550666100000912
通过执行一次前后向回代，即可获得方程的解。
实际中，通过调节截断误差εt可以调控叠层矩阵分解的精度。高精度(达到工程需求的精度)的叠层矩阵分解可用作直接解法，低精度(未达到工程需求的精度)的叠层矩阵分解可用作迭代解法的预条件。叠层矩阵分解精度越高，预条件的效果越好，但同时构造时间和内存消耗也越大。因此，选用何种精度的叠层矩阵分解作为预条件，以及选用直接解法还是预条件的迭代解法，都需要考虑多种因素进行权衡。一般地，直接解法更适合于多右边向量问题的求解，迭代解法更适合于单右边向量问题的求解。
求解出电流和电荷系数后，即可获得目标表面的电流及电荷分布，通过进一步计算即可提取出所需的电磁特性参数，可广泛用于芯片封装、集成电路、电磁散射、天线设计及电磁兼容等领域中的低频电磁特性仿真。
实施例1
采用本发明所述方法对一个具有复杂外形的飞机目标的低频电磁散射特性进行分析。该飞机目标的电尺寸为0.16λ×0.11λ×0.04λ，其中，λ表示入射电磁波的波长，如图7所示。目标受10MHz的入射平面波照射。采用三角形单元离散，生成40704条棱边和27136个面片，最小棱边长度仅为0.910×10-4λ。按照本发明所述步骤，构建基于叠层矩阵分解的预条件。采用低频多层快速多极子算法的层数取为5层，构造的叠层矩阵层数为18层，多极子个数P取为5，容许块重压缩精度εrec取为10-3。采用叠层矩阵分解预条件，来加速基于广义最小余量法(GMRES)的迭代解法的收敛，迭代收敛误差设置为10-6。叠层矩阵分解预条件中，分别测试截断精度εt为10-2和10-3两种情况。图8给出了迭代收敛曲线，并与常用的约束预条件进行了对比，可以看出，本发明所述方法能够显著提高迭代收敛速度。表1给出求解效率的对比，可以看出，约束预条件尽管构造时间短，但预条件效果不理想，收敛速度慢，因此求解时间长；而本发明所述方法则显著提高了求解效率。此外，从表1不难看出，对于单激励问题，低截断精度的预条件求解效率高；但是对于多激励问题(如：单站雷达散射截面分析等)，高截断精度的预条件甚至直接解法则更有效，因为叠层矩阵分解只需构造一次，而迭代过程则需进行多次。图9和图10表明，本发明所述叠层矩阵分解方法能够将计算复杂度和内存消耗分别降低到O(Nlog2N)和O(NlogN)。
表1几种预条件技术之间的性能比较
Figure BDA0001255066610000101
以上内容是结合具体的优选实施方式对本发明所作的进一步详细说明，不能认定本发明的具体实施只局限于这些说明。对于本发明所属技术领域的普通技术人员来说，在不脱离本发明构思的前提下，还可以做出若干简单推演或替代，都应当视为属于本发明的保护范围。

Claims (4)

1.一种基于叠层矩阵分解的低频电磁特性仿真方法，其特征在于，包括以下步骤：
第1步，针对目标表面建立电场积分方程，并引入电流连续性方程，采用三角形单元离散目标表面，用RWG基函数展开表面电流，用冲击基函数展开表面电荷，生成增量型电场积分方程的系统矩阵方程；
第2步，利用矩阵收缩方法或者缩减电荷自由度方法对该系统矩阵方程进行修正，消除电荷中性原理在极低频率时造成的矩阵奇异性；
第3步，建立增量型电场积分方程系统矩阵的叠层矩阵；
第4步，采用叠层矩阵格式的运算法则对该叠层矩阵进行上下三角分解；
第5步，将叠层矩阵的上下三角分解与前后向回代结合，采用预条件的迭代解法或者直接解法进行求解，获得目标表面的电流及电荷分布，通过计算提取所需的电磁特性参数；
第3步中，建立增量型电场积分方程系统矩阵的叠层矩阵的步骤如下：
第(1)步，分别针对定义在棱边上的RWG基函数簇和定义在面片上的冲击基函数簇，采用几何递归二分的方法，构造电流簇树和电荷簇树；
第(2)步，电流簇树和电荷簇树相互作用，构造电流-电流块树、电流-电荷块树、电荷-电流块树和电荷-电荷块树，并引入容许条件来划分容许块和非容许块，其中，非容许块以满阵形式表示，容许块以低秩分解矩阵形式表示；
第(3)步，将系统矩阵中各子矩阵的非零元素分别填入对应的块树中，其中，对于电流子矩阵和电荷子矩阵，采用矩量法构造非容许块，采用低频多层快速多极子技术构造容许块，并利用基于QR分解和奇异值分解的重压缩方法进一步对容许块进行压缩；对于其余子矩阵，将其非零元素对应填入非容许块或容许块，生成系统矩阵的叠层矩阵。
2.根据权利要求1所述的一种基于叠层矩阵分解的低频电磁特性仿真方法，其特征在于，第(3)步中，采用低频多层快速多极子技术构造容许块具体包括，提取低频多层快速多极子的聚合、转移和配置矩阵，将容许块表示成低秩分解矩阵的形式，并通过调节多极子数目来控制低秩分解的精度。
3.根据权利要求1所述的一种基于叠层矩阵分解的低频电磁特性仿真方法，其特征在于，第4步中，叠层矩阵格式的运算法则中涉及低秩分解矩阵运算，通过调节低秩分解矩阵运算中的截断精度来控制叠层矩阵上下三角分解的精度，高精度的分解用于构造直接解法，低精度的分解用于构造预条件来加速迭代解法的收敛，且不同精度的分解具有不同的加速效果。
4.根据权利要求3所述的一种基于叠层矩阵分解的低频电磁特性仿真方法，其特征在于，叠层矩阵上下三角分解的精度通过分解前后矩阵之差的二范数和分解前矩阵的二范数之比来定义，高精度指达到实际仿真精度要求的精度，低精度指低于实际仿真精度要求的精度。
CN201710186894.1A 2017-03-27 2017-03-27 一种基于叠层矩阵分解的低频电磁特性仿真方法 Active CN106991222B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201710186894.1A CN106991222B (zh) 2017-03-27 2017-03-27 一种基于叠层矩阵分解的低频电磁特性仿真方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201710186894.1A CN106991222B (zh) 2017-03-27 2017-03-27 一种基于叠层矩阵分解的低频电磁特性仿真方法

Publications (2)

Publication Number Publication Date
CN106991222A CN106991222A (zh) 2017-07-28
CN106991222B true CN106991222B (zh) 2020-04-21

Family

ID=59413377

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201710186894.1A Active CN106991222B (zh) 2017-03-27 2017-03-27 一种基于叠层矩阵分解的低频电磁特性仿真方法

Country Status (1)

Country Link
CN (1) CN106991222B (zh)

Families Citing this family (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN108170647B (zh) * 2017-12-13 2019-12-27 南京理工大学 一种低频电磁特性分析的自适应嵌套交叉近似方法
CN108629143A (zh) * 2018-05-16 2018-10-09 南京邮电大学 基于低秩分解的电磁有限元-边界元的直接求解方法
CN110826010B (zh) * 2019-10-10 2022-11-15 东南大学 介质表面积分方程多层快速多极子的效率倍增方法
CN113609646B (zh) * 2021-07-08 2022-04-12 中国人民解放军32215部队 一种复杂陆地环境与装备的耦合电磁散射特性建模仿真方法
CN113781332B (zh) * 2021-08-24 2023-09-29 中国电子产品可靠性与环境试验研究所（（工业和信息化部电子第五研究所）（中国赛宝实验室）） 电磁数据的降噪方法、装置、计算机设备和存储介质

Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102033985A (zh) * 2010-11-24 2011-04-27 南京理工大学 基于*-矩阵算法的高效时域电磁仿真方法
CN104731996A (zh) * 2013-12-24 2015-06-24 南京理工大学 一种快速提取电大尺寸金属腔体目标瞬态散射信号的仿真方法

Family Cites Families (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
FR2810136B1 (fr) * 2000-06-09 2002-10-18 Thomson Csf Algorithme de simulation electromagnetique, notamment des performances d'une antenne

Patent Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102033985A (zh) * 2010-11-24 2011-04-27 南京理工大学 基于*-矩阵算法的高效时域电磁仿真方法
CN104731996A (zh) * 2013-12-24 2015-06-24 南京理工大学 一种快速提取电大尺寸金属腔体目标瞬态散射信号的仿真方法

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
A matrix-free spectral rotation approach to the computation of electromagnetic fields generated by a surface current distribution;Tiberi G 等;《IEEE Antennas & Wireless Propagation Letters》;20050620;第4卷(第1期);第121-124页 *
复杂目标电磁散射算法与特性仿真研究;余定峰;《中国博士学位论文全文数据库 基础科学辑》;20140731;A005-22 *

Also Published As

Publication number Publication date
CN106991222A (zh) 2017-07-28

Similar Documents

Publication Publication Date Title
CN106991222B (zh) 一种基于叠层矩阵分解的低频电磁特性仿真方法
Chen et al. Discontinuous Galerkin time-domain methods for multiscale electromagnetic simulations: A review
CN111639447B (zh) 多级局部时间步进技术的任意高阶混合网格时域不连续伽辽金方法
US6185517B1 (en) Electromagnetic field intensity computing device
Chai et al. Direct matrix solution of linear complexity for surface integral-equation-based impedance extraction of complicated 3-D structures
Zemanian et al. Three-dimensional capacitance computations for VLSI/ULSI interconnections
CN111046603A (zh) 基于gpu并行加速特征基函数算法的电磁散射特性分析方法
Omar et al. A linear complexity direct volume integral equation solver for full-wave 3-D circuit extraction in inhomogeneous materials
Wu et al. A multiple-grid adaptive integral method for multi-region problems
Wan et al. Hierarchical matrix techniques based on matrix decomposition algorithm for the fast analysis of planar layered structures
CN102708229A (zh) 复杂分层媒质结构的矩阵分解结合新奇异值分解方法
CN102054094A (zh) 一种平面微带电路的快速多层方向性仿真方法
CN115018079A (zh) 量子电路、仿真方法、装置、设备及存储介质
CN114417769B (zh) 基于贝塞尔函数分段积分的集成电路电磁仿真方法及系统
CN116187258A (zh) 量子芯片版图的仿真方法、装置、计算设备及存储介质
US20090254873A1 (en) Circuit board analyzer and analysis method
EP3485404A1 (en) Eigen augmentation methods for electromagnetic modelling and simulation
Otto et al. Using the FEniCS package for FEM solutions in electromagnetics
Avula et al. Pole residue equivalent system solver (PRESS)
Sharma et al. A complete surface integral method for broadband modeling of 3D interconnects in stratified media
Liu et al. Layered ${\cal H} $-Matrix Based Inverse and LU Algorithms for Fast Direct Finite-Element-Based Computation of Electromagnetic Problems
Seo et al. Analyzing PEC Scattering Structure Using an IE-FFT Algorithm”
CN114239239A (zh) 用于弹目交会目标电磁特性快速模拟的直接稀疏求解方法
Li et al. MultiAIM: Fast electromagnetic analysis of multiscale structures using boundary element methods
CN115329973A (zh) 仿真方法、装置、设备及存储介质

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant