文章快速检索  
  高级检索
互相关系数自约束的重力三维反演与高效求解
梁生贤     
中国地质调查局成都地质调查中心, 成都 610081
摘要: 本文将拟合残差计算所得的互相关系数作为先验信息,与深度加权函数同时引入到重力正则化反演的模型约束中,以提升反演结果的可靠性。针对三维反演中的大型线性方程组问题,引入阻尼LSQR(最小二乘QR分解)算法,结合等效几何格架技术,将大矩阵按照模型单元划分为若干个子矩阵进行存储与运算。理论模型计算结果表明:同时利用互相关系数和深度加权的模型自约束反演,能较清晰地反映真实异常体;基于分块矩阵的阻尼LSQR算法求解线性方程组较直接法节省了几千甚至上万倍的存储量,且计算速率提高了数倍,可在普通计算机上实现较大规模的反演计算。将其应用于云南芦子园铁铅锌铜多金属矿床隐伏花岗岩体定位,取得了良好的效果。
关键词: 重力三维反演    互相关系数    自约束反演    阻尼LSQR算法    云南芦子园铁铅锌铜多金属矿床    花岗岩体    
A Self-Constrained 3D Inversion and Efficient Solution of Gravity Data Based on Cross-Correlation Coefficient
Liang Shengxian     
Chengdu Center, China Geological Survey, Chengdu 610081, China
Supported by National Key Research and Development Program of China(2016YFC060308) and Project of China Geological Survey(121201010000150007, 121201010000150014)
Abstract: In this study, we applied both cross-correlation coefficient from fitting residual and depth weight to the constrained model of regularization gravity inversion to improve the reliability. In our approach, for the solutions in the system linear equations of gravity 3D inversion, damping LSQR algorithm method was introduced. Combined with the equivalent geometric trellis technology, the large matrix was divided into several sub-matrixes for storage and computation according to the model unit. The results of the theoretical model indicate that through the application of the cross-correlation coefficient with the depth weight into the constraint inversion model at the same time, the true anomaly can be clearly reflected. The damping LSQR algorithm method is sufficient for improving thousands or even tens of thousands of times of storage, and for increasing several times of computation rate. Thus, a large scale inversion calculation can be realized on a normal computer. The developed model is applied to the buried granite positioning in Luziyuan area, and the result is consistent with the true situation.
Key words: gravity 3D inversion    cross-correlation coefficient    self-constrained inversion    damping LSQR algorithm    Luziyuan area of Yunnan Province    granite    

0 引言

重力反演就是根据观测场值求解场源的信息,是资料定量解释的重要环节之一。物性反演将模型空间离散化为若干个单元,只求解各单元相应的密度值,这种方法易于模拟复杂的地质体[1],逐渐成为重力三维反演的重要方向[2]。本文所讨论的就是物性反演。它是线性离散不适定问题,加之三维反演的数据量比待求解的模型量更少、系数矩阵(核矩阵)严重病态性和观测数据受误差影响等,多解性和不稳定性严重;因此,如何使得解模型更加符合实际情况是反演关注的主要问题之一。此外,当观测数据量很大时,系数矩阵大型稠密,反演涉及到大规模数据的存储与运算,在实际运算中须考虑到这一棘手问题。

在不适定问题求解方面,Tikhonov正则化方法是目前应用最为广泛的算法之一,它在一定意义上有效减少了解的多解性和不稳定性。已有研究表明,在正则化反演中尽可能地利用先验信息对场源施加约束,可使解模型更加符合实际情况,如:Li等[3-4]在重磁反演中通过加入深度加权函数来克服位场的“趋肤效应”;Paoletti等[5]提出位场自约束反演概念,指出位场本身包含了大量的潜在信息,利用这些潜在的信息对模型自约束可提高反演结果的可靠性;Sun等[6]利用磁场两个方向上的水平梯度数据计算互相关系数,得到空间加权函数,并将空间加权函数引入正则化反演的模型约束中,显示这种模型自约束反演结果对场源的边界刻画更清晰。

