一种基于横向平滑约束的瞬变电磁反演方法与流程

本发明涉及地球物理勘探领域，更具体的，涉及一种瞬变电磁反演方法。

背景技术：

瞬变电磁法作为一种地球物理勘探的有力手段，能够解决矿产勘查、能源、工程、水文、环境地质调查、考古探测等多种地球物理问题。它利用阶跃波或其他脉冲电流地下发射一次脉冲磁场，地下导电地质体在一次场激发下产生感应涡流，进而激发出二次磁场。当发射电流关断时，一次场消失，地下涡流随时间衰减，其衰减时间与导电地质体的电性参数(体积、结构、埋深、电阻率、电导率、介电常数)有关。通过观测二次场值，可以重建地下不均匀体的电性参数，这个过程称为反演。反之，利用已知电性参数的地下结构计算二次场变化的过程，称为正演。

在传统的一维瞬变电磁反演中，往往分别对一条测线上的每个测点单独反演，然后将这条测线上每个测点的反演结果结合作图，形成二维电性参数-深度断面成像。在实际工程中，某个测点的测量结果可能被外界因素(环境噪声、仪器噪声等)干扰，从而带来反演结果的突变。而实际大地相邻测点处的电性参数变化应该不大，特别是在测点较密集时，不应产生异常的电性参数突变。测点单独反演忽略了相邻测点间的约束关系，局部的电性参数突变会影响二维电性参数一深度断面的成像结果，给后期的解释工作带来困难。

2004年，Auken提出了基于Marquardt反演方法的横向约束反演。Marquardt单点反演算法的原理是直接对单个测点的测量数据和理论数据的拟合差求最小，得到此时对应的电阻率值。在此基础上，Auken将若干测点的测量数据组合起来构成输入数据矩阵，建立同时反演多个测点的Marquardt方程，并添加横向约束方程，使相邻测点间电阻率梯度变化最小。将这两个方程联立为方程组求解，得到基于Marquardt方法施加横向约束的反演结果。该方法通过相邻测点间电阻率梯度最小的约束关系，能抑制突变的反演结果，使成像结果更平滑可靠。

在实现本发明的过程中，申请人发现上述现有技术存在如下技术缺陷：

(1)Marquardt方法反演不稳定，极易陷入局部最优解，反演结果依赖初始模型，初始模型选取的不合理会导致反演结果的偏差甚至错误，因此基于Marquardt方法的横向约束平滑方法反演不稳定，适用性不强。

(2)采用正则化方法能寻求满足约束条件的全局最优解，得到稳定的反演结果。目前的正则化反演方法仅局限于单点反演中，无法在相邻测点间施加约束。

技术实现要素：

(一)要解决的技术问题

鉴于以上技术缺陷，本发明提供了一种能够抑制局部测点的跳变，使反演结果更加真实可靠性、同时具备正则化迭代稳定和全局最优的基于横向平滑约束的瞬变电磁反演方法。

(二)技术方案

本发明技术方案如下：

本发明提供了一种基于横向平滑约束的瞬变电磁反演方法，包括：

S1、向地下发射一次脉冲磁场，测量地下导电地质体在一次脉冲磁场激发下产生感应涡流而激发的二次磁场，所测量得到的二次磁场数据作为观测数据；

S2、构建反演初始模型；

S3、构建观测数据与正演理论数据之间的不匹配泛函；

S4、构建约束导电地质体的电性参数变化的纵向约束泛函；

S5、构建约束测点间导电地质体的电性参数变化的横向约束泛函；

S6、由不匹配泛函、纵向约束泛函和横向约束泛函构建反演目标函数；

S7、以初始反演模型为初始条件，求解反演目标函数取最小时对应的电性参数，根据该电性参数重建地下结构。

所述步骤S6构建的反演目标函数为：P＝P_t+αP_V+λP_H；式中：P_t为不匹配泛函；P_V为纵向约束泛函；P_H为横向约束泛函；α称为正则化因子，用于调整反演过程中纵向约束泛函占总目标函数的比重；λ称为横向约束因子，是大于0的正数，根据地层横向连续性取值。所述步骤S5构建的横向约束泛函是：其中代表横向电性参数梯度，|| ||_L2代表二范数。所述步骤S4构建的纵向约束泛函包括最小约束模型泛函、最平缓约束泛函、最光滑约束泛函、总变分约束泛函、最小支撑泛函，或者最小梯度支撑泛函。所述步骤S7求解反演目标函数的算法包括高斯牛顿法。所述电性参数包括导电地质体的电阻率、电导率或介电常数。

(三)有益效果

1、本发明通过施加横向约束泛函，增加了相邻测点间的约束关系，抑制了局部测点的跳变，使反演成像结果更真实可靠；

2、本发明基于正则化方法，具备正则化迭代稳定、全局最优的优点，通过选择不同的纵向电性参数约束条件，能得到纵向模型不同的反演结果。

附图说明

图1是根据本发明的一种基于横向平滑约束的瞬变电磁反演方法实施例流程图；

图2是根据本发明实施例的真实大地模型图；

图3是根据本发明实施例的测点1(10m)处的输入数据在瞬变电磁二次场响应图；

图4是无约束反演的Occam反演成像结果；

图5是根据本发明实施例的横向平滑约束的Occam反演成像结果。

具体实施方式

