CN114779365A - 一种离散函数拟合的重磁交叉梯度联合物性反演方法 - Google Patents

Abstract

本发明公开了一种离散函数拟合的重磁交叉梯度联合物性反演方法，所述方法包括：采集重力异常数据和磁力异常数据；确定研究区域及函数阶次并剖分；建立重力异常、磁力异常矩阵，得到重力异常灵敏度矩阵和磁力异常灵敏度矩阵；反演目标函数，计算得到函数参数A、B；依据函数参数正演重磁异常，根据拟合差判断是否为预期结果，不是即循环迭代，直至得到预期目标值，从而反演出地下物性。本发明通过构建地下半空间密度模型，引入物性函数，将运行时内存占用率从传统方法的31％降到7％，运行解算时间从传统方法的4993s降到647s，大幅提高了交叉梯度反演的反演效率。

Description

一种离散函数拟合的重磁交叉梯度联合物性反演方法

技术领域

本发明涉及一种离散函数拟合的重磁交叉梯度联合物性反演方法，属于物性反演技术领域。

背景技术

重力场及磁力场是地球的天然磁场，重力异常是重力场的扰动/磁力异常是磁力场的扰动，由地下物性差异产生。因此，通过重力异常/磁力异常即可以定量计算得到地下的物性差异，这种反推过程即为物性反演。

现有的物性反演中，通常包含两个步骤：1.将研究区域离散化，得到单位物性离散单元在观测点处响应的灵敏度矩阵；2.采用该灵敏度矩阵对观测数据进行反算得到各个离散单元的物性值。因此密度反演方法需要较少的先验信息就能够直观地得到地下密度的分布情况，被广泛应用于重力数据处理中。由于离散单元个数及观测数据点个数通常较多，会得到一个很大的灵敏度矩阵。100×100个观测数据100×100×100个剖分单元即需要建立一个37.25GB的核矩阵。同时，离散单元个数往往小于观测数据点个数，为欠定问题求解，因此需要迭代计算以保证计算精度，这导致密度反演需要多次大量的矩阵运算来实现。

这种将研究区域地下空间采用小网格单元进行剖分，然后计算每一个剖分单元物性的方式，剖分单元较少会影响反演结果的分辨率，而剖分单元过多会极大地降低计算效率，特别是在传统物性反演方法的基础上加入两种物性的交叉梯度，实现的交叉梯度反演方法，由于多种物性的引入，反演效率更为低下。

发明内容

本发明的目的是提供一种离散函数拟合的重磁交叉梯度联合物性反演方法，该反演方法将地下离散物性分布用函数形式表示，不仅解决了现有技术中对研究区域剖分不合理，导致分辨率与计算效率难以同时满足要求的问题，还解决了现有技术中交叉梯度反演方法，计算量大，计算效率低的问题。

为达到上述目的，本发明采用的技术方案是：一种离散函数拟合的重磁交叉梯度联合物性反演方法，所述方法包括：

S1：采集重力异常数据记为Δg，采集重力异常数据记为ΔT；

S2：将研究区域的密度及磁化率分布表示成随空间位置变化的函数，并确定研究区域及物性函数阶次s，物性函数的形式为：

其中ρ(x，z)为(x，z)点处密度，M(x，z)为(x，z)点处磁化强度，a_i，j，b_i，j均为函数系数；

S3：将研究区域剖分为n＝n_x×n_z个长方形单元，利用(1)式将剖分单元的物性用物性函数方式表示，其计算形式为：

其中，ρ＝[ρ₁ ρ₂ … ρ_n]^T为保存所有剖分单元密度的数组；

M＝[M₁ M₂ … M_n]^T为保存所有剖分单元磁化强度的数组；

A＝[a_0，0 … a_0，s a_1，0 … a_1，s-1 … a_i，j … a_s，0]^T为保存函数系数的数组

B＝[b_0，0 … b_0，s b_1，0 … b_1，s-1 … b_i，j … b_s，0]^T为保存函数系数的数组；

为保存物性函数各项的矩阵；(x_k，z_k)为第k个剖分单元中心点的坐标；k＝1，2，...，n；

S4：构建二维重磁反演求解目标函数：