在节省计算成本方面,主要有并行计算[7-8]、核矩阵压缩、重构[9-13]以及快速正演计算[1-2, 14]等方法,其中并行计算需要依靠计算机硬件设备,核矩阵压缩则会导致部分信号失真等。对于Tikhonov正则化反演,最终都归结为大型线性方程组问题的求解。Krylov子空间迭代法是求解大型线性方程组问题的一种有效途径,在位场反演中应用较多的为共轭梯度法[9-10, 15],近年来LSQR(最小二乘QR分解)算法[16]也越来越多地被应用于位场数据反演中[17-19];对于非病态线性方程组两者等效,而对于病态线性方程组后者求解更稳定。

本文从提升重力反演结果的可靠性和节省计算成本两个方面着手。其一,将快速扫描的互相关系数作为先验信息,通过处理引入到重力正则化反演的目标函数中,以提高反演结果的可靠性;其二,利用LSQR算法求解线性方程组问题,并对其进行相应的改进,与快速正演计算方法结合,以节省计算成本。最后通过理论模型和实际数据来展示上述方法的应用效果。

1 反演理论 1.1 正则化反演

重力正演公式可写为

(1)

式中:mM阶模型向量;dN阶数据向量;AN×M阶核矩阵,其元素A(i, j)为第j个模型单元在地表i处的重力响应核函数,在反演中保持不变。

物性反演的任务是根据观测数据向量反求模型向量。由于反演是不适定的,Tikhonov正则化目标函数构制如下:

(2)

式中:Δd为观测数据向量与模型响应向量之间的残差;Δm为模型修改向量;λ为正则化因子或拉格朗日因子,体现了数据拟合与模型约束之间的某种“折衷”;D为模型约束矩阵。此外,为了克服重磁观测幅值随场源深度增加而迅速衰减的“趋肤效应”,通常加入一个深度加权函数,在仅考虑深度方向的情况下可写为[3-4]

式中:zj为模型单元中心埋深,j∈[1, M];β取1.5~2.0。

根据极值条件,求目标函数Em, λ)关于模型修改向量Δm或ΔmT的偏导数,并令其等于0,可得模型空间迭代求解公式:

(3)

Siripunvaraporn等[20]基于Occam反演策略提出数据空间算法,在数据量N远小于模型量M的情况下,可大幅度提高计算效率,数据空间迭代求解公式为

(4)

式中,I为单位矩阵。

对于大规模重力数据三维正则化反演而言(数据量可能为几千甚至上万,模型量可能达几十万甚至百万),无论模型空间还是数据空间,直接法(奇异值分解)求解都需要耗费巨大的计算量和存储量(表 1)。

表 1 模型空间、数据空间与分块矩阵LSQR方法计算成本对比 Table 1 Comparison of computational cost of model space, data space and block matrix & LSQR method
方法 计算成本 时间复杂度 空间复杂度
矩阵与矩阵的乘积 矩阵与向量的乘积 矩阵求逆
模型空间 (M×N)×(N×M)和(M×M)×(M×N) (M×N)×(N×1) (M×M) O(2M2 N)+O(MN)+O(M3) O(M2+2MN)
数据空间 (N×M)×(M×N)和(M×N)×(N×N) (M×N)×(N×1) (N×N) O(2MN2)+O(MN)+O(N3) O(N2+2MN)
分块矩阵LSQR方法 (M×N)×(N×1)×it和(N×M)×(M×1)×it O(2iMN) O(2M+2N)
注:it为LSQR算法的迭代次数,一般而言it<<min(M, N);存储量不包含矩阵求逆所需空间;O表示时间或空间复杂度。
1.2 利用互相关系数与深度加权的约束反演

在反演中,尽可能地利用已知地质信息(包括地表地质、钻孔以及地震资料等)对场源施加约束是提高解的可靠性有效而实用的措施,但在有些情况下,地质信息并不充足。

互相关系数[21-22]表征了实测重力异常与模型单元核函数的相关程度,其绝对值越接近1,表示该单元对重力异常的贡献可能性越大,且具有计算简单而快速的特点[23-24]。因而,可根据互相关系数的绝对值判断模型单元的重要性,从而将相对重要的模型单元作为先验信息对模型进行加权约束。并鉴于利用剩余异常计算所得的互相关系数成像结果可靠性较高[23-26],利用每步反演迭代的拟合残差计算互相关系数,取绝对值和归一化后对模型进行约束:

(5)

式中:Dω为互相关系数约束矩阵;ωj为互相关系数,其数学表达式[23-24]

(6)

式中,Δdi为第i个观测数据与模型响应之间的残差。若初始模型为均匀半空间,则第一次迭代的Δdi=diobs(diobs为观测数据)。

此外,相邻单元的互相关系数往往是渐变的,整体具有平缓变化的特点,因而无需再加入模型粗糙度约束矩阵。考虑到位场“趋肤效应”,令D=Dω·Ddepth

我们根据文献[4]中的理论模型来展示不同算法的加权结果。该理论模型(图 1)由两个倾斜方向相反的脉状体异常体组成[4],其中:向西倾斜的长脉状体剩余密度为1.0 g/cm3,向东倾斜的短脉状体剩余密度为0.8 g/cm3。正演数据密度:东西方向点距50 m,南北方向点距100 m,共861个数据。将数据加入4%的随机噪声,模型水平方向剖分与数据网格一一对应,垂向剖分24层,共20 664个模型单元。

图 1 理论模型 Figure 1 Theory model

图 2为理论模型计算的单元权重结果。由图 2可见:仅利用深度加权时(图 2a),水平向各单元权重相同,加权值仅在垂向随深度增加而变大;仅利用互相关系数加权时(图 2b),模型单元权重依赖于互相关系数,而互相关系数对场源的准确位置反映较差,求解结果可能会出现偏差;同时利用深度与互相关系数加权时(图 2cde),由于一般情况下深度加权值的量级远大于互相关系数的量级,因而不至让反演结果过分地依赖于互相关系数,模型迭代求解从深部开始(图 2c),随着迭代进行,拟合残差逐渐变小,根据拟合残差计算所得的互相关系数逐渐趋于0附近,各单元权重对深度加权的依赖减小,相应的加权值也总体变小,其形态逐渐逼近理论模型(图 2de)。可见,同时利用互相关系数与深度加权的模型约束无疑使得先验信息更加丰富,其本质上属于自约束反演[5]的一种,有助于提高反演结果的可靠性。

a.深度加权;b.互相关系数加权;c.深度与初始互相关系数加权;d.反演迭代两次后的深度与互相关系数加权;e.反演迭代三次后的深度与互相关系数加权。黑色线框为理论模型。 图 2 理论模型不同加权结果 Figure 2 Synthetic model of different weighted results

在正则化反演过程中,拉格朗日因子(正则化因子)体现了模型约束与数据拟合之间的平衡:过大的拉格朗日因子往往偏重于模型约束中各单元体的重要性,而忽视了数据的拟合程度,会使得迭代收敛缓慢,甚至于出现不收敛的情况;过小的拉格朗日因子则偏重于数据拟合,而忽视了模型各单元应有权重的影响,会使得反演出现过拟合现象,无疑增加了解的非唯一性,使得求解结果不够稳定。目前常用的拉格朗日因子求取算法有广义偏差准则、广义交叉验证准则和L曲线法,3种计算方法都包含了大规模矩阵的多次计算,无疑增加了计算成本。本文拉格朗日因子求取采用在一定步长下逐次递减的方法[27],当出现迭代发散时,拉格朗日因子相应地增大并重新求解。该方法求解过程较为稳定,且避免了多次拉格朗日因子搜索所需的额外计算成本。

2 分块矩阵LSQR方法 2.1 阻尼LSQR算法

由于核矩阵为大型稠密的,直接法求解式(3)或式(4)的计算成本巨大。Krylov子空间迭代法仅出现矩阵与向量的乘积,具有收敛速度快、对计算机内存要求低等优势。CGLS(共轭梯度最小二乘法)和LSQR算法同属于Krylov子空间投影方法,两者的运算量与存储量相当,但当执行多次迭代或系数矩阵为病态时,LSQR的数值稳定性更好[16]

因此,我们引入LSQR算法求解反问题。令A*=AD-1,Δm*=DΔm;由于D为对角矩阵,有D-1D=I,则式(2)可写为

此时式(7)为阻尼最小二乘法,可引入阻尼LSQR算法求解线性方程组问题。在求得Δm*后,根据Δm=D-1Δm*反求模型修改量。

