CN112558164A - 一种基于偏差原理的大地电磁正则化反演方法及处理终端 - Google Patents
一种基于偏差原理的大地电磁正则化反演方法及处理终端 Download PDFInfo
- Publication number
- CN112558164A CN112558164A CN202011424836.6A CN202011424836A CN112558164A CN 112558164 A CN112558164 A CN 112558164A CN 202011424836 A CN202011424836 A CN 202011424836A CN 112558164 A CN112558164 A CN 112558164A
- Authority
- CN
- China
- Prior art keywords
- inversion
- magnetotelluric
- regularization
- principle
- observation data
- 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.)
- Pending
Links
Images
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V3/00—Electric or magnetic prospecting or detecting; Measuring magnetic field characteristics of the earth, e.g. declination, deviation
- G01V3/38—Processing data, e.g. for analysis, for interpretation, for correction
-
- 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
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Mathematical Physics (AREA)
- Data Mining & Analysis (AREA)
- Mathematical Optimization (AREA)
- Mathematical Analysis (AREA)
- Life Sciences & Earth Sciences (AREA)
- Pure & Applied Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Computational Mathematics (AREA)
- General Engineering & Computer Science (AREA)
- Software Systems (AREA)
- Databases & Information Systems (AREA)
- Algebra (AREA)
- Operations Research (AREA)
- Environmental & Geological Engineering (AREA)
- Geology (AREA)
- Remote Sensing (AREA)
- General Life Sciences & Earth Sciences (AREA)
- Geophysics (AREA)
- Measuring Magnetic Variables (AREA)
Abstract
本发明涉及一种基于偏差原理的大地电磁正则化反演方法及处理终端,所述方法包括,步骤1:获得大地电磁观测数据,并构建反演算子方程;步骤2:利用偏差原理求解出当前的反演算子方程中的正则化因子,得到新的反演算子方程;步骤3:求解出大地电磁观测数据中当前的灵敏度,利用所述当前的灵敏度和新的反演算子方程对大地电磁观测数据不断进行反演,最后一次的反演结果作为最终的反演结果。本发明基于偏差原理来计算出正则化因子,进而得到反演结果,使得不依赖于经验值,反演效果更好,并且通过收敛速度更快的偏差原理高阶收敛算法计算得到灵敏度,使得收敛速度得到快速提高,也即提高了反演速度。
Description
技术领域
本发明涉及地球科学中的地球物理探测技术领域,具体涉及一种基于偏差原理的大地电磁正则化反演方法及处理终端。
背景技术
大地电磁反演问题属于不适定问题,通常引入正则化反演方法来进行求解,通过加入先验的约束条件来增强反演过程的稳定性质,减少反演结果的非唯一性。常用的正则化反演算法有Occam反演、Rebocc反演、NLCG反演等,Occam反演能引入先验信息来去掉对已知构造边界处的光滑度或者加上对同类电性单元的电阻率差的限制,具有其他反演算法不具有的优点。但目前基于普通Occam反演算法进行大地电磁反演,仍然存在以下不足,主要包括仍然需要依赖反演这的经验和未考虑收敛速度的问题。例如,吴小平等认为,可以在主要拟合观测数据的同时得到一个光滑而不是最光滑的地电模型,提出了正则化因子递减方案,具体做法是给一个初始正则化因子和步长进行每次迭代,虽避免了线性搜索正则化因子,减少了冗余计算,但初始正则化因子和步长的选取仍然存有一定的随意性。若想得到一个满意结果,仍然需要依赖反演者的经验,未能摆脱经验选择方式。张罗磊等提出,初始正则化因子每次迭代随模型粗糙度的变化,在总体粗糙度没用减少时,以粗糙度变化率减少正则化因子。这种方法运算速度较快,但容易陷入局部极小的缺点。朴英哲等在搜索正则化因子的过程中加上了观测数据拟合差与其期望值的比较,能够排除多余的正演计算,使得反演速度加快,但未考虑收敛速度问题。
相关的参考文献如下:
吴小平,徐果明.大地电磁数据的Occam反演改进[J].地球物理学报,1998,41(4):547-554;
张罗磊,于鹏,王家林,等.光滑模型与尖锐边界结合的MT二维反演方法[J].地球物理学报,2009,52(6):1625-1632;
朴英哲,李桐林,刘永亮.在大地电磁二维Occam反演中求取拉格朗日乘子方法改进[J].吉林大学学报:地球科学版,2014,44(2):660-667。
发明内容
针对现有技术的不足,本发明的目的之一提供一种基于偏差原理的大地电磁正则化反演方法,其能够解决大地电磁反演的问题;
本发明的目的之二提供一种处理终端,其能够解决大地电磁反演的问题。
实现本发明的目的之一的技术方案为:一种基于偏差原理的大地电磁正则化反演方法,包括如下步骤:
步骤1:获得大地电磁观测数据,并构建反演算子方程;
步骤2:利用偏差原理求解出当前的反演算子方程中的正则化因子,得到新的反演算子方程;
步骤3:求解出大地电磁观测数据中当前的灵敏度,
步骤4:利用所述当前的灵敏度和新的反演算子方程对大地电磁观测数据进行反演,得到反演结果;
步骤5:重复步骤2-4,直至最后一次反演,最后一次反演的结果作为最终的反演结果。
进一步地,所述正则化因子按公式①得到:
其中,αk表示第k次反演的正则化因子。
进一步地,所述步骤2中,采用Morozov偏差原理求解出正则化因子。
实现本发明的目的之二的技术方案为:一种处理终端,其包括:
存储器,用于存储程序指令;
处理器,用于运行所述程序指令,以执行所述基于偏差原理的大地电磁正则化反演方法的步骤。
本发明的有益效果为:本发明基于偏差原理来计算出正则化因子,进而得到反演结果,使得不依赖于经验值,反演效果更好,并且通过收敛速度更快的偏差原理高阶收敛算法计算得到灵敏度,使得收敛速度得到快速提高,也即提高了反演速度。
附图说明
图1为较佳实施例的流程示意图;
图2是在低阻模型的示意图;
图3是在低阻模型下按现有的TM模型进行反演得到正则化因子的示意图;
图4是在低阻模型下按本发明进行反演得到正则化因子的示意图;
图5是在低阻模型下按TM模型和本发明进行反演得到均方根误差曲线图;
图6是在低阻模型的示意图;
图7是在低阻模型下按现有的TE模型进行反演得到正则化因子的示意图;
图8是在低阻模型下按本发明进行反演得到正则化因子的示意图;
图9是在低阻模型下按TE模型和本发明进行反演得到均方根误差曲线图;
图10是在高低阻模型的示意图;
图11是在高低阻模型下按现有的TM模型进行反演得到正则化因子的示意图;
图12是在高低阻模型下按本发明进行反演得到正则化因子的示意图;
图13是在高低阻模型下按TM模型和本发明进行反演得到均方根误差曲线图;
图14是在高低阻模型的示意图;
图15是在高低阻模型下按现有的TE模型进行反演得到正则化因子的示意图;
图16是在高低阻模型下按本发明进行反演得到正则化因子的示意图;
图17是在高低阻模型下按TE模型和本发明进行反演得到均方根误差曲线图;
图18为处理终端的示意图。
具体实施方式
下面,结合附图以及具体实施方案,对本发明做进一步描述。
如图1-图17所示,一种基于偏差原理的大地电磁正则化反演方法,包括如下步骤:
步骤1:获得大地电磁观测数据,并构建反演算子方程,以对大地电磁观测数据进行反演。
大地电磁反演算子方程如公式①:
F(m)=d------①
式中,d表示大地电测观测数据,为向量形式,向量中的值为一系列大地电磁实际观测得到的数据。m表示电阻率模型参数,为向量形式,F表示大地电磁观测数据正演响应算子。通常,野外测量得到的大地电磁观测数据带有一定的误差,即有dδ=d+δd,dδ表示带有误差后的大地电磁观测数据,δ表示误差系数,通常为常数。因此,公式①可以变成公式②:
dδ=Fα(m)------②
为了求解公式②的mα,现有方法中通常会通过正则化进行求解,其中Occam反演是一种较常用的正则化反演方法,可以用于求解公式②。在Occam反演中,Occam反演的展平泛函φ(m)如公式③:
因此,只需要求解出公式③中的正则化因子α,即可求解公式①,从而得到大地电磁的反演结果。目前传统的方法求解公式③中的正则化因子α,大多数是凭借经验直接预设一个固定值。这导致存在一个明显的缺陷和不足:需要依赖反演者的经验,未能摆脱经验选择方式。一旦反演者的经验不足,导致反演结果出现很大的偏差,反演结果质量差。
式中,dδ表示大地电磁观测数据,zα表示正则化解,A表示大地电磁正演算子。公式④可通过具有二阶收敛速度的Newton求解,即对公式④进行一阶泰勒展开,得到公式⑤:
令l(α)=0,便得到Newton迭代公式⑥:
式中,ξk表示介于α和αk之间的常数。公式⑦也即是二阶公式,联合公式⑥和⑦得到三阶收敛的迭代公式,也即得到公式⑧:
其中Ω′是稳定泛函,起着恢复解的稳定性的作用;F是公式②中的正演算子,F*是其伴随算子。方程组中的三个方程的系数矩阵是一样的,因此,只需要对系数矩阵做一次分解和回代后,便可以得到 和进而计算得到和从而可以求解出正则化因子α。
当然,以上只是本实施例给出求解基于偏差原理下的正则化因子,具体求解过程也可以通过其他数学手段进行求解得到,本质上属于可通过现有数学方法进行求解得到,因此,也不进行赘述。
所述αk表示第k次反演的正则化因子,αk+1表示第k+1次反演的正则化因子,也即每一次反演对应有一个正则化因子。因此,通过公式⑧,只需要给出初始值(也即给出α1的值),即可求解出任意一次反演的正则化因子,得到每次反演均更新得到的正则化因子,而不依赖于人为给出的经验值。
步骤3:在大地电磁观测数据反演中,灵敏度是一个重要的参数,灵敏度的计算最终可归结为电磁场对电导率的求导,在反演中,灵敏度有多个,并构成一个灵敏度矩阵,因此,灵敏度的计算也即是灵敏度矩阵的计算过程。通过基于偏差原理的高阶收敛算法,可以提高大地电磁观测数据反演过程的收敛速度,最终得到反演结果。灵敏度具体计算过程属于现有技术,在此不进行具体求解过程的介绍。
本实施例中,在每一次将正则化因子代入大地电磁观测数据构建的反演算子方程中进行反演过程中,均重新计算灵敏度,然后将灵敏度和对应的正则化因子代入到反演算子方程中进行反演,直至最后一次反演(也即最后一个k值),并将最后一次反演的结果作为最终的反演结果,也即得到大地的电性结构,至于反演次数可以通过人为预设,例如预设反演30次,也即αk中的k的最大值可以预设得到。也即,每一次反演,均会计算一次正则化因子和灵敏度并代入反演中,得到反演结果。
本发明基于偏差原理来计算出正则化因子,进而得到反演结果,使得不依赖于经验值,反演效果更好,并且通过收敛速度更快的偏差原理高阶收敛算法计算得到灵敏度,使得收敛速度得到快速提高,也即提高了反演速度。
将背景技术中由吴小平提出的“大地电磁数据的Occam反演改进”进行计算得到的正则化因子记为按方案一得到的正则化因子,按本发明基于偏差原理进行高阶收敛算法计算得到的正则化因子记为按方案二得到的正则化因子。
图2-图9是在低阻模型下分别按现有的TE和TM模型进行反演得到的反演结果图。图2和图6部分是低阻模型示意图,图3和图7部分是按方案一得到的正则化因子,图4和图8部分是按方案二得到的正则化因子,图5和图9部分是均方根误差曲线图。其中,研究区域中为电阻率100Ω·m的均匀半空间,存在电阻率为10Ω·m的低阻异常体,长宽为2.4km x1km。沿测线以300m的间距布置30个测点,正演频率以对数等间隔在0.01Hz-100Hz的频段中取25个频率,反演的初始模型取均匀半空间。将两种正则化因子的选取方案应用到TE、TM两种反演模式中。其中方案一的正则化因子的初始值为50,步长为2。
从图中可以看出,TM和TE反演结果均能显示异常体的顶部、底部埋深基本与理论模型一致。反演精度都较高。单独的模式下,将二种搜索正则化因子的方案加以对比,反演结果显示,方案二的反演效果最佳,精度最高。
从均方根误差曲线来看,两种反演模式下,方案一迭代初期下降较慢,方案二迭代初期下降迅速,且方案二比方案一收敛更快。TM模式下初始均方根误差均为0.4716,迭代8次后,方案一的均方根误差为0.0295,方案二的均方根误差为0.0278。TE模式下初始均方根误差均为0.4387,迭代8次后,方案一的均方根误差为0.1817,方案的二均方根误差为0.1417。这表明了方案二收敛迅速的优点,体现了算法的有效性和稳健性。
图2-图9均是低阻模型下的比较,下面给出高低阻模型的比较,如图10-图17。图10和图14高低阻模型示意图,图11和图15是按方案一得到的正则化因子,图12和图16部分是按方案二得到的正则化因子,图13和图17部分是均方根误差曲线图。其中,电阻率100Ω·m的均匀半空间中,存在高低组两个异常体,左低阻体长宽为1.5km x 1km,电阻率为10Ω·m。右高阻体长宽为1.5km1km,电阻率为1000Ω·m。频率、测点数目和距离和低阻模型的一致。反演初始模型为均匀半空间,将二种正则化因子的选取方案应用到TE、TM两种反演模式中。其中方案一的初始正则化因子取50,步长为2。
TM和TE两种模式的反演结果分别如图4和图5所示。从图中可以看出TE和TM两种反演模式基本都能反应出高低阻体的位置,但两种模式又有区别,TE模式对高低阻体的反映在纵向分辨率上比TM模式的低。将二种搜索正则化因子的方案加以对比,反演结果显示,方案二的反演效果最佳,精度最高。方案一迭代初期下降较慢,方案二迭代初期下降迅速,且方案二比方案一收敛更快。TM模式下初始均方根误差均为0.4940,迭代10次后,方案一均方根误差为0.0213,方案二均方根误差为0.0188。TE模式下初始均方根误差均为0.4784,迭代10次后,方案一的均方根误差为0.1226,方案二的均方根误差为0.1001。表明了方案二收敛迅速的优点,体现了算法的有效性和稳健性。
通过模型试算,均体现了本发明(也即方案二)的高效收敛。若将此算法进一步推广到三维大地电磁反演中,因此算法快速收敛的优点,有望进一步提高反演效率和缩短反演耗时。
如图18所示,本发明还涉及一种处理终端100,其包括:
存储器101,用于存储程序指令;
处理器102,用于运行所述程序指令,以执行所述基于偏差原理的大地电磁正则化反演方法的步骤。
本说明书所公开的实施例只是对本发明单方面特征的一个例证,本发明的保护范围不限于此实施例,其他任何功能等效的实施例均落入本发明的保护范围内。对于本领域的技术人员来说,可根据以上描述的技术方案以及构思,做出其它各种相应的改变以及变形,而所有的这些改变以及变形都应该属于本发明权利要求的保护范围之内。
Claims (4)
1.一种基于偏差原理的大地电磁正则化反演方法,其特征在于,包括如下步骤:
步骤1:获得大地电磁观测数据,并构建反演算子方程;
步骤2:利用偏差原理求解出当前的反演算子方程中的正则化因子,得到新的反演算子方程;
步骤3:求解出大地电磁观测数据中当前的灵敏度,
步骤4:利用所述当前的灵敏度和新的反演算子方程对大地电磁观测数据进行反演,得到反演结果;
步骤5:重复步骤2-4,直至最后一次反演,最后一次反演的结果作为最终的反演结果。
3.根据权利要求1所述的基于偏差原理的大地电磁正则化反演方法,其特征在于,所述步骤2中,采用Morozov偏差原理求解出正则化因子。
4.一种处理终端,其特征在于,其包括:
存储器,用于存储程序指令;
处理器,用于运行所述程序指令,以执行如权利要求1-3任一项所述基于偏差原理的大地电磁正则化反演方法的步骤。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202011424836.6A CN112558164A (zh) | 2020-12-08 | 2020-12-08 | 一种基于偏差原理的大地电磁正则化反演方法及处理终端 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202011424836.6A CN112558164A (zh) | 2020-12-08 | 2020-12-08 | 一种基于偏差原理的大地电磁正则化反演方法及处理终端 |
Publications (1)
Publication Number | Publication Date |
---|---|
CN112558164A true CN112558164A (zh) | 2021-03-26 |
Family
ID=75059804
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN202011424836.6A Pending CN112558164A (zh) | 2020-12-08 | 2020-12-08 | 一种基于偏差原理的大地电磁正则化反演方法及处理终端 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN112558164A (zh) |
Citations (3)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN110398782A (zh) * | 2019-07-17 | 2019-11-01 | 广州海洋地质调查局 | 一种重力数据和重力梯度数据联合正则化反演方法 |
CN111323830A (zh) * | 2020-01-14 | 2020-06-23 | 东华理工大学 | 一种基于大地电磁和直流电阻率数据的联合反演方法 |
CN111833412A (zh) * | 2020-07-16 | 2020-10-27 | 中北大学 | 一种基于分数滤波框架的Tikhonov正则化图像重建方法 |
-
2020
- 2020-12-08 CN CN202011424836.6A patent/CN112558164A/zh active Pending
Patent Citations (3)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN110398782A (zh) * | 2019-07-17 | 2019-11-01 | 广州海洋地质调查局 | 一种重力数据和重力梯度数据联合正则化反演方法 |
CN111323830A (zh) * | 2020-01-14 | 2020-06-23 | 东华理工大学 | 一种基于大地电磁和直流电阻率数据的联合反演方法 |
CN111833412A (zh) * | 2020-07-16 | 2020-10-27 | 中北大学 | 一种基于分数滤波框架的Tikhonov正则化图像重建方法 |
Non-Patent Citations (2)
Title |
---|
张昊: "大地电磁正则化反演方法研究", 《中国优秀博硕士学位论文全文数据库(博士) 基础科学辑》 * |
张衍敏: "基于正则化的多分散系纳米颗粒粒度反演优化方法研究", 《中国优秀博硕士学位论文全文数据库(硕士) 工程科技Ⅰ辑》 * |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN111323830B (zh) | 一种基于大地电磁和直流电阻率数据的联合反演方法 | |
CN109977585B (zh) | 一种高精度大地电磁正演模拟方法 | |
Haber et al. | On optimization techniques for solving nonlinear inverse problems | |
US8537158B2 (en) | Parallel triangle tessellation | |
CN113553748B (zh) | 一种三维大地电磁正演数值模拟方法 | |
Gutzmer et al. | Detection of discontinuities in scattered data approximation | |
CN104851080B (zh) | 一种基于tv的三维pet图像重建方法 | |
van den Doel et al. | Adaptive and stochastic algorithms for EIT and DC resistivity problems with piecewise constant solutions and many measurements | |
CN108287371A (zh) | 直流电阻率无单元法中的背景网格自适应剖分方法 | |
Harhanen et al. | Edge-enhancing reconstruction algorithm for three-dimensional electrical impedance tomography | |
CN109358379B (zh) | 修正总变分模型约束下基于泛函重构的地球物理反演方法 | |
CN112558164A (zh) | 一种基于偏差原理的大地电磁正则化反演方法及处理终端 | |
CN113516754B (zh) | 一种基于磁异常模量数据的三维可视化成像方法 | |
CN112666612B (zh) | 基于禁忌搜索的大地电磁二维反演方法 | |
CN116879964B (zh) | 一种时频电磁频率域数据自约束稳健电阻率反演方法 | |
CN115906559B (zh) | 一种基于混合网格的大地电磁自适应有限元正演方法 | |
CN117092702A (zh) | 孔-隧激发极化探水结构的施工方法及反演探水方法 | |
CN112084655A (zh) | 一种基于非单调线搜索的探地雷达参数反演方法 | |
CN111929348A (zh) | 水平多层结构大地环境下直流接地极地表电位计算方法 | |
CN110322415A (zh) | 基于点云的高精度表面三维重构方法 | |
CN113434942B (zh) | 一种多要素河槽尺度变化数据批量快速提取及计算方法 | |
Pei et al. | A modified L-curve method for choosing regularization parameter in electrical resistance tomography | |
CN114970289A (zh) | 三维大地电磁各向异性正演数值模拟方法、设备及介质 | |
CN114488314A (zh) | 一种基于陆地与水下直流电联合测量的地质反演方法 | |
Rozhenko | Comparison of radial basis functions |
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 | ||
RJ01 | Rejection of invention patent application after publication |
Application publication date: 20210326 |
|
RJ01 | Rejection of invention patent application after publication |