﻿ 非结构化网格与散度校正技术结合的三维CSEM矢量有限元正演
 石油地球物理勘探  2021, Vol. 56 Issue (4): 891-901  DOI: 10.13810/j.cnki.issn.1000-7210.2021.04.022 0
 文章快速检索 高级检索

### 引用本文

TANG Wenwu, DENG Juzhi, HUANG Qinghua. Three-dimensional CSEM forward modeling using edge-based finite element method based on unstructured meshes and divergence correction. Oil Geophysical Prospecting, 2021, 56(4): 891-901. DOI: 10.13810/j.cnki.issn.1000-7210.2021.04.022.

### 文章历史

① 北京大学地球与空间科学学院, 北京 100871;
② 东华理工大学核资源与环境国家重点实验室, 江西南昌 330013

Three-dimensional CSEM forward modeling using edge-based finite element method based on unstructured meshes and divergence correction
TANG Wenwu①② , DENG Juzhi , HUANG Qinghua
① School of Earth and Space Sciences, Peking University, Beijing 100871, China;
② State Key Laboratory of Nuclear Resources and Environment, East China University of Technology, Nanchang, Jiangxi 330013, China
Abstract: The iterative computation with a system of linear equations derived from the three-dimensional (3D) controlled-source electromagnetic(CSEM) forward modeling of the electric field equation suffers from slow convergence. Moreover, unstructured meshes can make the system of linear equations more ill-posed. In view of this, we propose an algorithm for finite element forward modeling based on unstructured tetrahedral meshes and divergence correction. Starting from the divergence equation of current density, we derive the corrected divergence equation of the potential on geo-electrical interfaces. Solving the system of linear equations is accelerated with the preconditioned quasi-minimal residual (QMR) method and the alternate divergence correction during the iteration. A three-layer medium model is subjected to the forward modeling under two conditions (with/without divergence correction) to verify the reliability of the proposed algorithm. The iterative convergence and the accuracy of numerical solution of the system of linear equations indicate that the divergence correction is effective to accelerate the iteration and improve the forward modeling accuracy. On this basis, a 3D geo-electric model is built, the electromagnetic response of which is employed for the comparison of numerical solutions between the proposed algorithm and the forward modeling based on the quadratic coupling potential equation. It further confirms the high accuracy of the algorithm in this study. The modeling of a complex oil and gas monitoring system demonstrates the application potential of the CSEM method in oil and gas monitoring.
Keywords: three-dimensional controlled-source EM    unstructured meshes    divergence correction    forward modeling
0 引言

