CN113505516A - 一种用于低频大地电磁法三维正演的边界截断层方法及装置 - Google Patents

一种用于低频大地电磁法三维正演的边界截断层方法及装置 Download PDF

Info

Publication number
CN113505516A
CN113505516A CN202110971573.9A CN202110971573A CN113505516A CN 113505516 A CN113505516 A CN 113505516A CN 202110971573 A CN202110971573 A CN 202110971573A CN 113505516 A CN113505516 A CN 113505516A
Authority
CN
China
Prior art keywords
boundary
layer
parameter
boundary truncation
truncation layer
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.)
Granted
Application number
CN202110971573.9A
Other languages
English (en)
Other versions
CN113505516B (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.)
Institute of Geology and Geophysics of CAS
Original Assignee
Institute of Geology and Geophysics of CAS
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 Institute of Geology and Geophysics of CAS filed Critical Institute of Geology and Geophysics of CAS
Publication of CN113505516A publication Critical patent/CN113505516A/zh
Application granted granted Critical
Publication of CN113505516B publication Critical patent/CN113505516B/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]
    • 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/11—Complex mathematical operations for solving equations, e.g. nonlinear equations, general mathematical optimization problems
    • G06F17/13—Differential equations
    • 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)
  • General Physics & Mathematics (AREA)
  • Mathematical Physics (AREA)
  • Theoretical Computer Science (AREA)
  • Mathematical Analysis (AREA)
  • Pure & Applied Mathematics (AREA)
  • Mathematical Optimization (AREA)
  • Computational Mathematics (AREA)
  • General Engineering & Computer Science (AREA)
  • Data Mining & Analysis (AREA)
  • Geometry (AREA)
  • Evolutionary Computation (AREA)
  • Computer Hardware Design (AREA)
  • Operations Research (AREA)
  • Algebra (AREA)
  • Databases & Information Systems (AREA)
  • Software Systems (AREA)
  • Hall/Mr Elements (AREA)
  • Complex Calculations (AREA)
  • Measuring Magnetic Variables (AREA)

Abstract

本申请公开了一种用于低频大地电磁法三维正演的方法和装置，所述方法包括：低频电磁波从目标区域中无反射地入射到各向异性边界截断层中，确定各向异性边界截断层的调制参数表达式；对所述边界截断层的调制参数进行渐变处理，获得透射波的吸收参数，确定有限厚度的边界截断层的最佳调制参数值；根据低频大地电磁正演采取的方法，将有限厚度的边界截断层的调制参数引入到电磁场麦克斯韦方程中，得到低频大地电磁法正演参数信息。可以对大地电磁法这类低频电磁方法的目标区域进行有效截断。

Description

