CN100429531C - 目标最小化的三维电磁快速反演方法 - Google Patents

Abstract

本发明是电磁勘探数据处理的目标最小化的三维电磁快速反演方法，步骤是根据已知地质资料，确定已知的目标区域，根据已知电导率参数，三维插值建立背景电导率组合；对组合用积分法计算背景场，根据生产要求反演未知区域，建立目标初始电导率组合计算异常场与背景场，合并得到总场，进而获得整个电导率参数组合的地面电磁场响应和视电导率、相位响应；对实测视电导率曲线及相位曲线进行三维插值，分布于网格节点；根据目标电导率参数组合的深度、尺寸选择反演频率，采用非线性共轭梯度法迭代反演；反演中只对目标电导率组合修改，拟合主要的反演频率得到目标体的地下电导率参数组合。

Description

目标最小化的三维电磁快速反演方法

技术领域

本发明涉及地球物理勘的大地电磁勘探数据处理技术，是一种目标最小化的三维电磁快速反演方法。

背景技术

大地电磁是一种地球物理勘探方法，一般通过在地面观测电磁场的变化来了解地下介质电导率的变化情况。上世纪九十年代，由大地电磁勘探方法(简称MT)发展出了一种新的电磁勘探方法，称为连续电磁阵列剖面法(简称CEMP)，由简单的大地电磁单点采集方法，发展为沿测线阵列武的多道采集方法；同时，由于计算机技术快速发展，二维有限元法、二维连续介质反演方法成熟起来，计算精度大大提高，可以比较精确地模拟实际地下电导率介质的情况，因而，使该方法由针对区域勘探为主，转换为针对目标勘探为主一类重要的油气勘探方法，在地震勘探困难区，如山前带、火成岩出露区、表层砾石区、碳酸盐岩地区、黄土塬覆盖区取得了较好的勘探效果，弥补了地震勘探资料的不足。

近几年随着油气勘探工作的不断深入，电磁勘探方法也面临更复杂的勘探问题，勘探工区地表地质条件更加复杂，地质目标更为隐秘，使电磁勘探方法面临着前所未有的挑战。对大地电磁以及连续电磁阵列剖面法反演精度要求越来越高，目前正在使用的二维有限元法或二维连续介质方法难以满足地质解释的需要。因而，电磁三维反演方法得到快速发展。

但是目前电磁三维反演方法理论研究较多，实际应用的例子并不多见，原因是：通常电磁三维反演方法是将地下介质统统作为未知的电导率参数处理，如果一个面积为10km×10km，勘探区域最大深度10km的三维勘探，为了精细地进行反演，可以将电导率参数组合剖分成20m×20m×20m单元，未知电导率参数就将会有500×500×500＝125,000,000个未知电导率参数，可见未知参数非常多。

未知参数越多，等值性越严重；在反演过程中组成的灵敏度矩阵就非常巨大，会占用大量的计算机内存，现在的计算机难以承受；在正演模拟过程中计算量非常巨大，要想得到理想的反演结果，可能会用一周或更长的时间。

如果增大剖分单元，未知参数减少了，但是很明显三维反演模拟精度会随之降低。这就是目前三维反演方法难以实用的根本原因。

开展大地电磁测深法或连续电磁阵列剖面法三维工作的地区一般都是开展过地震勘探工作的地区，甚至是开展过三维地震勘探工作的地区，已知资料多但对电磁法资料处理没有发挥作用。如测井资料、浅层地震资料，未得到充分利用。

发明内容

本发明目的是针对背景技术中存在的问题，提供一种利用地震或钻井已知信息作为非均匀电性介质背景，反演目标达到最小，能快速、高精度三维反演的目标最小化的三维电磁快速反演方法。

本发明采用以下技术方案实现，反演步骤如下：

1)在开展过地震、地质或钻井工作的地区，实测目标勘探工作地区视电阻率曲线及相位曲线，采用通常的方法进行三维插值，使之均匀分布于剖分单元的网格节点；

2)将电磁勘探工作区分为三部分：①地下介质已知的背景区域及电阻率值(σ_n)，②相邻的已知区域(D_b)及相应的电阻率参数(σ_b)，③由勘探任务确定的要进行反演的目标区域(D_a)，其相应的电阻率值参数(σ_a)由邻区已知电测井资料、地质露头资料确定；

3)采用通常的电磁一维解析公式计算已知的背景区域及电阻率值(σ_n)有关的电场(Eⁿ)和磁场(Hⁿ)，采用公式(7)(8)三维电磁积分公式计算已知区域(D_b)有关的电场E^Db和磁场H^Db：

