刘天佑，杨宇山，李媛媛

刘天佑，杨宇山，李媛媛_大人地球物理科学

刘天佑，杨宇山，李媛媛

（中国地质大学地球物理与空间信息学院，武汉，430074）

摘　要：介绍了20世纪80年代以来笔者追踪国际地球物理反演的进展，在重磁勘探反演方面所做的一些工作：①根据Backus-Gilbert反演理论与广义逆矩阵反演方法，实现了复杂磁异常的广义逆矩阵反演与居里等温面反演；②在反演的目标函数的改造优化方面，实现了频率域最优化反演与约束反演；利用小波多尺度分析方法实现了重磁异常多尺度反演；利用小波压缩技术求解大型积分方程，实现了重磁异常的曲化平；③研究了重磁、重震、井地磁测、井中磁测多分量的联合自动反演与交互反演方法；采用重震联合反演多层密度界面并应用于我国西部与南方海相油气勘探；采用2.5D、3D人机交互反演方法实现井地磁测联合反演、地面与井中磁测多分量的联合自动反演，并在我国危机矿山深部找矿中发挥了作用；④在非线性反演方面，根据BP人工神经网络原理，实现了密度界面反演方法，将群智能反演方法（粒子群、蚁群算法等）用于井地磁测的联合反演。这些工作扩大了重磁勘探方法的应用领域，并提高了重磁勘探方法的应用效果。

关键词：重磁勘探，Backus-Gilbert反演，广义逆矩阵反演，联合反演，群智能反演

前　言

反演问题是地球物理的核心问题。20世纪60～70年代，在西方Backus G E和Gilbert J F从连续模型出发，建立了地球物理反演理论；Wiggens R A，Jackson D和Franklin J N从离散模型出发，提出了广义线性反演方法，利用广义逆矩阵理论，使得Backus-Gilbert反演理论得到广泛的应用，成为近30年来地球物理反演卓有成效的理论成果。而在同期的前苏联，ТихоновАН提出了解不适定问题的正则化方法；СтраховВН，ЦирульскийАВ从复变函数理论出发，发展了重磁正反演理论；ГольцманФМ运用统计学方法，建立了重磁异常的统计解释理论。这些从不同角度，不同体系进行研究所获得的成果，极大地丰富了重磁反演的理论与方法。Backus-Gilbert反演理论和Тихонов正则化方法的贡献在于它们指出了地球物理观测数据是有限的和包含误差的，对于这类不适定问题应当用正则化方法解决，对获得的解也应当进行评价，这些理论与方法推动了我国20世纪80～90年代重磁反演的发展。

20世纪90年代以来，非线性反演方法，如小波分析、人工神经网络、模拟退火算法、遗传算法在重磁反演中得到较广泛的应用；而进入21世纪，一类基于模拟生物群体行为的智能优化算法（群智能算法）诞生了，如粒子群算法，源于对简单社会系统的模拟；蚁群算法，源于对蚂蚁群体觅食行为的模拟，它们通过个体与群体间信息交换与分工协作，调整群体及自身的运动方向及路径，体现了一种对环境的适应能力。这些方法也开始在重磁反演中应用。而随着我国经济发展对矿产资源的需求，为了适应危机矿山深部找矿、整装勘探与立体填图的需要，可视化反演、联合反演方法与物性反演得到了迅速的发展，并为我国固体矿产与油气勘探作出了贡献。

我们的主要工作是实现了重磁异常Backus-Gilbert反演与广义逆矩阵反演；实现了频率域最优化反演与约束反演；小波多尺度反演与小波压缩反演；实现了重磁、重震、井地磁测、井中磁测多分量的联合自动反演与交互反演方法；实现了群智能反演方法等。将这些方法应用于我国西部与南方海相油气勘探，以及危机矿山深部找矿，扩大与提高了重磁勘探方法的应用领域与应用效果。