一种用于低频大地电磁法三维正演的边界截断层方法及装置
技术领域
本申请涉及地球物理电磁探测领域，特别是涉及一种用于低频大地电磁法三维正演的边界截断层(Boundary Truncation Layer，BTL)方法及装置。
背景技术
目前，对电磁波的目标区域都需要人为地进行截断，由于当前计算机的计算能力有限，不可能进行无限区域的电磁场计算。在电磁波领域，以往模拟中对边界条件是采用大小渐增的网格来延伸几倍的趋肤深度，将二次场衰减为0，从而在模拟区域外边界处设置Dirichlet边界条件或Neumann边界条件，可以让外行波无反射地入内并使其快速衰减殆尽。但是此方式需要延拓距离足够远，而且针对于较宽的频率范围，延伸的网格大小需要缓慢增涨，故需要延拓的网格数目较多，需要较大的计算量以及消耗较多的计算时间。
发明内容
本申请的目的在于克服上述问题或者至少部分地解决或缓减解决上述问题。
本申请中提出了针对于扩散场的边界截断层公式，根据本申请的一个方面，提供了一种用于低频大地电磁法三维正演的边界截断层方法，可以对大地电磁法这类低频电磁方法的目标区域进行有效截断。所述方法包括：
低频电磁波从目标区域中无反射地入射到各向异性边界截断层中，确定各向异性边界截断层的调制参数表达式；
对所述边界截断层的调制参数进行渐变处理，获得透射波的吸收参数，确定有限厚度的边界截断层的最佳调制参数值；
根据低频大地电磁正演采取的方法，将有限厚度的边界截断层的调制参数引入到电磁场麦克斯韦方程中，得到低频大地电磁法正演参数信息。
优选地，确定各向异性边界截断层的调制参数表达式包括：
各向异性介质的电导率为σ·Λ，介电常数为ε·Λ，磁导率为μ·Λ，其中，Λ为调制关键参数，且满足如下的形式：
Figure BDA0003225905270000021
其中，sz为各向异性主轴参数；
设定时间因子为e-iωt，频域的麦克斯方程表达式为：
Figure BDA0003225905270000022
Figure BDA0003225905270000023
其中，E为电场矢量，H为磁场矢量，ω为角频率，σ为目标区域介质的电导率，μ为目标区域介质的磁导率，ε为目标区域介质的介电常数，
Figure BDA0003225905270000024
为虚数，其中σE表示扩散电流，iωεE表示位移电流；
推倒出：
Figure BDA0003225905270000025
令
Figure BDA0003225905270000026
计算出
Figure BDA0003225905270000027
设定边界截断层的调制参数si，i＝x,y,z，表达式为
Figure BDA0003225905270000028
其中，σi为边界截断层的关键参数，控制边界截断层中电磁波的衰减，得到
Figure BDA0003225905270000029
优选地，确定各向异性边界截断层的调制参数表达式之后还包括：对所述边界截断层的调制参数进行简化。
优选地，对所述边界截断层的调制参数进行简化包括：
在高频极限ωε＞＞σ时，si，i＝x,y,z的表达式：
Figure BDA00032259052700000210
当时间因子为eiωt时，得到
Figure BDA00032259052700000211
在低频极限ωε＜＜σ时，将
Figure BDA0003225905270000031
简化为
Figure BDA0003225905270000032
优选地，确定各向异性边界截断层的调制参数表达式之后还包括：
根据目标区域的电磁波表达式和所述调制参数确定边界截断层的透射波的表达式。
优选地，根据目标区域的电磁波表达式和所述调制参数确定的边界截断层透射波的表达式包括：
在扩散方程中，波数等于
Figure BDA0003225905270000033
根据相位调制原则，入射波在x方向上的波数kix等于透射波在x方向上波数ktx，在z方向上，入射波的波数kiz与透射波的波数ktz满足如下关系
ktz＝szkiz
在目标区域中，电磁波表示为
Figure BDA0003225905270000034
其中，ei表示电场的偏振方向，E0表示入射电场幅度；
得到边界截断层的透射波表示为
Figure BDA0003225905270000035
设定入射波的入射方向与分解面的夹角为θ，则z方向上的波数表达式为
Figure BDA0003225905270000036
得到
Figure BDA0003225905270000037
优选地，对所述边界截断层的调制参数进行渐变处理，获得透射波的吸收参数，确定有限厚度的边界截断层的最佳调制参数值包括：
将边界截断层的调制参数
Figure BDA0003225905270000038
修改为
Figure BDA0003225905270000039
的形式，对入射波吸收达到预设效果；
其中，κi和αi为两个控制参数，κi控制透射波的传播方向和对倏逝波的吸收，αi控制透射波在边界截断层中不同深度的衰减速度，σi、κi和αi共同协作，一同控制着透射波的衰减；
设置有限厚度的边界截断层，并让电磁波在边界截断层外边界处衰减为0，在边界截断层最外侧设置Direchlet边界条件；
对边界截断层参数si进行渐变处理，在目标区域与边界截断层界面处从1开始平缓上升，在远离界面处上升速率加快；对边界截断层参数si中的参数σi、κi和αi设置归纳如下：
Figure BDA0003225905270000041
其中，d为边界截断层的总厚度，x为边界截断层不同单元到目标区域与边界截断层分界面的距离；m为调节参数变化速率的一个常数；dz为边界截断层单元的厚度；κmax、σmax、αmax分别为控制参数σi、κi和αi大小的关键参数；κmax、σmax、αmax的最佳值κopt、σopt、αopt；通过数值试验选取。
优选地，根据低频大地电磁正演采取的方法，将有限厚度的边界截断层的调制参数引入到电磁场麦克斯韦方程中，得到低频大地电磁法正演参数信息包括：
对于二次场，采用边界截断层进行网格截断和吸收；
调制参数Λ作为介质的本构参数矩阵，将Λ引入到麦克斯韦方程中
Figure BDA0003225905270000042
Figure BDA0003225905270000043
电场的双旋度表达式为：
Figure BDA0003225905270000044
在目标区域中，调制参数Λ等于单位矩阵；在边界截断层中，调制参数Λ表示为：
Figure BDA0003225905270000045
地球背景模型的电磁波初始场和背景模型中包含异常体时的总场均符合电场的双旋度表达式，将总场公式减去初始场公式，可以得到大地电磁法二次场正演公式：
Figure BDA0003225905270000051
其中，σp为地球背景介质的电导率分布，σ为地球背景介质中加入异常体后的电导率分布，Ep为地球背景模型的初始场，采用一维解析解计算；
通过有限元方法求解二次场正演公式的偏微分方程，得到二次场后，再加上初始场Ep，得到电场总场场值；
通过电场和磁场计算得到阻抗信息，进一步计算出视电阻率和相位。
第二方面，本发明提供一种用于低频大地电磁法三维正演的边界截断层装置，包括：
参数确定模块，设置为低频电磁波从目标区域中无反射地入射到各向异性边界截断层中，确定各向异性边界截断层的调制参数表达式；
渐变处理模块，设置为对所述边界截断层的调制参数进行渐变处理，获得透射波的吸收参数，确定有限厚度的边界截断层的最佳调制参数值；
正演处理模块，设置为根据低频大地电磁正演采取的方法，将有限厚度的边界截断层的调制参数引入到电磁场麦克斯韦方程中，得到低频大地电磁法正演参数信息。
第三方面，本发明提供一种计算设备，包括存储器、处理器和存储在所述存储器内并能由所述处理器运行的计算机程序，其中，所述处理器执行所述计算机程序时实现上述的方法。
本申请采用基于扩散场的边界截断层，可以对大地电磁法这类低频电磁方法的目标区域进行有效截断。将该边界截断层应用到有限单元法大地电磁法正演中，求解目标区域内的场值分布，进而计算出正演所需的视电阻率和相位信息。
本发明实施例根据下文结合附图对本申请的具体实施例的详细描述，本领域技术人员将会更加明了本申请的上述以及其他目的、优点和特征。
附图说明
后文将参照附图以示例性而非限制性的方式详细描述本申请的一些具体实施例。附图中相同的附图标记标示了相同或类似的部件或部分。本领域技术人员应该理解，这些附图未必是按比例绘制的。附图中：
图1是本申请实施例的一种用于低频大地电磁法三维正演的边界截断层方法的示意性流程图；；
图2是本申请实施例的一种用于低频大地电磁法三维正演的边界截断层装置的结构示意图；
图3是本申请实施例的目标区域与边界截断层的示意图；
图4是本申请实施例的边界截断层截断的三维COMMEMI3D-1A模型图；
图5是本申请实施例的不同频率的趋肤深度示意图；
图6是本申请实施例的频率为0.1Hz时，电场场值对比图，其中，(a)在x＝125,y＝0m处，总场在z方向上的场值对比图，(b)在x＝125,z＝1250m处，二次场Ex在x方向上的场值对比图，(c)在x＝125,z＝1250m处，二次场在y方向上的场值对比，(d)在x＝125,y＝0m处，二次场在z方向上的场值对比图；
图7是本申请实施例的频率为10Hz时，电场场值对比图，其中，(a)在x＝125,y＝0m处，总场在z方向上的场值对比图，(b)在x＝125,z＝1250m处，二次场Ex在x方向上的场值对比图，(c)在x＝125,z＝1250m处，二次场在y方向上的场值对比，(d)在x＝125,y＝0m处，二次场在z方向上的场值对比图；；
图8是本申请实施例的频率为0.1Hz时，视电阻率和相位的正演结果图；
图9是本申请实施例的频率为10Hz时，视电阻率和相位的正演结果图；
图10是本申请实施例的一种计算设备的结构示意图；
图11是本申请实施例的一种计算机可读存储介质的结构示意图。
具体实施方式
图1是根据本申请一个实施例的一种用于低频大地电磁法三维正演的边界截断层方法的示意性流程图。本发明实施例的一种用于低频大地电磁法三维正演的边界截断层方法，包括：
S110、低频电磁波从目标区域中无反射地入射到各向异性边界截断层中，确定各向异性边界截断层的调制参数表达式；
S120、对所述边界截断层的调制参数进行渐变处理，获得透射波的吸收参数，确定有限厚度的边界截断层的最佳调制参数值；
S130、根据低频大地电磁正演采取的方法，将有限厚度的边界截断层的调制参数引入到电磁场麦克斯韦方程中，得到低频大地电磁法正演参数信息。
本发明实施例中，步骤S110确定各向异性边界截断层的调制参数表达式包括：
各向异性介质的电导率为σ·Λ，介电常数为ε·Λ，磁导率为μ·Λ，其中，Λ为调制关键参数，且满足如下的形式：
Figure BDA0003225905270000071
其中，sz为各向异性主轴参数；
设定时间因子为e-iωt，频域的麦克斯方程表达式为：
Figure BDA0003225905270000072
Figure BDA0003225905270000073
其中，E为电场矢量，H为磁场矢量，ω为角频率，σ为目标区域介质的电导率，μ为目标区域介质的磁导率，ε为目标区域介质的介电常数，
Figure BDA0003225905270000074
为虚数，其中σE表示扩散电流，iωεE表示位移电流；
推倒出：
Figure BDA0003225905270000075
令
Figure BDA0003225905270000076
计算出
Figure BDA0003225905270000077
设定边界截断层的调制参数si，i＝x,y,z，表达式为
Figure BDA0003225905270000078
其中，σi为边界截断层的关键参数，控制边界截断层中电磁波的衰减，得到
Figure BDA0003225905270000079
本发明实施例中，确定各向异性边界截断层的调制参数表达式之后还包括：对所述边界截断层的调制参数进行简化。
本发明实施例中，对所述边界截断层的调制参数进行简化包括：
在高频极限ωε＞＞σ时，si，i＝x,y,z的表达式：
Figure BDA00032259052700000710
当时间因子为eiωt时，得到
Figure BDA0003225905270000081
在低频极限ωε＜＜σ时，将
Figure BDA0003225905270000082
简化为
Figure BDA0003225905270000083
该公式可以应用于低频电磁场的截断。
本发明实施例中，sz为各向异性主轴参数，其公式的选择至关重要。确定各向异性边界截断层的调制参数表达式之后还包括：
根据目标区域的电磁波表达式和所述调制参数确定边界截断层的透射波的表达式。具体包括：
在扩散方程中，波数等于
Figure BDA0003225905270000084
根据相位调制原则，入射波在x方向上的波数kix等于透射波在x方向上波数ktx，在z方向上，入射波的波数kiz与透射波的波数ktz满足如下关系，
ktz＝szkiz
其中，ei表示电场的偏振方向，E0表示入射电场幅度，得到边界截断层中的透射波表示为
Figure BDA0003225905270000085
设定入射波的入射方向与分解面的夹角为θ，则z方向上的波数表达式为
Figure BDA0003225905270000086
得到
Figure BDA0003225905270000087
本发明实施例中，步骤S120，对所述边界截断层的调制参数进行渐变处理，获得透射波的吸收参数，确定有限厚度的边界截断层的最佳调制参数值包括：
为了提高边界截断层对低频电磁波的吸收效果，将边界截断层的调制参数
Figure BDA0003225905270000088
修改为
Figure BDA0003225905270000089
的形式，对入射波吸收达到预设效果；可以达到对入射波更好的吸收效果；
其中，κi和αi为两个控制参数，κi控制透射波的传播方向和对倏逝波的吸收，αi控制透射波在边界截断层中不同深度的衰减速度，σi、κi和αi共同协作，一同控制着透射波的衰减；
由于计算机的计算能力有限，需要设置有限厚度的边界截断层。让电磁波在边界截断层外边界处衰减为0，在边界截断层最外侧设置Direchlet边界条件；
对边界截断层参数si进行渐变处理，在目标区域与边界截断层界面处从1开始平缓上升，在远离界面处上升速率加快；对边界截断层参数si中的参数σi、κi和αi设置归纳如下：
Figure BDA0003225905270000091
其中，d为边界截断层的总厚度，x为边界截断层不同单元到目标区域与边界截断层分界面的距离；m为调节参数变化速率的一个常数；dz为边界截断层单元的厚度；κmax、σmax、αmax分别为控制参数σi、κi和αi大小的关键参数；κmax、σmax、αmax的最佳值κopt、σopt、αopt需要通过数值试验选取。
其中，si中的参数αi在三维大地电磁正演中，根据低频二次场强于频率高的二次场的要求。
本发明实施例中，步骤S130，根据低频大地电磁正演采取的方法，将有限厚度的边界截断层的调制参数引入到电磁场麦克斯韦方程中，得到低频大地电磁法正演参数信息包括：
对于二次场，采用边界截断层进行网格截断和吸收；
调制参数Λ作为介质的本构参数矩阵，将Λ引入到麦克斯韦方程中
Figure BDA0003225905270000092
Figure BDA0003225905270000093
电场的双旋度表达式为：
Figure BDA0003225905270000094
在目标区域中，调制参数Λ等于单位矩阵；在边界截断层中，调制参数Λ表示为：
Figure BDA0003225905270000101
地球背景模型的电磁波初始场和背景模型中包含异常体时的总场均符合电场的双旋度表达式，将总场公式减去初始场公式，可以得到大地电磁法二次场正演公式：
Figure BDA0003225905270000102
其中，σp为地球背景介质的电导率分布，σ为地球背景介质中加入异常体后的电导率分布，Ep为地球背景模型的初始场，采用一维解析解计算；
通过有限元方法求解二次场正演公式的偏微分方程，得到二次场后，再加上初始场Ep，得到电场总场场值；
通过电场和磁场计算得到阻抗信息，进一步计算出视电阻率和相位。
如图2所示，本发明实施例提供一种三维大地电磁法正演的装置，包括：
参数确定模块100，设置为低频电磁波从目标区域中无反射地入射到各向异性边界截断层中，确定各向异性边界截断层的调制参数表达式；
渐变处理模块200，设置为对所述边界截断层的调制参数进行渐变处理，获得透射波的吸收参数，确定有限厚度的边界截断层的最佳调制参数值；
正演处理模块300，设置为根据低频大地电磁正演采取的方法，将有限厚度的边界截断层的调制参数引入到电磁场麦克斯韦方程中，得到低频大地电磁法正演参数信息。
如图3所示，在x-z的平面上，目标区域在上侧，下侧为一各向异性介质，它们的分界面垂直于z轴。目标区域的电导率为，σ介电常数为，ε磁导率为。μ电磁波从目标区域入射到下侧的各向异性介质中，入射角为θ。倘若下侧的各向异性介质的电导率为σ·Λ，介电常数为ε·Λ，磁导率为μ·Λ，且Λ满足如下的形式
Figure BDA0003225905270000111
则电磁波可以从上侧的目标区域中无反射地入射到下侧各向异性介质中，上侧目标区域介质与下侧各向异性介质阻抗调制，因此下侧的各向异性介质可以称作边界截断调制层。sz理论上可以任意设定，并不会影响电磁波无反射地进入边界截断层。但是为了让进入边界截断层的电磁波快速衰减，本发明实施例还是需要对sz的公式进行精心定义。
首先，本发明实施例假设时间因子为e-iωt，频域的麦克斯方程可以写为
Figure BDA0003225905270000112
Figure BDA0003225905270000113
其中E为电场矢量，H为磁场矢量，ω为角频率，σ为介质的电导率，μ为介质的磁导率，ε为介质的介电常数，
Figure BDA0003225905270000114
为虚数。其中σE可以被视为扩散电流，iωεE可以被看做位移电流。结合公式(3)，公式(2)可以写为
Figure BDA0003225905270000115
定义
Figure BDA0003225905270000116
可以计算出
Figure BDA0003225905270000117
本发明实施例定义边界截断层的参数si(i＝x,y,z)为
Figure BDA0003225905270000118
σi为边界截断层的关键参数，可以控制边界截断调制层中电磁波的衰减。将公式(6)代入公式(7)，得到
Figure BDA0003225905270000119
公式(8)中分母中的σ是指目标区域介质的电导率，与分子中的σi具有完全不同的含义。此边界截断层形式适用于全频带的电磁波目标区域的截断。
在高频极限(ωε＞＞σ)时，公式(8)可以简化为
Figure BDA0003225905270000121
当时间因子为eiωt时，公式(9)应为
Figure BDA0003225905270000122
在低频极限(ωε＜＜σ)时，公式(8)可以简化为
Figure BDA0003225905270000123
在大地电磁中，由于研究的频率较低，往往扩散电流远大于位移电流，满足低频极限ωε＜＜σ。因此可以采用公式(11)的边界截断层进行大地电磁法正演模拟中的网格截断。可以看到本发明实施例提出的扩散场的边界截断层公式(11)与传统的用于高频电磁波的边界截断层公式(9)完全不同，它与目标区域的电导率有关，且与高频边界截断层相比具有不同的频率依赖关系。
在扩散方程中，波数等于
Figure BDA0003225905270000124
如图3所示，根据相位调制原则，入射波在x方向上的波数kix等于透射波在x方向上波数ktx。在z方向上，入射波的波数kiz与透射波的波数ktz满足如下关系
ktz＝szkiz (13)
在目标区域中，电磁波可以表示为
Figure BDA0003225905270000125
所以边界截断层中的透射波可以表示为
Figure BDA0003225905270000126
假设入射波的入射方向与分解面的夹角为θ，则z方向上的波数可以写为
Figure BDA0003225905270000127
将公式(11)代入到公式(15)中，可得到
Figure BDA0003225905270000128
与目标区域的电磁波公式(14)相比，可以看到边界截断层中电磁波公式(16)的前两个指数项为电磁波的自然衰减项，最后一个指数项为边界截断层的额外衰减项。可以看出，边界截断层的衰减能力与参数σz有关，且与入射角度θ有关。电磁波垂直入射时，衰减能力最强。电磁波斜入射时，衰减能力较弱。
考虑到x、y、z三个不同的方向的边界截断层时，得到更一般的公式
Figure BDA0003225905270000131
其中sx、sy与sz为x、y、z方向上的边界截断层参数，均可写为公式(11)。当边界截断层的分界面垂直于x轴时，sx采取公式(11)的形式，公式(17)中的sy和sz看作为1。当边界截断层的分界面垂直于y或z轴时,采用类似的方式。
为了提高边界截断层的参数si对透射波的吸收能力，本发明实施例将公式(11)修改为
Figure BDA0003225905270000132
其中κi和αi为新增添的两个参数，无具体物理含义。κi主要控制透射波的传播方向和对倏逝波的吸收，αi主要控制透射波在边界截断层中不同深度的衰减速度。σi、κi和αi共同协作，一同控制着透射波的衰减。
由于计算机计算能力有限，本发明实施例只能设置有限厚度的边界截断层，并让电磁波快速衰减为0，在边界截断层最外侧设置Direchlet边界条件。为了减小因目标区域与边界截断层在界面处物性骤变所引起的发射，本发明实施例还需要对边界截断层参数si进行一个渐变处理，在目标区域与边界截断层界面处从1开始平缓上升，在远离界面处上升速率加快。对边界截断层参数si中的参数σi、κi和αi设置归纳如下：
Figure BDA0003225905270000133
其中d为边界截断层的总厚度，x为边界截断层不同单元到目标区域与边界截断层分界面的距离。m为调节参数变化速率的一个常数。dz为边界截断层单元的厚度。κmax、σmax、αmax分别为控制参数σi、κi和αi大小的关键参数。κmax、σmax、αmax的最佳值κopt、σopt、αopt可通过数值试验选取。
大地电磁法正演中，通常包含两种算法：总场法和二次场法。总场法同时考虑地球背景介质和异常体介质，进行总场计算。二次场法先计算地球背景介质的初始场值，然后再计算异常体所产生的二次场值，最后将初始场值与二次场值相加得到总场场值。对于二次场，可以看作是外行波，从而采用边界截断层进行网格截断和吸收。矩阵Λ可以看做是介质的本构参数矩阵，将Λ引入到麦克斯韦方程中
Figure BDA0003225905270000141
Figure BDA0003225905270000142
结合公式(20)和(21)，写成电场的双旋度形式
Figure BDA0003225905270000143
在目标区域中，矩阵Λ等于单位矩阵。在边界截断层中，矩阵Λ采用公式(17)的形式。地球背景模型的电磁波初始场和背景模型中包含异常体时的总场都可以写成公式(22)形式，将总场公式减去初始场公式，可以得到大地电磁法二次场正演公式
Figure BDA0003225905270000144
其中σp为地球背景介质的电导率分布，σ为地球背景介质中加入异常体后的电导率分布。Ep为地球背景模型的初始场，采用一维解析解计算。通过有限元方法求解公式(23)的偏微分方程，得到二次场Es后，再加上初始场Ep，得到电场总场场值。通过电场和磁场计算得到阻抗信息，进一步计算出视电阻率和相位。
在大地电磁法正演中，需要考虑地面上方的空气层。空气层的电阻率往往大于106Ω·m，所以空气层的电阻率通常要远大于常见的地球介质电阻率。将边界截断层应用到大地电磁正演中时，空气层的存在是一个很大的障碍。由于边界截断层公式中的分母上存在一个其所包围介质的电导率，通常所包围的介质电导率一致或者悬殊较小时才能取得较好的电磁波吸收效果。在传统的边界截断层的应用中，关键参数是介电常数和磁导率，不同介质的介电常数和磁导率一般悬殊不大，所以边界截断层包围的介质物性差异对边界截断层使用效果的影响有限。但是在大地电磁法正演中，空气层与地球的电导率往往悬殊几千倍甚至更大。所以这给边界截断层的应用带来了很大的难度。幸运的是，经过大量的数值试验发现，如果把所有方向上的边界截断层参数分母中的介质电导率都设置成地球背景介质的电导率后，通过仔细寻找κopt、σopt、αopt参数后，即使存在空气层，也可以取得较好的网格截断效果。选取得到的κopt、σopt、αopt参数，可适用于大地电磁法的宽频带和常用的地球背景介质电阻率范围，改变频率或模型背景介质电阻率后，无需再次寻找最佳参数。
为了测试边界截断层在三维大地电磁法正演中的网格截断效果，本发明实施例建立图4所示的COMMEMI3D-1A模型。COMMEMI3D-1A模型在1997年提出，征集了不同研究者的正演结果，给出了正演结果的大致范围。该模型的空气电阻率设置为108Ω·m，地球背景介质的电阻率设置为100Ω·m。在背景介质中存在一个顶面埋深为250m的低阻棱柱体，其电阻率为0.5Ω·m。低阻棱柱体在x、y、z方向上的长度分别为1000、2000、2000m。本发明实施例将模拟区域离散成立方体单元，立方体单元在x、y方向上均为250m。空气层在z方向上剖分4层，每层厚度为250m，所以空气层的厚度为1000m。在地下介质中，为了计算出可靠的视电阻率和相位信息，在地表下的3层网格的z方向厚度设置为5m，地下第4层和第5层网格的厚度为117.5m。从地下第6层网格起，z方向上的网格厚度均为250m。地下介质z方向上的网格数为15，x和y方向上有效模拟区域的网格数为16，所以在地表的有效观测范围为：x方向上-2000m到2000m,y方向上-2000m到2000m。本发明实施例在模拟区域所有外侧上加载6层网格的边界截断调制层，每层厚度为50m，所以每侧边界截断层的总厚度为300m，远小于该介质中低频电磁波的趋肤深度。所以，本模型包含模拟区域和边界截断调制层区域的网格剖分总数为282831，有限元方程组中的未知数数目为78039。在三维模拟中，本发明实施例采用的边界截断调制层参数为m＝1，κopt＝38.57，σopt＝0.2357，αopt＝100。本发明实施例计算7个频点，分别为0.001、0.01、0.1、1、10、100、1000Hz。
为了与传统方法对比，本发明实施例采用传统的网格延拓方法计算出一个参考解。当地球背景电阻率为100Ω·m时，频率范围0.001～1000Hz内的7个频点的趋肤深度如图5所示。可以看到，不同频率的电磁波在该介质中的趋肤深度悬殊巨大，0.001Hz时趋肤深度可达159km，而1000Hz时仅为0.16km。在大地电磁模拟中，传统的网格截断方案是在模拟区域外围加载大小逐渐增加的一定数量的网格。考虑到不同频率趋肤深度不一，在每个频率的趋肤深度内最少需要2～3个采样点，因此在大地电磁法这样频带极宽的电磁方法模拟中，延拓网格的大小需要逐渐增加才能取得较好的精度。对于二次场法，在延拓的网格的最外侧边界上通常设置Dirichlet边界条件，也就是场值设置为0。所以延拓的网格的总长度需要2～3倍趋肤深度才能满足需求。在该模型中，本发明实施例让每侧的网格延拓数量为16，延拓网格的大小从250m起逐渐以1.55倍的速率递增。每侧网格延拓的总长度为780km，远大于频率为0.001Hz时在该介质中的趋肤深度159km。所以该传统方法作为参考解是可靠的且精确的。网格延拓方法的网格总数量为484851，有限元方程组中的未知数数目为367059，远大于边界截断调制层方法的未知数。在求解计算时间方面，边界截断调制层方法耗时为13分钟，传统的网络延拓方法耗时为4小时36分钟。因此传统的网格延拓方法对延拓的网格数目要求较多，占用内存较大，而且计算时间较长。
本发明实施例接下来仅展示0.1和10Hz时的正演结果。图6为频率为0.1Hz时，有限元求解出的场值对比。图6(a)为边界截断层法与传统网格延拓方法的总场场值对比，可以看出，两者高度吻合，其中紫色虚线为初始场值。图6(b)、(c)和(d)为电场的二次场值在x、y、z方向上的对比，蓝色虚线且带‘o’标记的为边界截断调制层方法的数据，红色实线且带实心点标记的为传统网格延拓方法的数据。从图中可以看出，边界截断层方法求解出来的二次电场与传统大距离网格延拓方法的二次电场拟合的很好，说明本发明实施例提出的边界截断调制层是有效的且具有很高的计算精度的截断边界，而且边界截断调制层法可以有效节省计算量和计算时间。图7(a)、(b)、(c)和(d)是频率为10Hz时的场值对比，与图6(a)、(b)、(c)和(d)具有相同的含义，仅频率不一致。
图8和图9分别为0.1和10Hz时的视电阻率和相位正演结果。图8(a)和8(c)分别为视电阻率和相位的平面图，从图中可以看出，在异常体上方对应的地面观测位置处存在明显的低阻异常。图8(b)和8(d)为y＝0处的x方向上的视电阻率和相位曲线。绿色实线为任政勇(2013)文中的结果，洋红色误差棒为1997年文献中不同研究者给出的视电阻率和相位范围。蓝色虚线且带‘o’标记的为边界截断调制层方法的数据，红色实线且带实心点标记的为传统网格延拓方法的数据。可以看到，边界截断调制层方法的正演结果与传统网格延拓方法的正演结果高度一致，说明边界截断调制层方法应用到大地电磁法正演网格的截断中是可以取得非常精确的正演结果的。采用边界截断调制层方法和网格延拓方法计算出的结果与任政勇的数据存在微小差异，这可能是因为采用了不同的阻抗计算方式造成的，这个差异产生的原因是可以理解和接受的。
本发明提出的边界截断调制层，它可以应用于大地电磁法的正演中。传统的边界截断调制层主要适用于位移电流为主的高频电磁波，不再适用于以扩散电流为主的低频电磁波正演中。本发明提出的边界截断调制层可以很好地截断低频扩散场电磁波的模拟区域，对外行波的吸收性能稳定，且正演的结果精度较高。传统的边界截断调制层参数在边界截断调制层中往往以多项式方式离散，本发明中提出的边界截断调制层以指数方式离散可以取得更佳的效果。空气层对于边界截断调制层应用到大地电磁法正演中是一个很大的障碍，本文在设置空气层外侧边界截断调制层参数时，采用地球背景电导率进行近似设置，通过数值发现，这可能是最好的妥协方式之一。在三维大地电磁法正演中，采用6层边界截断调制层就可以取得较好的网格截断效果。当然边界截断调制层层数越多，截断效果越好，但是计算量也会随之提升。不同的边界截断调制层层数需要不同的最佳边界截断调制层参数，但是对于某层数的边界截断调制层，一旦最佳的边界截断调制层参数已经找好，则可以应用于大地电磁法中常用的频率范围和背景电阻率的范围。正演中频率参数和背景电阻率发生改变时，无需再次寻找最优的边界截断调制层参数。所以本发明实施例在今后大地电磁法正演中可以采用固定的边界截断调制层层数，一直采用相同的最优的边界截断调制层参数，从而可以节省工作量且可以取得满意的网格截断效果。通常，传统的网格延拓方法针对宽频带的大地电磁法频率范围，需要非常长的网格延拓距离，延拓距离甚至可以数百千米。而且还要求延拓的网格大小是渐变的，从而会导致网格延拓方法对延拓的网格数量较大。所以取得一个可靠的结果，可能会占用非常大的内存和计算时间。而把边界截断调制层应用到大地电磁法正演中，6层边界截断调制层就可以取得精度非常高的正演结果，计算占用内存小，计算用时少。所以边界截断调制层方法在大地电磁法正演中是一个十分值得推荐的方法。
本申请实施例还提供了一种计算设备，参照图10，该计算设备包括存储器1120、处理器1110和存储在所述存储器1120内并能由所述处理器1110运行的计算机程序，该计算机程序存储于存储器1120中的用于程序代码的空间1130，该计算机程序在由处理器1110执行时实现用于执行任一项根据本发明的方法步骤1131。
本申请实施例还提供了一种计算机可读存储介质。参照图11，该计算机可读存储介质包括用于程序代码的存储单元，该存储单元设置有用于执行根据本发明的方法步骤的程序1131′，该程序被处理器执行。
本申请实施例还提供了一种包含指令的计算机程序产品。当该计算机程序产品在计算机上运行时，使得计算机执行根据本发明的方法步骤。
在上述实施例中，可以全部或部分地通过软件、硬件、固件或者其任意组合来实现。当使用软件实现时，可以全部或部分地以计算机程序产品的形式实现。所述计算机程序产品包括一个或多个计算机指令。在计算机加载和执行所述计算机程序指令时，全部或部分地产生按照本申请实施例所述的流程或功能。所述计算机可以是通用计算机、专用计算机、计算机网络、获取其他可编程装置。所述计算机指令可以存储在计算机可读存储介质中，或者从一个计算机可读存储介质向另一个计算机可读存储介质传输，例如，所述计算机指令可以从一个网站站点、计算机、服务器或数据中心通过有线(例如同轴电缆、光纤、数字用户线(DSL))或无线(例如红外、无线、微波等)方式向另一个网站站点、计算机、服务器或数据中心进行传输。所述计算机可读存储介质可以是计算机能够存取的任何可用介质或者是包含一个或多个可用介质集成的服务器、数据中心等数据存储设备。所述可用介质可以是磁性介质，(例如，软盘、硬盘、磁带)、光介质(例如，DVD)、或者半导体介质(例如固态硬盘Solid State Disk(SSD))等。
专业人员应该还可以进一步意识到，结合本文中所公开的实施例描述的各示例的单元及算法步骤，能够以电子硬件、计算机软件或者二者的结合来实现，为了清楚地说明硬件和软件的可互换性，在上述说明中已经按照功能一般性地描述了各示例的组成及步骤。这些功能究竟以硬件还是软件方式来执行，取决于技术方案的特定应用和设计约束条件。专业技术人员可以对每个特定的应用来使用不同方法来实现所描述的功能，但是这种实现不应认为超出本申请的范围。
本领域普通技术人员可以理解实现上述实施例方法中的全部或部分步骤是可以通过程序来指令处理器完成，所述的程序可以存储于计算机可读存储介质中，所述存储介质是非短暂性(英文：non-transitory)介质，例如随机存取存储器，只读存储器，快闪存储器，硬盘，固态硬盘，磁带(英文：magnetic tape)，软盘(英文：floppy disk)，光盘(英文：optical disc)及其任意组合。
以上所述，仅为本申请较佳的具体实施方式，但本申请的保护范围并不局限于此，任何熟悉本技术领域的技术人员在本申请揭露的技术范围内，可轻易想到的变化或替换，都应涵盖在本申请的保护范围之内。因此，本申请的保护范围应该以权利要求的保护范围为准。

