CN111103628A - 大地电磁te极化模式对电场数据二维反演方法和装置 - Google Patents

大地电磁te极化模式对电场数据二维反演方法和装置 Download PDF

Info

Publication number
CN111103628A
CN111103628A CN202010035274.XA CN202010035274A CN111103628A CN 111103628 A CN111103628 A CN 111103628A CN 202010035274 A CN202010035274 A CN 202010035274A CN 111103628 A CN111103628 A CN 111103628A
Authority
CN
China
Prior art keywords
inversion
electric field
field data
jacobian matrix
magnetotelluric
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
CN202010035274.XA
Other languages
English (en)
Other versions
CN111103628B (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.)
Guilin University of Technology
Original Assignee
Guilin 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 Guilin University of Technology filed Critical Guilin University of Technology
Priority to CN202010035274.XA priority Critical patent/CN111103628B/zh
Publication of CN111103628A publication Critical patent/CN111103628A/zh
Application granted granted Critical
Publication of CN111103628B publication Critical patent/CN111103628B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G01MEASURING; TESTING
    • G01VGEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V3/00Electric or magnetic prospecting or detecting; Measuring magnetic field characteristics of the earth, e.g. declination, deviation
    • G01V3/40Electric or magnetic prospecting or detecting; Measuring magnetic field characteristics of the earth, e.g. declination, deviation specially adapted for measuring magnetic field characteristics of the earth