1 Backus-Gilbert反演理论与广义逆矩阵反演

20世纪70年代末以来，笔者在国内最早开始Backus-Gilbert反演理论与Тихонов正则化方法在重磁异常反演方面的研究，利用广义逆矩阵的选择法反演复杂模型磁异常，并应用于沉积变质铁矿的解释^［1～3］，利用Backus-Gilbert反演理论与Тихонов正则化方法反演磁性体上、下界面^［4～6］。

求解重磁反问题，归结为解以下的非线性矛盾方程组：

在（x）^（0）点展为泰勒级数并略去高次项，最后可以表示为：

式（2）的广义逆解为：

式中：A＝UΛV^T称为A的奇异值分解。给出奇异值截断标准，使对角矩阵Λ^-1中大于λp的取为1/λi，小于λp的取为零，上述方法称奇异值截断法。除了采用该方法以外，还可以采用在对角矩阵加阻尼因子的方法。

图1是对一个向斜模型采用阻尼最小二乘法、变尺度法与广义逆矩阵法反演结果的对比，广义逆矩阵法完全收敛，而其余两种方法不收敛。

2　对目标函数改造的反演方法（小波多尺度反演与频率域最优化反演）^{［1，7～10］}

待求解的未知参数多、目标函数维数大，及目标函数多极值是影响最优化反演方法不收敛的原因。20世纪70年代末以来，笔者通过对目标函数性质的分析，改造目标函数，实现了频率域磁异常振幅谱与相位谱的反演，该方法降低了参数空间的维数，并且与磁化方向无关；利用线性规划中多阶段决策理论，通过对目标函数的小波多尺度分解，实现了多尺度反演，克服了广义线性反演方法局部寻优的困难^［1，7^，9］。

图1　向斜模型阻尼最小二乘法、变尺度法与广义逆矩阵法反演结果对比

高斯法、最速下降法、阻尼最小二乘法等是反演二度板状体、三度棱柱体等地质体重磁异常常用的最优化方法。由于目标函数的复杂性，这类方法往往不收敛。经过分析可知，阻尼最小二乘法等用来求解复杂重磁异常不收敛的原因是由于目标函数不是简单的凸函数；目标函数的复杂性不仅与参数空间的维数有关，而且与参与计算的采样点有关。

20世纪90年代以来，小波分析方法在地球物理数据处理方面得到广泛应用。小波的多尺度分析是Mallat根据多尺度分析理论发展起来的一种方法。通常，图像的细节以不同的分辨率描述了图像的不同物理结构。在粗的分辨率下，这些细节对应了较大的结构，它们提供图像粗的轮廓。因此，可以采用由粗到细的策略来寻求目标函数真正的极小值来解决了直接用复杂目标函数求极小不容易收敛的问题。

应用非线性规划中的多阶段决策思想，把求解复杂重磁异常最优化问题分解为多个阶段决策问题。为了保证寻找极小值的成功，通过小波多尺度分解把目标函数分解成不同阶的目标函数，由于从低阶到高阶目标函数由简单变复杂，先用低阶小波分解的目标函数粗略找出其极小值的范围，再逐次改用高阶小波分解的目标函数缩小极小值的范围，最终找到目标函数极小值。

采用小波多尺度分析方法来实现重磁异常的最优化反演的方法原理如下，重磁异常反演的最优化选择法可以归结为对如下目标函数求极小值问题。

式中：ΔT_k（k＝1，2，…，m）是观测值；f（x，b₁，b₂，…，b_n）是n个未知参数b₁，b₂，…，b_n产生的模型理论值。

多级系统的分步最优化常用的方法是Bellman（1957、1962）提出的“动态规划”，它采用多阶段决策的方法来实现序列系统的决策最优化。我们这里的最优化反演过程的多阶段决策与Bellman的动态规划有所区别，我们利用多级系统的分步最优化思想，把复杂的重磁异常最优化问题分解为多阶段决策问题，每一阶段目标函数是通过小波多尺度分解来实现的。其具体方法如下：