$E^{D_b}\left(r_j\right)=\underset{D_b}{\∫}{\hat{G}}_E\left(r_j\right|r)\·\left(\mathrm{\Δ}{\mathrm{\σ}}_{b}E\left(r\right)\right)\mathrm{dv}---\left(7\right)$

$H^{D_b}\left(r_j\right)=\underset{D_b}{\∫}{\hat{G}}_H\left(r_j\right|r)\·\left(\mathrm{\Δ}{\mathrm{\σ}}_{b}E\left(r\right)\right)\mathrm{dv}---\left(8\right)$

E^Db和H^Db：分别表示与已知区域有关的电场和磁场其中：

和

分别为电场和磁场的格林算子，Δσ_b＝σ_b-σ_n已知区域电阻率值与相应的背景电阻率值之差，E(r)为任意一点的电场；

4)采用如下(9)(10)三维电磁积分公式计算目标区域(D_a)的E^Da和磁场H^Da

$E^{D_a}\left(r_j\right)=\underset{D_a}{\∫}{\hat{G}}_E\left(r_j\right|r)\·\left(\mathrm{\Δ}{\mathrm{\σ}}_{a}E\left(r\right)\right)\mathrm{dv}---\left(9\right)$

$H^{D_a}\left(r_j\right)=\underset{D_a}{\∫}{\hat{G}}_H\left(r_j\right|r)\·\left(\mathrm{\Δ}{\mathrm{\σ}}_{a}E\left(r\right)\right)\mathrm{dv}---\left(10\right)$

E^Da和H^Da：分别表示与目标区域有关的电场和磁场，

和

分别为电场和磁场的格林算子，Δσ_a＝σ_a-σ_n目标区域电阻率值与背景电阻率值之差，E(r)为任意一点的电场；

5)根据上述3)、4)的结果计算地面上的视电阻率及相位相应；

6)根据目标区域大小、埋深、规模及由5)计算的目标区域对应的视电阻率曲线和频率相位曲线，确定5～20个频点的视电阻率曲线及相位曲线进行反演；

7)采用常规的非线性共轭梯度法，针对目标区域电阻率值参数进行三维迭代反演，得到电阻率参数；

8)根据背景电阻率参数(σ_n)、已知区域电阻率(σ_b)以及通过三维反演得到电阻率参数(σ_a)得到整个地下介质的电阻率分布情况，建立电磁三维反演数据体，作为地质解释的依据。

本发明还采用如下技术方案实现：

所述的(σ_a)(σ_b)代表区域内每一个单元的电阻率值，根据目标地区反演的深度范围内已有钻井资料确定其值。

所述的所述的(σ_a)(σ_b)代表区域内每一个单元的电阻率值用一个任意电阻率值充填。

所述的5～20个频点的视电阻率曲线及相位曲线频率范围一定大于该目标区域对应的频率范围。

所述的非线性共轭梯度反演过程是，第一步修改未知区域的电阻率参数，计算异常场；第二步利用已计算过的背景场及未知区域的异常场，计算视电阻率曲线和相位曲线；第三步计算视电阻率曲线及相位曲线拟合误差；重复上述三步直到拟合误差满足要求，得到未知区域的电阻率参数。

本发明可以充分利用诸如地震或钻井已知信息，作为已知非均匀电性介质背景，在整个反演过程前只需正演计算一次，由于已知电导率参数不参与反演，反演目标达到最小，从而达到快速、高精度三维反演的目的。

附图说明

图1为本发明三维地下电导率参数组合示意。

具体实施方式

以下结合附图详细说明本发明。

本发明第一步根据地震、钻井等已知资料，确定已知的区域，同时对非研究区根据已知资料确定电导率参数，通过三维插值建立三维背景电导率参数组合；对三维背景电导率参数组合采用积分方程法计算背景场，整个反演过程只需计算一次；根据实际生产要求确定三维反演未知区域，并建立研究目标的初始电导率参数组合；对初始电导率参数组合采用积分方程法计算异常场与背景场，合并得到总场，进而获得整个电导率参数组合的地面电磁场响应和视电导率、相位响应；对实测视电导率曲线及相位曲线进行三维插值，使之均匀分布于电导率参数组合网格节点；根据目标电导率参数组合的深度、尺寸选择主要的反演频率，使计算量进一步减少；采用非线性共轭梯度法迭代反演进行三维反演；反演中只对要研究初始电导率参数组合进行修改，只需拟合主要的反演频率；最终得到要研究目标体的地下电导率参数组合。

本发明仅对研究目标进行反演，将研究目标作为异常区域，同时给定初始地下电导率参数组合，利用三维积分方程的方法计算异常场。