首先得到四边形重磁正演公式：

其中，Δg_t，ΔT_t分别为第k个剖分单元在第t个观测点处所产生的重力异常及磁异常，γ为万有引力常数6.672×10^-11m³/(kg·s²)，μ₀为真空磁导率4π×10^-7H/m，ρ_k为第k个剖分单元的密度，M_k为第k个剖分单元的磁化强度，X_i＝x_t-x_k-(i-0.5)dx，Z_i＝z_t-z_k-(j-0.5)dz，(x_t，z_t)为第t个观测点的坐标，t＝1，2，...，R，R为观测点个数，(x_k，z_k)为第k个剖分单元中心点的坐标，dx为一个剖分单元的x方向的长度，dz为一个剖分单元的z方向的长度；

将(3)式写成以下形式：

其中，Δg＝[Δg₁ Δg₂ … Δg_R]^T为保存所有重力异常数据的数组；

ΔT＝[ΔT₁ ΔT₂ … ΔT_R]^T为保存所有重力异常数据的数组，G_g，G_T分别为重磁异常正演的核矩阵，其中，第t行第k个元素为：

将(2)式带入(4)式得：

令G_g·P＝V_g、G_T·P＝V_T，其中，Vg，VT为重磁反演的灵敏度矩阵，构建出交叉梯度反演目标求解目标函数：

其中，

表示二范数的平方；

S5：已知V_g，V_T，P，Δg，ΔT，对于(7)式所示的目标函数求得其极小值即可得到物性函数系数A、B，将其代入(2)式得到地下物性分布，采用迭代求解方法计算，依据函数参数正演重磁异常，根据拟合差判断是否为预期结果，不是即循环迭代，直至得到预期目标值。

进一步地，所述物性函数阶次s取5。

进一步地，设置n_x＝50，n_z＝20。

进一步地，设置n_x＝60，n_z＝20。

由于上述技术方案的运用，本发明与现有技术相比具有下列优点：

本发明一种离散函数拟合的重磁交叉梯度联合物性反演方法，通过构建地下半空间密度模型，引入物性函数，将运行时内存占用率从传统方法的31％降到7％，运行解算时间从传统方法的4993s降到647s，大幅提高了交叉梯度反演的反演效率。

附图说明

附图1为本发明一种离散函数拟合的重磁交叉梯度联合物性反演方法的流程示意图；

附图2为地下空间离散化剖分示意图；

附图3为实施例1的重力异常数据及反演结果示意图；

附图4为实施例1的磁力异常数据及反演结果示意图；

附图5为实施例2的重力异常数据及反演结果示意图；

附图6为实施例2的磁力异常数据及反演结果示意图。

具体实施方式

实施例1：

一种离散函数拟合的重磁交叉梯度联合物性反演方法，参照附图1-4，所述方法包括：

S1：采集重力异常数据记为Δg，采集重力异常数据记为ΔT；

S2：将研究区域的密度及磁化率分布表示成随空间位置变化的函数，并确定研究区域及物性函数阶次s＝5，物性函数的形式为：

其中ρ(x，z)为(x，z)点处密度，M(x，z)为(x，z)点处磁化强度，a_i，j，b_i，j均为函数系数；

S3：将研究区域剖分为n＝n_x×n_z个长方形单元，这里，n_x＝50、n_z＝20，利用(1)式将剖分单元的物性用物性函数方式表示，其计算形式为：

其中，ρ＝[ρ₁ ρ₂ … ρ_n]^T为保存所有剖分单元密度的数组；

M＝[M₁ M₂ … M_n]^T为保存所有剖分单元磁化强度的数组；

A＝[a_0，0 … a_0，s a_1，0 … a_1，s-1 … a_i，j … a_s，0]^T为保存函数系数的数组

B＝[b_0，0 … b_0，s b_1，0 … b_1，s-1 … b_i，j … b_s，0]^T为保存函数系数的数组；

为保存物性函数各项的矩阵；(x_k，z_k)为第k个剖分单元中心点的坐标；k＝1，2，...，n；

S4：构建二维重磁反演求解目标函数：

首先得到四边形重磁正演公式：