（1）对目标函数的非线性问题线性化，求Jacobi矩阵及残差。

（2）对Jacobi矩阵及残差向量进行多尺度分解，并形成新的核函数矩阵和残差向量。

（3）用奇异值分解法求取模型修正量。

（4）判断所得结果是否满足这一分解尺度精度要求。若满足则完成这一尺度反演，否则返回（1）。

（5）判断是否完成反演；否则以第四步得到的模型作为初始模型转第一步。

下面是理论模型的试算结果：水平方向上相互叠加的两个下延有限倾斜板状体模型，板状体大小、位置、磁性参数、初值及反演结果如图2所示。采用Harr小波，最大尺度为4，迭代精度为10^-4，待反演未知参数有14个，即两个模型的参数：板状体上顶中心位置坐标（X_i，Z_i），板宽度2b_i，下延深度2l_i，有效磁化倾角i_si，有效磁化强度Msi，板状体倾角α_i（i＝1，2）。

图2　板状体模型不同方法反演结果

（1）理论模型和多尺度反演结果；（2）初始模型；（3）阻尼最小二乘法反演结果；（4）模型理论值及多尺度反演拟合值；（5）初始模型产生的异常曲线

由图2可以看出，尽管理论模型的磁异常曲线不具明显双峰值异常特征，但多尺度反演结果很好地收敛到理论模型的位置。

3　联合反演与约束反演^{［1，2，11～18］}

多种地球物理资料的联合反演是减少多解性的重要措施。20世纪80年代以来，我们研究了重磁、重震、井地磁测、井中磁测多分量联合自动反演与交互反演方法，并采用重震联合反演多层密度界面应用于我国西部与南方海相油气勘探；采用惩罚函数法对反演参数进行约束，充分利用了已知信息；采用2.5D、3D人机交互反演方法实现井地磁测联合反演；地面与井中磁测多分量的联合自动反演；在我国危机矿山深部找矿中发挥了作用。

3.1　重磁联合反演

在地球物理资料处理与解释中，我们常常会遇到所谓“重磁同源”的情况，如盆地的结晶基底，既是高密度，又是强磁性；又如火成岩侵入体，通常是强磁性，也可能是高密度的；固体矿产中的铁矿床，更是典型的高密度、强磁性。对于这类重磁同源问题，以往人们用对应分析、相关分析方法计算磁源重力异常方法进行定性分析与解释，重磁联合反演为我们定量解释这类地质问题提供了新途径。

我们在“八五”国家科技攻关项目“新疆塔里木盆地北部综合物探评价圈闭的方法技术研究”项目中，曾采用重磁联合反演和交互反演方法研究塔里木盆地北部含油气构造圈闭。

重磁联合反演的目标函数表示为：

式中：W^G、W^T表示重力与磁法的权重系数；下标1、2分别表示观测值和模型理论值；带＊号表示用平均值进行归一化。

利用加阻尼因子的广义逆矩阵方法就可以求解式（5）的联合反演问题。

3.1.1　界面模型

如图3所示模型，密度界面剩余密度Δσ＝0.1g/cm³，磁性界面磁化强度M＝1 000×10^-3 A/m，磁化倾角I＝90°，测线磁方位角A′＝0°，界面平均深度3.2km，剖面点数m＝25，初始模型参数n＝15，每个模型的初始深度均取z＝3.2km，反演中取重力与磁法的权系数相同，经过5次迭代，得到如图结果。图下部的虚线为理论模型的深度，实线为联合反演的结果，重力资料与磁法资料的联合反演，减少了反演的多解性，提高了反演结果的可靠性。

3.1.2　孤立地质体模型

采用任意多边形截面二度水平棱柱体作为地质模型，此类模型适合于反演孤立的地质体，如火成岩体、矿体等等。