本发明对实际资料进行网格差值建立对应于电导率参数组合的实际数据体，根据研究目标的深度确定视电导率曲线及阻抗相位曲线频率段，计算拟合误差。

如在一个目标地区施工，一般地震浅层反射资料很好，而深层地质资料不好，就可以将浅层资料作为已知资料处理，浅层4km资料已知，精细反演时的地下电导率参数组合中未知电导率参数个数就变为原来的50％；同理，如在山地施工，一般平原区地震反射资料较好，而山区地震资料不好，也可将平原区作为已知资料处理，来减少未知参数。本发明的最大特点是充分利用已知资料，缩小反演区域，省时省力，减少了多解性，精度较高。

本发明地下电导率参数组合如图1，背景地下电导率参数组合为三层，电阻率分别为σ₁，σ₂，σ₃，以及已知的三维异常区域D_b，电阻率为σ_b，令σ_b＝σ_n+Δσ_b

其中：σ_n为背景电阻率值，目标区域D_a，电阻率为σ_a，

令σ_a＝σ_n+Δσ_a

麦克斯韦方程可写成：

$\▿\×H={\mathrm{\σ}}_{n}E+j={\mathrm{\σ}}_{n}E+j^{{\mathrm{\Δ\σ}}_b}+j^{{\mathrm{\Δ\σ}}_a}+j^e---\left(1\right)$

$\▿\×E=i{\mathrm{\ω\μ}}_{0}H---\left(2\right)$

$j^{{\mathrm{\Δ\σ}}_a}=\left\{\begin{array}{cc}{\mathrm{\Δ\σ}}_{a}E,& r\∈D_a\\ 0,& r\∉D_a\end{array}\right.---\left(3\right)$

$j^{{\mathrm{\Δ\σ}}_b}=\left\{\begin{array}{cc}{\mathrm{\Δ\σ}}_{b}E,& r\∈D_b\\ 0,& r\∉D_b\end{array}\right.---\left(4\right)$

电磁三维积分积分方程为

$E\left(r_j\right)=E^n+\underset{D_b}{\∫}{\hat{G}}_E\left(r_j\right|r)\·\left({\mathrm{\Δ\σ}}_{b}E\left(r\right)\right)\mathrm{dv}---\left(5\right)$

$+\underset{D_a}{\∫}{\hat{G}}_E\left(r_j\right|r)\·\left({\mathrm{\Δ\σ}}_{a}E\left(r\right)\right)\mathrm{dv}$

$H\left(r_j\right)=H^n+\underset{D_b}{\∫}{\hat{G}}_H\left(r_j\right|r)\·\left({\mathrm{\Δ\σ}}_{b}E\left(r\right)\right)\mathrm{dv}---\left(6\right)$

$+\underset{D_a}{\∫}{\hat{G}}_H\left(r_j\right|r)\·\left({\mathrm{\Δ\σ}}_{a}E\left(r\right)\right)\mathrm{dv}$

公式(5)、(6)中的Eⁿ、Hⁿ分别表示背景电场、磁场，第二项分别表示已知区域的电场、磁场，第三项分别表示异常区域的电场、磁场。

本发明具体实施步骤是：

1)在开展过地震、地质或钻井工作的地区，实测目标勘探工作地区视电阻率曲线及相位曲线，采用通常的方法进行三维插值，使之均匀分布于剖分单元的网格节点；

在野外实际大地电磁或连续阵列剖面法观测时，往往受到地形起伏、河流、胡泊、村庄、电站、煤矿、高压线、高速路影响，测点疏密不均，为了能够有效地反映和控制地下异常区域电导率参数组合，对实测的视电导率曲线及相位曲线进行三维插值，使之在地面上均匀规则分布。

2)将电磁勘探工作区分为三部分：①地下介质已知的背景区域及电阻率值(σ_n)，②相邻的已知区域(D_b)及相应的电阻率参数(σ_b)，③由勘探任务确定的要进行反演的目标区域(D_a)，其相应的电阻率值参数(σ_a)由邻区电测井资料、地质露头资料确定。

MT或CEMP开展目标勘探工作的地区，一般都开展过地震、地质或钻井工作，都有或多或少的已知的物探地质资料。为了充分利用这些已知资料，将这些已知资料引入到三维反演过程中。

(D_b)、(D_a)在反演前，都被剖分成了多个小的正方体，就是所说的单元，每个小的单元被认为是均匀的，可以用一个电阻率值表示，上面所说的(σ_a)(σ_b)代表一系列的电阻率值，包括相应区域内每一个单元的电阻率值。如要反演的深度范围为4-8千米，如果目标区域附近有钻井资料，可以参照给出大概的值。如果根本不知道任何参考信息，可以直接用一个电阻率值充填。在反演过程中修改目标区域的电阻率值，使之尽量符合客观实际。