方程A*×Δm*d关于阻尼因子λ最小二乘问题的LSQR算法见表 2[28]

表 2 阻尼LSQR算法 Table 2 Algorithm damping LSQR
初始化:
  
  
for i=1, 2, 3,…
  
  
  
  
  
  
  
   达到迭代终止条件后退出
注:μ1N维向量;ω1ν1、Δm0*M维向量;β1α1ρ1*φ1*为实数。具体含义见文献[28]。
2.2 系数矩阵分块运算

由上述阻尼LSQR算法可见,它仅涉及到矩阵与向量的乘积。若将A*按列划分为k个子矩阵,且向量b(表 2中的μiνi)及子向量bi的阶数与对应的矩阵及子矩阵列数一致,则有:

假设面积数据是网格化水平分布的,并且模型剖分单元与数据网格采取一一对应的关系,则根据平移等效性和互换对称性,在同一层模型各单元之间,核矩阵A的元素具有特定的规律,即等效几何格架;实际中只需计算和存储每一层的第一个单元在所有观测点处的重力响应核函数,同一层的其他元素可根据下式进行索引查找[2]

(8)

式中:ij为任意模型单元在x方向和y方向上的排列号;kl为任意数据观测点在x方向和y方向上的排列号。

若子矩阵按照模型单元划分,根据式(8)可快速获取某个模型单元在所有观测点处的重力响应核函数。大型稠密系数矩阵不再被完整地表示出来,且由于D为对角矩阵,不会明显增加额外的计算量,实现分块矩阵LSQR方法与等效几何格架技术的结合。

在具体计算成本方面,分块矩阵LSQR方法的空间复杂度仅为O(2M+2N),时间复杂度为O(2itMN),由于一般情况下迭代次数i远小于数据量N和模型量M,相比直接法求解式(3)、(4),分块矩阵LSQR方法节省了大量存储空间和计算时间(表 1)。假设测区网格数据为100×100的规模,模型为100×100×50的三维网格,LSQR迭代次数为5 000次,则分块矩阵LSQR方法的计算量不到数据空间的一半,反演速率至少提高了2倍;假设数据以双精度存储,则模型空间求解需1.89 TB存储量,数据空间求解需75 GB存储量,分块矩阵LSQR方法则仅需约8 MB存储量,在普通计算机上就能实现较大规模的三维反演计算。

3 反演实例 3.1 理论模型

理论模型及数据密度、网格剖分等见2.2节。我们分别利用DdepthDωDdepth对模型进行约束反演,附加密度约束范围为0~1 g/cm3,最终迭代反演的均方误差分别为0.031、0.112 mGal。

图 3ab分别为利用DdepthDωDdepth进行模型约束的反演。两种模型约束的反演异常在浅部差别不大,但在深度为400 m以下则差别逐渐变大:利用Ddepth进行模型约束的反演异常呈较正的“v”字型(图 3a),这与真实模型不一致;而利用DωDdepth进行模型约束的反演异常呈“y”字型(图 3b),较清晰地反映了真实异常体的基本轮廓。

a.利用深度加权函数进行模型约束的反演结果;b.利用互相关系数与深度加权进行模型约束的反演结果。 图 3 理论模型两种约束方法的反演结果对比图 Figure 3 Comparison of the inversion results of the synthetic model with two model constraints
3.2 实例

芦子园铁铅锌铜多金属矿位于保山—镇康地块的南端,成矿类型为岩浆热液型,已有地质、物探及化探工作均认为该区成矿作用与隐伏中酸性岩体密切相关[29]。但由于矿区地表未见岩体出露,且在矿区南东部出露的印支期木场花岗岩体与成矿作用无明显关系,因而针对岩体的研究总体较少,譬如岩体的埋深、规模及其空间形态等均鲜有研究。区内主要出露寒武系、奥陶系、志留系、泥盆系、石炭系、二叠系的碳酸盐岩、碎屑岩,石炭系、三叠系的火山岩以及第四系。物性资料(图 4)表明:木厂出露的花岗岩体密度常见值为2.59 g/cm3,明信坝出露石英闪长玢岩密度常见值为2.66 g/cm3,两者相对于区内广泛出露的碳酸盐岩(密度通常在2.73~2.75 g/cm3)围岩表现为低密度特征。在不同类型岩石中,泥岩、粉砂岩密度最小,其次为中酸性岩体,矽卡岩化、矿石以及火山岩密度最高,碳酸盐岩密度居中;区内第四系以及新近系、古近系、白垩系和侏罗系规模有限。综合密度测试结果可知,重力测量在区内寻找隐伏中酸性岩体具备较好的物性条件。我们以云南地质调查院提供的芦子园地区1:5布格重力数据为例,进行三维反演以推测区内隐伏中酸性岩体的空间分布特征。