Landscapes

  • Engineering & Computer Science (AREA)
  • Remote Sensing (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Environmental & Geological Engineering (AREA)
  • Geology (AREA)
  • Physics & Mathematics (AREA)
  • General Life Sciences & Earth Sciences (AREA)
  • General Physics & Mathematics (AREA)
  • Geophysics (AREA)
  • Measuring Magnetic Variables (AREA)

Abstract

本申请涉及一种大地电磁TE极化模式对电场数据二维反演方法和装置。所述方法包括:根据测线中多个观测点沿X轴的不同频率观测电场数据分量以及观测点沿X轴的不同频率计算电场数据分量,构建大地电磁TE极化模式下的反演目标函数;通过反演目标函数,得到反演目标函数对应的反演迭代方程组;通过引入归一化电场数据分量误差和用计算得到的磁场Hy分量代替实测磁场Hy分量,使得实际拟合的是张量阻抗数据;根据归一化电场数据分量误差的数学表达式,计算雅克比矩阵,根据雅克比矩阵和归一化电场数据分量误差计算反演迭代方程组,根据所述反演迭代方程组对电场数据进行二维反演。采用本方法能够提高二维反演的计算效率和反演精度。

Description

大地电磁TE极化模式对电场数据二维反演方法和装置
技术领域
本申请涉及地磁场二维反演技术领域,特别是涉及一种大地电磁TE极化模式对电场数据二维反演方法和装置。
背景技术
大地电磁测深法(MT)是前苏联学者Tikhonov(1950)和Cagniard(1953)于20世纪50年代提出来的,该方法是通过天然交变电磁场研究地球模型电性性质及其分布特征的一种地球物理探测方法。它通过研究大地对天然电磁场的频率响应,获得地下不同深度介质电阻率分布信息。大地电磁测深方法具有极大的穿透深度和分辨能力,在地壳于上地幔深部构造、矿产于地热资源、油气勘探、环境监测等领域得到了广泛的应用。
大地电磁二维结构中,共包括TE(横电波)和TM(横磁波)两种极化模式,TE极化模式表示沿构造走向方向只有电场分量。目前,大地电磁TE极化模式反演算法中,主要采用视电阻率、阻抗相位及张量阻抗数据进行反演。其中,视电阻率表征了电磁场实部信息,阻抗相位表征了电磁场虚部信息,为了充分利用电磁场全部信息,通常采用视电阻率和阻抗相位联合反演的方式。大地电磁正演计算结果是电磁场分量,从理论上讲,直接对电磁场数据进行拟合是最直接和最有效的方法。由于天然电磁场源具有随机性和不稳定性,因此无法获得稳定的电场数据,常规做法是将观测得到的电磁场分量数据换算到只与地下介质相关的张量阻抗数据或视电阻率与阻抗相位数据。
以大地电磁TE极化模式视电阻率和阻抗相位数据联合反演为例,首先,需要联合反演才能完全利用电磁场全部信息,而联合反演不仅仅增加了一倍的计算量,更是增加了反演的计算难度;其次,与直接对电场数据反演相比,由于视电阻率和阻抗相位在计算过程中增加大量的空间导数和倒数计算,从而增加了反演本身的非线性程度和偏导数矩阵的计算难度。我们知道,非线性程度越高,采用线性(非线性)正则化反演算法得到的反演结果分辨率越低。张量阻抗数据反演与此相同,也同样具有增加反演本身非线性程度和降低反演分辨率的特点。
发明内容
基于此,有必要针对上述技术问题,提供一种能够解决大地地磁二维反演计算复杂以及反演分辨率低问题的大地电磁TE极化模式对电场数据二维反演方法和装置。
一种大地电磁TE极化模式对电场数据二维反演方法,所述方法包括:
根据测线中多个观测点沿X轴方向的不同频率观测电场数据分量以及观测点沿X轴的不同频率计算电场数据分量,构建大地电磁TE极化模式下的反演目标函数;
通过所述反演目标函数,得到所述反演目标函数对应的反演迭代方程组;
根据张量阻抗数据,求解所述反演迭代方程组对应的归一化电场数据分量误差,根据所述归一化电场数据分量误差的数学表达式,计算雅克比矩阵;
根据所述雅克比矩阵和所述归一化电场数据分量误差迭代求解所述反演迭代方程组,从而实现大地电磁TE极化模式对电场数据二维反演。
在其中一个实施例中,还包括:根据测线中多个观测点沿X轴方向的不同频率观测电场数据分量以及观测点沿X轴的不同频率计算电场数据分量,构建大地电磁TE极化模式下的反演目标函数为:
Figure BDA0002365780200000021
其中,ΦTE(m)表示反演目标函数,m=[m1,m2,…,mM]T表示二维地电模型向量,M表示反演节点的数量,mj表示第j个反演节点的电导率,上标T表示转置,
Figure BDA0002365780200000022
表示观测点沿X轴的观测电场数据分量,
Figure BDA0002365780200000023
表示观测点沿X轴的计算电场数据分量,*表示共轭算子,λ表示正则化因子,L表示二阶差分算子。
在其中一个实施例中,还包括:将所述反演目标函数ΦTE(m)在二维地电模型向量m0进行二阶泰勒-拉格朗日展开,去除展开结果中的三阶及以上的高阶项,得到展开表达式;对所述展开表达式的两端对二维地电模型向量m计算一阶导数,通过求解
Figure BDA0002365780200000031
得到所述反演目标函数对应的反演迭代方程组为:
Figure BDA0002365780200000032
其中,G表示雅克比矩阵。
在其中一个实施例中,还包括:根据所述张量阻抗数据,得到相互正交的大地电场和磁场分量的关系为:
Ei=ZijHj i,j∈y,x
根据大地电场和磁场分量的关系式,得到所述反演迭代方程组对应的归一化电场数据分量误差为:
Figure BDA0002365780200000033
设置野外观测时所述观测点磁场沿Y方向分量的初始值
Figure BDA0002365780200000034
对上述归一化电场数据分量误差进行简化,得到:
Figure BDA0002365780200000035
在其中一个实施例中,还包括:根据所述归一化电场属于分量误差的数学表达式,计算雅克比矩阵为:
Figure BDA0002365780200000036
其中,所述雅克比矩阵为N×M的矩阵。
在其中一个实施例中,还包括:根据所述雅克比矩阵,计算所述雅克比矩阵与向量的乘积;根据所述雅克比矩阵,构建海森矩阵,计算所述海森矩阵与向量的乘积;根据所述雅克比矩阵与向量的乘积、所述海森矩阵与向量的乘积以及所述归一化电场数据分量误差,计算所述反演迭代方程组。
在其中一个实施例中,还包括:根据RODI法求解角频率为ω,极化模式为TE时所述雅克比矩阵中的局部因子为:
Figure BDA0002365780200000037
其中,i=1,2,...,N表示观测点序列,j=1,...,M表示反演节点序列,
Figure BDA0002365780200000038
ip表示第i个观测点对应的正演网格节点编号,向量中第ip个元素为1,其余为零,K-1表示有限元正演生成的对称、稀疏复系数矩阵,E表示正演生成的节点电场向量,
Figure BDA0002365780200000041
表示L×L大小的稀疏矩阵;根据所述局部因子,计算所述雅克比矩阵与向量的乘积为:
Figure BDA0002365780200000042
Figure BDA0002365780200000043
其中,向量x=(x1,x2,...,xM)T,y=(y1,y2,...,yN)T
Figure BDA0002365780200000044
Figure BDA0002365780200000045
表示每次迭代中正演计算得到的测点磁场分量。
在其中一个实施例中,还包括:根据所述雅克比矩阵,构建海森矩阵为:
H≈GTG*
其中,H表示海森矩阵;根据所述雅克比矩阵与向量的乘积,计算得到所述海森矩阵与向量的乘积。
一种大地电磁TE极化模式对电场数据二维反演装置,所述装置包括:
目标函数构建模块,用于根据测线中多个观测点沿X轴方向的不同频率观测电场数据分量以及观测点沿X轴的不同频率计算电场数据分量,构建大地电磁TE极化模式下的反演目标函数;
迭代方程构建模块,用于通过所述反演目标函数,得到所述反演目标函数对应的反演迭代方程组;
雅克比矩阵计算模块,用于根据张量阻抗数据,求解所述反演迭代方程组对应的归一化电场数据分量误差,根据所述归一化电场数据分量误差的数学表达式,计算雅克比矩阵;
二维反演模块,用于根据所述雅克比矩阵和所述归一化电场数据分量误差迭代求解所述反演迭代方程组,从而实现大地电磁TE极化模式对电场数据二维反演。
一种计算机设备,包括存储器和处理器,所述存储器存储有计算机程序,所述处理器执行所述计算机程序时实现以下步骤:
根据测线中多个观测点沿X轴方向的不同频率观测电场数据分量以及观测点沿X轴的不同频率计算电场数据分量,构建大地电磁TE极化模式下的反演目标函数;
通过所述反演目标函数,得到所述反演目标函数对应的反演迭代方程组;
根据张量阻抗数据,求解所述反演迭代方程组对应的归一化电场数据分量误差,根据所述归一化电场数据分量误差的数学表达式,计算雅克比矩阵;
根据所述雅克比矩阵和所述归一化电场数据分量误差迭代求解所述反演迭代方程组,从而实现大地电磁TE极化模式对电场数据二维反演。
上述大地电磁TE极化模式对电场数据二维反演方法、装置和计算机设备,通过构建对电场数据的反演目标函数以及根据反演目标函数计算得到反演迭代方程组,通过对反演目标函数归一化处理,使得反演中实际拟合的为阻抗相位数据,与直接对张量阻抗数据反演相比,该方法不会额外增加反演目标函数的非线性程度,保证了反演的稳定性和精度;与对视电阻率和相位数据联合反演相比,该反演方法除了不会额外增加反演目标函数的非线性程度,且一次反演就能包含了电磁场的全部信息,简化了反演计算的难度并提高了反演的效率。
附图说明
图1为一个实施例中大地电磁TE极化模式对电场数据二维反演方法的流程示意图;
图2为另一个实施例中大地电磁TE极化模式对电场数据二维反演方法的流程示意图;
图3为一个实施例中二维理论模型示意图;
图4为一个实施例中大地电磁TE极化模式正演图。
图5为一个实施例中大地电磁TE极化模式迭代收敛图;
图6为一个实施例中大地电磁TE极化模式二维反演剖面图(初始模型100Ωm);
图7为一个实施例中大地电磁TE极化模式二维反演剖面图(初始模型400Ωm);
图8为一个实施例中大地电磁TE极化模式二维反演剖面图(初始模型20Ωm);
图9为一个实施例中大地电磁TE极化模式对电场数据二维反演装置的结构框图;
图10为一个实施例中计算机设备的内部结构图。
具体实施方式
为了使本申请的目的、技术方案及优点更加清楚明白,以下结合附图及实施例,对本申请进行进一步详细说明。应当理解,此处描述的具体实施例仅仅用以解释本申请,并不用于限定本申请。
在一个实施例中,如图1所示,提供了一种大地电磁TE极化模式对电场数据二维反演方法,包括以下步骤:
步骤102,根测线中多个观测点沿X轴方向的不同频率观测电场数据分量以及观测点沿X轴的不同频率计算电场数据分量,构建大地电磁TE极化模式下的反演目标函数。
反演节点指的是反演时所采用的节点,一般的,进行反演时,需要通过线性系统生成反演网格模型,反演节点为反演网格模型中设置的节点。
观测点指的是在进行正演时,正演网格模型中当前观测的点,本实施例中,采用正演模拟计算,得到观测点的张量阻抗数据。具体的,三维电磁场分量包括有:Ex、Ey、Ez、Hx、Hy、Hz共6个分量,TE极化模式下,包含:沿X方向的电场Ex,沿Y轴方向的磁场Hy和沿Z轴方向的磁场Hz
步骤104,通过反演目标函数,得到反演目标函数对应的反演迭代方程组。
通过反演迭代方程组,进行迭代,可以进行电场数据的二维反演。
步骤106,根据张量阻抗数据,求解反演迭代方程组对应的归一化电场数据分量误差,根据归一化电场数据分量误差的数学表达式,计算雅克比矩阵。
张量阻抗数据是将大地电磁TE模式下的观测数据转换得到,观测数据包括视电阻率和阻抗相位,视电阻率和阻抗相位即可以转化为张量阻抗。
本步骤中,采用张量阻抗数据可以求解反演迭代方程组对应的归一化电场数据分量误差。
步骤108,根据雅克比矩阵和归一化电场数据分量误差迭代求解反演迭代方程组,从而实现大地电磁TE极化模式对电场数据二维反演。
本步骤中,通过计算反演迭代方程组中的未知参数,可以得到用于电场数据二维反演的迭代方程,从而通过迭代方程的迭代,进行二维反演。
上述大地电磁TE极化模式对电场数据二维反演方法中,通过构建对电场数据的反演目标函数以及根据反演目标函数计算得到反演迭代方程组,通过对反演目标函数归一化处理,使得反演中实际拟合的为阻抗相位数据,与直接对张量阻抗数据反演相比,该方法不会额外增加反演目标函数的非线性程度,保证了反演的稳定性和精度;与对视电阻率和相位数据联合反演相比,该反演方法除了不会额外增加反演目标函数的非线性程度,且一次反演就能包含了电磁场的全部信息,简化了反演计算的难度并提高了反演的效率。
在其中一个实施例中,构建反演目标函数的方式包括:根据测线中多个观测点沿X轴方向的不同频率观测电场数据分量以及观测点沿X轴的不同频率计算电场数据分量,构建大地电磁TE极化模式下的反演目标函数为:
Figure BDA0002365780200000071
其中,ΦTE(m)表示反演目标函数,m=[m1,m2,...,mM]T表示二维地电模型向量,M表示反演节点的数量,mj表示第j个反演节点的电导率,上标T表示转置,
Figure BDA0002365780200000072
表示观测点沿X轴的观测电场数据分量,
Figure BDA0002365780200000073
表示观测点沿X轴的计算电场数据分量,*表示共轭算子,λ表示正则化因子,L表示二阶差分算子。
本实施例中,分别取归一化电场数据的转置和共轭算子,可以确保得到的数据为实数,另一方面,也利用要了数据的虚部和实部,保证了计算的准确性。
具体的,二阶差分算子可以选择拉普拉斯算子。
在其中一个实施例中,将反演目标函数ΦTE(m)在二维地电模型向量m0进行二阶泰勒-拉格朗日展开,去除展开结果中的三阶及以上的高阶项,得到展开表达式;对展开表达式的两端对二维地电模型向量m计算一阶导数,通过求解
Figure BDA0002365780200000074
得到反演目标函数对应的反演迭代方程组为:
Figure BDA0002365780200000075
其中,G表示雅克比矩阵。
本实施例中,去除高阶项之后,并保留了二阶以下的项,便于进行数学计算,另外,本反演迭代方程组无需对视电阻率和阻抗相位数据反演。
在其中一个实施例中,计算归一化电场数据分量误差的步骤包括:根据张量阻抗数据,得到相互正交的大地电场和磁场分量的关系为:
Ei=ZijHj i,j∈y,x
根据大地电场和磁场分量的关系式,得到反演迭代方程组对应的归一化电场数据分量误差为:
Figure BDA0002365780200000081
设置野外观测时所述观测点沿Y轴方向磁场分量的初始值
Figure BDA0002365780200000082
对上述归一化电场数据分量误差进行简化,得到:
Figure BDA0002365780200000083
本实施例中,在TE极化模式下,利用地电场和磁场分量的关系,将归一化电场数据分量误差进行转化,以及根据野外观测时所述观测点沿Y轴磁场分量的初始值
Figure BDA0002365780200000084
省去归一化电场数据分量误差中的磁场分量,便于进行反演计算。
在其中一个实施例中,根据归一化电场数据分量误差的数学表达式,计算雅克比矩阵为:
Figure BDA0002365780200000085
其中,雅克比矩阵为N×M的矩阵。
在另一个实施例中,根据雅克比矩阵,计算雅克比矩阵与向量的乘积,根据雅克比矩阵,构建海森矩阵,计算海森矩阵与向量的乘积;根据雅克比矩阵与向量的乘积、海森矩阵与向量的乘积以及归一化电场数据分量误差,计算反演迭代方程组。至此,计算出反演迭代方程组中的未知项,可以利用反演迭代方程组进行二维反演。
具体的,首先需要计算角频率为ω,极化模式为TE时,根据Rodi法可以得到:
Figure BDA0002365780200000091
其中,i=1,2,...,N表示观测点序列,j=1,...,M表示反演节点序列,
Figure BDA0002365780200000092
ip表示第i个观测点对应的正演网格节点编号,向量中第ip个元素为1,其余为零,K-1表示有限元正演生成的对称、稀疏复系数矩阵,E表示正演生成的节点电场向量,
Figure BDA0002365780200000093
表示L×L大小的稀疏矩阵,L为正演网格节点个数。
然后,根据局部因子,计算雅克比矩阵与向量的乘积为:
Figure BDA0002365780200000094
Figure BDA0002365780200000095
其中,向量x=(x1,x2,...,xM)T,y=(y1,y2,...,yN)T
Figure BDA0002365780200000096
Figure BDA0002365780200000097
表示每次迭代中正演计算得到的测点磁场分量。
在另一个具体实施例中,根据雅克比矩阵,构建海森矩阵为:
H≈GTG*
其中,H表示海森矩阵;根据雅克比矩阵与向量的乘积,计算得到海森矩阵与向量的乘积。
具体的,计算海森矩阵与向量的乘积可以转化为Gx*(G*x=Gx*)和GTy的计算过程,可以令:
Figure BDA0002365780200000098
每个频率生成的系数矩阵K在正演计算中都采用Pardiso_64位求解器分解成了LU矩阵存储在内存中,上式的求解只需要回代即可。然后再根据向量b可以得到大小为的Gx*向量。
然后,令:
Figure BDA0002365780200000099
由于K-1为对称矩阵,则有(KT)-1=K-1,可以转化为:
Figure BDA0002365780200000101
对上式采用回代求解,根据求解所得的向量b可以得到大小为M×1的GTy向量。从而得到海森矩阵与向量的乘积。
具体的,可以采用共轭梯度法求解大地电磁TE极化模式下的反演迭代方程组。
至于如何采用反演迭代方程组进行电场数据二维反演,以下以一个具体实施例进行说明。
在其中一个实施例中,提供一种大地电磁TE极化模式对电场数据二维反演,包括以下步骤:
步骤202,根据野外观测得到的视电阻率和阻抗相位数据换算成张量阻抗相位数据。
本步骤中,具体的换算表达式为:
Figure BDA0002365780200000102
其中,ω表示角频率,i为虚数单位,μ0为磁导率,ρTE为视电阻率参数,φTE为阻抗相位数据,π为圆周率。
步骤204,给定初始模型,大地电磁TE极化模式进行有限元正演计算,得到观测点的沿X轴方向电场分量和张量阻抗数据。
步骤206,将野外观测得到的张量阻抗相位数据、观测点沿X轴方向的电场分量以及张量阻抗数据进行归一化数据拟合,得到拟合误差。
步骤208,将拟合误差作为迭代终止条件,当迭代终止时,得到二维反演结果。
本步骤中,若迭代未终止,则利用反演迭代方程组计算得到初始模型的修改量,继续进行迭代,从而得到二维反演结果。
本实施例中,采用正演的方式进行数值模拟,然后通过归一化计算结果作为迭代终止条件,得到反演结果。
以下一个具体实施例对本发明的应用场景进行说明。
本应用场景下,对起伏地形下二维隆起构造模型大地电磁TE极化模式二维反演。
如图3所示为理论模型示意图,背景电阻率为100Ωm,隆起构造地形下浅层低阻体电阻率为20Ωm,深部高阻层电阻率为1000Ωm。采用TE极化模式进行有限元数值模拟,在0.1Hz~1000Hz对数区间上均匀选取40个频率进行计算,如图4所示为正演得到的视电阻率和阻抗相位与频率组合剖面图。在合成数据反演计算中以此为野外观测数据。
反演过程如下:
1)根据测线观测系统生成反演网格模型;
2)由合成数据换算成野外观测阻抗数据。
对电场数据反演中,通过取近似磁场分量消元后,数据误差转变成计算张量阻抗归一化数据。
3)给定初始模型(100Ωm的均匀半空间模型)或反演更新模型,进行正演数值模拟计算,得到测点计算张量阻抗数据,然后计算数据拟合误差,当满足迭代要求,跳出循环,否则向下迭代。
4)计算反演迭代方程组右端相目标函数的梯度。
5)采用共轭梯度算法求解反演迭代方程组,获得模型修改量,更新电阻率模型,返回步骤3继续迭代。
如图5所示为反演迭代收敛曲线图,如图6所示为最终反演剖面图。如图7所示为初始模型为采用400Ωm的均匀半空间模型反演剖面图,如图8为初始模型为采用20Ωm的均匀半空间模型时反演剖面图。从图中可以看出,本发明提出的大地电磁TE极化模式对电场数据二维反演方法具有收敛稳定,收敛精度高的特点,反演结果分辨率高,静态效应影响非常小,起伏地形和浅层低阻体不影响深部构造的成像,此外,本方法对于初始模型的依赖性非常小。
应该理解的是,虽然图1和2的流程图中的各个步骤按照箭头的指示依次显示,但是这些步骤并不是必然按照箭头指示的顺序依次执行。除非本文中有明确的说明,这些步骤的执行并没有严格的顺序限制,这些步骤可以以其它的顺序执行。而且,图1和2中的至少一部分步骤可以包括多个子步骤或者多个阶段,这些子步骤或者阶段并不必然是在同一时刻执行完成,而是可以在不同的时刻执行,这些子步骤或者阶段的执行顺序也不必然是依次进行,而是可以与其它步骤或者其它步骤的子步骤或者阶段的至少一部分轮流或者交替地执行。
在一个实施例中,如图9所示,提供了一种大地电磁TE极化模式对电场数据二维反演装置,包括:目标函数构建模块902、迭代方程构建模块904、雅克比矩阵计算模块906和二维反演模块908,其中:
目标函数构建模块902,用于根据测线中多个观测点沿X轴方向的不同频率观测电场数据分量以及观测点沿X轴的不同频率计算电场数据分量,构建大地电磁TE极化模式下的反演目标函数;
迭代方程构建模块904,用于通过所述反演目标函数,得到所述反演目标函数对应的反演迭代方程组;
雅克比矩阵计算模块906,用于根据张量阻抗数据,求解所述反演迭代方程组对应的归一化电场数据分量误差,根据所述归一化电场数据分量误差的数学表达式,计算雅克比矩阵;
二维反演模块908,用于根据所述雅克比矩阵和所述归一化电场数据分量误差迭代求解所述反演迭代方程组,从而实现大地电磁TE极化模式对电场数据二维反演。
在其中一个实施例中,目标函数构建模块902还用于根据测线中多个观测点沿X轴方向的不同频率观测电场数据分量以及观测点沿X轴的不同频率计算电场数据分量,构建大地电磁TE极化模式下的反演目标函数为:
Figure BDA0002365780200000121
其中,ΦTE(m)表示反演目标函数,m=[m1,m2,…,mM]T表示二维地电模型向量,M表示反演节点的数量,mj表示第j个反演节点的电导率,上标T表示转置,
Figure BDA0002365780200000122
表示观测点沿X轴的观测电场数据分量,
Figure BDA0002365780200000123
表示观测点沿X轴的计算电场数据分量,*表示共轭算子,λ表示正则化因子,L表示二阶差分算子。
在其中一个实施例中,迭代方程构建模块904还用于将所述反演目标函数ΦTM(m)在二维地电模型向量m0进行二阶泰勒-拉格朗日展开,去除展开结果中的三阶及以上的高阶项,得到展开表达式;对所述展开表达式的两端对二维地电模型向量m计算一阶导数,通过求解
Figure BDA0002365780200000131
得到所述反演目标函数对应的反演迭代方程组为:
Figure BDA0002365780200000132
其中,G表示雅克比矩阵。
在其中一个实施例中,雅克比矩阵计算模块906还用于根据所述张量阻抗数据,得到相互正交的大地电场和磁场分量的关系为:
Ei=ZijHj i,j∈y,x
根据大地电场和磁场分量的关系式,得到所述反演迭代方程组对应的归一化电场数据分量误差为:
Figure BDA0002365780200000133
设置野外观测时所述观测点观测电场数据分量的初始值
Figure BDA0002365780200000134
对上述归一化电场数据分量误差进行简化,得到:
Figure BDA0002365780200000135
在其中一个实施例中,二维反演模块908还用于根据所述归一化电场数据分量误差的数学表达式,计算雅克比矩阵为:
Figure BDA0002365780200000136
其中,所述雅克比矩阵为N×M的矩阵。
在其中一个实施例中,二维反演模块908还用于根据所述雅克比矩阵,计算所述雅克比矩阵与向量的乘积;根据所述雅克比矩阵,构建海森矩阵,计算所述海森矩阵与向量的乘积;根据所述雅克比矩阵与向量的乘积、所述海森矩阵与向量的乘积以及所述归一化电场数据分量误差,计算所述反演迭代方程组。
在其中一个实施例中,二维反演模块908还用于根据RODI法求解角频率为ω,极化模式为TE时所述雅克比矩阵中的局部因子为:
Figure BDA0002365780200000137
其中,i=1,2,...,N表示观测点序列,j=1,...,M表示反演节点序列,
Figure BDA0002365780200000141
ip表示第i个观测点对应的正演网格节点编号,向量中第ip个元素为1,其余为零,K-1表示有限元正演生成的对称、稀疏复系数矩阵,E表示正演生成的节点电场向量,
Figure BDA0002365780200000142
表示L×L大小的稀疏矩阵;根据所述局部因子,计算所述雅克比矩阵与向量的乘积为:
Figure BDA0002365780200000143
Figure BDA0002365780200000144
其中,向量x=(x1,x2,...,xM)T,y=(y1,y2,...,yN)T
Figure BDA0002365780200000145
Figure BDA0002365780200000146
表示每次迭代中正演计算得到的测点磁场分量。
在其中一个实施例中,二维反演模块908还用于根据所述雅克比矩阵,构建海森矩阵为:
H≈GTG*
其中,H表示海森矩阵;根据所述雅克比矩阵与向量的乘积,计算得到所述海森矩阵与向量的乘积。
关于大地电磁TE极化模式对电场数据二维反演装置的具体限定可以参见上文中对于大地电磁TE极化模式对电场数据二维反演方法的限定,在此不再赘述。上述大地电磁TE极化模式对电场数据二维反演装置中的各个模块可全部或部分通过软件、硬件及其组合来实现。上述各模块可以硬件形式内嵌于或独立于计算机设备中的处理器中,也可以以软件形式存储于计算机设备中的存储器中,以便于处理器调用执行以上各个模块对应的操作。
在一个实施例中,提供了一种计算机设备,该计算机设备可以是终端,其内部结构图可以如图10所示。该计算机设备包括通过系统总线连接的处理器、存储器、网络接口、显示屏和输入装置。其中,该计算机设备的处理器用于提供计算和控制能力。该计算机设备的存储器包括非易失性存储介质、内存储器。该非易失性存储介质存储有操作系统和计算机程序。该内存储器为非易失性存储介质中的操作系统和计算机程序的运行提供环境。该计算机设备的网络接口用于与外部的终端通过网络连接通信。该计算机程序被处理器执行时以实现一种大地电磁TE极化模式对电场数据二维反演方法。该计算机设备的显示屏可以是液晶显示屏或者电子墨水显示屏,该计算机设备的输入装置可以是显示屏上覆盖的触摸层,也可以是计算机设备外壳上设置的按键、轨迹球或触控板,还可以是外接的键盘、触控板或鼠标等。
本领域技术人员可以理解,图10中示出的结构,仅仅是与本申请方案相关的部分结构的框图,并不构成对本申请方案所应用于其上的计算机设备的限定,具体的计算机设备可以包括比图中所示更多或更少的部件,或者组合某些部件,或者具有不同的部件布置。
在一个实施例中,提供了一种计算机设备,包括存储器和处理器,该存储器存储有计算机程序,该处理器执行计算机程序时实现上述实施例中方法的步骤。
在一个实施例中,提供了一种计算机可读存储介质,其上存储有计算机程序,计算机程序被处理器执行时实现上述实施例中方法的步骤。
本领域普通技术人员可以理解实现上述实施例方法中的全部或部分流程,是可以通过计算机程序来指令相关的硬件来完成,所述的计算机程序可存储于一非易失性计算机可读取存储介质中,该计算机程序在执行时,可包括如上述各方法的实施例的流程。其中,本申请所提供的各实施例中所使用的对存储器、存储、数据库或其它介质的任何引用,均可包括非易失性和/或易失性存储器。非易失性存储器可包括只读存储器(ROM)、可编程ROM(PROM)、电可编程ROM(EPROM)、电可擦除可编程ROM(EEPROM)或闪存。易失性存储器可包括随机存取存储器(RAM)或者外部高速缓冲存储器。作为说明而非局限,RAM以多种形式可得,诸如静态RAM(SRAM)、动态RAM(DRAM)、同步DRAM(SDRAM)、双数据率SDRAM(DDRSDRAM)、增强型SDRAM(ESDRAM)、同步链路(Synchlink)DRAM(SLDRAM)、存储器总线(Rambus)直接RAM(RDRAM)、直接存储器总线动态RAM(DRDRAM)、以及存储器总线动态RAM(RDRAM)等。
以上实施例的各技术特征可以进行任意的组合,为使描述简洁,未对上述实施例中的各个技术特征所有可能的组合都进行描述,然而,只要这些技术特征的组合不存在矛盾,都应当认为是本说明书记载的范围。
以上所述实施例仅表达了本申请的几种实施方式,其描述较为具体和详细,但并不能因此而理解为对发明专利范围的限制。应当指出的是,对于本领域的普通技术人员来说,在不脱离本申请构思的前提下,还可以做出若干变形和改进,这些都属于本申请的保护范围。因此,本申请专利的保护范围应以所附权利要求为准。