Claims (10)

1.一种用于低频大地电磁法三维正演的边界截断层方法，包括：
低频电磁波从目标区域中无反射地入射到各向异性边界截断层中，确定各向异性边界截断层的调制参数表达式；
对所述边界截断层的调制参数进行渐变处理，获得透射波的吸收参数，确定有限厚度的边界截断层的最佳调制参数值；
根据低频大地电磁正演采取的方法，将有限厚度的边界截断层的调制参数引入到电磁场麦克斯韦方程中，得到低频大地电磁法正演参数信息。
2.根据权利要求1所述的方法，其特征在于，确定各向异性边界截断层的调制参数表达式包括：
各向异性介质的电导率为σ·Λ，介电常数为ε·Λ，磁导率为μ·Λ，其中，Λ为调制关键参数，且满足如下的形式：
Figure FDA0003225905260000011
其中，sz为各向异性主轴参数；
设定时间因子为e-iωt，频域的麦克斯方程表达式为：
Figure FDA0003225905260000012
Figure FDA0003225905260000013
其中，E为电场矢量，H为磁场矢量，ω为角频率，σ为目标区域介质的电导率，μ为目标区域介质的磁导率，ε为目标区域介质的介电常数，
Figure FDA0003225905260000014
为虚数，其中σE表示扩散电流，iωεE表示位移电流；
推倒出：
Figure FDA0003225905260000015
令
Figure FDA0003225905260000016
计算出
Figure FDA0003225905260000017
设定边界截断层的调制参数si，i＝x,y,z，表达式为
Figure FDA0003225905260000021
其中，σi为边界截断层的关键参数，控制边界截断层中电磁波的衰减，得到
Figure FDA0003225905260000022
3.根据权利要求2所述的方法，其特征在于，确定各向异性边界截断层的调制参数表达式之后还包括：对所述边界截断层的调制参数进行简化。
4.根据权利要求3所述的方法，其特征在于，对所述边界截断层的调制参数进行简化包括：
在高频极限
Figure FDA0003225905260000023
时，si，i＝x,y,z的表达式：
Figure FDA0003225905260000024
当时间因子为eiωt时，得到
Figure FDA0003225905260000025
在低频极限
Figure FDA0003225905260000026
时，将
Figure FDA0003225905260000027
简化为
Figure FDA0003225905260000028
5.根据权利要求2所述的方法，其特征在于，确定各向异性边界截断层的调制参数表达式之后还包括：
根据目标区域的电磁波表达式和所述调制参数确定边界截断层的透射波的表达式。
6.根据权利要求5所述的方法，其特征在于，根据目标区域的电磁波表达式和所述调制参数确定的边界截断层透射波的表达式包括：
在扩散方程中，波数等于
Figure FDA0003225905260000029
根据相位调制原则，入射波在x方向上的波数kix等于透射波在x方向上波数ktx，在z方向上，入射波的波数kiz与透射波的波数ktz满足如下关系
ktz＝szkiz
在目标区域中，电磁波表示为
Figure FDA00032259052600000210
其中，ei表示电场的偏振方向，E0表示入射电场幅度；
得到边界截断层的透射波表示为
Figure FDA0003225905260000031
设定入射波的入射方向与分解面的夹角为θ，则z方向上的波数表达式为
Figure FDA0003225905260000032
得到
Figure FDA0003225905260000033
7.根据权利要求4所述的方法，其特征在于，对所述边界截断层的调制参数进行渐变处理，获得透射波的吸收参数，确定有限厚度的边界截断层的最佳调制参数值包括：
将边界截断层的调制参数
Figure FDA0003225905260000034
修改为
Figure FDA0003225905260000035
的形式，对入射波吸收达到预设效果；
其中，κi和αi为两个控制参数，κi控制透射波的传播方向和对倏逝波的吸收，αi控制透射波在边界截断层中不同深度的衰减速度，σi、κi和αi共同协作，一同控制着透射波的衰减；
设置有限厚度的边界截断层，并让电磁波在边界截断层外边界处衰减为0，在边界截断层最外侧设置Direchlet边界条件；
对边界截断层参数si进行渐变处理，在目标区域与边界截断层界面处从1开始平缓上升，在远离界面处上升速率加快；对边界截断层参数si中的参数σi、κi和αi设置归纳如下：
Figure FDA0003225905260000036
其中，d为边界截断层的总厚度，x为边界截断层不同单元到目标区域与边界截断层分界面的距离；m为调节参数变化速率的一个常数；dz为边界截断层单元的厚度；κmax、σmax、αmax分别为控制参数σi、κi和αi大小的关键参数；κmax、σmax、αmax的最佳值κopt、σopt、αopt；通过数值试验选取。
8.根据权利要求7所述的方法，其特征在于，根据低频大地电磁正演采取的方法，将有限厚度的边界截断层的调制参数引入到电磁场麦克斯韦方程中，得到低频大地电磁法正演参数信息包括：
对于二次场，采用边界截断层进行网格截断和吸收；
调制参数Λ作为介质的本构参数矩阵，将Λ引入到麦克斯韦方程中
Figure FDA0003225905260000041
Figure FDA0003225905260000042
电场的双旋度表达式为：
Figure FDA0003225905260000043
在目标区域中，调制参数Λ等于单位矩阵；在边界截断层中，调制参数Λ表示为：
Figure FDA0003225905260000044
地球背景模型的电磁波初始场和背景模型中包含异常体时的总场均符合电场的双旋度表达式，将总场公式减去初始场公式，可以得到大地电磁法二次场正演公式：
Figure FDA0003225905260000045
其中，σp为地球背景介质的电导率分布，σ为地球背景介质中加入异常体后的电导率分布，Ep为地球背景模型的初始场，采用一维解析解计算；
通过有限元方法求解二次场正演公式的偏微分方程，得到二次场Es后，再加上初始场Ep，得到电场总场场值；
通过电场和磁场计算得到阻抗信息，进一步计算出视电阻率和相位。
9.一种用于低频大地电磁法三维正演的边界截断层装置，包括：
参数确定模块，设置为低频电磁波从目标区域中无反射地入射到各向异性边界截断层中，确定各向异性边界截断层的调制参数表达式；
渐变处理模块，设置为对所述边界截断层的调制参数进行渐变处理，获得透射波的吸收参数，确定有限厚度的边界截断层的最佳调制参数值；
正演处理模块，设置为根据低频大地电磁正演采取的方法，将有限厚度的边界截断层的调制参数引入到电磁场麦克斯韦方程中，得到低频大地电磁法正演参数信息。
10.一种计算设备，包括存储器、处理器和存储在所述存储器内并能由所述处理器运行的计算机程序，其中，所述处理器执行所述计算机程序时实现如权利要求1-8中任一项所述的方法。
CN202110971573.9A 2021-02-25 2021-08-24 一种用于低频大地电磁法三维正演的边界截断层方法及装置 Active CN113505516B (zh)