图 4 云南芦子园地区岩石密度测试统计结果 Figure 4 Statistical results of rock density in Luziyuan, Yunnan

反演前采用矩形窗口滑动平均法求取剩余重力异常(图 5)。利用网格化的剩余重力异常数据进行三维反演,共4 148个数据,模型水平方向剖分与数据一一对应,垂向剖分34层,共141 032个网格单元。在普通Thinkpad台式机电脑上(处理器:Intel i5-2400CPU,3.10 GH;内存:4 GB)反演迭代3次,共耗时51 h 54 min,最终反演的均方误差为0.104 mGal。

图 5 云南芦子园地区剩余重力异常图 Figure 5 Residual gravity anomaly contour of Luziyuan, Yunnan

由于中酸性岩体表现为低密度,为方便起见,我们只提取三维反演结果的低密度体(根据物性测试统计结果,并结合MT测量与密度成像结果的形态特征[30-31],低密度体取剩余密度小于0.1 g/cm3),并将其套合在地质图上(图 6)。由于区内砾岩、砂岩、粉砂岩等低密度体规模有限,且表现为低电阻率特征,而中酸性岩体表现为高电阻率特征,结合MT测量结果[30-31]以及其他地质、物探、化探[29-32]证据,推测这种大规模低密度体主要反映了中酸性岩体。

1.第四系砂、砾、黏土;2.新近系南林组砾岩、砂岩、粉砂岩;3.古近系勐腊组砾岩、砂岩、泥岩;4.白垩系南新组砂砾岩;5.侏罗系中统粉砂岩、页岩;6.三叠系粉砂岩、灰岩、玄武岩;7.二叠系灰岩;8.石炭系玄武岩、灰岩;9.泥盆系灰岩;10.志留系栗柴坝组灰岩;11.奥陶系灰岩;12.寒武系上统碳酸盐岩;13.碱长花岗岩;14.铁矿;15.铅锌矿;16.铅锌铜矿;17.锡矿;18.推测中酸性岩体。 图 6 云南芦子园地区地质及重力三维反演结果图 Figure 6 Results of 3D gravity inversion and geology of Luziyuan, Yunnan

图 6可见:研究区外围南东部出露的印支期木场花岗岩体与反演结果的南东部低密度体在平面位置上一致,总体倾向南东,显示与成矿作用关系不大,这与地质、化探认识一致;而研究区北东的低密度体在总体趋势上与镇康复背斜一致(走向北东),显示复背斜对岩浆侵位的控制;已知矿床在平面位置上均位于岩体形成的凹入部位或转折部位,为岩浆热液成矿的最有利地段,其中,芦子园、天生桥一带的侵入岩枝隆起中心为天生桥,并向芦子园矿段侧伏,这与天生桥到芦子园的“背形”隆起条带状矽卡岩型磁铁矿(铁矿体自天生桥向芦子园呈隐伏延伸,埋深逐渐加大)形态一致。根据上述推断的岩体空间分布特征与已知地质信息[29-32],显示了本算法的实用性。

4 结果与结论

1) 将拟合残差计算所得的互相关系数作为先验信息,与深度加权函数同时引入正则化反演模型约束中,既体现了互相关系数在模型约束中的作用,又不至让反演结果过分依赖于互相关系数。理论模型反演结果表明,这种自约束反演方法对真实异常体的轮廓反映较为清晰,解模型更加符合实际情况。