1 方法原理 1.1 三维电磁正演边值问题

 $\nabla \times {\boldsymbol{E}} = - \frac{{\partial {\boldsymbol{B}}}}{{\partial t}} = {{\rm{i}}_{\omega \mu }}_0{\boldsymbol{H}}$ (1)
 $\nabla \times {\boldsymbol{H}} = \sigma {\boldsymbol{E}} + {{\boldsymbol{J}}_外}$ (2)

 $\nabla \times \nabla \times {\boldsymbol{E}} = {\rm{i}}\omega {\mu _0}\nabla \times {\boldsymbol{H}}$ (3)

 $\nabla \times \nabla \times {\boldsymbol{E}} - {\rm{i}}\omega {\mu _0}\sigma {\boldsymbol{E}} = {{\rm{i}}_\omega }{\mu _0}{{\boldsymbol{J}}_外}$ (4)

 ${\boldsymbol{E}} = {{\boldsymbol{E}}_{\rm{p}}} + {{\boldsymbol{E}}_{\rm{s}}}$ (5)

 $\nabla \times \nabla \times {{\boldsymbol{E}}_{\rm{p}}} - {\rm{i}}\omega {\mu _0}{\sigma _{\rm{p}}}{{\boldsymbol{E}}_{\rm{p}}} = {\rm{i}}\omega {\mu _0}{{\boldsymbol{J}}_外}$ (6)

 $\nabla \times \nabla \times {{\boldsymbol{E}}_{\rm{s}}} - {\rm{i}}\omega {\mu _0}\sigma {{\boldsymbol{E}}_{\rm{s}}} = {\rm{i}}\omega {\mu _0}\Delta \sigma {{\boldsymbol{E}}_{\rm{p}}}$ (7)

 $\left\{ {\begin{array}{*{20}{l}} {\nabla \times \nabla \times {{\boldsymbol{E}}_s} - {\mathop{\rm i}\nolimits} \omega {\mu _0}\sigma {{\boldsymbol{E}}_s} = {\mathop{\rm i}\nolimits} \omega {\mu _0}\Delta \sigma {{\boldsymbol{E}}_{\rm{p}}}\;\;\;\;\; \in \mathit{\Omega }}\\ {{\boldsymbol{n}} \times {{\boldsymbol{E}}_{\rm{s}}} = 0\;\;\;\;\;\; \in \mathit{\Gamma }} \end{array}} \right.$ (8)

1.2 基于伽辽金有限元法的边值问题求解

 图 1 研究区域示意图 Idl代表接地电流源，黑点代表测点

 图 2 四面体剖分单元示意图 数字1~4表示节点编号，N1~N6表示对应各条棱边的矢量形函数

 ${{\boldsymbol{E}}_{\rm{p}}} = \sum\limits_{i = 1}^6 {{{\boldsymbol{N}}_i}} {E_{{\rm{p}}i}}$ (9)
 ${{\boldsymbol{E}}_{\rm{s}}} = \sum\limits_{i = 1}^6 {{{\boldsymbol{N}}_i}} {E_{{\rm{s}}i}}$ (10)

 ${{\boldsymbol{N}}_i} = \left( {{{\boldsymbol{L}}_{{i_1}}}\nabla {{\boldsymbol{L}}_{{i_2}}} - {{\boldsymbol{L}}_{{i_2}}}\nabla {{\boldsymbol{L}}_{{i_1}}}} \right){l_i}$ (11)

 $\begin{array}{l} \int_\mathit{\Omega} {\left( {\nabla \times \nabla \times {{\boldsymbol{E}}_{\rm{s}}} - {\rm{i}}\omega {\mu _0}\sigma {{\boldsymbol{E}}_{\rm{s}}}} \right)} \cdot {{\boldsymbol{N}}_i}{\rm{d}}\mathit{\Omega }\\ \;\;\;\;\;\;\;\; = \int_\mathit{\Omega} {\rm{i}} \omega {\mu _0}\Delta \sigma {{\boldsymbol{E}}_{\rm{p}}} \cdot {{\boldsymbol{N}}_i}{\rm{d}}\mathit{\Omega } \end{array}$ (12)

 $\begin{array}{l} {\rm{LHS}} = \int_\mathit{\Omega} {{{\boldsymbol{N}}_i}} \cdot \left( {\nabla \times \nabla \times {{\boldsymbol{E}}_{\rm{s}}}} \right){\rm{d}}\mathit{\Omega } - \\ \;\;\;\;\;\;\;\;\;\;\;\sum\limits_{j = 1}^6 {\int_\mathit{\Omega} {\rm{i}} } \omega {\mu _0}\sigma {{\boldsymbol{N}}_i} \cdot {{\boldsymbol{N}}_j}{{\boldsymbol{E}}_{{\rm{s}}j}}{\rm{d}}\mathit{\Omega } \end{array}$ (13)

 $- {\boldsymbol{A}} \cdot (\nabla \times {\boldsymbol{B}}) + {\boldsymbol{B}} \cdot (\nabla \times {\boldsymbol{A}}) = \nabla \cdot ({\boldsymbol{A}} \times {\boldsymbol{B}})$

 (14)

 (15)

 $\begin{array}{l} \sum\limits_{j = 1}^6 {\left[ {\int_\mathit{\Omega} {\left( {\nabla \times {{\boldsymbol{N}}_i}} \right)} \cdot \left( {\nabla \times {{\boldsymbol{N}}_j}} \right){\rm{d}}\mathit{\Omega} - } \right.} \\ \;\;\;\;\left. {\int_\mathit{\Omega} {\rm{i}} \omega {\mu _0}\sigma {{\boldsymbol{N}}_i} \times {{\boldsymbol{N}}_j}{\rm{d}}\mathit{\Omega} } \right]{E_{{\rm{s}}j}}\\ \;\;\;\; = {E_{{\rm{p}}j}}\left( {\sum\limits_{j = 1}^6 {\rm{i}} \omega {\mu _0}\int_\mathit{\Omega} \Delta \sigma {{\boldsymbol{N}}_i} \cdot {{\boldsymbol{N}}_j}{\rm{d}}\mathit{\Omega} } \right) \end{array}$ (16)

 ${\boldsymbol{K}}{{\boldsymbol{E}}_{\rm{s}}} = {\boldsymbol{M}}{{\boldsymbol{E}}_{\rm{p}}}$ (17)

 $\left\{ {\begin{array}{*{20}{l}} {{K_{ij}} = \int_\mathit{\Omega} {\left[ {\left( {\nabla \times {{\boldsymbol{N}}_i}} \right) \cdot \left( {\nabla \times {{\boldsymbol{N}}_j}} \right) - {\mathop{\rm i}\nolimits} \omega {\mu _0}\sigma {{\boldsymbol{N}}_i} \cdot {{\boldsymbol{N}}_j}} \right]} {\rm{d}}\mathit{\Omega} }\\ {{M_{ij}} = {\rm{i}}\omega {\mu _0}\Delta \sigma \int_\mathit{\Omega} {{{\boldsymbol{N}}_i}} \cdot {{\boldsymbol{N}}_j}{\rm{d}}\mathit{\Omega} } \end{array}} \right.$ (18)

1.3 散度校正技术

 $\nabla \cdot \left( {\sigma {{\boldsymbol{E}}_{\rm{s}}} + \Delta \sigma {{\boldsymbol{E}}_{\rm{p}}}} \right) = 0$ (19)

 $\nabla \cdot (\sigma \nabla \phi ) = 0\;\;\;\; \in \Omega$ (20)

 $\begin{array}{l} \;\;\;\;\;\;\;\;\;\;{\boldsymbol{n}} \cdot \left( {{\sigma _i}\nabla {\phi _i} - {\sigma _j}\nabla {\phi _j}} \right)\left| {{s_{ij}} = } \right.\\ {\boldsymbol{n}} \cdot \left( {{\sigma _i}{{\boldsymbol{E}}_{si}} - {\sigma _j}{{\boldsymbol{E}}_{sj}} + \Delta {\sigma _i}{{\boldsymbol{E}}_{pi}} - \Delta {\sigma _j}{{\boldsymbol{E}}_{pj}}} \right)\;\;\; \in \Gamma \end{array}$ (21)
 图 3 电导率分界面上的散度校正示意图

 $\phi = \sum\limits_{k = 1}^M {{\phi _k}} {v_k}$ (22)

 ${{\boldsymbol{E}}_{\rm{p}}} = \sum\limits_{j = 1}^N {{E_{{\rm{p}}j}}} {{\boldsymbol{\omega }}_j}$ (23)
 ${{\boldsymbol{E}}_{\rm{s}}} = \sum\limits_{j = 1}^N {{E_{{\rm{s}}j}}} {{\boldsymbol{\omega }}_j}$ (24)

 ${\boldsymbol{Dp}} = {\boldsymbol{t}}$ (25)

 $\left\{ \begin{array}{l} {D_{hk}} = \int_\mathit{\Omega} \sigma \nabla {v_h} \cdot \nabla {v_k}{\rm{d}}\mathit{\Omega} \\ {P_k} = {\phi _k}\\ {t_h} = \sum\limits_{j = 1}^N {\left( {{E_{sj}}\int_S {{\sigma _h}} {{\boldsymbol{\omega }}_j} \cdot {\boldsymbol{n}}{\rm{d}}{\boldsymbol{S}} + } \right.} \\ \;\;\;\;\;\;\left. {{E_{{\rm{p}}j}}\int_S \sigma {v_h}{{\boldsymbol{\omega }}_j} \cdot {\boldsymbol{n}}{\rm{d}}{\boldsymbol{S}}} \right)\;\;\;\;\;h, k = 1, \cdots , M \end{array} \right.$

 ${\boldsymbol{E}}_{\rm{s}}^{\rm{c}} = {{\boldsymbol{E}}_{\rm{s}}} - \nabla \phi$ (26)

2 数值模拟实验

2.1 水平层状介质模型

 图 4 水平三层介质模型yOz面示意图

 图 5 水平层状介质模型线性方程组迭代求解收敛曲线对比 (a)无散度校正；(b)施加散度校正

 图 6 水平层状介质模型二次电场分量数值解与解析解对比(1Hz) (a)Exs实部；(b)Exs虚部；(c)Eys实部；(d)Eys虚部
2.2 三维低阻体模型

 图 7 三维地电模型示意图

 图 8 三维地电模型网格剖分切面图 (a)全区域yOz切面；(b)全区域xOy切面；(c)图a中蓝色方框区域放大图；(d)图b中蓝色方框区域放大图紫红色区域代表空气，绿色区域代表地下围岩，蓝色区域代表低阻体，白色线条代表地面测线位置

 图 9 三维地电模型不同频率时的迭代收敛曲线 (a)无散度校正；(b)施加散度校正。图中标注的时间t为模拟耗时

 图 10 基于二次耦合势(a)及二次电场(b)计算的8Hz时Exs数值解实部(上)及虚部(下)

 图 11 8Hz时y=3km剖面分别基于二次电场和二次耦合势计算得到的Exs数值解对比 (a)实部分量；(b)虚部分量；(c)实部分量相对误差；(d)虚部分量相对误差

 图 12 32Hz时y=3km剖面分别基于二次电场和二次耦合势计算得到的Exs数值解对比 (a)实部分量；(b)虚部分量

 图 13 16Hz时y=5km剖面分别基于二次电场和二次耦合势计算得到的Exs数值解对比 (a)实部分量；(b)虚部分量；(c)实部分量相对误差；(d)虚部分量相对误差

 图 14 8Hz时y=3km剖面分别基于二次电场和二次耦合势计算得到的Eys数值解对比 (a)实部分量；(b)虚部分量

 图 15 32Hz时y=3km剖面分别基于二次电场和二次耦合势计算得到的Eys数值解对比 (a)实部分量；(b)虚部分量
2.3 三维油气监测模型

 图 16 三维水驱油气监测模型yOz断面图

 图 17 三维油气监测模型电场水平分量Ex(上)和Ey(下)振幅 (a)油气占比100%；(b)油气占比50%；(c)无油气

 图 18 三维油气监测模型x=-1.75km(上)和x=-1.00km(下)剖面电场分量Ex(a)和Ey(b)振幅曲线
3 结论

 [1] Nabighian M, Corbett J. Electromagnetic Methods in Applied Geophysics: Theory[M]. SEG, 1988. [2] He Z, Liu X, Qiu W, et al. Mapping reservoir boundary by borehole-surface TFEM: Two case studies[J]. The Leading Edge, 2005, 24(9): 896-900. DOI:10.1190/1.2056379 [3] 邓居智, 陈辉, 殷长春, 等. 九瑞矿集区三维电性结构研究及找矿意义[J]. 地球物理学报, 2015, 58(12): 4465-4477. DENG Juzhi, CHEN Hui, YIN Changchun, et al. Three-dimensional electrical structures and significance for mineral exploration in the Jiujiang-Ruichang District[J]. Chinese Journal of Geophysics, 2015, 58(12): 4465-4477. DOI:10.6038/cjg20151211 [4] 汤井田, 何继善. 可控源音频大地电磁法及其应用[M]. 湖南长沙: 中南大学出版社, 2005. [5] Wannamaker P E, Hohmann G W, and San Filipo W A. Electromagnetic modeling of three-dimensional bodies in layered earths using integral equations[J]. Geophysics, 1984, 49(1): 60-74. DOI:10.1190/1.1441562 [6] Wannamaker P E. Advances in three-dimensional magnetotelluric modeling using integral equations[J]. Geophysics, 1991, 56(11): 1716-1728. DOI:10.1190/1.1442984 [7] 李静和, 何展翔, 孟淑君, 等. 三维地形频率域井筒电磁场区域积分方程法模拟[J]. 物理学报, 2019, 68(14): 202-211. LI Jinghe, HE Zhanxiang, MENG Shujun, et al. Domain decomposition based integral equation modeling of 3D topography in frequency domain for well electromagnetic[J]. Acta Physica Sinica, 2019, 68(14): 202-211. [8] Newman G A, Alumbaugh D L. Frequency-domain modelling of airborne electromagnetic responses using staggered finite differences[J]. Geophysical Prospecting, 1995, 43(8): 1021-1042. DOI:10.1111/j.1365-2478.1995.tb00294.x [9] 沈金松. 用交错网格有限差分法计算三维频率域电磁响应[J]. 地球物理学报, 2003, 46(2): 281-288. SHEN Jinsong. Modeling of 3D electromagnetic responses in frequency domain by using staggered-grid finite difference method[J]. Chinese Journal of Geophysics, 2003, 46(2): 281-288. DOI:10.3321/j.issn:0001-5733.2003.02.024 [10] Streich R. 3D finite-difference frequency-domain mode-ling of controlled-source electromagnetic data: Direct solution and optimization for high accuracy[J]. Geophysics, 2009, 74(5): F95-F105. DOI:10.1190/1.3196241 [11] 邓居智, 谭捍东, 陈辉, 等. CSAMT三维交错采样有限差分数值模拟[J]. 地球物理学进展, 2011, 26(6): 2026-2032. DENG Juzhi, TAN Handong, CHEN Hui, et al. CSAMT 3D modeling using staggered grid finite difference method[J]. Progress in Geophysics, 2011, 26(6): 2026-2032. DOI:10.3969/j.issn.1004-2903.2011.06.017 [12] Haber E, Ascher U, Aruliah D, et al. Fast simulation of 3D electromagnetics problems using potentials[J]. Journal of Computational Physics, 2000, 163: 150-171. DOI:10.1006/jcph.2000.6545 [13] Jahandari H and Farquharson C. Finite-volume model-ling of geophysical electromagnetic data on unstructured grids using potentials[J]. Geophysical Journal International, 2015, 202: 1859-1876. DOI:10.1093/gji/ggv257 [14] 陈辉, 殷长春, 邓居智. 基于Lorenz规范条件下磁矢势和标势耦合方程的频率域电磁法三维正演[J]. 地球物理学报, 2016, 59(8): 3087-3097. CHEN Hui, YIN Changchun, DENG Juzhi. A finite volume solution to 3D frequency-domain electromagnetic modelling using Lorenz-gauged magnetic vector and scalar potentials[J]. Chinese Journal of Geophy-sics, 2016, 59(8): 3087-3097. [15] Badea E A, Everett M E, Newman G A, et al. Finite-element analysis of controlled-source electromagnetic induction using Coulomb-gauged potentials[J]. Geophysics, 2001, 66(3): 786-799. DOI:10.1190/1.1444968 [16] 徐志锋, 吴小平. 可控源电磁三维频率域有限元模拟[J]. 地球物理学报, 2010, 53(8): 1931-1939. XU Zhifeng, WU Xiaoping. Controlled source electromagnetic 3-D modeling in frequency domain for finite element method[J]. Chinese Journal of Geophysics, 2010, 53(8): 1931-1939. DOI:10.3969/j.issn.0001-5733.2010.08.019 [17] Puzyrev V, Koldan J, Puente J, et al. A parallel finite-element method for three-dimensional controlled-source electromagnetic forward modelling[J]. Geophysical Journal International, 2013, 193(2): 678-693. DOI:10.1093/gji/ggt027 [18] Ansari S, Farquharson C G. 3D finite-element forward modeling of electromagnetic data using vector and scalar potentials and unstructured grids[J]. Geophy-sics, 2014, 79(4): E149-E165. [19] Cai H, Hu X, Li J, et al. Parallelized 3D CSEM mode-ling using edge-based finite element with total field formulation and unstructured mesh[J]. Computers & Geosciences, 2017, 99: 125-134. [20] 刘长生, 汤井田, 任政勇, 等. 基于非结构化网格的三维大地电磁自适应矢量有限元模拟[J]. 中南大学学报(自然科学版), 2010, 41(5): 1855-1859. LIU Changsheng, TANG Jingtian, REN Zhengyong, et al. Three-dimension magnetotellurics modeling by adaptive edge finite-element using unstructured meshes[J]. Journal of Central South University (Science and Technology), 2010, 41(5): 1855-1859. [21] 杨军, 刘颖, 吴小平. 海洋可控源电磁三维非结构矢量有限元数值模拟[J]. 地球物理学报, 2015, 58(8): 2827-2838. YANG Jun, LIU Ying, WU Xiaoping. 3D simulation of marine CSEM using vector finite element method on unstructured grids[J]. Chinese Journal of Geophysics, 2015, 58(8): 2827-2838. [22] 殷长春, 张博, 刘云鹤, 等. 面向目标自适应三维大地电磁正演模拟[J]. 地球物理学报, 2017, 60(1): 327-336. YIN Changchun, ZHANG Bo, LIU Yunhe, et al. A goal-oriented adaptive algorithm for 3D magnetotelluric forward modeling[J]. Chinese Journal of Geophy-sics, 2017, 60(1): 327-336. [23] 刘颖, 李予国, 韩波. 可控源电磁场三维自适应矢量有限元正演模拟[J]. 地球物理学报, 2017, 60(12): 4874-4886. LIU Ying, LI Yuguo, HAN Bo. Adaptive edge finite element modeling of the CSEM field on unstructured grids[J]. Chinese Journal of Geophysics, 2017, 60(12): 4874-4886. DOI:10.6038/cjg20171227 [24] 叶益信, 李予国, 刘颖, 等. 基于局部加密非结构网格的海洋可控源电磁法三维有限元正演[J]. 地球物理学报, 2016, 59(12): 4747-4758. YE Yixin, LI Yuguo, LIU Ying, et al. 3D finite element modeling of marine controlled-source electromagnetic fields using locally refined unstructured meshes[J]. Chinese Journal of Geophysics, 2016, 59(12): 4747-4758. DOI:10.6038/cjg20161233 [25] 赵宁, 王绪本, 余刚, 等. 面向目标自适应海洋可控源电磁三维矢量有限元正演[J]. 地球物理学报, 2019, 62(2): 779-788. ZHAO Ning, WANG Xuben, YU Gang, et al. 3D MCSEM parallel goal-oriented adaptive vector finite element modeling[J]. Chinese Journal of Geophysics, 2019, 62(2): 779-788. [26] Smith J T. Conservative modeling of 3D electromagnetic fields, part Ⅱ: Biconjugate gradient solution and an accelerator[J]. Geophysics, 1996, 61(5): 1319-1324. DOI:10.1190/1.1444055 [27] 陈辉, 邓居智, 谭捍东, 等. 大地电磁三维交错网格有限差分数值模拟中的散度校正方法研究[J]. 地球物理学报, 2011, 54(6): 1649-1659. CHEN Hui, DENG Juzhi, TAN Handong, et al. Study on divergence correction method in three dimensional magnetotelluric modeling with staggered-grid finite difference method[J]. Chinese Journal of Geophysics, 2011, 54(6): 1649-1659. DOI:10.3969/j.issn.0001-5733.2011.06.025 [28] Farquharson C G, Miensopust M P. Three-dimensional finite-element modelling of magnetotelluric data with a divergence correction[J]. Journal of Applied Geophysics, 2011, 75: 699-710. DOI:10.1016/j.jappgeo.2011.09.025 [29] Liu Y H, Yin C C. Electromagnetic divergence correction for 3D anisotropic EM modeling[J]. Journal of Applied Geophysics, 2013, 77: 19-27. [30] 秦策, 王绪本, 赵宁, 等. 频率域电磁法三维有限元正演线性方程组迭代算法[J]. 地球物理学报, 2020, 63(8): 3180-3190. QIN Ce, WANG Xuben, ZHAO Ning, et al. Research on the iterative solver of linear equations in three-dimensional finite element forward modeling for frequency domain electromagnetic method[J]. Chinese Journal of Geophysics, 2020, 63(8): 3180-3190. [31] 汤文武, 李耀国, 柳建新, 等. 基于二次电场的可控源电磁法三维矢量有限元正演模拟[J]. 石油物探, 2015, 54(6): 665-673. TANG Wenwu, LI Yaoguo, LIU Jianxin, et al. Three-dimensional controlled-source electromagnetic forward modeling by edge-based finite element using secondary electrical field[J]. Geophysical Prospecting for Petroleum, 2015, 54(6): 665-673. [32] Jin J. The Finite Element Method in Electromagne-tics[M]. Wiley Interscience, 2002. [33] Saad Y. Iterative Methods for Sparse Linear Systems[M]. SIAM, 2003. [34] Hang Si. TetGen: a delaunay-based quality tetrahedral mesh generator[J]. ACM Transactions on Mathematical Software, 2015, 41(2): 1-36. [35] Tang W W, Li Y G, Swidinsky A, et al. Three-dimensional controlled-source electromagnetic modelling with a well casing as a grounded source: a hybrid method of moments and finite element scheme[J]. Geo-physical Prospecting, 2015, 63(4): 1491-1507.