Applications Claiming Priority (2)

Application Number Priority Date Filing Date Title
CN2021102074283 2021-02-25
CN202110207428 2021-02-25

Publications (2)

Publication Number Publication Date
CN113505516A true CN113505516A (zh) 2021-10-15
CN113505516B CN113505516B (zh) 2022-05-17

Family

ID=75911251

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202110971573.9A Active CN113505516B (zh) 2021-02-25 2021-08-24 一种用于低频大地电磁法三维正演的边界截断层方法及装置

Country Status (2)

Country Link
CN (1) CN113505516B (zh)
AU (2) AU2021101629A4 (zh)

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN114970289B (zh) * 2022-07-25 2022-10-25 中南大学 三维大地电磁各向异性正演数值模拟方法、设备及介质

Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20110114856A1 (en) * 2008-07-15 2011-05-19 Danmarks Tekniske Universitet All-optical control of thz radiation in parallel plate waveguides
CN103969627A (zh) * 2014-05-26 2014-08-06 苏州市数字城市工程研究中心有限公司 基于fdtd的探地雷达大规模三维正演模拟方法
US9538292B1 (en) * 2012-03-23 2017-01-03 Coleridge Design Associates Llc Speaker with voice coil and field coil
US20180143941A1 (en) * 2016-11-23 2018-05-24 Hubei University Of Technology Method of solving stress based on force boundary and balance condition
CN109407160A (zh) * 2018-12-19 2019-03-01 中国科学院地质与地球物理研究所 一种基于wem矢量有限元法进行三维正演的方法
CN110058315A (zh) * 2019-05-29 2019-07-26 中南大学 一种三维各向异性射频大地电磁自适应有限元正演方法
CN110119586A (zh) * 2019-05-21 2019-08-13 中煤科工集团西安研究院有限公司 轴向电导率各向异性瞬变电磁三分量三维fdtd正演方法