2) 引入阻尼LSQR算法求解反问题,由于一般情况下迭代次数远小于数据量和模型量,反演速率可提高数倍。算法中只涉及到矩阵与向量的乘积,便于实现分块运算,结合等效几何格架技术,将原矩阵按照模型单元划分为若干个子矩阵进行存储与运算,极大地节省了反演对计算机存储空间的需求。

3) 将本文的重力三维反演方法应用于云南芦子园隐伏花岗岩体的定位中,共4 148个数据,141 032个模型单元,在普通计算机上运算仅需约52 h。反演结果显示已知矿床均位于推测岩体形成的凹入部位内侧,为岩浆热液成矿的最有利地段,验证了方法的有实用性。

致谢: 王永华教授在程序编写及论文撰写期间提供了帮助,在此表示感谢。
参考文献
[1]
侯遵泽, 杨文采, 刘家琦. 中国大陆地壳密度差异多尺度反演[J]. 地球物理学报, 1998, 41(5): 642-651.
Hou Zunze, Yang Wencai, Liu Jiaqi. Multiscale Inversion of the Density Contrast Within the Crust of China[J]. Chinese Journal Geophysics, 1998, 41(5): 642-651. DOI:10.3321/j.issn:0001-5733.1998.05.007
[2]
姚长利, 郝天珧, 管志宁, 等. 重磁遗传算法三维反演中高速计算及有效存储方法技术[J]. 地球物理学报, 2003, 46(2): 252-258.
Yao Changli, Hao Tianyao, Guan Zhining, et al. High Speed Computation and Efficient Storage in 3D Gravity and Magnetic Inversion Based on Genetic Algorithms[J]. Chinese Journal Geophysics, 2003, 46(2): 252-258. DOI:10.3321/j.issn:0001-5733.2003.02.020
[3]
Li Y G, Oldenburg D W. 3-D Inversion of Magnetic Data[J]. Geophysics, 1996, 61(2): 394-408. DOI:10.1190/1.1443968
[4]
Li Y G, Oldenburg D W. 3-D Inversion of Gravity Data[J]. Geophysics, 1998, 63(1): 109-119. DOI:10.1190/1.1444302
[5]
Paoletti V, Ialongo S, Florio G, et al. Self-Constrain-ed Inversion of Potential Fields[J]. Geophysical Journal International, 2013, 195(2): 854-869. DOI:10.1093/gji/ggt313
[6]
Sun S D, Chen C. A Self-Constrained Inversion of Magnetic Data Based on Correlation Method[J]. Journal of Applied Geophysics, 2016, 135: 8-16. DOI:10.1016/j.jappgeo.2016.09.022
[7]
Moorkamp M, Jegen M, Roberts A, et al. Massively Parallel Forward Modeling of Scalar and Tensor Gravimetry Data[J]. Computers & Geosciences, 2010, 36(5): 680-686.
[8]
陈召曦, 孟小红, 郭良辉, 等. 基于GPU并行的重力、重力梯度三维正演快速计算及反演策略[J]. 地球物理学报, 2012, 55(12): 4069-4077.
Chen Zhaoxi, Meng Xiaohong, Guo Lianghui, et al. Three-Dimensional Fast Forward Modeling and Inversion Strategy for Large Scale Gravity and Gravimetry Data Based on GPU[J]. Chinese Journal Geophysics, 2012, 55(12): 4069-4077. DOI:10.6038/j.issn.0001-5733.2012.12.019
[9]
Portniaguine O, Zhdanov M S. 3-D Magnetic Inversion with Data Compression and Image Focusing[J]. Geophysics, 2002, 67(5): 1532-1541. DOI:10.1190/1.1512749
[10]
Li Y G, Oldenburg D W. Fast Inversion of Large-Scale Magnetic Data Using Wavelet Transforms and Alogarithmic Barrier Method[J]. Geophysical Journal International, 2003, 152(2): 251-265. DOI:10.1046/j.1365-246X.2003.01766.x
[11]
刘天佑. 位场勘探数据处理新方法[M]. 北京: 科学出版社, 2007: 47-49.
Liu Tianyou. New Data Processing Methods for Potential Field Exploration[M]. Beijing: Science Press, 2007: 47-49.
[12]
Davis K, Li Y G. Fast Solution of Geophysical Inversion Using Adaptive Mesh, Space-Filling Curves and Wavelet Compression[J]. Geophysical Journal International, 2011, 185(1): 157-166. DOI:10.1111/gji.2011.185.issue-1
[13]
Davis K, Li Y G. Efficient 3D Inversion of Magnetic Data via Octree-Mesh Discretization, Space-Filling Curves, and Wavelets[J]. Geophysics, 2013, 78(5): J61-J73. DOI:10.1190/geo2012-0192.1
[14]
姚长利, 郑元满, 张聿文. 重磁异常三维物性反演随机子域方法技术[J]. 地球物理学报, 2007, 50(5): 1576-1583.
Yao Changli, Zheng Yuanman, Zhang Yuwen. 3-D Gravity and Magnetic Inversion for Physical Properties Using Stochastic Subspaces[J]. Chinese Journal Geophysics, 2007, 50(5): 1576-1583. DOI:10.3321/j.issn:0001-5733.2007.05.035
[15]
Pilkington M. 3-D Magnetic Imaging Using Conjugate Gradients[J]. Geophysics, 1997, 62(4): 1132-1142. DOI:10.1190/1.1444214
[16]
Paige C C, Saunders M A. LSQR:An Algorithm for Aparse Linear Equations and Sparse Least Squares[J]. Acm Transactions on Mathematical Software, 1982, 8(1): 43-71. DOI:10.1145/355984.355989
[17]
Abedi M, Gholami A, Norouzi G H, et al. Fast Inversion of Magnetic Data Using Lanczos Bidiagonalization Method[J]. Journal of Applied Geophysics, 2013, 90: 126-137. DOI:10.1016/j.jappgeo.2013.01.008
[18]
Meng Z H, Li F T, Zhang D L, et al. Fast 3D Inversion of Airborne Gravity-Gradiometry Data Using Lanczos Bidiagonalization Method[J]. Journal of Applied Geophysics, 2016, 132: 211-228. DOI:10.1016/j.jappgeo.2016.07.013
[19]
Rezaie M, Moradzadeh A, Kalateh A N. Fast 3D Inversion of Gravity Data Using Solution Space Priorconditioned Lanczos Bidiagonalization[J]. Journal of Applied Geophysics, 2017, 136: 42-50. DOI:10.1016/j.jappgeo.2016.10.019
[20]
Siripunvaraporn W, Egbert G. An Efficient Data-Subspace Inversion Method for 2-D Magnetotelluric Data[J]. Geophysics, 2000, 65(3): 791-803. DOI:10.1190/1.1444778
[21]
Patella D. Introduction to Ground Surface Self-Potential Tomography[J]. Geophysical Prospecting, 1997, 45(4): 653-681. DOI:10.1046/j.1365-2478.1997.430277.x
[22]
Mauriello P, Patella D. Localization of Maximum-Depth Gravity Anomaly Sources by a Distribution of Equivalent Point masses[J]. Geophysics, 2001, 66(5): 1431-1437. DOI:10.1190/1.1487088
[23]
郭良辉, 孟小红, 石磊, 等. 重力和重力梯度数据三维相关成像[J]. 地球物理学报, 2009, 52(4): 1098-1106.
Guo Lianghui, Meng Xiaohong, Shi Lei, et al. 3-D Correlation Imaging for Gravity Gradiometry Data[J]. Chinese Journal Geophysics, 2009, 52(4): 1098-1106. DOI:10.3969/j.issn.0001-5733.2009.04.027
[24]
孟小红, 刘国峰, 陈召曦, 等. 基于剩余异常成像的相关重磁物性反演方法[J]. 地球物理学报, 2012, 55(1): 304-309.
Meng Xiaohong, Liu Guofeng, Chen Zhaoxi, et al. 3-D Gravity and Magnetic Inversion for Physical Properties Based on Residual Anomaly Correlation[J]. Chinese Journal Geophysics, 2012, 55(1): 304-309. DOI:10.6038/j.issn.0001-5733.2012.01.030
[25]
阮帅, 王绪本, 高永才, 等. 位场正则化下延在自然电位概率成像中的应用[J]. 物化探计算技术, 2006, 28(4): 332-336.
Ruan Shuai, Wang Xuben, Gao Yongcai, et al. Application of Potential Field Normalized Downward Extrapolation in SP Potential Tomography[J]. Computing Techniques for Geophysical and Geochemical Exploration, 2006, 28(4): 332-336.
[26]
王绪本, 高永才.基于正则化向下延拓的重磁概率成像方法研究[C]//《地球物理学报》编辑部.重磁数据处理解释应用研讨会论文集.杭州: 浙江大学出版社, 2008: 49-54.
Wang Xuben, Gao Yongcai. Probability Tomography of Gravity and Magnetic Anomalies Based on Regularized Downward Continuation[C]//"Chinese Journal Geophysics" Editorial Department. Collection of Papers on the Workshop on Gravity and Magnetic Interpretation and Applications. Hangzhou: Zhejiang University Press, 2008: 49-54.
[27]
吴小平, 徐果明. 大地电磁数据的Occam反演改进[J]. 地球物理学报, 1998, 41(4): 547-553.
Wu Xiaoping, Xu Guoming. Improvement of Occam's Inversion for MT Data[J]. Chinese Journal Geophysics, 1998, 41(4): 547-553. DOI:10.3321/j.issn:0001-5733.1998.04.013
[28]
杨文采, 杜剑渊. 层析成像新算法及其在工程检测上的应用[J]. 地球物理学报, 1994, 37(2): 239-244.
Yang Wencai, Du Jianyuan. A New Algorithm of Seismic Tomography with Application to Engineering Detsctions[J]. Chinese Journal Geophysics, 1994, 37(2): 239-244. DOI:10.3321/j.issn:0001-5733.1994.02.012
[29]
蒋成兴, 卢映祥, 陈永清, 等. 滇西南芦子园超大型铅锌多金属矿床成矿模式与综合找矿模型[J]. 地质通报, 2013, 32(11): 1832-1844.
Jiang Chengxing, Lu Yingxiang, Chen Yongqing, et al. Metallogenic Model and Integrated Prospecting Pattern of the Luziyuan Pb-Zn Polymetallic Deposit, Southwestern Yunnan Province[J]. Geological Bulletin of China, 2013, 32(11): 1832-1844. DOI:10.3969/j.issn.1671-2552.2013.11.016
[30]
Liang S X, Jiao Y J, Guo J. Prediction of Hidden Granites in the Luziyuan Area of Yunnan Province and Prospecting Direction[J]. Acta Geologica Sinica (English Edition), 2015, 89(5): 1781-1782. DOI:10.1111/acgs.2015.89.issue-5
[31]
吾守艾力·肉孜, 梁生贤, 邹光富, 等. AMT与重力方法在云南芦子园地区隐伏岩体勘查中的应用[J]. 物探与化探, 2015, 39(3): 525-529.
Wushouaili Rouzi, Liang Shengxian, Zou Guangfu, et al. The Application of AMT and 3D Gravity Methods to the Prospecting for Concealed Rock Body in the Luziyuan Area[J]. Geophysical and Geochemical Exploration, 2015, 39(3): 525-529.
[32]
董文伟, 陈少玲. 镇康芦子园铅锌矿床特征及成因[J]. 云南地质, 2007, 26(4): 404-410.
Dong Wenwei, Chen Shaoling. The Characteristics & Genesis of Luziyuan Pb-Zn Deposit, Zhenkang[J]. Yunnan Geology, 2007, 26(4): 404-410. DOI:10.3969/j.issn.1004-1885.2007.04.004
http://dx.doi.org/10.13278/j.cnki.jjuese.20170167
吉林大学主办、教育部主管的以地学为特色的综合性学术期刊
0

文章信息

梁生贤
Liang Shengxian
互相关系数自约束的重力三维反演与高效求解
A Self-Constrained 3D Inversion and Efficient Solution of Gravity Data Based on Cross-Correlation Coefficient
吉林大学学报(地球科学版), 2018, 48(5): 1473-1482
Journal of Jilin University(Earth Science Edition), 2018, 48(5): 1473-1482.
http://dx.doi.org/10.13278/j.cnki.jjuese.20170167

文章历史

收稿日期: 2017-12-31

相关文章

工作空间