起伏地形下大地电磁二维正演数值模拟方法、装置、设备及介质制造方法及图纸

起伏地形下大地电磁二维正演数值模拟方法、装置、设备及介质，包括：构建内部包含勘探目标的二维矩形模型；对二维矩形模型进行网格剖分并进行电导率赋值，得到刻画任意电导率分布的二维异常体模型；根据频率参数和二维异常体模型，构建关于观测点上待求解的场的积分表达式；根据积分表达式，加载一维边界条件以及二维大地电磁场源，得到场的控制方程；将二维异常体模型沿Z方向划分为多层子区域，在各层子区域分别沿水平方向构建不同稀疏度的插值算子，所有插值算子合成后得到一个总的稀疏采样算子，基于场的控制方程求解得到修正后的场，进而计算TE模式和TM模式下的阻抗、视电阻率和相位。该方法能够满足大规模电磁数据精细、快速正演成像的需求。快速正演成像的需求。快速正演成像的需求。

【技术实现步骤摘要】
起伏地形下大地电磁二维正演数值模拟方法、装置、设备及介质


[0001]本专利技术属于强磁体数值模拟
，特别涉及一种起伏地形下大地电磁二维正演数值模拟方法、装置、设备及介质。

技术介绍

[0002]大地电磁测深法是通过天然交变电磁场研究地下电性结构的一种勘查方法，具有穿透深度大和分辨能力高等特点，被广泛应用于矿产与地热资源勘查、油气勘探和环境监测等领域。
[0003]大地电磁测深法不仅可进行平坦地形的探测，也可观测山区丘陵地带。在实际生产中，走向方向无限延伸的地质体可看作二维构造，以便研究起伏地形畸变带来的影响。起伏地形对二维数据影响很大，如果不考虑起伏地形引起的大地电磁场畸变，将使大地电磁资料产生显著误差，其解释结果必然偏离实际构造。因此，研究二维起伏地形对山区大地电磁勘探工作的影响具有现实意义。
[0004]对于起伏地形的模拟，目前本领域研究人员的工作多集中于有限单元法、有限差分法和有限积分法。其中有限单元法采用非结构网格能够更好的模拟复杂地形，其精度较高，然而其模拟需要计算和存储的未知量远高于其他两类方法。有限差分法计算效率较高，但是其计算精度相对其他两类方法较差。有限积分法作为上述两种方法的平衡，其基本思想是将求解区域剖分成一系列不交叉的控制积分单元，使得各个积分单元满足守恒原理，然后对控制方程在积分单元内进行积分离散。
[0005]然而，无论选用何种方法，起伏地形的精细离散都需要大量的网格单元，必然使得计算所需的时间和内存急剧增加，此时计算效率和计算精度难以维持平衡。
专利
技术实现思路