Claims (10)

1.一种大地电磁TE极化模式对电场数据二维反演方法,所述方法包括:
根据测线中测点观测到的沿X轴方向不同频率实测电场数据分量以及理论模型计算得到的不同频率计算电场数据分量,构建大地电磁TE极化模式下的最小二范数反演目标函数;
通过所述反演目标函数,得到所述反演目标函数对应的反演迭代方程组;
根据张量阻抗数据,求解所述反演迭代方程组对应的归一化电场数据分量误差,根据所述归一化电场数据分量误差的数学表达式,计算雅克比矩阵;
根据所述雅克比矩阵和所述归一化电场数据分量误差迭代求解所述反演迭代方程组,从而实现大地电磁TE极化模式对电场数据二维反演。
2.根据权利要求1所述的方法,其特征在于,所述根据测线中测点观测到的沿X轴方向不同频率实测电场数据分量以及理论模型计算得到的不同频率计算电场数据分量,构建大地电磁TE极化模式下的最小二范数反演目标函数,包括:
根据测线中测点观测到的沿X轴方向不同频率实测电场数据分量以及理论模型计算得到的不同频率计算电场数据分量,构建大地电磁TE极化模式下的最小二范数反演目标函数:
Figure FDA0002365780190000011
其中,ΦTE(m)表示反演目标函数,m=[m1,m2,…,mM]T表示二维地电模型向量,M表示反演节点的数量,mj表示第j个反演节点的电导率,上标T表示转置,
Figure FDA0002365780190000012
表示观测点沿X轴的观测电场数据分量,
Figure FDA0002365780190000013
表示观测点沿X轴的计算电场数据分量,*表示共轭算子,λ表示正则化因子,L表示二阶差分算子。
3.根据权利要求2所述的方法,其特征在于,通过所述反演目标函数,得到所述反演目标函数对应的反演迭代方程组,包括:
将所述反演目标函数ΦTE(m)在二维地电模型向量m0进行二阶泰勒-拉格朗日展开,去除展开结果中的三阶及以上的高阶项,得到展开表达式;
对所述展开表达式的两端对二维地电模型向量m计算一阶导数,通过求解
Figure FDA0002365780190000014
得到所述反演目标函数对应的反演迭代方程组为:
Figure FDA0002365780190000021
其中,G表示雅克比矩阵。
4.根据权利要求3所述的方法,其特征在于,所述根据张量阻抗数据,求解所述反演迭代方程组对应的归一化电场数据分量误差,包括:
根据所述张量阻抗数据,得到相互正交的大地电场和磁场分量的关系为:
Ei=ZijHj i,j∈y,x
根据大地电场和磁场分量的关系式,得到所述反演迭代方程组对应的归一化电场数据分量误差为:
Figure FDA0002365780190000022
设置野外观测时所述观测点观沿Y轴磁场分量的初始值
Figure FDA0002365780190000023
对上述归一化电场数据分量误差进行简化,得到:
Figure FDA0002365780190000024
5.根据权利要求4所述的方法,其特征在于,根据所述归一化电场数据分量误差的数学表达式,计算雅克比矩阵,包括:
根据所述归一化电场数据分量误差的数学表达式,计算雅克比矩阵为:
Figure FDA0002365780190000025
其中,所述雅克比矩阵为N×M的矩阵。
6.根据权利要求5所述的方法,其特征在于,根据所述雅克比矩阵和所述归一化电场数据分量误差计算所述反演迭代方程组,包括:
根据所述雅克比矩阵,计算所述雅克比矩阵与向量的乘积;
根据所述雅克比矩阵,构建海森矩阵,计算所述海森矩阵与向量的乘积;
根据所述雅克比矩阵与向量的乘积、所述海森矩阵与向量的乘积以及所述归一化电场数据分量误差,计算所述反演迭代方程组。
7.根据权利要求6所述的方法,其特征在于,根据所述雅克比矩阵,计算所述雅克比矩阵与向量的乘积,包括:
根据RODI法求解角频率为ω,极化模式为TE时所述雅克比矩阵中的局部因子为:
Figure FDA0002365780190000031
其中,i=1,2,...,N表示观测点序列,j=1,...,M表示反演节点序列,
Figure FDA0002365780190000032
ip表示第i个观测点对应的正演网格节点编号,向量中第ip个元素为1,其余为零,K-1表示有限元正演生成的对称、稀疏复系数矩阵,E表示正演生成的节点电场向量,
Figure FDA0002365780190000033
表示L×L大小的稀疏矩阵;
根据所述局部因子,计算所述雅克比矩阵与向量的乘积为:
Figure FDA0002365780190000034
Figure FDA0002365780190000035
其中,向量x=(x1,x2,...,xM)T,y=(y1,y2,...,yN)T
Figure FDA0002365780190000036
Figure FDA0002365780190000037
表示每次迭代中正演计算得到的测点磁场分量。
8.根据权利要求7所述的方法,其特征在于,根据所述雅克比矩阵,构建海森矩阵,计算所述海森矩阵与向量的乘积,包括:
根据所述雅克比矩阵,构建海森矩阵为:
H≈GTG*
其中,H表示海森矩阵;
根据所述雅克比矩阵与向量的乘积,计算得到所述海森矩阵与向量的乘积。
9.一种大地电磁TE极化模式对电场数据二维反演装置,其特征在于,所述装置包括:
目标函数构建模块,用于根据测线中多个观测点沿X轴方向的不同频率观测电场数据分量以及观测点沿X轴的不同频率计算电场数据分量,构建大地电磁TE极化模式下的反演目标函数;
迭代方程构建模块,用于通过所述反演目标函数,得到所述反演目标函数对应的反演迭代方程组;
雅克比矩阵计算模块,用于根据张量阻抗数据,求解所述反演迭代方程组对应的归一化电场数据分量误差,根据所述归一化电场数据分量误差的数学表达式,计算雅克比矩阵;
二维反演模块,用于根据所述雅克比矩阵和所述归一化电场数据分量误差迭代求解所述反演迭代方程组,从而实现大地电磁TE极化模式对电场数据二维反演。
10.一种计算机设备,包括存储器和处理器,所述存储器存储有计算机程序,其特征在于,所述处理器执行所述计算机程序时实现权利要求1至8中任一项所述方法的步骤。
CN202010035274.XA 2020-01-14 2020-01-14 大地电磁te极化模式对电场数据二维反演方法和装置 Active CN111103628B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202010035274.XA CN111103628B (zh) 2020-01-14 2020-01-14 大地电磁te极化模式对电场数据二维反演方法和装置

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202010035274.XA CN111103628B (zh) 2020-01-14 2020-01-14 大地电磁te极化模式对电场数据二维反演方法和装置