图4是截面为任意六边形的水平棱柱体不同资料的反演结果对比。磁化强度M＝2 000×10^-3 A/m，剖面内有效磁化倾角i_s＝70°，剩余密度Δσ＝1.0g/cm³。由反演结果可知，联合反演迭代19次，基本上收敛于真值。而单一磁法资料反演迭代27次，下部有两个角点收敛不太好。单一重力资料反演则需迭代37次，收敛也不好，由此可以看出联合反演的效果。可以分析，对于有误差的重力和磁法数据，联合反演能够互相弥补数据中的缺陷，预计会获得更好的效果。

图3　重磁联合反演单一密度与磁性界面图

图4　不同资料反演结果对比图

3.2　重震联合反演

重震联合广义反演的模型为N个共同具有密度差异与波阻抗差异的多层界面。设G为由m_g个数据组成的重力观测值，S为m_s个地震旅行时数据，G^f为模型正演重力值，S^f为模型正演地震数据，由最小二乘原理知，求模型参数向量，使

与

最小，则重震联合反演可表示为

式中：W_G、W_S为相对加权矩阵，由于重力观测值和地震旅行时是两类不同量纲的数据，因此对G、G^f，S、S^f数据必须作归一化处理。根据非线性问题的广义线性反演方法，在给定的初始模型m°附近作线性展开，将非线性问题变为迭代求解线性反问题：

式中：X为模型参数修改向量；A^G和A^S分别为模型响应对参数的偏导数矩阵；BG和BS为重力残差向量和地震残差向量。

式（9）中参数修改量X为：

重力模型为二度半任意多边形截面水平棱柱体。

图5中（a）是给定的地下地质模型，（b）是进行重震联合广义反演的结果，反演结果与实际界面的较一致。可以看出，对于多层界面的模型，层数越多，反演的多解性越严重。越是浅部的地层的界面，反演的精度越高；越是深部的地层界面的反演精度就越低。可以说，层数越多对初始模型的要求越严，或者说初始模型的选取对反演结果的影响越大。我们这里是采用零偏移距地震反射旅行时与重力资料进行联合反演，这种特征尤为突出，若改用非零偏移距反射旅行时与重力进行联合反演，这种特性将相对变弱。

图5　重震联合反演结果

3.3　井地联合反演^{［19～21］}

3.3.1　人机交互反演

老矿山人文干扰严重，它们“屏蔽”了深部的有用信息，制约了地面高精度磁测方法深部找矿的效果。而井中磁测能够避开地面及浅部磁性体的干扰，但是井中磁测受钻井的限制，控制的范围有限。如果将地面高精度磁测资料和井中磁测资料结合起来，进行联合反演解释，势必能够发挥两种方法的优点，达到优势的互补。

我们采用任意形状三度体模型，数值积分法近似计算，在Windows环境下，用Visual C语言，OpenGL函数编制了SWMI3D软件，实现了井地磁测资料人机交互反演，并以大冶铁矿的实例说明了方法的有效性。

对于任意形态三度体磁场的计算可以采用数值积分方法。Talwani M（1965）用面元、线元与点元法实现了不规则形体数值积分方法。我们用两组相互垂直的截面把任意三度体分割成许多小棱柱体，每个棱柱体相当于一个直立线元。沿Z轴用解析方法实现一重积分，求出各线元的作用值，然后在垂直线元的X、Y方向分别作数值积分，即可得出整个地质体的近似磁场作用值。

所谓“可视化”反演是指实际观测磁异常曲线、拟合的地质模型及其正演计算的磁异常曲线在计算机屏幕上始终以图形方式显示，我们可以直接对地质模型进行操作（修改、反演），实时地计算修改后的地质模型所产生的磁异常曲线并与实际观测的磁异常曲线进行比较，通过地质模型形态、物性参数的不断修改，使得正演计算磁异常曲线与实际观测值的磁异常曲线逐步拟合以达到对地下地质体的反演解释的目的。