为使本发明的目的、技术方案和优点更加清楚明白，以下结合具体实施例，并参照附图，对本发明作进一步的详细说明。

本发明采用的技术方案是，一种基于横向平滑约束的瞬变电磁反演方法，包括以下步骤：

S1、向地下发射一次脉冲磁场，测量地下导电地质体在一次脉冲磁场激发下产生感应涡流而激发的二次磁场，所测量得到的二次磁场数据作为观测数据。

该步骤由瞬变电磁法观测二次场数据，包括采集时间T及其对应的场值d_obs。其中，场值d_obs＝[d_obs1，d_obs2，...，d_obsM]^T由多个测点的观测数据拼接而成，d_obsi＝[d₁，d₂，...，d_L]，M为参与反演的测点数，L为时间测道总数，d_i为某个采样时间处的采样数据。每个采样点的数据要满足衰减特性，需剔除早期包含一次场的部分和晚期由于噪声产生震荡的部分。

S2、构建反演初始模型。

构建反演初始模型需要输入反演参数。反演参数包括发射电流I，发射线圈半径R，发射线圈匝数N_t，观测偏移距r，接收线圈有效面积S。反演初始模型包括测点数M，反演层数N，各层厚度h＝[h₁，h₂，...，h_N-1]，反演初始电性参数，m＝[m_1，1，m_1，2，...，m_1，N，m_2，1，m_2，2，...，m_2，N，...，m_M，1，m_M，2，...，m_M，N]^T。反演初始电性参数可以是导电地质体的电阻率、电导率、介电常数。由于随深度变化反演分辨率降低，通常令h_i+1/h_i＞1。

S3、构建观测数据与正演理论数据之间的不匹配泛函

不匹配泛函可写为

P_t＝||W_d(F(m)-d_obs)||² (1)

式中，W_d为描述数据比重的权值矩阵，F(m)是通过反演参数正演得到的理论数据，由M个正演数据拼接而成。W_d可选择为

σ_i，j为第i个测点第j个时间测道的噪声。

S4、构建约束电性参数变化的纵向约束泛函。

纵向约束泛函可选择但不仅限于以下几种：

最小约束模型泛函

最平缓约束泛函

最光滑约束泛函

总变分约束泛函

最小支撑泛函

最小梯度支撑泛函

其中m为需要进行约束的电性参数矩阵，m_ref为参考模型，代表纵向电性参数梯度，β称为聚焦因子，是一个远小于1的正数。

S5、构建约束测点间电性参数变化的横向约束泛函。

为使相邻测点间横向变化尽量平滑，泛函可写为

S6、由不匹配泛函、纵向约束泛函和横向约束泛函构建反演目标函数。

反演目标函数可写为

P＝P_t+αP_V+λP_H (10)

式中α称为正则化因子，用于调整反演过程中纵向约束泛函占总目标函数的比重。λ称为横向约束因子，是大于0的正数。当已知地层横向连续性较好时，可以选择稍大的λ值。

S7、以初始反演模型为初始条件，求解目标函数取最小时对应的电性参数，根据该电性参数重建地下结构。

该步骤可采用高斯牛顿法求解，将正演函数F(m)近似写为

F(m)＝F(m₀)+J₀(m₁-m₀) (11)

式中m₀为某步迭代的初始条件，m₁为本次迭代的待求电性参数，J₀为某步迭代时正演函数对初始条件的导数矩阵，其每个元素为

将公式(13)代入目标函数中，令得迭代公式

每步迭代后用新求得的结果更新初始条件，当拟合差小于期望拟合差时，迭代终止。每步迭代过程需在一定区间内动态选择α值，使目标函数达到最小。选择使用的最优化算法包括但不仅限于高斯牛顿法。

实例

本实施例中，通过所阐述的横向平滑约束算法，对瞬变电磁二次场数据进行反演。

(1)输入由瞬变电磁法测得的二次场数据。设置大地模型如图1所示，电性参数使用导电地质体的电阻率，分别为300Ωm，100Ωm，300Ωm。每个测点处的观测数据应具备衰减特性，如图2所示。

(2)输入反演参数并构建初始模型。发射线圈匝数为1，发射电流为5A，发射半径R＝50m，观测点偏移距r＝0m，接收线圈有效面积为2000m²。设测点数M＝19，反演模型为30层，厚度增长比例h_i+1/h_i＝1.04，第一层厚度h₁＝13m，每层的初始电阻率为100Ωm。

(3)构建观测数据与正演理论数据之间的不匹配泛函。在此假设各时间测道数据所占比重一致，即权值矩阵为单位对角阵。

(4)构建约束电阻率变化的纵向约束泛函。在此使用Occam反演中应用的最平缓约束泛函其中

(5)构建约束测点间电阻率变化的横向约束泛函其中

每一行中1和-1相隔N-1个0。

(6)由不匹配泛函、纵向约束泛函和横向约束泛函构建反演目标数。设置α的搜索范围为10^-2-10¹，λ＝1.1。

(7)以初始模型为初始条件，采用高斯牛顿法求解目标函数取最小时对应的电阻率，得到如图4所示的反演结果。图3为未施加横向平滑约束的反演结果。

以上所述的具体实施例，对本发明的目的、技术方案和有益效果进行了进一步详细说明，应理解的是，以上所述仅为本发明的具体实施例而已，并不用于限制本发明，凡在本发明的精神和原则之内，所做的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。