Publications (2)

Publication Number Publication Date
CN111103628A true CN111103628A (zh) 2020-05-05
CN111103628B CN111103628B (zh) 2022-03-11

Family

ID=70426242

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202010035274.XA Active CN111103628B (zh) 2020-01-14 2020-01-14 大地电磁te极化模式对电场数据二维反演方法和装置

Country Status (1)

Country Link
CN (1) CN111103628B (zh)

Cited By (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111856596A (zh) * 2020-08-05 2020-10-30 中国海洋大学 层状介质电阻率各向异性海洋可控源电磁快速反演方法
CN111856597A (zh) * 2020-08-05 2020-10-30 中国海洋大学 拖曳式海洋电磁地层电阻率与接收站位置联合反演方法
CN113484920A (zh) * 2021-08-17 2021-10-08 成都理工大学 一种频域电磁测深资料二维结构化反演方法

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20120123683A1 (en) * 2009-01-20 2012-05-17 Statoil Asa Csem survey method
CN102798898A (zh) * 2012-08-20 2012-11-28 中国地质科学院矿产资源研究所 大地电磁场非线性共轭梯度三维反演方法
CN106643744A (zh) * 2016-12-29 2017-05-10 武汉大学 一种基于四程中继跟踪模式的远月面着陆器精密定位方法
CN110187398A (zh) * 2019-07-11 2019-08-30 中南大学 一种寻找井间目标体的多电极系探测方法

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20120123683A1 (en) * 2009-01-20 2012-05-17 Statoil Asa Csem survey method
CN102798898A (zh) * 2012-08-20 2012-11-28 中国地质科学院矿产资源研究所 大地电磁场非线性共轭梯度三维反演方法
CN106643744A (zh) * 2016-12-29 2017-05-10 武汉大学 一种基于四程中继跟踪模式的远月面着陆器精密定位方法
CN110187398A (zh) * 2019-07-11 2019-08-30 中南大学 一种寻找井间目标体的多电极系探测方法

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
TIFFANY TJONG 等: "Two Dimensional Finite Element Based Magnetotelluric Inversion using Singular Value Decomposition Method on Transverse Electric Mode", 《THE INTERNATIONAL CONFERENCE ON THEORETICAL AND APPLIED PHYSICS》 *
康敏 等: "大地电磁二维反演方法分析对比", 《地球物理学进展》 *

Cited By (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111856596A (zh) * 2020-08-05 2020-10-30 中国海洋大学 层状介质电阻率各向异性海洋可控源电磁快速反演方法
CN111856597A (zh) * 2020-08-05 2020-10-30 中国海洋大学 拖曳式海洋电磁地层电阻率与接收站位置联合反演方法
CN111856597B (zh) * 2020-08-05 2023-03-21 中国海洋大学 拖曳式海洋电磁地层电阻率与接收站位置联合反演方法
CN113484920A (zh) * 2021-08-17 2021-10-08 成都理工大学 一种频域电磁测深资料二维结构化反演方法
CN113484920B (zh) * 2021-08-17 2023-05-19 成都理工大学 一种频域电磁测深资料二维结构化反演方法

Also Published As

Publication number Publication date
CN111103628B (zh) 2022-03-11

Similar Documents

Publication Publication Date Title
CN111103627B (zh) 大地电磁tm极化模式对电场数据二维反演方法和装置
CN111103628B (zh) 大地电磁te极化模式对电场数据二维反演方法和装置
Rücker et al. pyGIMLi: An open-source library for modelling and inversion in geophysics
CN112287534B (zh) 基于nufft的二维磁异常快速正演模拟方法和装置
Ford et al. A nonhydrostatic finite-element model for three-dimensional stratified oceanic flows. Part I: Model formulation
Gavete et al. Generalized finite differences for solving 3D elliptic and parabolic equations
CN110346834B (zh) 三维频率域可控源电磁的正演方法、系统
CN113051779B (zh) 一种三维直流电阻率法数值模拟方法
CN111638556B (zh) 基于地空分解策略的大地电磁正演方法及装置、存储介质
Bera et al. A MATLAB-based boundary data simulator for studying the resistivity reconstruction using neighbouring current pattern
Alvarez-Aramberri et al. Dimensionally adaptive hp-finite element simulation and inversion of 2D magnetotelluric measurements
CN114065511A (zh) 起伏地形下大地电磁二维正演数值模拟方法、装置、设备及介质
CN116842813A (zh) 一种三维大地电磁正演数值模拟方法
CN114970289B (zh) 三维大地电磁各向异性正演数值模拟方法、设备及介质
Egbert Hybrid conjugate gradient-Occam algorithms for inversion of multifrequency and multitransmitter EM data
Dumont et al. The best of two worlds: The expedite boundary element method
Boonchaisuk et al. Two-dimensional direct current (DC) resistivity inversion: Data space Occam's approach
CN113076678B (zh) 一种频率域二度体重力异常快速数值模拟方法和装置
Chen et al. New post-processing method for interpretation of through casing resistivity (TCR) measurements
CN113779818B (zh) 三维地质体其电磁场数值模拟方法、装置、设备及介质
CN115755199A (zh) 一种实用的非结构网格三维电磁反演平滑正则化方法
Wubs et al. The performance of implicit ocean models on B-and C-grids
Chen et al. Method for imposing boundary conditions on Reissner–Mindlin plates for analysis using structured background mesh
Guo et al. Microwave inversion for sparse data using descent learning technique
Klees et al. Calculation of strongly singular and hypersingular surface integrals

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
EE01 Entry into force of recordation of patent licensing contract

Application publication date: 20200505

Assignee: Guangxi Aigui Intelligent Technology Co.,Ltd.

Assignor: GUILIN University OF TECHNOLOGY

Contract record no.: X2022450000132

Denomination of invention: Two dimensional inversion method and device of magnetotelluric TE polarization model for electric field data

Granted publication date: 20220311

License type: Common License

Record date: 20221123

EE01 Entry into force of recordation of patent licensing contract