可视化人机交互反演包括如下几方面：

（1）采用任意三度体模型，其正演计算用数值积分法，即用辛普生积分和梯形积分实现三度体磁场三重积分的近似计算。

（2）为了便于修改模型，修改的过程是在剖面内完成，对x、y不同方向剖面逐条修改拟合，每一次修改拟合的正演计算都采用数值积分对任意三度模型进行计算。

（3）初始模型由已知的勘探线所控制的矿体、围岩形成，也可以根据先验的信息解释人员直接在屏幕前用鼠标建立。

（4）在Windows环境下，用Visual C语言，OpenGL函数实现立体模型与平面组合模型的旋转、移动、放大、缩小，以及任意选择剖面、断面进行精细反演解释。

可视化交互反演界面，如图6所示。

图6　系统可视化交互反演界面

大冶铁矿属接触交代型铁矿，接触带整体呈“S”形展布，矿体沿接触带分布。前人对ZK19-1-15井中三分量磁测资料采用传统的半定量方法推断在其标高-220～-320m之间和-400～-600m之间分别有一旁侧矿体（Fe4和Fe3），但没有确定其空间位置，并认为Fe3矿体在19-1剖面的旁侧。旁侧矿体（Fe4和Fe3）是否存在？Fe3矿体向东延伸，是否到19-1剖面终止？为了深入解释，我们采用3D任意形状地质体井-地磁测资料联合反演方法，对该矿段的6条测线（17线、18线、19-0线、19-1线、21线和22线）和3个钻孔（ZK18-12、ZK19-1-15和ZK22-16）的资料进行解释，如图7所示。通过对相邻的18、21、22线的分析发现，在与Fe3对应的标高位置，各测线均有一铁矿存在，而且铁矿的形状、产状均与Fe3非常相似。考虑到从18线到22线，接触带发生“S”形扭曲，对应标高位置的铁矿也随之发生扭曲。标高-400～-600m之间的矿体在走向上是连接在一起的，也就是说Fe3矿体实际上是穿过了19-1剖面，而19-1-15孔只不过是从Fe3的旁边穿过，并没有打到Fe3。2008年布置钻孔19-1-17孔，在井深638.32～689.40m处打到一厚层铁矿，视厚度约40m，证实了上述推断。图8是ZK19-1-15孔井地联合反演井中曲线拟合结果，可见，添加了Fe3矿体，并考虑它的空间位置的井中正演曲线较好地拟合了实测曲线。

图7　大冶铁矿某矿段测线及钻孔分布图

（底图：ΔZ为磁异常；实线为2.5D反演解释剖面）

图8　大冶铁矿ZK19-1-15孔井地联合反演井中曲线拟合结果

3.3.2　井地联合成像^［22］

根据磁异常反演磁性体的问题转为求解磁化率在地下三维空间的分布，这种反演方法称为视磁化率成像法。它可以将地下三维空间划分为诸多小的长方体磁性模型，通过最优化反演方法，在空间域对所建立的目标函数进行求解，最终得到地下三维空间的磁性分布情况。

图9（a）、（b）、（c）分别是地面磁测资料单独反演、井中磁测资料单独反演与3D井地联合反演的结果的对比。从图中的结果说明，3D井地磁测资料的联合反演，充分利用了地面磁测资料横向分辨率高和井中磁测资料纵向分辨率高的优点，反演出的结果最接近理论模型，要优于单一的地面或井中磁测资料的反演结果。

图9　井地联合反演结果对比

（图中的粉色线框为理论模型的实际位置）

山东某铁矿为接触交代型铁矿，矿体埋深超过1 000m，单独利用地面磁测资料无法解释深部矿体。图10是利用地面磁测资料与ZK18-11、ZK18-9和ZK18-17的井中磁测资料联合反演的结果。可以看到，磁化率成像的结果与钻遇的条带状磁铁矿体产状一致。