其中，Δg_t，ΔT_t分别为第k个剖分单元在第t个观测点处所产生的重力异常及磁异常，γ为万有引力常数6.672×10^-11m³/(kg·s²)，μ₀为真空磁导率4π×10^-7H/m，ρ_k为第k个剖分单元的密度，M_k为第k个剖分单元的磁化强度，X_i＝x_t-x_k-(i-0.5)dx，Z_i＝z_t-z_k-(j-0.5)dz，(x_t，z_t)为第t个观测点的坐标，t＝1，2，...，R，R为观测点个数，(x_k，z_k)为第k个剖分单元中心点的坐标，dx为一个剖分单元的x方向的长度，dz为一个剖分单元的z方向的长度；

将(3)式写成以下形式：

其中，Δg＝[ΔΔg₁ Δg₂ … Δg_R]^T为保存所有重力异常数据的数组；

ΔT＝[ΔT₁ ΔT₂ … ΔT_R]^T为保存所有重力异常数据的数组，G_g，G_T分别为重磁异常正演的核矩阵，其中，第t行第k个元素为：

将(2)式带入(4)式得：

令G_g·P＝V_g、G_T·P＝V_T，其中，Vg，VT为重磁反演的灵敏度矩阵，构建出交叉梯度反演目标求解目标函数：

其中，

表示二范数的平方；

S5：已知V_g，V_T，P，Δg，ΔT，对于(7)式所示的目标函数求得其极小值即可得到物性函数系数A、B，将其代入(2)式得到地下物性分布，采用迭代求解方法计算，依据函数参数正演重磁异常，根据拟合差判断是否为预期结果，不是即循环迭代，直至误差小于真实值的1％。

实施例2：

一种离散函数拟合的重磁交叉梯度联合物性反演方法，参照附图1-2、5-6，所述方法包括：

S1：采集重力异常数据记为Δg，采集重力异常数据记为ΔT；

S2：将研究区域的密度及磁化率分布表示成随空间位置变化的函数，并确定研究区域及物性函数阶次s＝5，物性函数的形式为：

其中ρ(x，z)为(x，z)点处密度，M(x，z)为(x，z)点处磁化强度，a_i，j，b_i，j均为函数系数；

S3：将研究区域剖分为n＝n_x×n_z个长方形单元，这里，n_x＝60、n_z＝20，利用(1)式将剖分单元的物性用物性函数方式表示，其计算形式为：

其中，ρ＝[ρ₁ ρ₂ … ρ_n]^T为保存所有剖分单元密度的数组；

M＝[M₁ M₂ M₂ … M_n]^T为保存所有剖分单元磁化强度的数组；

A＝[a_0，0 … a_0，s a_1，0 … a_1，s-1 … a_i，j … a_s，0]^T为保存函数系数的数组

B＝[b_0，0 … b_0，s b_1，0 … b_1，s-1 … b_i，j … b_s，0]^T为保存函数系数的数组；

为保存物性函数各项的矩阵；(x_k，z_k)为第k个剖分单元中心点的坐标；k＝1，2，...，n；

S4：构建二维重磁反演求解目标函数：

首先得到四边形重磁正演公式：

其中，Δg_t，ΔT_t分别为第k个剖分单元在第t个观测点处所产生的重力异常及磁异常，γ为万有引力常数6.672×10^-11m³/(kg·s²)，μ₀为真空磁导率4π×10^-7H/m，ρ_k为第k个剖分单元的密度，M_k为第k个剖分单元的磁化强度，X_i＝x_t-x_k-(i-0.5)dx，Z_i＝z_t-z_k-(j-O.5)dz，(x_t，z_t)为第t个观测点的坐标，t＝1，2，...，R，R为观测点个数，(x_k，z_k)为第k个剖分单元中心点的坐标，dx为一个剖分单元的x方向的长度，dz为一个剖分单元的z方向的长度；

将(3)式写成以下形式：

其中，Δg＝[Δg₁ Δg₂ … Δg_R]^T为保存所有重力异常数据的数组；

ΔT＝[ΔT₁ ΔT₂ … ΔT_R]^T为保存所有重力异常数据的数组，G_g，G_T分别为重磁异常正演的核矩阵，其中，第t行第k个元素为：