然后，要确定研究对象，即要进行反演的区域D_a，比如要反演深度大于3km含油气构造，根据已知资料确定其规模，延伸深度，从而得到要反演的研究区域，根据相邻区域的电测井资料、地质露头资料，确定异常区域电导率值参数；

3)对背景、已知区域的三维地下电导率参数组合采用三维电磁积分方程法计算相对应的背景电磁场；

根据2)确定的背景电导率参数，根据公式(5)、(6)分别计算Eⁿ、Hⁿ，根据2)确定的已知区域及电导率参数分别计算第一个积分项，得到相应的背景电磁场。

电场E^Db和磁场H^Db由三维电磁积分公式计算：

$E^{D_b}\left(r_j\right)=\underset{D_b}{\∫}{\hat{G}}_E\left(r_j\right|r)\·\left(\mathrm{\Δ}{\mathrm{\σ}}_{b}E\left(r\right)\right)\mathrm{dv}---\left(7\right)$

$H^{D_b}\left(r_j\right)=\underset{D_b}{\∫}{\hat{G}}_H\left(r_j\right|r)\·\left(\mathrm{\Δ}{\mathrm{\σ}}_{b}E\left(r\right)\right)\mathrm{dv}---\left(8\right)$

其中：

和

分别为电场和磁场的格林算子，Δσ_b＝σ_b-σ_n已知区域电阻率值与相应的背景电阻率值之差，E(r)为任意一点的电场。

4)采用积分方程法计算异常区域对应的电磁场，与已经计算过的背景场，合并得到总的电磁场；

根据2)确定的未知区域及电导率参数分别计算公式(5)、(6)中的第二个积分项，得到相应的未知区域的电磁场。与3)计算得到的背景电磁场合并得到电场(公式(5)的结果)和总的磁场(公式(6)的结果)。

由三维电磁积分公式计算目标区域(D_a)的E^Da和磁场H^Da：

$E^{D_a}\left(r_j\right)=\underset{D_a}{\∫}{\hat{G}}_E\left(r_j\right|r)\·\left(\mathrm{\Δ}{\mathrm{\σ}}_{a}E\left(r\right)\right)\mathrm{dv}---\left(9\right)$

$H^{D_a}\left(r_j\right)=\underset{D_a}{\∫}{\hat{G}}_H\left(r_j\right|r)\·\left(\mathrm{\Δ}{\mathrm{\σ}}_{a}E\left(r\right)\right)\mathrm{dv}---\left(10\right)$

和

分别为电场和磁场的格林算子，Δσ_a＝σ_a-σ_n目标区域电阻率值与背景电阻率值之差，E(r)为任意一点的电场。

5)根据上述3)、4)的计算结果计算地面上的视电阻率及相位相应；

6)根据目标区域大小、埋深、规模及由5)计算的目标区域对应的视电阻率曲线和频率相位曲线，确定5～20个频点的视电阻率曲线及相位曲线进行反演；

实际MT及CEMP方法应用于石油勘探的频率范围一般为320～0.00055Hz，频带范围较宽，当未知区域规模较小时，未知区域电导率参数只与几个或十几个频点的视电导率曲线或相位曲线有关，因此只需选择与未知区域有关的频段范围进行反演。

如果目标区域深度为4000-8000m，则根据深度经验公式 $D=356\sqrt{\frac{1}{\mathrm{\σf}}}$ (其中D表示勘探深度，σ为地下介质平均电阻率，f表示频率)反推频率范围。如果假设σ值为0.1，则对应的频率范围0.079-0.02Hz，(而实际的频率范围一般为320-0.00055Hz)，一般在此频率范围内选择5-10个频点，小于或大于此频率范围分别延长3-5个频点，作为反演频率范围。

5～20个频点的视电阻率曲线及相位曲线频率范围一定大于该目标区域对应的频率范围。

7)采用非线性共轭梯度法，针对异常区域电导率值参数进行三维迭代反演，得到未知区域的电导率参数；

非线性共轭梯度法，是目前比较流行的一种非线性的反演方法。在反演过程中，第一步修改未知区域的电导率参数，计算异常场；第二步利用已计算过的背景场及未知区域的异常场，计算视电导率曲线和相位曲线；第三步计算视电导率曲线及相位曲线拟合误差；重复上述三步直到拟合误差满足要求，得到未知区域的电导率参数。