图10　山东某矿18线井地联合成像

3.4　井中磁测三分量联合反演^［23］

井中磁测是在钻孔中测量磁场的3个分量，近年来，井中三分量磁测井仪器的定位精度与分辨力有了大幅的提高，为井中三分量的反演解释提供了条件。

三分量磁测数据H_ax、H_ay、Z_a联合反演的目标函数：

式中：分别为三分量磁测H_ax、H_ay、Z_a的权重系数矩阵。求目标函数E的极小值，得到方程（11）的最小方差解

式中：，J为对称阵。

由于系数矩阵J存在小的特征值或者特征值为零，即方程是病态或者奇异的，这造成解的不稳定或者解的不存在。改善解的稳定性的有效办法是在矩阵J上加一阻尼因子λ²I，加阻尼因子后的预测解为：

图11为井底异常，曲线不完整。理论模型为直立长方体，深度500m，对数据均加50nT的高斯白噪声。分别采用三分量H_ax、H_ay、Z_a联合反演与Z_a单分量反演对比如表1所示，三分量反演结果如图11所示。在测井曲线不完整且加噪声的情况下，Z_a单独反演结果发散。而三分量联合反演结果可以大致确定模型体的位置和大小。

图11　三分量联合反演结果（剖面位置y＝550m）

其中：红、蓝、绿3条实线分别为H_ax、H_ay、Z_a的“观测值”；红、蓝、绿3条虚线分别为H_ax、H_ay、Z_a的反演结果；黑色线框为理论模型，虚线黑线框为反演结果，粉线为钻孔

表1　井底异常、一个钻孔、有噪声情况反演结果

4　非线性反演方法^{［20，22～27］}

由于非线性反演具有全局寻优的特点，20世纪90年代以来，开始在地球物理反演中得到应用。我们根据BP人工神经网络原理，实现了密度界面反演方法^［8］，将群智能反演方法（粒子群、蚁群算法等）用于井地磁测的联合反演。

粒子群算法（Particle　Swarm　Optimization，简称PSO）是继遗传算法、蚁群算法之后的又一种新的群体智能算法^［1］，它是以模拟鸟群集智能为特征，以求解连续变量优化问题为背景的一种优化算法。

PSO源于对鸟群捕食行为的研究，一群鸟在随机搜寻食物，如果这个区域里只有一块食物，那么找到食物的最简单有效的策略就是搜寻目前离食物最近的鸟的周围区域。PSO求解优化问题时，问题的解对应于搜索空间中鸟的位置，称这些鸟为“粒子”（particle）或“主体”（agent）。每个粒子都有自己的位置和速度（决定飞行的方向和距离），还有一个由被优化函数决定的适应值。各个粒子记忆、追随当前的最优粒子，在解空间中搜索，在每一次迭代中，粒子通过跟踪两个“极值”来更新自己：第一个就是粒子本身所找到的最好解，叫做个体极值点（用pbest表示其位置），另一个极值点是整个种群目前找到的最好解，称为全局极值点（用gbest表示其位置）。在找到这两个最好解后，粒子根据一定的方程来更新自己的速度和位置。

上述思想可抽象为，在一个D维的目标搜索空间中，由m个粒子组成一个群落，其中第i个粒子表示为一个D维的向量x^→_i＝（x_i1，x_i2，…，x_iD），i＝1，2，…，m，即第i个粒子在D维空间的位置。将^-x_i带入一个目标函数就可以计算出其适应值，根据适应值的大小衡量^-x_i的优劣。第i个粒子的“飞翔”速度也是一个D维的向量，记为v^→_i＝（v_i1，v_i2，…，v_iD）。

记第i个粒子迄今为止搜索到的最优位置（pbest）为p^→_i＝（p_i1，p_i2，…，p_iD），整个粒子群迄今为止搜索到的最优位置（gbest）为_p^→g＝（p_g1，p_g2，…，p_gD），PSO算法采用下面的公式对粒子的速度和位置更新，方程为

