2018, Vol. 61
三维反演(Li and Oldenburg, 1996，2003；Pilkington, 1997, 2009；Portniaguine and Zhdanov, 2002；姚长利等, 2003, 2007)是磁测数据定量解释的重要方法，在金属矿勘探中扮演着重要的角色.但是在实际矿区的应用中，传统的磁总场异常反演方法依然存在两个问题：一是地面磁异常反演的深度分辨率较低，由于磁异常随距离的三次方衰减，因此对于同一个场源体，随埋深的增加，其磁异常迅速衰减，而对应反演结果中，则表现为深部场源体的成像效果差.二是金属矿(及含矿岩石)中常常包含强剩磁，此时矿体的磁化方向可能与当地地磁场方向偏差很大，此时若简单将地磁场方向当作总磁化方向来进行反演，反演结果会出现较大偏差甚至是完全错误的.针对上述两个问题，国内外学者已分别进行了大量研究.
在改善反演深度分辨率方面，主要的途径是引入额外的地质或地球物理信息.常用的手段有三种.第一种是通过添加物性信息来改善反演结果.已知的物性信息可以作为参考模型或者等式约束(Aster et al., 2005；Menke，2012)加入到反演的目标函数中.此外，物性信息也可以作为边界(即物性上下界)约束(Last and Kubik, 1983；Guillen and Menichetti, 1984；Li and Oldenburg, 1996，2003；Portniaguine and Zhdanov, 1999，2002)或者离散值约束(Camacho et al., 2000；Krahenbuhl and Li, 2006，2009；Uieda and Barbosa, 2012)来改善反演结果.第二种则是通过添加构造方向信息来改善反演结果.Guillen和Menichetti(1984)在反演中极小化场源体关于构造方向的转动惯量.Li和Oldenburg(2000a)利用旋转矩阵将模型目标函数中的空间导数旋转至任意方向，然后通过调整旋转后的导数的权重实现了简单地质体的构造倾向信息的加入.Lelièvre和Oldenburg(2009a)将Li和Oldenburg(2000a)的方法扩展至了可适用任意形态的地质构造的情况.Zhou等(2014)则实现了一种图像引导的反演，引导的图像可以是地质图或者高分辨率的地震数据，该方法通过从引导图中提取构造信息加入到反演中，来获得与引导图一致的反演结果.上述两种方法均可认为是引入了额外的地质信息，第三种则是引入额外的地球物理信息，通过多地球物理数据联合反演来改善反演结果.Lines等(1988)提出了一种序贯反演方法，该方法对每种数据独立反演，但在迭代的过程中通过共享信息来改善每种方法的反演结果.Li和Oldenburg(2000b)则利用广义加权函数将地面磁异常和井中三分量异常联系在一起，结合了地面磁异常横向分辨率高以及井中三分量数据纵向分辨率高的优点，实现了两者的联合反演.Gallardo和Meju(2004)通过交叉梯度方法将不同类型的地球物理数据联系起来，实现了不同类型数据的联合反演.上述三类手段也可以结合使用，来进一步改善反演结果(Guillen and Menichetti, 1984；Barbosa and Silva, 1994；Lelièvre et al., 2012).
在减少剩磁影响方面，主要的方法有两类.第一类方法仅利用磁异常中的振幅信息.该类方法首先将磁异常转换为受磁化方向影响小的非负转换量，然后反演这个转换量.Paine等(2001)将磁异常转换为两种转换量，然后将其假设为化极磁异常进行反演.Shearer(2005)、Li等(2010)、李泽林等(2015)利用磁异常模量的三维反演来减少反演中的剩磁影响.Pilkington和Beiki(2013)则反演归一化磁源强度来压制剩磁的影响.这类方法的优势是其适用于多场源复杂磁异常的反演，不足之处是只能反演出场源体的大概位置，而不能反演出其产状，主要原因是场源体的产状信息主要包含在磁异常的相位中，但转换量很大程度的消除了磁异常的相位信息(Li et al., 2010).第二类方法是磁化强度矢量反演，该方法试图最大程度上同时利用磁异常中的振幅和相位信息来恢复场源磁化强度的三分量.对于孤立场源，可以假设其磁化方向单一，首先估计场源的磁化方向(Phillips，2005；Dannemiller and Li, 2006；Gerovska et al., 2009)，然后将得到的磁化方向用于标准的三维反演中.对于多场源的情况，通常需要同时反演磁化强度的大小和方向.王妙月等(2004)提出了简单层状模型的二维磁化强度矢量反演.刘双等(2013)利用磁异常模量实现了井中二维磁化强度矢量反演.Liu等(2013, 2015)提出并改进了二维磁异常的磁化强度矢量反演方法.Lelièvre和Oldenburg(2009b)提出了直角坐标系和球坐标系下直接反演磁化强度矢量的三维算法，但是在不添加剩磁信息的情况下，这种方法反演结果中的磁化方向通常与真实值偏差较大.Li和Sun(2016)采用模糊C均值聚类来改善三维磁化强度矢量反演的效果，但该方法需要给定聚类的个数，若聚类个数不对，则可能产生错误的结果.总而言之，对于磁异常矢量反演而言，二维反演参数少(仅包含一个标量和一个方向)，效果较好，三维反演参数多(一个标量和两个方向)，需要添加特定的约束来减少反演参数增多带来的更加明显的非唯一性.
尽管前人已经对上述两个问题分别进行了比较深入的研究，但尚未尝试同时解决这两个问题.本文在前人研究的基础上，通过井地磁异常模量联合反演的方法，在减少剩磁影响的同时，也改善了深部成像效果.比较而言，该方法需要的控制参数少，无需加入额外的地质信息，且可用于多场源复杂磁异常的反演，具有较强的适用性.
本文首先回顾了常用的受磁化方向影响小的转换量，然后给出井地磁异常模量联合反演的目标函数及优化方法，最后通过模型试验对本方法进行了验证.
1 井地磁异常模量联合反演方法 1.1 受磁化方向影响小的转换量重力异常极大值往往对应场源的中心点，且形状规则的场源产生的重力异常往往也具有对称性.但是，对于磁异常而言，由于地磁场和磁性体磁化强度均存在方向性，这种方向性会改变磁异常的相位，进而改变磁异常的形态和幅值，使得磁异常的解释复杂化.当总磁化方向已知时，可以采用化极(Baranov，1957)来简化磁异常的解释，使得化极后的磁异常具有跟重力异常类似的特性.然而，当强剩磁存在时，总磁化方向未知，此时常规化极无法实现.
针对上述问题，前人已经做了一些研究.其中，解决该问题的一类方法是将磁异常转化为一个受磁化方向影响小的转换量，来达到类似化极的效果.常用的转换量有解析信号振幅(Nabighian，1972；Roest et al., 1992)、归一化磁源强度(Wilson，1985)和磁异常模量(侯重初，1979；Stavrev and Gerovska, 2000)等.
上述转换量在二维情况下均不受磁化方向的影响，在三维情况下则受磁化方向影响较小，与场源平面位置对应关系良好.理论上讲，这些转换量均可用于联合反演.本文选择磁异常模量进行联合反演，主要的原因是磁异常模量与磁异常的量级一致，在转换的过程中不会放大噪声，保持了长波长信息.
对于地面磁测而言，测量数据通常为总场异常ΔT，因此需要利用ΔT换算磁异常模量.对于中高纬度的平面网格数据，可通过频率域转换的方法将总场异常转化为异常三分量，然后投影合成得到模量数据.利用频率域方法换算三分量仅需已知地磁场方向，而不对磁异常做任何假设(即是否包含剩磁)，因此换算过程不会引入明显误差.对于其他情况下的观测数据，常规分量转换不易实现，可以采用空间域的等效源方法(Dampney，1969)得到模量数据.利用等效源方法换算三分量需要假设等效源的磁化方向，我们通常采用地磁场方向作为假设的磁化方向，由于利用等效源计算三分量时并不改变这个假设的磁化方向，因此用等效源方法换算三分量也不会引入明显误差.对于井中磁测而言，测量数据目前主要还是三分量数据，因此可以直接投影合成模量数据.
1.2 正演问题在反演中，我们将地下模型空间划分为三维长方体单元，反演的最终目标即为确定每个单元的物性.与常规反演不同，我们首先定义有效磁化率κ=μ0M/T0，其中，μ0=4π×10-7H/m为真空中的磁导率，M为有效磁化强度，单位为A/m，T0为地磁场强度，单位为T.如果假设观测的模量数据(包括地面和井中数据)为d，有效磁化率向量为m，正演问题可以表述如下：
|
(1) |
其中F为磁异常模量关于有效磁化率向量m的非线性正演算子.理论上讲，磁异常模量的正演需要正确的磁化方向，但真实的磁化方向是未知的，由于磁异常模量受磁化方向的影响小，因此在正演时，可以采用地磁场方向作为假设的磁化方向，这是模量反演中需要的(Li et al., 2010).显然，这个磁化方向的选择已经不再像常规磁异常分量反演那么重要了.此外，通常每一个井中磁测数据都会位于某一个物性单元的内部或边界位置，此时，常规的解析正演公式是不适用的.因此，需要对这种情况下的核函数进行特殊处理.我们采用的方法是将该数据的点位移动一小段距离，使得其位于对应物性单元的外部，来避免奇异性.对于每一个井中磁测数据，这种处理方法仅改变核函数中一个元素的值，对于最终的反演结果是几乎没有影响的.
1.3 反演问题磁测反演问题通常表述为约束优化问题，这里，我们参考Shearer(2005)和Li等(2010)在地面磁异常模量反演中选取的目标函数和约束项，将井地磁异常模量反演问题表述如下：
|
(2) |
其中μ为正则化参数.ϕd(m)为数据拟合差，其矩阵形式如下：
|
(3) |
其中Wd为对角数据加权矩阵，其第i个对角元素为第i个观测数据的噪声标准差的倒数.ϕm(m)为模型目标函数项，其体积分形式如下：
|
(4) |
其中，αs、αx、αy、αz为对应项的权重，w(r)为广义深度加权函数，将在后面进行详细讨论.从反演角度来讲，(4)式可以认为是最小模型和最平坦模型的组合.从正则化角度来讲，(4)式则可以认为是0阶和1阶吉洪诺夫正则化的组合.对(4)式进行有限差分近似，即可得到其矩阵形式：
|
(5) |
其中Ws为最小模型项的加权矩阵，Wx、Wy、Wz为最平坦模型项在x、y、z方向的加权矩阵，Wm为组合模型加权矩阵.对于磁异常三维反演，Li和Oldenburg(1996)建议模型目标函数中最小模型项和最平坦模型项应保持近似相同的权重.但对于井地模量联合反演，我们建议将最小模型项的系数设置为较小的权重，甚至为0，其他三项的系数则设置为1，即αs=0、αx=αy=αz=1，来获得一个最光滑解.这也是磁异常联合反演和模量联合反演的一个主要差别.对此，我们给出如下分析：磁异常(地面总场或井中三分量)反演假设不包含剩磁，即磁化方向是完全正确的，地下真实场源和井、地磁测数据之间均是严格对应的，此时不管选择最小模型、最平坦模型、还是两者的组合，都可以获得较好的反演结果；但对于剩磁影响下的磁异常模量反演而言，此时我们并不知道真实的磁化方向，而模量在一定程度上依然受到磁化方向的影响且地面和井中模量所受到的影响可能不同，此时地面模量和井中模量的反演结果会存在不一致性.数值模拟表明：如果在井地磁异常模量联合反演中仅包含最小模型项，则反演结果往往包含大量的冗余构造和人工痕迹；如果反演中仅采用最平坦模型项，可以很大程度地减少和平滑冗余构造，改善反演结果.也就是说，最小模型约束项对磁化方向的敏感度要大于最平坦模型约束项的敏感度.两种选择的具体反演效果将在后续章节进行展示.
Z为对角广义深度加权矩阵.常规的深度加权矩阵(Li and Oldenburg, 1996)仅与模型的深度有关，其仅适用于地面数据的反演，而不能用于井地联合反演.这里，我们参考Li和Oldenburg(2000b)在地面磁异常和井中三分量联合反演中的策略，给出磁异常模量井地联合反演的广义深度加权函数，分别为基于核函数的广义深度加权函数
|
(6) |
以及基于距离的广义深度加权函数
|
(7) |
其中wj为Z的第j个对角元素，N为观测数据(包括地面和井中模量数据)的个数，Fij为单位有效磁化率的第j个模型单元对第i个观测点的磁异常模量正演响应，rij为第j个模型单元与第i个观测点之间的距离.β为与磁异常模量衰减速率相关的权重系数，通常取1即可.这两种广义深度加权函数均与所有观测数据有关，因此可以用于井地联合反演.将(5)式和(3)式代入(2)式中，得到(2)式的矩阵形式：
|
(8) |
(8) 式是不等式约束非线性反演问题.传统的磁总场异常反演属于不等式约束线性反演问题，对于该问题的求解，当前已有多种算法成功应用到了重磁领域，如非线性投影法(Li and Oldenburg, 1996；Lelièvre and Oldenburg, 2006)，内点法(Li and Oldenburg, 2003)，以及梯度投影法(Lelièvre et al., 2009).这些方法理论上也可以扩展到不等式约束非线性反演问题的求解，本文选择内点法求解该问题，详细的算法见附录A.
2 模型试验为了验证井地磁异常模量反演算法的有效性，我们进行了针对性的模型试验.在所有模型试验中，定义x轴指向北，y轴指向东，z轴垂直向上.
2.1 倾斜薄板模型第一个模型为倾斜薄板模型，该模型的空间位置和几何形态与Li和Oldenburg(1998)所使用的组合板状体模型中的右侧板状体相同.其南北向的延伸距离为500至1500 m.垂向延伸则从-100到-1000 m.板状体的倾角为45°.地磁场强度为50000 nT，地磁场倾角和偏角分别为65°和25°.模型的总磁化倾角和偏角分别为90°和25°，有效磁化率为0.05SI.地磁场方向和总磁化方向的偏差为25°，因此，本例主要用于模拟剩磁影响不太大时的情况.地表观测数据为21×21的网格数据，网格间距为100 m.地面磁测数据中加入了标准差为2%数据绝对值加上2 nT的高斯噪声.地面磁异常及对应的模量数据见图 1.同时，计算了三口垂直井的三分量数据，每个井中有19个点，垂向间隔为50 m，三个钻井的参数见表 1.井中磁测数据见图 2，其中包含了2%的高斯噪声.需要说明的是，随着深度的变化，测量数据离场源的距离急剧变化，磁异常曲线变化很大，为了减少显示上的混叠效应，图 2的绘制是采用的间隔5 m的井中数据，但反演依旧是采用间隔为50 m的数据.
|
图 1 (a) 薄板模型所产生的磁异常数据，(b)对应的模量数据 三口垂直井的平面投影见十字符号. Fig. 1 (a) The total-field anomaly produced by the thin dipping slab model. (b) The amplitude data computed from the total field anomaly in Fig. 1a The positions of three vertical boreholes are indicated by the cross symbols. |
|
表 1 倾斜板状体模型三个钻井的参数 Table 1 Parameters of three boreholes of the dipping slab model |
|
图 2 薄板模型产生的井中三分量数据及对应的模量数据 Fig. 2 The three-component borehole data and corresponding amplitude data produced by the thin dipping slab model |
对于该模型，本文分别利用井中三分量和地面磁异常数据、地面磁异常模量数据、井中和地面磁异常模量数据这三种数据进行反演，将地下剖分为20×20×10的立方体，反演结果的中心切片如图 3所示.图 3a为井中三分量和地面磁异常数据联合反演的结果，由于没有使用正确的有效磁化方向，反演结果较差.图 3b为地面模量数据的反演结果，可以看出，采用模量减少了剩磁影响，改善了反演结果.但随着深度的增加，反演结果的分辨率越来越低.图 3c为仅包含最小模型项的井中和地面磁异常模量数据的反演结果，可以看出，尽管深部分辨率得到了改善，但又产生了大量虚假的冗余构造.图 3d则为仅包含最平坦模型项的井中和地面磁异常模量数据的反演结果，该结果既改善了深部分辨率，又没有冗余构造，效果较好.因此，在下文第二个模型中，我们在井地磁异常模量联合反演时依然选择仅包含最平坦模型项.
|
图 3 倾斜薄板模型不同方法的反演结果中心切片图 (a)井中三分量和地面磁异常数据联合反演结果; (b)地面磁异常模量数据反演结果; (c)仅包含最小模型项的井中和地面磁异常模量数据的反演结果; (d)仅包含最平坦模型项的井中和地面磁异常模量数据的反演结果.黑框代表了真实模型的位置. Fig. 3 Inversion results of anomalies produced by the thin dipping slab model by using different methods (a) The result of joint inversion of surface and three-component borehole magnetic data; (b) The result of inversion of surface magnetic amplitude data; (c) The result produced by joint inversion of surface and borehole magnetic amplitude data by only incorporating the smallest term; (d) The result produced by joint inversion of surface and borehole magnetic amplitude data by only incorporating the flattest term. The position of the true model is indicated by the black boxes. |
第二个模型采用的是两个相距75 m长方体模型，其大小、埋深、磁性参数各不相同，详见表 2.该模型的空间位置和几何形态与Li和Oldenburg(2000b)所使用的组合模型完全一致，但总磁化方向不相同.地磁场强度为50000 nT，地磁场倾角和偏角分别为65°和25°.由表 2可见，模型的磁化方向均与地磁场方向不同，用来模拟剩磁影响较大的情况.观测数据的范围为x：-300到300 m，y：-300到300 m，x方向数据间隔为100 m，y方向数据间隔为25 m，这样，地面观测数据为175个，地面磁测数据及换算得到的模量见图 4.同时，计算了三口垂直井的三分量数据，每个井中有16个点，垂向间隔为25 m，三个钻井的参数见表 3，井中磁测数据见图 5.所有的数据均加入了2 nT的高斯噪声.
|
|
表 2 两个长方体的顶点坐标和磁性参数 Table 2 Vertex coordinates and magnetic parameters of the two rectangular prisms |
|
图 4 (a) 组合长方体模型所产生的磁异常数据，(b)对应的模量数据 三口垂直井的平面投影见十字符号.黑框为模型的水平投影位置. Fig. 4 (a) The total-field anomaly produced by the model composed of two rectangular prisms. (b) The amplitude data computed from the total-field anomaly in Fig. 4a The positions of three vertical boreholes are indicated by the cross symbols. The horizontal projection of the true model is indicated by the black boxes. |
|
|
表 3 组合长方体模型三个钻井的参数 Table 3 Parameters of three boreholes of the model composed of two rectangular prisms |
|
图 5 组合长方体模型产生的井中三分量数据及对应的模量数据 Fig. 5 The three-component borehole data and corresponding amplitude data produced by the model composed of two rectangular prisms |
对于该模型，本文分别利用地面磁异常数据、井中三分量磁测数据、井中三分量和地面磁异常数据、地面模量数据、井中模量数据以及井中和地面模量数据这六种数据进行反演，将地下剖分为24×24×16的立方体，前三种数据的反演结果如图 6所示，后三种数据的反演结果如图 7所示.
|
图 6 组合长方体模型磁异常数据不同方法的反演结果切片图 (a)地面磁异常数据的反演结果; (b)井中三分量的反演结果; (c)井中三分量和地面磁异常数据联合反演结果.黑框代表了模型的真实位置. Fig. 6 Inversion results of magnetic anomalies produced by the model composed of two rectangular prisms by using different methods (a) Results of inversion of surface total-field data; (b) Results of inversion of three-component borehole data; (c) Results produced by joint inversion of surface and three-component borehole data. The position of the true model is indicated by the black boxes. |
|
图 7 组合长方体模型磁异常模量数据不同方法的反演结果切片图 (a)地面磁异常模量的反演结果; (b)井中磁异常模量的反演结果; (c)井中和地面磁异常模量联合反演结果.黑框代表了模型的真实位置. Fig. 7 Inversion results of magnetic amplitude data produced by the model composed of two rectangular prisms by using different methods (a) Results of inversion of surface amplitude data; (b) Results of inversion of borehole amplitude data; (c) Results produced by joint inversion of surface and borehole amplitude data. The position of the true model is indicated by the black boxes. |
从前三种数据的反演结果来看，地面磁异常的反演结果显示为一个场源，与实际情况不符；井中三分量反演对深部棱柱体的成像效果较好，浅部棱柱体不足，且存在部分虚假构造；井中三分量和地面磁异常联合反演的结果与井中三分量的反演效果类似，说明此时井中三分量数据对反演结果的贡献较大，其对浅部棱柱体的成像效果有所改善，但总体而言，浅部棱柱体的成像依然很差.
从后三种数据的反演结果来看，地面磁异常模量反演的结果和地面磁异常反演结果类似，均显示为单个场源，与实际情况不符；井中磁异常模量的反演结果则完全错误，原因主要是由磁异常换算模量数据是从矢量数据到标量数据的转换，原始三分量数据中的相位信息丢失，因此模量在减少剩磁效应的同时也增大了反演的非唯一性，此时仅凭借一维井中模量数据难以控制整个下半空间，导致反演结果错误；井地磁异常模量联合反演的效果较好，由于模量依然受磁化方向的影响，因此，反演结果的位置略有偏移，但从该反演结果中可以清晰地识别深部和浅部的两个棱柱体，整体成像效果较好.
3 结论金属矿区磁异常的反演和解释面临深度分辨率低和受剩磁影响大的问题，针对上述问题，我们将磁异常模量反演算法扩展到了井地联合反演中，模型试验表明，井地磁异常模量联合反演算法能够在有效的减少剩磁影响的同时提高反演的深度分辨率，是一种有效的联合反演算法.此外，我们强调磁异常模量在减少剩磁影响的同时也会增加反演解释的非唯一性，这是在模量数据解释中需要额外注意的.
附录A 井地磁异常模量反演的内点法内点法是近年来求解大规模不等式约束优化问题中最流行的方法之一，其基本思想是将不等式约束优化问题转化为一系列等式约束优化问题进行求解.对于(8)式，我们直接给出采用内点法优化的目标函数
|
(A1) |
其中，λ为障碍参数，ϕλ(m)为障碍函数，其形式如下：
|
(A2) |
通过添加对数障碍项，不等式约束就隐含到了目标函数中，障碍项的作用是沿着模型可行域的边界形成一个障碍并防止模型在目标函数优化的过程中越界.对于固定的μ而言，该方法求解一系列非线性优化问题，在迭代过程中，λ不断减小，随着λ趋近于0，障碍项占总目标函数的比重不断减小，最终解也就趋近了(8)式的解.将(A2)式代入(A1)式中，得到(A1)式的矩阵形式如下：
|
(A3) |
该式的优化包含两重迭代：其中外部迭代过程是根据不同的正则化参数μ和对应数据拟合差之间的关系来确定与已知数据噪声水平的相一致的μ，合适的μ值可以通过一维线性搜索来确定，我们不再进行详细说明；而内部迭代则是对于固定的μ值，求解对应的m，这里假设已知合适的μ值，而重点阐述内部迭代过程.采用高斯牛顿法对(A3)式进行极小化，可以得到第k次迭代时目标函数下降方向的求解公式：
|
(A4) |
式中，Δm为目标函数的下降方向.Jk为第k次迭代时的雅可比矩阵，其详细形式及推导过程见Li等(2010)和李泽林等(2015)，在此不再详述.Xk为对角矩阵，其主对角线上的元素为mk-1.1表示所有元素均为1的列向量.利用(A4)式，结合合适的迭代步长选择以及障碍参数更新策略，即可迭代求解最终结果.这里，依旧采用Li和Oldenburg(2003)的方法.对于固定的μ值，详细的迭代算法如下：
初始化：m0=0.001，λ0=ϕ2(m0)/ϕλ(m0)，k=1.
重复：
(1) 由mk-1生成Xk和Jk.
(2) 采用预处理共轭梯度法求解(A4)式得到Δm.
(3) mk←mk-1+γβΔm，其中γ=0.925，
(4) λk←(1-min[β, γ])λk-1.
(5) k←k+1.
直到收敛.
对上述算法，我们给出部分解释.由于算法迭代的要求，模型初值m0不能为0，因此给了一个比较接近0的初值.模型初值给定后，障碍参数初值λ0的选择是用来平衡(A1)式中障碍函数与另一项的大小.β的选择是为了保证更新后的所有模型依然为正值.γ设置为0.925的原因则是使得更新后的模型应该距离边界(即0)有一定的距离，来增强迭代的稳定性.λ的更新则依赖于迭代步长：如果步长很大(非常接近1)，则说明模型整体离边界(即0)还有很远，那么可以尽可能的减小障碍项的作用，即λ可以减小很多；如果步长很小，则说明模型中有部分值已经离边界很近，那么需要尽可能保持目前障碍项的作用，即λ应减小的尽可能少.更详细的解释及方法的原理可以参考Li和Oldenburg(2003)和Nocedal和Wright(2006).
当障碍函数充分小，障碍项对总目标函数的贡献很小的时候，迭代停止.我们在实际编程中通常设置三个附加准则来保证迭代成功，一是目标函数的改变量小于百分之一，二是障碍项与总目标函数的比值小于百分之一，三是达到预先设置的最大迭代次数.
Aster R C, Borchers B, Thurber C H. 2005. Parameter Estimation and Inverse Problems. Aster R C, Borchers B, Thurber C H. 2005. Parameter Estimation and Inverse Problems. New York: Academic Press..
|
Baranov V. 1957. A new method for interpretation of aeromagnetic maps:pseudo-gravimetric anomalies. Geophysics, 22(2): 359-382. DOI:10.1190/1.1438369 |
Barbosa V C F, Silva J B C. 1994. Generalized compact gravity inversion. Geophysics, 59(1): 57-68. |
Camacho A G, Montesinos F G, Vieira R. 2000. Gravity inversion by means of growing bodies. Geophysics, 65(1): 95-101. DOI:10.1190/1.1444729 |
Dampney C N G. 1969. The equivalent source technique. Geophysics, 34(1): 39-53. |
Dannemiller N, Li Y G. 2006. A new method for determination of magnetization direction. Geophysics, 71(6): L69-L73. DOI:10.1190/1.2356116 |
Gallardo L A, Meju M A. 2004. Joint two-dimensional DC resistivity and seismic travel time inversion with cross-gradients constraints. Journal of Geophysical Research:Solid Earth, 109(B3): :B03311. DOI:10.1029/2003JB002716 |
Gerovska D, Aráuzo-Bravo M J, Stavrev P. 2009. Estimating the magnetization direction of sources from southeast Bulgaria through correlation between reduced-to-the-pole and total magnitude anomalies. Geophysical Prospecting, 57(4): 491-505. DOI:10.1111/gpr.2009.57.issue-4 |
Guillen A, Menichetti V. 1984. Gravity and magnetic inversion with minimization of a specific functional. Geophysics, 49(8): 1354-1360. DOI:10.1190/1.1441761 |
Hou Z C. 1979. Using potential field transformation to build an interpretation system of gravity and magnetic anomalies. Geophysical and Geochemical Exploration (in Chinese), 3(3): 1-10. |
Krahenbuhl R A, Li Y G. 2006. Inversion of gravity data using a binary formulation. Geophysical Journal International, 167(2): 543-556. DOI:10.1111/gji.2006.167.issue-2 |
Krahenbuhl R A, Li Y G. 2009. Hybrid optimization for lithologic inversion and time-lapse monitoring using a binary formulation. Geophysics, 74(6): I55-I65. DOI:10.1190/1.3242271 |
Last B J, Kubik K. 1983. Compact gravity inversion. Geophysics, 48(6): 713-721. DOI:10.1190/1.1441501 |
Lelièvre P G, Oldenburg D W. 2006. Magnetic forward modelling and inversion for high susceptibility. Geophysical Journal International, 166(1): 76-90. DOI:10.1111/gji.2006.166.issue-1 |
Lelièvre P G, Oldenburg D W. 2009a. A comprehensive study of including structural orientation information in geophysical inversions. Geophysical Journal International, 178(2): 623-637. DOI:10.1111/gji.2009.178.issue-2 |
Lelièvre P G, Oldenburg D W. 2009b. A 3D total magnetization inversion applicable when significant, complicated remanence is present. Geophysics, 74(3): L21-L30. |
Lelièvre P G, Oldenburg D W, Williams N C. 2009. Integrating geological and geophysical data through advanced constrained inversions. Exploration Geophysics, 40(4): 334-341. |
Lelièvre P G, Farquharson C G, Hurich C A. 2012. Joint inversion of seismic traveltimes and gravity data on unstructured grids with application to mineral exploration. Geophysics, 77(1): K1-K15. |
Li Y G, Oldenburg D W. 1996. 3-D inversion of magnetic data. Geophysics, 61(2): 394-408. DOI:10.1190/1.1443968 |
Li Y G, Oldenburg D W. 1998. 3-D inversion of gravity data. Geophysics, 63(1): 109-119. |
Li Y G, Oldenburg D W. 2000a. Incorporating geological dip information into geophysical inversions. Geophysics, 65(1): 148-157. |
Li Y G, Oldenburg D W. 2000b. Joint inversion of surface and three-component borehole magnetic data. Geophysics, 65(2): 540-552. |
Li Y G, Oldenburg D W. 2003. Fast inversion of large-scale magnetic data using wavelet transforms and a logarithmic barrier method. Geophysical Journal International, 152(2): 251-265. DOI:10.1046/j.1365-246X.2003.01766.x |
Li Y G, Shearer S E, Haney M M, et al. 2010. Comprehensive approaches to 3D inversion of magnetic data affected by remanent magnetization. Geophysics, 75(1): L1-L11. |
Li Y G, Sun J J. 2016. 3D magnetization inversion using fuzzy c-means clustering with application to geology differentiation. Geophysics, 81(5): J61-J78. DOI:10.1190/geo2015-0636.1 |
Li Z L, Yao C L, Zheng Y M, et al. 2015. 3D data-space inversion of magnetic amplitude data. Chinese Journal of Geophysics (in Chinese), 58(10): 3804-3814. DOI:10.6038/cjg20151030 |
Lines L R, Schultz A K, Treitel S. 1988. Cooperative inversion of geophysical data. Geophysics, 53(1): 8-20. |
Liu S, Feng J, Gao W L, et al. 2013. 2D inversion for borehole magnetic data in the presence of significant remanence and demagnetization. Chinese Journal of Geophysics (in Chinese), 56(12): 4297-4309. DOI:10.6038/cjg20131232 |
Liu S, Hu X Y, Liu T Y, et al. 2013. Magnetization vector imaging for borehole magnetic data based on magnitude magnetic anomaly. Geophysics, 78(6): D429-D444. DOI:10.1190/geo2012-0454.1 |
Liu S, Hu X Y, Xi Y F, et al. 2015. 2D sequential inversion of total magnitude and total magnetic anomaly data affected by remanent magnetization. Geophysics, 80(3): K1-K12. |
Menke W. 2012. Geophysical Data Analysis:Discrete Inverse Theory:MATLAB Edition. London: Academic Press.
|
Nabighian M N. 1972. The analytic signal of two-dimensional magnetic bodies with polygonal cross-section; its properties and use for automated anomaly interpretation. Geophysics, 37(3): 507-517. DOI:10.1190/1.1440276 |
Nocedal J, Wright S. 2006. Numerical Optimization. 2nd ed. New York: Springer.
|
Paine J, Haederle M, Flis M. 2001. Using transformed TMI data to invert for remanently magnetised bodies. Exploration Geophysics, 32(4): 238-242. |
Phillips J D. 2005. Can we estimate total magnetization directions from aeromagnetic data using Helbig's integrals?. Earth, Planets and Space, 57(8): 681-689. DOI:10.1186/BF03351848 |
Pilkington M. 1997. 3D magnetic imaging using conjugate gradients. Geophysics, 62(4): 1132-1142. DOI:10.1190/1.1444214 |
Pilkington M. 2009. 3-D magnetic data-space inversion with sparseness constraints. Geophysics, 74(1): L7-L15. |
Pilkington M, Beiki M. 2013. Mitigating remanent magnetization effects in magnetic data using the normalized source strength. Geophysics, 78(3): J25-J32. DOI:10.1190/geo2012-0225.1 |
Portniaguine O, Zhdanov M S. 1999. Focusing geophysical inversion images. Geophysics, 64(3): 874-887. DOI:10.1190/1.1444596 |
Portniaguine O, Zhdanov M S. 2002. 3D magnetic inversion with data compression and image focusing. Geophysics, 67(5): 1532-1541. DOI:10.1190/1.1512749 |
Roest W R, Verhoef J, Pilkington M. 1992. Magnetic interpretation using the 3-D analytic signal. Geophysics, 57(1): 116-125. |
Shearer S E. 2005. Three-dimensional inversion of magnetic data in the presence of remanent magnetization[Master's thesis]. Colorado: Colorado School of Mines. https://www.researchgate.net/publication/33999135_Three-dimensional_inversion_of_magnetic_data_in_the_presence_of_remanent_magnetization
|
Stavrev P, Gerovska D. 2000. Magnetic field transforms with low sensitivity to the direction of source magnetization and high centricity. Geophysical Prospecting, 48(2): 317-340. DOI:10.1046/j.1365-2478.2000.00188.x |
Uieda L, Barbosa V C F. 2012. Robust 3D gravity gradient inversion by planting anomalous densities. Geophysics, 77(4): G55-G66. DOI:10.1190/geo2011-0388.1 |
Wang M Y, Di Q Y, Xu K, et al. 2004. Magnetization vector inversion equations and 2D forward and inversed model study. Chinese Journal of Geophysics (in Chinese), 47(3): 528-534. |
Wilson H S. 1985. Analysis of the magnetic gradient tensor:Defense Research Establishment Pacific. Technical Memorandum: 5-13. |
Yao C L, Hao T Y, Guan Z N, et al. 2003. High-speed computation and efficient storage in 3-D gravity and magnetic inversion based on genetic algorithms. Chinese Journal of Geophysics (in Chinese), 46(2): 252-258. |
Yao C L, Zheng Y M, Zhang Y W. 2007. 3-D gravity and magnetic inversion for physical properties using stochastic subspaces. Chinese Journal of Geophysics (in Chinese), 50(5): 1576-1583. |
Zhou J Y, Revil A, Karaoulis M, et al. 2014. Image-guided inversion of electrical resistivity data. Geophysical Journal International, 197(1): 292-309. DOI:10.1093/gji/ggu001 |
侯重初. 1979. 利用位场转换建立一个重磁异常解释系统. 物探与化探, 3(3): 1-10. |
李泽林, 姚长利, 郑元满, 等. 2015. 数据空间磁异常模量三维反演. 地球物理学报, 58(10): 3804-3814. DOI:10.6038/cjg20151030 |
刘双, 冯杰, 高文利, 等. 2013. 强剩磁强退磁条件下的二维井中磁测反演. 地球物理学报, 56(12): 4297-4309. DOI:10.6038/cjg20131232 |
王妙月, 底青云, 许琨, 等. 2004. 磁化强度矢量反演方程及二维模型正反演研究. 地球物理学报, 47(3): 528-534. DOI:10.3321/j.issn:0001-5733.2004.03.025 |
姚长利, 郝天珧, 管志宁, 等. 2003. 重磁遗传算法三维反演中高速计算及有效存储方法技术. 地球物理学报, 46(2): 252-258. DOI:10.3321/j.issn:0001-5733.2003.02.020 |
姚长利, 郑元满, 张聿文. 2007. 重磁异常三维物性反演随机子域法方法技术. 地球物理学报, 50(5): 1576-1583. DOI:10.3321/j.issn:0001-5733.2007.05.035 |