8)根据背景、已知区域及通过三维反演得到的异常区域的地下电导率参数组合建立电磁三维电磁反演数据体，作为地质解释的依据。

根据背景电导率参数(σ_n)、已知区域电导率(σ_b)以及通过三维反演得到电导率参数(σ_a)得到整个地下介质的电导率分布情况，建立电磁三维反演数据体，作为地质解释的依据。

Claims (4)

1、一种目标最小化的三维电磁快速反演方法，其特征采用如下步骤：

1)在开展过地震、地质或钻井工作的地区，实测目标勘探工作地区视导电率曲线及相位曲线，采用通常的方法进行三维插值，使之均匀分布于剖分单元的网格节点；

2)将电磁勘探工作区分为三部分：①地下介质已知的背景区域及导电率值(σ_n)，②相邻的已知区域(D_b)及相应的导电率参数(σ_b)，③由勘探任务确定的要进行反演的目标区域(D_a)，其相应的导电率值参数(σ_a)由邻区已知电测井资料、地质露头资料确定；

3)采用通常的电磁一维解析公式计算已知的背景区域及导电率值(σ_n)有关的电场(Eⁿ)和磁场(Hⁿ)，采用公式(7)(8)三维电磁积分公式计算已知区域(D_b)有关的电场E^Db和磁场H^Db：

$E^{D_b}\left(r_j\right)=\underset{D_b}{\∫}{\hat{G}}_E\left(r_j\right|r)\·\left(\mathrm{\Δ}{\mathrm{\σ}}_{b}E\left(r\right)\right)\mathrm{dv}---\left(7\right)$

$H^{D_b}\left(r_j\right)=\underset{D_b}{\∫}{\hat{G}}_H\left(r_j\right|r)\·\left(\mathrm{\Δ}{\mathrm{\σ}}_{b}E\left(r\right)\right)\mathrm{dv}---\left(8\right)$

E^Db和H^Db：分别表示与已知区域有关的电场和磁场其中：

和

分别为电场和磁场的格林算子，Δσ_b＝σ_b-σ_n已知区域导电率值与相应的背景导电率值之差，E(r)为任意一点的电场；

4)采用如下三维电磁积分公式计算目标区域(D_a)的E^Da和磁场H^Da

$E^{D_a}\left(r_j\right)=\underset{D_a}{\∫}{\hat{G}}_E\left(r_j\right|r)\·\left(\mathrm{\Δ}{\mathrm{\σ}}_{a}E\left(r\right)\right)\mathrm{dv}---\left(9\right)$

$H^{D_a}\left(r_j\right)=\underset{D_a}{\∫}{\hat{G}}_H\left(r_j\right|r)\·\left(\mathrm{\Δ}{\mathrm{\σ}}_{a}E\left(r\right)\right)\mathrm{dv}---\left(10\right)$

E^Da和H^Da：分别表示与目标区域有关的电场和磁场，和

分别为电场和磁场的格林算子，Δσ_a＝σ_a-σ_n目标区域导电率值与背景导电率值之差，E(r)为任意一点的电场；

5)根据上述3)、4)的计算结果计算地面上的视导电率及相位相应；

6)根据目标区域大小、埋深、规模及计算的目标区域对应的视导电率曲线和频率相位曲线，确定5～20个频点的视导电率曲线及相位曲线进行反演，所述的5～20个频点的视导电率曲线及相位曲线频率范围一定大于该目标区域对应的频率范围；

7)采用常规的非线性共轭梯度法，针对目标区域导电率值参数进行三维迭代反演，得到导电率参数；

8)根据背景导电率参数(σ_n)、已知区域导电率(σ_b)以及通过三维反演得到导电率参数(σ_a)得到整个地下介质的导电率分布情况，建立电磁三维反演数据体，形成地质解释的依据。

2、根据权利要求1所述的目标最小化的三维电磁快速反演方法，其特征在于：所述的(σ_a)(σ_b)代表区域内每一个单元的导电率值，根据目标地区反演的深度范围内已有钻井资料确定其值。

3、根据权利要求1或2所述的目标最小化的三维电磁快速反演方法，其特征在于：所述的所述的(σ_a)(σ_b)代表区域内每一个单元的导电率值用一个任意导电率值充填。

4、根据权利要求1所述的目标最小化的三维电磁快速反演方法，其特征在于：所述的非线性共轭梯度反演过程是，第一步修改未知区域的导电率参数，计算异常场；第二步利用已计算过的背景场及未知区域的异常场，计算视导电率曲线和相位曲线；第三步计算视导电率曲线及相位曲线拟合误差；重复上述三步直到拟合误差满足要求，得到未知区域的导电率参数。

Images