Patent Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20110114856A1 (en) * 2008-07-15 2011-05-19 Danmarks Tekniske Universitet All-optical control of thz radiation in parallel plate waveguides
US9538292B1 (en) * 2012-03-23 2017-01-03 Coleridge Design Associates Llc Speaker with voice coil and field coil
CN103969627A (zh) * 2014-05-26 2014-08-06 苏州市数字城市工程研究中心有限公司 基于fdtd的探地雷达大规模三维正演模拟方法
US20180143941A1 (en) * 2016-11-23 2018-05-24 Hubei University Of Technology Method of solving stress based on force boundary and balance condition
CN109407160A (zh) * 2018-12-19 2019-03-01 中国科学院地质与地球物理研究所 一种基于wem矢量有限元法进行三维正演的方法
CN110119586A (zh) * 2019-05-21 2019-08-13 中煤科工集团西安研究院有限公司 轴向电导率各向异性瞬变电磁三分量三维fdtd正演方法
CN110058315A (zh) * 2019-05-29 2019-07-26 中南大学 一种三维各向异性射频大地电磁自适应有限元正演方法

Non-Patent Citations (6)

* Cited by examiner, † Cited by third party
Title
LAURA MARIA SCHMIDT 等: "3D Boundary Conditions in 3D Finite-Element Electromagnetic Forward Modelling", 《GEM 2019 XI"AN: INTERNATIONAL WORKSHOP AND GRAVITY, ELECTRICAL & MAGNETIC METHODS AND THEIR APPLICATIONS》 *
ZHENGYONG REN 等: "A goal-oriented adaptive finite-element approach for plane wave 3-D electromagnetic modelling", 《GEOPHYSICAL JOURNAL INTERNATIONAL》 *
余翔: "时域瞬变电磁三维有限差分正演及广义逆矩阵反演研究", 《中国博士学位论文全文数据库 基础科学辑》 *
喻国: "大地电磁三维任意电各向异性正演模拟及应用研究", 《中国优秀硕士学位论文全文数据库 基础科学辑》 *
杨良勇 等: "基于棱边基有限元的三维WEM各向异性正演", 《2019年中国地球科学联合学术年会论文集（三十一）——专题84：超深层（油气）重磁电震勘探技术、专题85：中国钾盐矿产基地成矿规律与深部探测技术示范》 *
薛帅 等: "耦合PML边界条件的三维大地电磁二次场有限差分法", 《地球物理学报》 *