式中：i＝1，2，…，m；d＝1，2，…，D；惯性权重w^{［3，4，7］}是用来控制粒子以前速度对当前速度的影响，它将影响粒子的全局和局部搜索能力。较大的w值有利于全局搜索，较小的w值有利于局部搜索。选择一个合适的w值可以平衡全局和局搜索能力，这样可以以最少的叠代次数找到最优解。对PSO全局和局部之搜索能力的平衡主要被惯性权重控制；学习因子（又称加速因子）c₁和c₂是非负常数，分别调节向全局最好粒子和个体最好粒子方向飞行的最大步长，若太小，粒子可能远离目标区域，若太大则会导致突然向目标区域飞去，或飞过目标区域。合适的c₁，c₂可以加快收敛且不易陷入局部最优，通常c₁＝c₂＝2；r₁和r₂是［0，1］之间的随机数，v_id∈［-v_dmax，＋v_dmax］，即粒子速度被一个最大速度所限制。假设搜索空间第d维定义为［-x_dmax，＋x_dmax］，则通常v_dmax＝kx_dmax，0.1≤k≤1.0，粒子规模m不用取太大，几十个就足够用了，有时候还要视具体情况和问题来定，粒子的维度D，是由优化问题所决定的，也就是解空间的维度。

迭代中止条件根据具体问题一般选为达到最大迭代次数或（和）粒子群搜索到的全局最优解所满足误差的要求。

式（14）的第一部分为粒子先前的速度；第二部分为“认知”部分，表示粒子本身的思考；第三部分为“社会”部分，表示粒子间的信息共享与相互合作。“认知”部分可解释为一个得到加强的随机行为在将来更有可能出现。这里的行为即“认知”，并假设获得正确的知识是得到加强的，这样一个模型假定粒子被激励着去减小误差。“社会”部分可解释为粒子本身是被其他粒子所模仿，即作为一个种群中粒子运动的参考。

图12是两个板状体模型产生的双峰值的叠加异常，通常的最优化反演方法，如阻尼最小二乘法无法收敛。利用粒子群算法，种群规模，当对反演精度要求不高时取40～80，当精度要求较高时取N＝120，迭代500次，得到2个板状体的参数与理论模型的对比见表2和图12，可见，理论模型与反演结果较为一致，说明非线性反演方法由于全局寻优，在多峰异常反演中的效果。

表2　PSO算法反演板状体模型参数结果

图12　PSO算法磁异常反演结果

5　结　论

重磁勘探是一种利用天然场的方法，其纵向分辨率低，反演解释具多解性，深入研究重磁异常反演方法具有重要意义。笔者追踪了20世纪80年代以来国际重磁异常反演的进展，实现了Backus-Gilbert反演与广义逆矩阵反演、频率域最优化反演与约束反演、小波多尺度反演与压缩反演；实现了重磁、重震、井地磁测、井中磁测多分量的联合自动反演与交互反演^{［25～27］}。并将这些方法应用于我国油气勘探与危机矿山深部找矿，提高了重磁勘探方法的应用效果。

致谢：参加本文重磁反演方法与应用研究的有博士、硕士研究生宣霖康、方晓梅、张世晖、耿喜哲、吴招才、魏伟、冯杰、李曼、习宇飞、曾琴琴、刘双、张大莲等。

参考文献

［1］刘天佑.重磁异常反演理论与方法［M］.武汉：中国地质大学出版社，1992.

［2］刘天佑.位场勘探数据处理新方法［M］.北京：科学出版社，2007.

［3］刘天佑，管志宁.对二维磁异常利用广义逆矩阵的选择法［J］.地球物理学报，1986，29（4）：390-398.