[0006]针对目前起伏地形下大地电磁二维正数值模拟方法大多存在计算效率和计算精度难以平衡的问题，本专利技术旨在提出一种起伏地形下大地电磁二维正演数值模拟方法、装置、设备及介质，以满足起伏地形精细化正演计算的需求。
[0007]为实现上述技术目的，本专利技术提出的技术方案为：
[0008]一方面，本专利技术提供一种起伏地形下大地电磁二维正演数值模拟方法，包括：
[0009]勘探目标为二维异常体，构建内部包含所述勘探目标的二维矩形模型；
[0010]对所述二维矩形模型沿x、z方向进行网格剖分，剖分成若干个小矩形单元，得到二维长方体模型的网格剖分参数，根据二维异常体的电导率分布，给每个小矩形单元的电导率赋值，每一个小矩形单元的电导率为常值，不同矩形单元的电导率值不同，得到刻画任意电导率分布的二维异常体模型；
[0011]根据频率参数和二维异常体模型，构建关于观测点上待求解的场的积分表达式
其中F代表观测点上待求解的场，表示梯度算子，表示散度算子，ω表示角频率，通过ω＝2πf求取，f表示给定的频率，μ表示磁导率，其值为4π
×
10
‑7，σ表示相邻两个小矩形单元的加权平均电导率，Ω表示小矩形单元的面积；
[0012]基于一维边界条件以及二维大地电磁场源，得到场的控制方程CF＝b，其中C表示系数矩阵，表示为b由二维大地电磁场源和一维边界条件构成的向量，为已知量；
[0013]将所述二维异常体模型沿Z方向划分为多层子区域，在各层子区域分别沿水平方向构建不同稀疏度的插值算子；
[0014]将各层子区域对应的插值算子合成后得到一个总的稀疏采样算子S，系数矩阵C乘上总的稀疏采样算子得到新的系数矩阵C
new
，由C
new
F
new
＝b求解得到修正后的场F
new
。
[0015]二维矩形模型的网格剖分参数包括二维矩形模型x、z方向网格剖分得到的小矩形单元的数量N
x
、N
z
，各个小矩形x、z方向的棱边长度，以及对各小矩形单元进行编号后得到的各小矩形单元的编号和坐标位置。
[0016]另一方面，本专利技术提供一种起伏地形下大地电磁二维正演数值模拟装置，包括：
[0017]第一模块，用于构建内部包含勘探目标的二维矩形模型，所述勘探目标为二维异常体；
[0018]第二模块，用于对所述二维矩形模型沿x、z方向进行网格剖分，剖分成若干个小矩形单元，得到二维长方体模型的网格剖分参数，根据二维异常体的电导率分布，给每个小矩形单元的电导率赋值，每一个小矩形单元的电导率为常值，不同矩形单元的电导率值不同，得到刻画任意电导率分布的二维异常体模型；
[0019]第三模块，用于根据频率参数和二维异常体模型，构建关于观测点上待求解的场的积分表达式其中F代表观测点上待求解的场，表示梯度算子，表示散度算子，ω表示角频率，通过ω＝2πf求取，f表示给定的频率，μ表示磁导率，其值为4π
×
10
‑7，σ表示相邻两个小矩形单元的加权平均电导率，Ω表示小矩形单元的面积；
[0020]第四模块，用于基于一维边界条件以及二维大地电磁场源，得到场的控制方程CF＝b，其中C表示系数矩阵，表示为b由二维大地电磁场源和一维边界条件构成的向量，为已知量；
[0021]第五模块，用于将所述二维异常体模型沿Z方向划分为多层子区域，在各层子区域分别沿水平方向构建不同稀疏度的插值算子；
[0022]第六模块，将各层子区域对应的插值算子合成后得到一个总的稀疏采样算子S，系数矩阵C乘上总的稀疏采样算子得到新的系数矩阵C
new
，由C
new
F
new
＝b求解得到修正后的场F
new
。
[0023]进一步地，还包括第七模块，用于根据修正后的场F
new
，计算阻抗、视电阻率和相
位。
[0024]另一方面，本专利技术提供一种计算机设备，包括存储器和处理器，存储器存储有计算机程序，处理器执行计算机程序时实现所述起伏地形下大地电磁二维正演数值模拟方法中的步骤。
[0025]再一方面，本专利技术还提供一种计算机可读存储介质，其上存储有计算机程序，计算机程序被处理器执行时实现所述起伏地形下大地电磁二维正演数值模拟方法中的步骤。
[0026]与现有技术相比，本专利技术的优点在于：
[0027]本专利技术可以根据实际需要，将整个研究区域灵活的分为多个子区域，在每个子区域水平方向施加不同“稀疏度”的插值算子，z方向计算点保持不变，通过分区域稀疏采样直接减少了离散电磁场的数量。因此，本专利技术大大提高了电磁正演的计算效率、减少了计算所需的内存；另外该方法不改变原有模型电导率单元的数量，转变为降低离散计算点数量，与已有主流的方法相比，这更容易用于二维反演算法。
附图说明
[0028]图1是本专利技术一实施例中流程图；
[0029]图2是本专利技术一实施例中复杂山峰地形模型示意图；
[0030]图3为本专利技术一实施例得到的TE模式下计算的视电阻率、相位的参考解、数值解以及相对误差等值线图，其中(a)为视电阻率的参考解，(b)为视电阻率的数值解，(c)为视电阻率的相对误差；(d)为相位的参考解，(e为相位的数值解，本文档来自技高网...

【技术保护点】

【技术特征摘要】
1.起伏地形下大地电磁二维正演数值模拟方法，其特征在于，包括：勘探目标为二维异常体，构建内部包含所述勘探目标的二维矩形模型；对所述二维矩形模型沿x、z方向进行网格剖分，剖分成若干个小矩形单元，得到二维长方体模型的网格剖分参数，根据二维异常体的电导率分布，给每个小矩形单元的电导率赋值，每一个小矩形单元的电导率为常值，不同矩形单元的电导率值不同，得到刻画任意电导率分布的二维异常体模型；根据频率参数和二维异常体模型，构建关于观测点上待求解的场的积分表达式其中F代表观测点上待求解的场，表示梯度算子，表示散度算子，ω表示角频率，通过ω＝2πf求取，f表示给定的频率，μ表示磁导率，其值为4π
×
10
‑7，σ表示相邻两个小矩形单元的加权平均电导率，Ω表示小矩形单元的面积；基于一维边界条件以及二维大地电磁场源，得到场的控制方程CF＝b，其中C表示系数矩阵，表示为b由二维大地电磁场源和一维边界条件构成的向量，为已知量；将所述二维异常体模型沿Z方向划分为多层子区域，在各层子区域分别沿水平方向构建不同稀疏度的插值算子；将各层子区域对应的插值算子合成后得到一个总的稀疏采样算子S，控制方程的系数矩阵C乘上总的稀疏采样算子得到新的系数矩阵C
new
，由C
new
F
new
＝b求解得到修正后的场F
new
。2.根据权利要求2所述的起伏地形下大地电磁二维正演数值模拟方法，其特征在于，根据修正后的场F
new
，计算阻抗、视电阻率和相位：TE模式时，修正后的场F
new
代表电场，TE模式下的阻抗Z
TE
、视电阻率ρ
TE
和相位φ
TE
为：TE模式TM模式时，修正后的场F
new
代表磁场，TM模式下的阻抗Z
TM
、视电阻率ρ
TM
和相位φ
TM
为：
TM模式其中：表示z方向长度,Im和Re表示虚部和实部，arctan表示反三角函数。3.根据权利要求1或2所述的起伏地形下大地电磁二维正演数值模拟方法，其特征在于，系数矩阵C通过以下方法获得：二维异常体模型中各个小矩形单元x、z方向的棱边长度，构成长度元矩阵L；二维异常体模型中各个小矩形单元的表面积，构成面积元矩阵S；系数矩阵C可表示为其中G为梯度算子，实质为二阶导数，G
T
为梯度算子转置，表示散度算子，S
‑1表示面积元矩阵S的倒数，L
‑1表示长度元矩阵L的倒数，表示相邻长度元矩阵L的平均。4.根据权利要求1或2所述的起伏地形下大地电磁二维正演数值模拟方法，其特征在于，将所述二维异常体模型沿Z方向划分为五层子区域，由上至下分别为高空层子区域、空气
‑
地表层子区域、浅地层子区域、中地层子区域和深地层子区域；对于高空层子区域和浅地层子区域，其插值算子为1阶稀疏采样算子，表示为S1＝s1，s1为一阶插值矩阵；对于空气
‑
地表层子区域其插值算子为单位矩阵I，其矩阵大小为(n,n)，主对角线的元素为1；对于中地层子区域其插值算子为2阶稀疏采样算子，表示为S2＝s1s2，s2为二阶插值矩阵；对于深地层子区域其插值算子为3阶稀疏采样算子，表示为S3＝s1s2s3，s3为三阶插值矩阵。5.根据权利要求4所述的起伏地形下大地电磁二维正演数值模拟方法，其特征在于，构建一阶插值矩阵s1，方法是：先构建一个一阶初始矩阵，该矩阵大小为(n,n)，主对角线的元素为1；对子区域中的各小矩形单元，构建x方向一阶插值算子系数为：其中分别表示子区域中x方向上第2j
‑
1、2j+1个小矩形单元上计算点对应的插值算子系数，j用于计数，j为正整数，j＝1,2,3...；而Δx
2j
‑1,Δx
2j
,Δx
2j+1
分别代表第
2j、2j
‑
1、2j+1个小矩形单元沿x方向的棱边长度；将所构建的一阶插值算子系数分别施加到对应子区域中的相应小矩形单元上，有将所构建的一阶插值算子系数分别施加到对应子区域中的相应小矩形单元上，有和表示沿x方向上第2j
‑
1、2j+1个小矩形单元上计算点对应的场，该计算点仍然会保留，记录下其位置(2j
‑
1，2j)和(2j+1，2j)，然后在一阶初始矩阵里面找到对应的位置并分别赋值插值系数位置并分别赋值插值系数表示x方向第2j个小矩形单元上计算点被插值的场，该计算点将会被去除，记录下其位置(2j，2j)，然后在一阶初始矩阵里面找到对应的位置并赋值0，以此规律，直至一阶初始矩阵中所有位置均赋值完毕，得到最终的一阶插值矩阵s1。6.根据权利要求5所述的起伏地形下大地电磁二维正演数值模拟方法，其特征在于，构建二阶插值矩阵s2，方法是：先构建一个二阶初始矩阵，其矩阵大小为(n,n)，主对角线元素为1；对中地层子区域中的各小矩形单元，构建x方向二阶插值系数为：其中分别表示中地层子区域中x方向上第4j
‑
3、4j+1个小矩形单元上计算点对应的插值系数；Δx
4j
‑3,Δx
4...

【专利技术属性】
技术研发人员：李健，柳建新，郭荣文，王永斐，王旭龙，袁洋，
申请(专利权)人：中南大学，
类型：发明
国别省市：