CN106777472B - 基于拉盖尔多项式的减少分裂误差的完全匹配层实现方法 - Google Patents

基于拉盖尔多项式的减少分裂误差的完全匹配层实现方法 Download PDF

Info

Publication number
CN106777472B
CN106777472B CN201611009431.XA CN201611009431A CN106777472B CN 106777472 B CN106777472 B CN 106777472B CN 201611009431 A CN201611009431 A CN 201611009431A CN 106777472 B CN106777472 B CN 106777472B
Authority
CN
China
Prior art keywords
field component
ζmax
calculating
updating
coefficient
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
CN201611009431.XA
Other languages
English (en)
Other versions
CN106777472A (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.)
Xian University of Technology
Original Assignee
Xian University of Technology
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 Xian University of Technology filed Critical Xian University of Technology
Priority to CN201611009431.XA priority Critical patent/CN106777472B/zh
Publication of CN106777472A publication Critical patent/CN106777472A/zh
Application granted granted Critical
Publication of CN106777472B publication Critical patent/CN106777472B/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
    • G06F30/23—Design optimisation, verification or simulation using finite element methods [FEM] or finite difference methods [FDM]
    • Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02E—REDUCTION OF GREENHOUSE GAS [GHG] EMISSIONS, RELATED TO ENERGY GENERATION, TRANSMISSION OR DISTRIBUTION
    • Y02E60/00—Enabling technologies; Technologies with a potential or indirect contribution to GHG emissions mitigation

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Computer Hardware Design (AREA)
  • Evolutionary Computation (AREA)
  • Geometry (AREA)
  • General Engineering & Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Complex Calculations (AREA)

Abstract

本发明公开了一种基于拉盖尔多项式的减少分裂误差的完全匹配层实现方法，具体按照以下步骤实施：输入模型文件；初始化参数和设置参数；添加场源到y方向上的电场分量系数中，使用因式分裂的WLP‑FDTD方法计算电场分量系数
Figure DDA0001153990550000011
记为初始场值
Figure DDA0001153990550000012
更新计算整个计算区域的y方向上电场分量系数
Figure DDA0001153990550000013
更新计算整个计算区域的x方向上电场分量系数
Figure DDA0001153990550000014
判断迭代次数k是否达到预设值；更新计算整个计算区域的磁场分量系数
Figure DDA0001153990550000015
更新计算整个计算区域的电磁场分量系数的辅助变量；更新计算观测点处的电磁场分量；判断拉盖尔多项式的阶数q是否达到预设值。本发明的基于拉盖尔多项式的减少分裂误差的完全匹配层实现方法，计算速度快、精度高，且对低频与凋落波的吸收更加有效。

Description