Also Published As

Publication number Publication date
AU2021101629A4 (en) 2021-05-20
CN113505516B (zh) 2022-05-17
AU2021204490A1 (en) 2022-09-08

Similar Documents

Publication Publication Date Title
CN113221393B (zh) 一种基于非结构有限元法的三维大地电磁各向异性反演方法
CN108875211B (zh) 一种基于fetd与fdtd耦合的二维模型的探地雷达正演方法
CN112949134B (zh) 基于非结构有限元方法的地-井瞬变电磁反演方法
CN113505516B (zh) 一种用于低频大地电磁法三维正演的边界截断层方法及装置
Tabatabaeenejad et al. Coherent scattering of electromagnetic waves from two-layer rough surfaces within the Kirchhoff regime
CN110210129B (zh) 自适应有限元gpr频率域正演方法
CN110502785B (zh) 一种三维时域计算波导s参数的电磁数值方法
Iqbal et al. A novel wavelet-Galerkin method for modeling radio wave propagation in tropospheric ducts
Brooke et al. PECan: A Canadian parabolic equation model for underwater sound propagation
Peeters et al. Calculation of MoM interaction integrals in highly conductive media
Bahar Propagation of VLF radio waves in a model earth‐ionosphere waveguide of arbitrary height and finite surface impedance boundary: Theory and experiment
Wang et al. An analysis of narrow-angle and wide-angle shift-map parabolic equations
Dey et al. Electromagnetic response of two‐dimensional inhomogeneities in a dissipative half‐space for Turam interpretation
Xu et al. The mode‐conversion coefficient of ELF wave propagation due to land‐to‐sea transition
Meng et al. Numerical analysis of multicomponent responses of surface-hole transient electromagnetic method
Warecka et al. Hybrid method analysis of unshielded guiding structures
Feliziani et al. Finite-difference time-domain modeling of thin shields
Meier Date: May 17, 2021 File: SearchPath_Documentation_EN_US_Letter. docx
Nazari et al. On the use of combined surface integral equations for the analysis of high contrast penetrable objects
Meier Date: May 17, 2021 File: SearchPath_Documentation_EN_DIN_A4. docx
CN117094198A (zh) 一种利用ade-pml法快速计算目标电磁场分布的方法
Yang et al. Comparison of response characteristics between the electromagnetic method of “Earth-ionosphere” mode and traditional magnetotellurics
Batista et al. Moving-window propagation model based on an unconditionally stable FDTD method
Wu et al. A conformal symplectic multi-resolution time-domain method with ADE-PML
Wait Influence of lower ionosphere on propagation of VLF waves to great distances

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