［4］刘天佑.利用广义逆矩阵的窗口系数法反演重磁界面［J］.地球科学——中国地质大学学报，1989，14（6）：645-656.

［5］刘天佑，宣霖康.利用巴克斯-吉尔伯特反演理论计算二维磁性体下界面［J］.现代地质，1989，4（3）：454-463.

［6］刘天佑，方晓梅.一种连续模型的居里面反演方法［J］.地球科学——中国地质大学学报，1987，12（6）：647-656.

［7］管志宁，刘天佑.二维磁异常频率域最优化反演方法［J］.物化探计算技术，1985，（4）：300-309.

［8］耿喜哲，刘天佑，丁艳红.小波神经网络在重磁资料反演中的应用前景［J］.物探与化探，2001，25（2）：102-108.

［9］刘天佑，朱立新，邓荣来.复杂重磁异常多阶段决策的最优化反演［J］.石油地球物理勘探，2002，37（4）：349-353.

［10］Yang Yushan，Li Yuanyuan，Liu Tianyou.Continuous wavelet transform，theoretical aspects and application to aeromagnetic data at the Huanghua Depression，Dagang Oilfield，China［J］.Geophysical Prospecting，2010，58（4）：669-684.

［11］刘天佑，杨宇山，李媛媛.重力、地震联合反演方法研究深层油气藏的分布［J］.天然气工业，2005，25（5）：34-37.

［12］魏伟，刘天佑.梯度法解释复杂二维断裂重力异常［J］.物探与化探，2005，29（4）：347-351.

［13］魏伟，刘天佑.台阶重力异常的梯度解释［J］.石油地球物理勘探，2005，40（2）：238-244.

［14］刘天佑，刘大为.磁测资料处理新方法及在危机矿山挖潜中的应用［J］.物探与化探，2006，30（5）：377-381.

［15］杨宇山，刘天佑.任意形状地质体数值积分法重磁场三维可视化［J］.地质与勘探，2006，42（5）：79-83.

［16］魏伟，吴招才.基于AutoCAD平台三维可视化规则几何形体磁场反演［J］.工程地球物理学报，2006，3（1）：54-59.

［17］中国地质大学“八五”科技攻关项目“新疆塔里木盆地北部综合物探评价圈闭的技术方法研究”.

［18］刘天佑，朱铉.综合地球物理数据处理新方法在西部油气勘探中的应用［J］.勘探地球物理进展，2006，29（2）：104-108.

［19］刘双，张大莲，刘天佑，等.井地磁测资料联合反演及应用［J］.地质与勘探，2008，44（6）：69-72.

［20］曾琴琴.位场资料群智能反演方法研究［D］.武汉：中国地质大学博士学位论文，2011.

［21］Yang Yushan，Li Yuanyuan，Liu Tianyou，et al.Interactive 3Dforward modeling of total field surface and three-component borehole magnetic data for the Daye ironore deposit（Central China）［J］.Journal of Applied Geophysics，2011，75：254-263.

［22］冯杰.井地磁测联合反演研究［D］.武汉：中国地质大学博士学位论文，2010.

［23］习宇飞.磁场三分量联合反演研究［D］.武汉：中国地质大学博士学位论文，2011.

［24］张大莲，刘天佑，陈石羡，詹应林.粒子群算法在磁测资料井地联合反演中的应用［J］.物探与化探，2009，33（5）：571-575.

［25］张大莲，刘双，陶德益，杨坤彪.井中磁测与地面磁测资料联合反演［J］.工程地球物理学报，2008，5（1）：60-64.

［26］张大莲，刘天佑.混合最优化算法在磁测资料井地联合反演中的应用［J］.第四届环境与工程地球物理国际会议论文集，2008，825-828.

［27］冯杰，刘天佑，杨宇山.3D井地磁测联合反演技术及其在危机矿山找矿中的应用［J］.武汉大学学报，信息科学版，2010，35（12）：1 436-1 439.