基于拉盖尔多项式的减少分裂误差的完全匹配层实现方法
技术领域
本发明属于计算电磁学技术领域，具体涉及一种基于拉盖尔多项式的减少分裂误差的完全匹配层实现方法。
背景技术
众所周知，时域有限差分(Finite-difference time-domain，FDTD)方法的时间步长受柯西稳定性条件的限制，这限制了FDTD方法在精细结构模型中的应用。为了消除柯西稳定性条件的限制，人们提出了无条件稳定时域有限差分方法，比如：交替方向隐式(Alternating-Direction-Implicit，ADI)的时域有限差分(ADI-FDTD)方法、局部一维(locally one dimensional，LOD)时域有限差分(LOD-FDTD)方法和基于拉盖尔多项式的时域有限差分(WLP-FDTD)方法。在这些方法中，WLP-FDTD方法既能消除柯西稳定性条件的限制，而且又能解决ADI-FDTD方法在使用较大的时间步长时产生很大的色散误差这个难题，因此WLP-FDTD方法在求解精细结构模型下的电磁场问题时，具有一定的优越性。然而，这种传统的WLP-FDTD方法在求解精细结构的电磁场问题时，会产生一个大型的稀疏矩阵方程，直接求解此方程会使得计算较复杂，内存消耗较大。之前提出了一种因式分裂的WLP-FDTD方法，这种方法计算速度得到了一定的提升，但是存在分裂误差。
另外，由于计算机容量的限制，电磁场的计算只能在有限区域进行。为了能模拟开域电磁波传播过程，在计算区域的截断边界处必须给出吸收边界条件。现有的直角坐标系下吸收边界主要有：Mur吸收边界，分裂场的PML(Split-field PML)，单轴各向异性PML(Uniaxial Anisotropic PML，UPML)。以上三种吸收边界均可以应用于WLP-FDTD方法的电磁场计算中，但上述吸收边界，对低频以及凋落波的吸收效果都不理想。
发明内容
本发明的目的是提供一种基于拉盖尔多项式的减少分裂误差的完全匹配层实现方法，计算速度较快、精度高、而且对于低频和凋落波具有很好的吸收效果。
本发明所采用的技术方案是，基于拉盖尔多项式的减少分裂误差的完全匹配层实现方法，具体按照以下步骤实施：
步骤1：输入模型文件；
步骤2：初始化参数和设置参数；
步骤3：添加场源到y方向上的电场分量系数中，使用因式分裂的WLP-FDTD方法计算电场分量系数
Figure BDA0001153990530000021
记为初始场值
Figure BDA0001153990530000022
步骤4：更新计算整个计算区域的y方向上电场分量系数
Figure BDA0001153990530000023
步骤5：更新计算整个计算区域的x方向上电场分量系数
Figure BDA0001153990530000024
步骤6：将k+1赋值给k，并判断迭代次数k是否达到预设值，若未达到预设值，则返回步骤4，若达到预设值，则执行步骤7；
步骤7：更新计算整个计算区域的磁场分量系数
Figure BDA0001153990530000025
步骤8：更新计算整个计算区域的电磁场分量系数的辅助变量；
步骤9：更新计算观测点处的电磁场分量；
步骤10：将q+1赋值给q，并判断拉盖尔多项式的阶数q是否达到预设值，若未达到预设值，则返回步骤3，若达到预设值，则结束。
本发明的特点还在于：
步骤1输入模型文件，具体为：
计算区域大小Nx×Ny，其中Nx为x方向的网格数，Ny为y方向的网格数；空间步长Δζ，ζ＝x，y，x为横坐标，y为纵坐标；时间步长Δt；真空中的电导率σ，磁导率μ0，介电常数ε0；吸收边界层数NPML与相关参数κζmax，σζmax，αζmax；其中，κζmax取整数，κζmax取值范围为[1，60]；αζmax取值范围为[0，1)；σζmax/σopt取值范围为(0，12]；仿真计算时长Tf；迭代次数k，k≥0且为整数；加权拉盖尔多项式的阶数q，q≥0且为整数；时间尺度因子s，s取值范围为[109，1013]；观测点；场源参数。
步骤2初始化参数和设置参数具体为：
初始化的参数包括：
将整个计算区域的电磁场分量系数
Figure BDA0001153990530000031
整个计算区域的电磁场分量系数的和
Figure BDA0001153990530000032
整个计算区域的辅助变量
Figure BDA0001153990530000033
拉盖尔多项式
Figure BDA0001153990530000034
均初始化为零，其中Fη＝Ex,Ey,Hz，ζ＝x，y，
Figure BDA0001153990530000035
初始化PML系数(C1ζ,C2ζ,C3ζ)，具体为：
C1ζ＝1/(1+0.5ε0s)
C2ζ＝1
C3ζ＝ε0/μ0
其中，ζ＝x，y，ε0是空气中的介电常数，s为时间尺度因子，s取值范围为[109，1013]；
设置的参数具体为：
设置CFS-PML吸收边界的参数，具体为：
σζ＝σζmax|ζ-ζ0|m/dm
κζ＝1+(κζmax-1)|ζ-ζ0|m/dm
αζ＝αζmaxζ0/d
式中ζ＝x，y，ζ0为PML层与非PML截面位置，d是PML吸收边界的厚度，κζmax取整数，κζmax取值范围为[1，60]；αζmax取值范围为[0，1)；σζmax根据σopt来设置，σζmax/σopt取值范围为(0，12]；σopt＝(m+1)/150πΔζ，m取值范围为[1，20]，其中m取值为4时边界的吸收效果最好，Δζ取值范围
Figure BDA0001153990530000041
λ为源的波长；
设置PML系数，具体为：
Figure BDA0001153990530000042
步骤3中所添加的场源的表达式为：
Figure BDA0001153990530000043
其中，Tc，Td为场源参数。
步骤4具体为：
步骤4.1：电场分量系数
Figure BDA0001153990530000044
在计算区域的方程为：
Figure BDA0001153990530000051
其中，i表示横坐标上的第i个计算网格，j表示纵坐标上的第j个计算网格；
步骤4.2：使用追赶法求解步骤4.1的方程，得到整个计算区域的电场分量系数
Figure BDA0001153990530000052
步骤5具体为：
步骤5.1：电场分量系数
Figure BDA0001153990530000053
在计算区域的方程为：
Figure BDA0001153990530000061
步骤5.2：使用追赶法求解步骤5.1的方程，得到整个计算区域的x方向上电场分量系数
Figure BDA0001153990530000062
步骤7更新计算整个计算区域的磁场分量系数
Figure BDA0001153990530000063
具体更新公式为：
Figure BDA0001153990530000064
步骤8更新计算整个计算区域的电磁场分量系数的辅助变量，具体更新公式为：
Figure BDA0001153990530000065
Figure BDA0001153990530000066
其中，Fη＝Ex,Ey,Hz；ζ＝x，y。
步骤9更新计算观测点处的电磁场分量，具体更新公式为：
Figure BDA0001153990530000071
其中，U表示电磁场分量Ex,Ey,Hz，Uq表示q阶电磁场分量系数，
Figure BDA0001153990530000072
是q阶加权拉盖尔多项式，
Figure BDA0001153990530000073
是带有时间尺度因子s＞0的扩展时间，
Figure BDA0001153990530000074
是q阶拉盖尔多项式。
本发明的有益效果是：
①本发明基于拉盖尔多项式的减少分裂误差的完全匹配层实现方法，在直角坐标系下，通过用加权拉盖尔多项式表示电磁场分量，来解时域麦克斯韦方程，使得在更新计算整个计算区域的电磁场分量系数时不涉及到时间步长，只是在最后计算观测点处的电磁场分量时用到时间步长，因此计算过程中时间步长可以取得比柯西稳定性条件限制的时间步长更大；
②本发明基于拉盖尔多项式的减少分裂误差的完全匹配层实现方法，在求解电磁场分量系数时，将大型稀疏矩阵方程分裂成两个三对角矩阵方程，同时使用迭代的方案，使得它在计算时比传统的WLP-FDTD方法更简单、计算速度更快、内存消耗更少、计算精度高，而且可以对大区域的电磁场问题进行求解；
③本发明基于拉盖尔多项式的减少分裂误差的完全匹配层实现方法，在设置PML系数时，由于采用了CFS因子，并且通过调整CFS因子中的参数，可以使得该吸收边界对低频与凋落波的吸收更加有效；
④本发明基于拉盖尔多项式的减少分裂误差的完全匹配层实现方法，由于采用了复扩展坐标系，使得PML在实现时避免了场的分裂且与媒质无关。
附图说明
图1是本发明完全匹配层实现方法的流程图；
图2是本发明实施例中的计算模型的结构示意图；
图3是本发明的方法与传统FDTD方法和因式分裂的WLP-FDTD方法在观测点处时域波形对比图；
图4是本发明实施例中观测点的不同吸收边界相对反射误差。
具体实施方式
下面结合附图和具体实施方式对本发明进行详细说明。
本发明基于拉盖尔多项式的减少分裂误差的完全匹配层实现方法，所依据的原理为：首先导出复扩展直角坐标系下，PML中电磁场所满足的麦克斯韦方程；然后利用一种迭代的因式分裂的WLP-FDTD方法推导出电磁场分量系数的更新方程；最后求解观测点处的电磁场分量。
在求电磁场分量系数更新方程时，首先需要推导出复扩展直角坐标系下，PML中电磁场满足的麦克斯韦方程；
复扩展直角坐标系下，PML中电磁场满足的麦克斯韦方程为：
Figure BDA0001153990530000081
其中，
Figure BDA0001153990530000082
表示电场矢量，
Figure BDA0001153990530000083
表示磁场矢量，j为虚数单位，ω为角频率，μ为介质的磁导率，ε为介质的介电常数，
Figure BDA0001153990530000084
为修正后的微分算子，可以写成：
Figure BDA0001153990530000085
sx、sy和sz是坐标扩展变量，可以表示成：
sζ＝kζ+σζ/jωε (3)
其中ζ表示x、y、z，kζ、σζ和αζ为PML的有关的参数。
本发明仅考虑简单无耗媒质中二维横电波情况，于是复扩展直角坐标下的麦克斯韦方程可以写成：
Figure BDA0001153990530000091
Figure BDA0001153990530000092
Figure BDA0001153990530000093
其中Ex,Ey分别表示x,y方向的电场，Hz分别表示z方向的磁场。
然后，求出电磁场分量系数的更新方程；
为了计算方便，引入下面几个辅助变量：
Figure BDA0001153990530000094
将(3)代入(7)，然后利用jω→t的变换，可以得到四组方程，这里仅给出第一个方程：
Figure BDA0001153990530000095
由于电磁场分量和其对时间的一阶偏导可以展开成一系列的电磁场分量系数与加权拉盖尔多项式的函数之和，公式如下：
Figure BDA0001153990530000096
上式中U表示电磁场分量Ex,Ey,Hz，Uq表示q阶电磁场分量系数，
Figure BDA0001153990530000097
是q阶加权拉盖尔多项式，
Figure BDA0001153990530000098
是带有时间尺度因子s＞0的扩展时间，
Figure BDA0001153990530000099
是q阶拉盖尔多项式。将(9)代入(4)、(5)和(6)中，然后让方程两边同乘以
Figure BDA0001153990530000101
可以得到：
Figure BDA0001153990530000102
Figure BDA0001153990530000103
Figure BDA0001153990530000104
上面三式中，b＝2ε0/μ0，q是加权拉盖尔多项式的阶数，Dx和Dy分别是沿x和y方向上的微分算子，
Figure BDA00011539905300001013
和
Figure BDA00011539905300001014
是q阶电磁场分量系数，Ciζ,i＝1,2,3；ζ＝x,y是与坐标网格有关的PML系数，计算式为：
Figure BDA0001153990530000105
Figure BDA0001153990530000106
和
Figure BDA0001153990530000107
是拉盖尔域中电磁场分量和辅助变量的低阶和，公式如下：
Figure BDA0001153990530000108
上式中的辅助变量
Figure BDA0001153990530000109
的计算式为：
Figure BDA00011539905300001010
将(10)、(11)、(12)写成一个矩阵形式的方程，如下：
Figure BDA00011539905300001011
式中
Figure BDA00011539905300001012
Figure BDA0001153990530000111
Figure BDA0001153990530000112
令
Figure BDA0001153990530000113
则(16)式变为：
Figure BDA0001153990530000114
添加一个微扰项
Figure BDA0001153990530000115
到上式，于是得到
Figure BDA0001153990530000116
上式可以分裂为下面两式
Figure BDA0001153990530000117
式中
Figure BDA0001153990530000118
是一个非物理中间量，为了解(22)式，使上式中I为单位对角矩阵，A和B的等式如下
Figure BDA0001153990530000119
式中
Figure BDA00011539905300001110
将(23)式代入(22)式化简后得到
Figure BDA00011539905300001111
将上式展开后得到
Figure BDA0001153990530000121
Figure BDA0001153990530000122
Figure BDA0001153990530000123
Figure BDA0001153990530000124
Figure BDA0001153990530000125
将(29)式代入(27)式和(30)式，(26)式代入(28)式得：
Figure BDA0001153990530000126
Figure BDA0001153990530000127
Figure BDA0001153990530000128
对上面三式进行中心差分，离散化后，得到：
Figure BDA0001153990530000131
Figure BDA0001153990530000132
Figure BDA0001153990530000141
上面三式中，i表示横坐标上的第i个计算网格，j表示纵坐标上的第j个计算网格，k表示电磁场系数迭代次数；在整个计算区域上，(34)式和(35)式可以写成三对角矩阵差分方程，与传统WLP-FDTD方法相比，这种迭代的WLP-FDTD方法将大型稀疏矩阵方程的求解转变成两个三对角矩阵方程的求解，于是可以使用追赶法，非常简单的解得整个计算区域电磁场分量系数，最后通过(9)式解得观测点的电磁场分量。
本发明基于拉盖尔多项式的减少分裂误差的完全匹配层实现方法，如图1所示，具体按照以下步骤实施：
步骤1：输入模型文件；
输入的模型文件具体为：计算区域大小Nx×Ny，其中Nx为x方向的网格数，Ny为y方向的网格数；空间步长Δζ，ζ＝x，y，x为横坐标，y为纵坐标；时间步长Δt；真空中的电导率σ，磁导率μ0，介电常数ε0；吸收边界层数NPML与相关参数κζmax，σζmax，αζmax；其中，κζmax取整数，κζmax取值范围为[1，60]；αζmax取值范围为[0，1)；σζmax/σopt取值范围为(0，12]；仿真计算时长Tf；迭代次数k，k≥0且为整数；加权拉盖尔多项式的阶数q，q≥0且为整数；时间尺度因子s，s取值范围为[109，1013]；观测点；场源参数。
步骤2：初始化参数和设置参数；
初始化的参数包括：
将整个计算区域的电磁场分量系数
Figure BDA0001153990530000142
整个计算区域的电磁场分量系数的和
Figure BDA0001153990530000151
整个计算区域的辅助变量
Figure BDA0001153990530000152
拉盖尔多项式
Figure BDA0001153990530000153
均初始化为零，其中Fη＝Ex,Ey,Hz，ζ＝x，y，
Figure BDA0001153990530000154
初始化PML系数(C1ζ,C2ζ,C3ζ)，具体为：
C1ζ＝1/(1+0.5ε0s)
C2ζ＝1
C3ζ＝ε0/μ0
其中，ζ＝x，y，ε0是空气中的介电常数，s为时间尺度因子，s取值范围为[109，1013]；
设置的参数具体为：
设置CFS-PML吸收边界的参数，具体为：
σζ＝σζmax|ζ-ζ0|m/dm
κζ＝1+(κζmax-1)|ζ-ζ0|m/dm
αζ＝αζmaxζ0/d
式中ζ＝x，y，ζ0为PML层与非PML截面位置，d是PML吸收边界的厚度，κζmax取整数，κζmax取值范围为[1，60]；αζmax取值范围为[0，1)；σζmax根据σopt来设置，σζmax/σopt取值范围为(0，12]；σopt＝(m+1)/150πΔζ，m取值范围为[1，20]，其中m取值为4时边界的吸收效果最好，Δζ取值范围
Figure BDA0001153990530000155
λ为源的波长；
设置PML系数，具体为：
Figure BDA0001153990530000156
步骤3：添加场源到y方向上的电场分量系数中，使用因式分裂的WLP-FDTD方法计算电场分量系数
Figure BDA0001153990530000161
记为初始场值
Figure BDA0001153990530000162
其中，所添加场源的表达式为：
Figure BDA0001153990530000163
其中，Tc，Td为场源参数。
步骤4：更新计算整个计算区域的y方向上电场分量系数
Figure BDA0001153990530000164
步骤4.1：电场分量系数
Figure BDA0001153990530000165
在计算区域的方程为：
Figure BDA0001153990530000166
其中，i表示横坐标上的第i个计算网格，j表示纵坐标上的第j个计算网格；
步骤4.2：使用追赶法求解步骤4.1的方程，得到整个计算区域的电场分量系数
Figure BDA0001153990530000167
步骤5：更新计算整个计算区域的x方向上电场分量系数
Figure BDA0001153990530000171
步骤5.1：电场分量系数
Figure BDA0001153990530000172
在计算区域的方程为：
Figure BDA0001153990530000173
步骤5.2：使用追赶法求解步骤5.1的方程，得到整个计算区域的x方向上电场分量系数
Figure BDA0001153990530000174
步骤6：将k+1赋值给k，并判断迭代次数k是否达到预设值，若未达到预设值，则返回步骤4，若达到预设值，则执行步骤7；
步骤7：更新计算整个计算区域的磁场分量系数
Figure BDA0001153990530000175
具体更新公式为：
Figure BDA0001153990530000176
步骤8：更新计算整个计算区域的电磁场分量系数的辅助变量，具体更新公式为：
Figure BDA0001153990530000181
Figure BDA0001153990530000182
其中，Fη＝Ex,Ey,Hz；ζ＝x，y。
步骤9：更新计算观测点处的电磁场分量，具体更新公式为：
Figure BDA0001153990530000183
其中，U表示电磁场分量Ex,Ey,Hz，Uq表示q阶电磁场分量系数，
Figure BDA0001153990530000184
是q阶加权拉盖尔多项式，
Figure BDA0001153990530000185
是带有时间尺度因子s＞0的扩展时间，
Figure BDA0001153990530000186
是q阶拉盖尔多项式。
步骤10：将q+1赋值给q，并判断拉盖尔多项式的阶数q是否达到预设值，若未达到预设值，则返回步骤3，若达到预设值，则结束。
实施例
点源辐射的计算
按照本发明的方法步骤进行实施，如图2所示，实验中整个计算区域为50×50网格，网格大小为1cm×1cm，即Δx＝Δy＝1cm。四个边界采用8层网格的PML吸收边界，计算中所加的源位于网格(25，25)，所加场源的表达式如下：
Figure BDA0001153990530000187
其中，Tc＝5ns,Td＝1ns。观测点位于(40，40)网格处。时间步长Δt＝117.85p,s迭代次数为k＝2，加权拉盖尔多项式的阶数q＝120，时间扩展因子s＝2.5×1010，整个仿真时间为Tf＝16.5ns，PML吸收边界参数κζmax＝19，σζmax＝0.8×σopt，αζmax＝0.0161，m＝3。采用本发明方法计算的观测点处的电场分量Ey与采用传统FDFD方法及因式分裂的WLP-FDTD方法计算的结果参见图3。从图3中可见，本发明方法与传统FDTD方法计算结果一致，验证了本发明方法的正确性，且本发明方法的计算精度比因式分裂的WLP-FDTD方法的计算精度更高。图4为观测点的不同吸收边界相对反射误差，其计算公式可以表示为：
Figure BDA0001153990530000191
其中，Epml为当存在SC-PML吸收边界时，观测点的时域波形，Eref(t)为参考波形，max|Eref(t)|为参考波形绝对值的最大值。由图4可知，带有CFS因子的SC-PML吸收边界最大的反射误差为-70dB，它比无CFS因子的PML吸收边界的吸收效果更好，说明CFS因子可以改善吸收边界的性能。