将(2)式带入(4)式得：

令G_g·P＝V_g、G_T·P＝V_T，其中，Vg，VT为重磁反演的灵敏度矩阵，构建出交叉梯度反演目标求解目标函数：

其中，

表示二范数的平方；

S5：已知V_g，V_T，P，Δg，ΔT，对于(7)式所示的目标函数求得其极小值即可得到物性函数系数A、B，将其代入(2)式得到地下物性分布，采用迭代求解方法计算，依据函数参数正演重磁异常，根据拟合差判断是否为预期结果，不是即循环迭代，直至误差小于真实值的1％。

Claims (4)

1.一种离散函数拟合的重磁交叉梯度联合物性反演方法，其特征在于，所述方法包括：

S1：采集重力异常数据记为Δg，采集重力异常数据记为ΔT；

S2：将研究区域的密度及磁化率分布表示成随空间位置变化的函数，并确定研究区域及物性函数阶次s，物性函数的形式为：

其中，ρ(x，z)为(x，z)点处密度，M(x，z)为(x，z)点处磁化强度，a_i，j，b_i，j均为函数系数；

S3：将研究区域剖分为n＝n_x×n_z个长方形单元，利用(1)式将剖分单元的物性用物性函数方式表示，其计算形式为：

其中，ρ＝[ρ₁ ρ₂ … ρ_n]^T为保存所有剖分单元密度的数组；

M＝[M₁ M₂ … M_n]^T为保存所有剖分单元磁化强度的数组；

A＝[a_0，0 … a_0，s a_1，0 … a_1，s-1 … a_i，j … a_s，0]^T为保存函数系数的数组

B＝[b_0，0 … b_0，s b_1，0 … b_1，s-1 … b_i，j … b_s，0]^T为保存函数系数的数组；

为保存物性函数各项的矩阵；(x_k，z_k)为第k个剖分单元中心点的坐标；k＝1，2，...，n；

S4：构建二维重磁反演求解目标函数：

首先得到四边形重磁正演公式：

其中，Δg_t，ΔT_t分别为第k个剖分单元在第t个观测点处所产生的重力异常及磁异常，γ为万有引力常数6.672×10^-11m³/(kg·s²)，μ₀为真空磁导率4π×10^-7H/m，ρ_k为第k个剖分单元的密度，M_k为第k个剖分单元的磁化强度，X_i＝x_t-x_k-(i-0.5)dx，Z_i＝z_t-z_k-(j-0.5)dz，(x_t，z_t)为第t个观测点的坐标，t＝1，2，...，R，R为观测点个数，(x_k，z_k)为第k个剖分单元中心点的坐标，dx为一个剖分单元的x方向的长度，dz为一个剖分单元的z方向的长度；

将(3)式写成以下形式：

其中，Δg＝[Δg₁ Δg₂ … Δg_R]^T为保存所有重力异常数据的数组；

ΔT＝[ΔT₁ ΔT₂ … ΔT_R]^T为保存所有重力异常数据的数组，G_g，G_T分别为重磁异常正演的核矩阵，其中，第t行第k个元素为：

将(2)式带入(4)式得：

令G_g·P＝V_g、G_T·P＝V_T，其中，Vg，VT为重磁反演的灵敏度矩阵，构建出交叉梯度反演目标求解目标函数：

其中，

表示二范数的平方；

S5：已知V_g，V_T，P，Δg，ΔT，对于(7)式所示的目标函数求得其极小值即可得到物性函数系数A、B，将其代入(2)式得到地下物性分布，采用迭代求解方法计算，依据函数参数正演重磁异常，根据拟合差判断是否为预期结果，不是即循环迭代，直至得到预期目标值。

2.根据权利要求1所述的一种离散函数拟合的重磁交叉梯度联合物性反演方法，其特征在于，所述物性函数阶次s取5。

3.根据权利要求2所述的一种离散函数拟合的重磁交叉梯度联合物性反演方法，其特征在于，设置n_x＝50，n_z＝20。

4.根据权利要求2所述的一种离散函数拟合的重磁交叉梯度联合物性反演方法，其特征在于，设置n_x＝60，n_z＝20。

Images