Claims (6)

1.基于拉盖尔多项式的减少分裂误差的完全匹配层实现方法，其特征在于，具体按照以下步骤实施：
步骤1：输入模型文件；
步骤2：初始化参数和设置参数
初始化的参数包括：
将整个计算区域的电磁场分量系数
Figure FDA0002259299250000011
整个计算区域的电磁场分量系数的和
Figure FDA0002259299250000012
整个计算区域的辅助变量
Figure FDA0002259299250000013
拉盖尔多项式
Figure FDA0002259299250000014
均初始化为零，其中Fη＝Ex,Ey,Hz，ζ＝x，y，
Figure FDA0002259299250000015
初始化PML系数(C1ζ,C2ζ,C3ζ)，具体为：
C1ζ＝1/(1+0.5ε0s)
C2ζ＝1
C3ζ＝ε0/μ0
其中，ζ＝x，y，ε0是空气中的介电常数，s为时间尺度因子，s取值范围为[109，1013]；
设置的参数具体为：
设置CFS-PML吸收边界的参数，具体为：
σζ＝σζmax|ζ-ζ0|m/dm
κζ＝1+(κζmax-1)|ζ-ζ0|m/dm
αζ＝αζmaxζ0/d
式中ζ＝x，y，ζ0为PML层与非PML截面位置，d是PML吸收边界的厚度，κζmax取整数，κζmax取值范围为[1，60]；αζmax取值范围为[0，1)；σζmax根据σopt来设置，σζmax/σopt取值范围为(0，12]；σopt＝(m+1)/150πΔζ，m取值范围为[1，20]，其中m取值为4时边界的吸收效果最好，Δζ取值范围
Figure FDA0002259299250000021
λ为源的波长；
设置PML系数，具体为：
C1ζ＝1/(κζαζ+σζ+0.5κζε0s)
C2ζ＝1+2αζ/(ε0s)
C3ζ＝ε0/μ0+2αζ/(μ0s)；
步骤3：添加场源到y方向上的电场分量系数中，使用因式分裂的WLP-FDTD方法计算电场分量系数
Figure FDA0002259299250000022
记为初始场值
Figure FDA0002259299250000023
所述添加的场源的表达式为：
Figure FDA0002259299250000024
其中，Tc，Td为场源参数；
步骤4：更新计算整个计算区域的y方向上电场分量系数
Figure FDA0002259299250000025
步骤4.1：电场分量系数
Figure FDA0002259299250000026
在计算区域的方程为：
Figure FDA0002259299250000027
其中，i表示横坐标上的第i个计算网格，j表示纵坐标上的第j个计算网格；
步骤4.2：使用追赶法求解步骤4.1的方程，得到整个计算区域的电场分量系数
Figure FDA0002259299250000031
步骤5：更新计算整个计算区域的x方向上电场分量系数
Figure FDA0002259299250000032
步骤6：将k+1赋值给k，并判断迭代次数k是否达到预设值，若未达到预设值，则返回步骤4，若达到预设值，则执行步骤7；
步骤7：更新计算整个计算区域的磁场分量系数
Figure FDA0002259299250000033
步骤8：更新计算整个计算区域的电磁场分量系数的辅助变量；
步骤9：更新计算观测点处的电磁场分量；
步骤10：将q+1赋值给q，并判断拉盖尔多项式的阶数q是否达到预设值，若未达到预设值，则返回步骤3，若达到预设值，则结束。
2.根据权利要求1所述的基于拉盖尔多项式的减少分裂误差的完全匹配层实现方法，其特征在于，所述步骤1输入模型文件，具体为：
计算区域大小Nx×Ny，其中Nx为x方向的网格数，Ny为y方向的网格数；空间步长Δζ，ζ＝x，y，x为横坐标，y为纵坐标；时间步长Δt；真空中的电导率σ，磁导率μ0，介电常数ε0；吸收边界层数NPML与相关参数κζmax，σζmax，αζmax；其中，κζmax取整数，κζmax取值范围为[1，60]；αζmax取值范围为[0，1)；σζmax/σopt取值范围为(0，12]；仿真计算时长Tf；迭代次数k，k≥0且为整数；加权拉盖尔多项式的阶数q，q≥0且为整数；时间尺度因子s，s取值范围为[109，1013]；观测点；场源参数。
3.根据权利要求1所述的基于拉盖尔多项式的减少分裂误差的完全匹配层实现方法，其特征在于，所述步骤5具体为：
步骤5.1：电场分量系数
Figure FDA0002259299250000041
在计算区域的方程为：
Figure FDA0002259299250000042
步骤5.2：使用追赶法求解步骤5.1的方程，得到整个计算区域的x方向上电场分量系数
Figure FDA0002259299250000043
4.根据权利要求3所述的基于拉盖尔多项式的减少分裂误差的完全匹配层实现方法，其特征在于，所述步骤7更新计算整个计算区域的磁场分量系数
Figure FDA0002259299250000044
具体更新公式为：
Figure FDA0002259299250000045
5.根据权利要求1所述的基于拉盖尔多项式的减少分裂误差的完全匹配层实现方法，其特征在于，所述步骤8更新计算整个计算区域的电磁场分量系数的辅助变量，具体更新公式为：
Figure FDA0002259299250000051
Figure FDA0002259299250000052
其中，Fη＝Ex，Ey，Hz；ζ＝x，y。
6.根据权利要求5所述的基于拉盖尔多项式的减少分裂误差的完全匹配层实现方法，其特征在于，所述步骤9更新计算观测点处的电磁场分量，具体更新公式为：
Figure FDA0002259299250000053
其中，U表示电磁场分量Ex，Ey，Hz，Uq表示q阶电磁场分量系数，
Figure FDA0002259299250000054
是q阶加权拉盖尔多项式，
Figure FDA0002259299250000055
是带有时间尺度因子s＞0的扩展时间，
Figure FDA0002259299250000056
是q阶拉盖尔多项式。
CN201611009431.XA 2016-11-16 2016-11-16 基于拉盖尔多项式的减少分裂误差的完全匹配层实现方法 Active CN106777472B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201611009431.XA CN106777472B (zh) 2016-11-16 2016-11-16 基于拉盖尔多项式的减少分裂误差的完全匹配层实现方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201611009431.XA CN106777472B (zh) 2016-11-16 2016-11-16 基于拉盖尔多项式的减少分裂误差的完全匹配层实现方法

Publications (2)

Publication Number Publication Date
CN106777472A CN106777472A (zh) 2017-05-31
CN106777472B true CN106777472B (zh) 2020-05-22

Family

ID=58969418

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201611009431.XA Active CN106777472B (zh) 2016-11-16 2016-11-16 基于拉盖尔多项式的减少分裂误差的完全匹配层实现方法

Country Status (1)

Country Link
CN (1) CN106777472B (zh)

Families Citing this family (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN110598225B (zh) * 2018-06-13 2023-03-24 南京理工大学 一种改进的分裂场时域有限差分算法
CN111259587A (zh) * 2020-01-16 2020-06-09 中国人民解放军陆军工程大学 一种三维交替迭代无条件稳定fdtd算法
CN113917526B (zh) * 2020-07-10 2023-07-28 中国石油化工股份有限公司 基于非分裂完全匹配层吸收边界的正演模拟方法
CN114417667A (zh) * 2022-01-17 2022-04-29 厦门大学 基于时域有限差分的双曲超材料的完美匹配层方法

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN104794289A (zh) * 2015-04-23 2015-07-22 西安理工大学 一种扩展直角坐标系下完全匹配吸收边界的实现方法
CN104809286A (zh) * 2015-04-23 2015-07-29 西安理工大学 一种等离子体中扩展坐标的完全匹配吸收边界的实现方法
CN104809343A (zh) * 2015-04-23 2015-07-29 西安理工大学 一种等离子体中使用电流密度卷积完全匹配层的实现方法

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN104794289A (zh) * 2015-04-23 2015-07-22 西安理工大学 一种扩展直角坐标系下完全匹配吸收边界的实现方法
CN104809286A (zh) * 2015-04-23 2015-07-29 西安理工大学 一种等离子体中扩展坐标的完全匹配吸收边界的实现方法
CN104809343A (zh) * 2015-04-23 2015-07-29 西安理工大学 一种等离子体中使用电流密度卷积完全匹配层的实现方法

Also Published As

Publication number Publication date
CN106777472A (zh) 2017-05-31

Similar Documents

Publication Publication Date Title
CN106777472B (zh) 基于拉盖尔多项式的减少分裂误差的完全匹配层实现方法
Özgün et al. MATLAB-based finite element programming in electromagnetic modeling
Luo et al. An accurate, fast, matrix-free implicit method for computing unsteady flows on unstructured grids
Vipiana et al. EFIE modeling of high-definition multiscale structures
CN113158527B (zh) 一种基于隐式fvfd计算频域电磁场的方法
Rallabhandi et al. Sonic-boom mitigation through aircraft design and adjoint methodology
Li et al. A low dissipation numerical scheme for implicit large eddy simulation
CN104794289B (zh) 一种扩展直角坐标系下完全匹配吸收边界的实现方法
CN111079278B (zh) 三维时域杂交间断伽辽金方法外加电磁源项的处理方法
Li et al. A new fast multipole boundary element method for two dimensional acoustic problems
Misaka et al. Numerical study on jet-wake vortex interaction of aircraft configuration
Li et al. POD-based model order reduction with an adaptive snapshot selection for a discontinuous Galerkin approximation of the time-domain Maxwell's equations
CN107016184B (zh) 一种二维高精度迭代的非磁化等离子体中的实现方法
Mao et al. Calderón preconditioned spectral-element spectral-integral method for doubly periodic structures in layered media
Druskin et al. An extended Krylov subspace model-order reduction technique to simulate wave propagation in unbounded domains
Reckinger et al. Adaptive volume penalization for ocean modeling
Hicken et al. The role of dual consistency in functional accuracy: error estimation and superconvergence
CN108090296B (zh) 基于高阶辛紧致格式的波导全波分析方法
CN112307639B (zh) 一种基于高品质算法的贝伦格完全匹配层仿真方法
Sattarzadeh et al. 3D implicit mesh-less method for compressible flow calculations
Prinn Efficient finite element methods for aircraft engine noise prediction
Acosta et al. Finite difference on grids with nearly uniform cell area and line spacing for the wave equation on complex domains
CN109190140B (zh) 基于浸入式边界方法的连续元音生成方法
Jin et al. Investigation of high-order cell-centered finite difference method for aeroacoustics
CN116401921B (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