地球物理学进展  2017, Vol. 32 Issue (2): 781-785   PDF    
带地形三维复电阻率正演研究
薛园兵1, 李京谋2, 李桐林1, 朱成1, 王少博1     
1. 吉林大学地球探测科学与技术学院, 长春 130026
2. 吉林大学地球科学学院, 长春 130061
摘要:在复电阻率法找矿过程中,地形对勘测结果的影响不可忽略.本文利用有限差分法二次场模拟研究了带地形的三维复电阻率正演,对求解区域采用不规则正六面体网格剖分,利用不同块体在地表之上的累加来模拟地形起伏的变化.其正演结果与积分方程法得到的结果相比较,所得结果基本一致,证明本文正演算法的可靠性和准确性.并设计了不同正演响应模型,通过对比分析正演响应结果可以得到:在低频时凹地形对视电阻率曲线影响要比高频时大,而地形对电场相位的影响在高频时更明显.
关键词复电阻率    三维正演    带地形    
Research on forward modeling of three-dimensional complex electromagnetic method with topography
XUE Yuan-bing1 , LI Jing-mou2 , LI Tong-lin1 , ZHU Cheng1 , WANG Shao-bo1     
1. College of Geoexploration Science and Technology, Jilin University, Changchun 130026, China
2. College of Earth Sciences, Jilin University, Changchun 130061, China
Abstract: In the complex resistivity prospecting process, the terrain plays a non negligible effect on the survey. In this paper, the terrain is regarded as abnormal body in air and uses finite difference with the secondary field simulation of the terrain 3D complex resistivity forward. Split the solution area with irregular hexahedral mesh, and the variation of topographic relief is simulated by the cumulative superposition of different blocks on the surface of the earth. Firstly, this method is tested to be right in the volume integral equation forward. Then, we use the design of different parameters of the theoretical model to get the results of forward modeling. The final design of the valley model under different frequency get forward response, comparing these forward response results, it is known that at low frequency the influence of concave terrain on the apparent resistivity curve is larger than that of high frequency, and the effect of topography on the phase of the electric field is more obvious at high frequency.
Key words: complex resistivity     three-dimensional forward modeling     topography    
0 引言

复电阻率法,又称频谱激电法 (SIP) 自从20世纪70年代发展起来后,很多学者都对其一、二、三维的正演做了研究,而对于三维复电阻率的正演,早些时候的研究主要集中于均匀半空间的情况,而对于非均匀半空间则忽略了电磁效应只考虑激电效应.罗延忠等 (1986) 用积分方程法研究了非各向同性二层大地中三维极化体的正演响应,并且分析了激电与电磁效应同时存在对响应影响.张辉等 (2006)在同时考虑激电效应和电磁效应情况下用积分方程法研究了均匀围岩中三维极化体的电磁场响应;蔡军涛等 (2007)采用基于三角单元剖分的有限单元法进行了复电阻率二维数值模拟研究;徐凯军 (2008)用有限元实现了2.5维复电阻率正演模拟;赵楠 (2008) 利用有限元法研究了人工源2.5维复电阻率正演;许仕良 (2009) 利用积分方程法研究了三维复电阻率正演;杨毅 (2010)利用积分方程法实现了基于磁性源的三维复电阻率正演;梁盛军 (2011)用有限差分实现了平地形下的三维复电阻率正演;范翠松 (2013)研究了基于有限元的2.5维复电阻率正演.

上述这些研究基本全是在平地形上进行的,然而起伏地形更符合实际情况,对带地形的复电阻率正演研究更具有实用价值,地形起伏会影响电流分布,从而影响到观测数据,使数据产生畸变,难以反映地下真实地质体分布.我国又是一个多山地国家,山区占国土总面积三分之二,要想在这种复杂的地形条件下获取地下地质体分布情况就急需研究起伏地形下的复电阻率数值模拟及异常特征.在这一方面的研究少之又少,其中:李建平 (2008)利用积分方程法研究了带地形条件下均匀围岩中三维极化体复电祖率电磁场响应;赵广茂 (2009)利用有限元法研究了带地形条件下复电阻率2.5维电磁场响应;尚晓等 (2013)研究了起伏地表条件下复电阻率法2.5维有限元数值模拟的一些规律;范翠松等 (2012)也利用有限元法实现了2.5维复电阻率正演研究;欧鸥 (2015)提出了在起伏地形条件下、三角单元剖分、复电导率分块连续变化、2.5维复电阻率有限元数值模拟方法.他们用积分方程法和有限元法分别实现了三维和2.5维带地形的正演,但是在求解复杂地形的速度以及精度都需要提高.

基于这些需求,作者采用交错网格有限差分法来进行带地形的三维复电阻率正演模拟,通过与三维复电阻率积分方程法正演结果比较,验证了本文算法的正确性.并且对不同正演模型做了不同频率下的正演计算,分析了正演响应特征,总结了带地形三维复电阻率正演模拟规律.

1 正演 1.1 正演理论

从电磁勘探的基础方程麦克斯韦方程组出发来推导复电阻率正演的控制方程.复电阻率法的工作频率相对来说比较低,所以可以忽略位移电流影响.假设地下介质磁导率和空气一样,电磁场的时谐因子为eiωt,其中i=,于是符合三维复电阻率的麦克斯韦方程组可以表示为

(1)
(2)

在此式中,EH分别是电场和磁场强度,μ为磁导率,本文中等同于真空中的磁导率 (μ0=4π×10-7 H/m),σ为不同介质中的电导率,Jp是源的电流密度.

本文正演中没有直接求电磁总场,而是分为求解一次场和二次场,即背景场和散射场,所以上面两式可以写成二次场形式的麦克斯韦方程组:

(3)
(4)

式中,EsHs分别是散射电、磁场强度,其中,Js(r)=(σ(r)-σp(r))Ep(r),Ep为背景电场,σp为背景电导率,Js为异常体处的电流密度.

将 (4) 式代入取旋度后的 (3) 式中得到一个基于二次场的亥姆霍兹方程为

(5)

上式即为三维复电阻率正演的控制方程,对此式进行交错采样的形式离散,可以得到正演线性方程组Ax=b,x是待求解向量,是不同方向的二次电场列向量,b是发射源的一次场形成的等效源向量,A是与模型参数有关的大型线性稀疏复对称矩阵.求解此线性方程时采用狄利克雷边界条件,采用雅克比预处理,最后用QMR迭代得到计算节点二次场.再通过观测数据插值算子和电场旋度方程可以得到任意位置的二次电磁场,并将其与背景场相加得到总电磁场.

1.2 地形

本文采用交错网格有限差分法进行带地形的正演模拟,对起伏地形而言,更精细的地表网格剖分可以更准确地模拟地形对电磁场的影响,虽然不规则单元有限元法在模拟地形方面比有限差分法准确,但是计算速度却比较慢.此外为了减少地形对正演收敛的影响,本文还加入了基于电流密度的散度校正,以实现加快收敛的目的.

对于地形的具体处理如下:将地形视为异常体处理,所以在正演中有两种类型的异常体,分别是地下异常体和地形异常体.这就需要在构建亥姆霍兹方程 (5) 的右端项时,电流密度项中不仅含有地下异常体产生的电流密度,还要加上地形产生的散射电流密度.将两部分电流密度相加代入到右端项中,便生成考虑地形的正演控制方程.

1.3 Cole-Cole模型

以上正演流程是以实电阻率为基础进行推导演算的,考虑到岩矿石的激电效应,由复电阻率理论出发,可以用Cole-Cole模型来计算岩矿石的复电阻率.公式为

其中,ρ(iω) 表示复电阻率,ρ0表示零频电阻率,m表示极化率,τ表示时间常数,c表示频率相关系数.本文中复电阻率就是用上式求得的.将正演方程组中的实电阻率换成上式表示的复电阻率,便得到复电阻率线性方程组.

2 正演结果验证及分析 2.1 纯地形正演响应验证

我们通过设计一个凹地形来验证本文带地形正演的准确性.凹地形如下图所示 (图 1),此地堑的长、宽、高分别为200 m、200 m、100 m,中心点坐标为 (-200, 0, 50),假设均匀大地的电阻率为100 Ω,空气电阻率为108 Ω·m,极化率,频率相关系数,时间常数设置为0,发射源坐标为 (-2000, 0, 0) 发射频率为1 Hz,偶极源长度为10 m,电流大小为10 A,接收测线沿y=0方向布置,x范围从-1000 m到1000 m.

图 1 凹地形验证模型 Figure 1 A concave terrain model to verify

对此凹地形进行正演模拟得到积分方程和有限差分的正演结果,并进行比较得到下列视电阻率以及电场相位结果曲线图.

从图中不难看出,在凹地形处,两种正演方法的电场相位都有稍微的变化,并且有限差分方法跟积分方程法得到的结果基本一致.而图 2中的视电阻率更能体现出这一点,两种方法都表现出一种从平地形到凹地形中点时候视电阻率先降低再升高的趋势,并且曲线形态一致,说明本文的有限差分算法在起伏地形下的正演结果也是正确的.

图 2 视电阻率曲线图 Figure 2 The figure of apparent resistivity curve

图 3 电场相位曲线图 Figure 3 The Phase curve of electric field
2.2 正演结果特征分析

首先,我们先研究一下平地形下异常体的复电阻率响应.为此,我们设计如图 4的模型一:地下含有一个低阻复电阻率模型,复电阻率参数为:零频电阻率10 Ω·m;极化率0.3;时间常数100;频率相关系数0.3,低阻体长、宽均为300 m,高180 m,背景场电阻率为100 Ω·m,发射源长度为100 m,电流强度10 A,中心点坐标位于 (-1500, 0, 0),接收测线过中心原点O,沿着y=0的方向放置,其范围是x=-1000 m到1000 m.

图 4 平地形下异常体模型 Figure 4 The abnormal body model in the flat ground

接着我们讨论在地形起伏情况下正演响应会有如何变化,我们设计如下正演模型二.图 5是我们所设计的凹地形模型示意图.凹地形的长宽高为200 m、200 m、90 m.发射源长度为100 m,电流强度10 A,中心点坐标位于 (-1500, 0, 0),背景电阻率为100 Ω·m.在模型二基础上在距凹地形底部60 m处放置一个异常体得到模型三,凹地形下低阻异常体的复电阻率参数为:零频电阻率10 Ω·m;极化率0.3;时间常数100;频率相关系数0.3,低阻体长、宽均为300 m,高180 m,接收测线过中心原点O,沿着y=0的方向放置,其范围是x=-1000 m到1000 m.

图 5 凹地形模型 Figure 5 The abnormal body model under concave topography

设置发射频率为1 Hz,上述两个模型正演得到的视电阻率曲线如图 6所示.从图 6不难看出,相较于平地形而言,在低阻体上面放置一个凹地形会使异常体处的正演视电阻率偏高,所以凹地形会在正演模拟中起到一个高阻的作用.接下来我们讨论在凹地形情况下,不同参数对正演响应结果的影响,依然采用模型二的凹地形.

图 6 视电阻率曲线比较 Figure 6 The comparison of apparent resistivity curves

在不考虑含有地下异常体,且其他收发参数不变情况下将背景复电阻率参数设置为:零频电阻率100 Ω·m;极化率0.3;时间常数100;频率相关系数0.3,得到平地形和凹地形下的电场相位曲线图如上图 7所示,可以看出地形对相位影响很小;其他参数不变情况下,改变发射频率得到上图 8结果,可以看出随着频率增加,相位值有变小的趋势,并且在地形处异常也比较明显.

图 7 电场相位曲线图比较 Figure 7 The comparison of electric field phase curves

图 8 不同频率下电场相位曲线图 Figure 8 The electric phase curves under different frequencies

使模型三在固定其他参数不变的情况下,将发射频率 (Hz) 分别设置为0.01、1、100,得到正演响应如图 9所示.同样的参数设定下,平地形正演模拟得到结果如图 10所示.将平地形结果视为背景场得到的结果,用这个结果抵消掉凹地形中的背景场,便得到不同频率下凹地形对结果的影响,如图 11所示.从图 11来看,在频率较低时,地形对视电阻率的影响更大一些.

图 9 凹地形不同发射频率下正演视电阻率曲线 Figure 9 The apparent resistivity curves at different frequencies in concave terrain

图 10 平地形下不同发射频率下正演视电阻率曲线 Figure 10 The apparent resistivity curves at different frequencies in the flat shape

图 11 不同频率下地形对正演模拟影响 Figure 11 The influence of topography on the apparent resistivity curve at different frequencies
3 结论

本文通过有限差分法得到了带地形三维复电阻率正演结果,并用积分方程法验证了本算法的正确性.通过设计不同模型,得出如下结论:

(1) 低频时凹地形对相位的影响不明显,但是随着发射频率的增加,凹地形对相位影响开始变大.

(2) 对视电阻率来说,在低频时地形对其影响要大于高频时的影响.本文通过不同正演模型,总结了带地形下正演响应规律,为野外实际工作提供了依据.

致谢 感谢审稿专家提出的修改意见和编辑部的大力支持!
参考文献
[] CAI Jun-Tao, RUAN Bai-Yao, ZHAO Guo-Ze, et al. 2007. Two-dimensional modeling of complex resistivity using finite element method[J]. Chinese Journal of Geophysics (in Chinese), 50(6): 1869–1876. DOI:10.3321/j.issn:0001-5733.2007.06.030
[] CHEN Hui, DENG Ju-Zhi, TAN Han-Dong, et al. 2011. Study on divergence correction method in three-dimensional magnetotelluric modeling with staggered-grid finite difference method[J]. Chinese Journal of Geophysics (in Chinese), 54(6): 1649–1659. DOI:10.3969/j.issn.0001-5733.2011.06.025
[] FAN Cui-Song. 2013. Research on complex resistivity forward and inversion with finite element method and its application (in Chinese)[Ph. D. thesis]. Changchun:Jilin University.
[] FAN Cui-Song, LI Tong-Lin, YAN Jia-Yong. 2012. Research and application on 2.5D SIP inversion[J]. Chinese Journal of Geophysics (in Chinese), 55(12): 4044–4050. DOI:10.6038/j.issn.0001-5733.2012.12.016
[] LI Jian-Ping. 2008. Research on electromagnetic modeling of 3D complex resistivity with topography (in Chinese)[Ph.D.thesis]. Changchun:Jilin University.
[] LIANG Sheng-Jun. 2011. Research on complex resistivity 3D modeling and inversing problem (in Chinese)[Ph. D. thesis]. Beijing:China University of Geosciences (Beijing).
[] LIU Yun-He. 2011. Research on 3-D controlled source electromagnetic method inversion using nonlinear conjugate gradients (in Chinese)[Ph. D. thesis]. Changchun:Jilin University.
[] OU Ou. 2015. 2.5 dimensional numerica modeling and inversion of the frequency domain resistivity method under the topography condition (in Chinese)[Ph. D. thesis]. Chengdu:Chengdu University of Technology.
[] SHANG Xiao, LI Yong, LIN Pin-Rong, et al. 2013. A study on finite element modeling and anomaly features of complex resistivity method with surface topography[J]. Computing Techniques for Geophysical and Geochemical Exploration (in Chinese), 35(1): 1–9.
[] SHEN Jin-Song. 2003. Modeling of 3-D electromagnetic responses in frequency domain by using staggered grid finite difference method[J]. Chinese Journal of Geophysics (in Chinese), 46(2): 281–288. DOI:10.3321/j.issn:0001-5733.2003.02.024
[] XU Kai-Jun. 2007. Study on 2.5D complex resistivity electromagnetic forward and inversion (in Chinese)[Ph. D. thesis]. Changchun:Jilin University.
[] XU Kai-Jun, LI Tong-Lin. 2006. The forward modeling of three-dimensional geoelectric field of vertical finite line source by finite-difference method[J]. Journal of Jilin University (Earth Science Edition) (in Chinese), 36(1): 137–141, 147.
[] YANG Yi. 2010. 3D electromagnetic simulation of complex resistivity of magnetic source (in Chinese)[MSc thesis]. Changchun:Jilin University.
[] ZHANG Hui. 2006. Research of complex resistivity 3D electromagnetic forward and inversion (in Chinese)[Ph. D. thesis]. Changchun:Jilin University.
[] ZHANG Hui, LI Tong-Lin, DONG Rui-Xia. 2006. Modeling electromagnetic responses of complex resistivity 3-D body using volume integral equation method[J]. Coal Geology & Exploration (in Chinese), 34(1): 70–74.
[] ZHAO Guang-Mao. 2009. Research of complex resistivity 2.5D electromagnetic forward and inversion with topography (in Chinese)[Ph. D. thesis]. Changchun:Jilin University.
[] ZHU Cheng. 2016. Research on modeling and inversion of three-dimensional controlled source electromagnetic method with topography in frequency-domain (in Chinese)[MSc thesis]. Changchun:Jilin University.
[] ZHU Cheng, LI Tong-lin, FAN Cui-song. 2014. Research on divergence correction method in 3D numerical modeling of 3D controlled source electromagnetic fields[J]. Global Geology (in Chinese), 33(4): 916–924.
[] 蔡军涛, 阮百尧, 赵国泽, 等. 2007. 复电阻率法二维有限元数值模拟[J]. 地球物理学报, 50(6): 1869–1876. DOI:10.3321/j.issn:0001-5733.2007.06.030
[] 陈辉, 邓居智, 谭捍东, 等. 2011. 大地电磁三维交错网格有限差分数值模拟中的散度校正方法研究[J]. 地球物理学报, 54(6): 1649–1659. DOI:10.3969/j.issn.0001-5733.2011.06.025
[] 范翠松. 2013. 基于有限元法的复电阻率正反演研究及应用[博士论文]. 长春: 吉林大学. http: //cdmd. cnki. com. cn/Article/CDMD-10183-1013193112. htm
[] 范翠松, 李桐林, 严加永. 2012. 2.5维复电阻率反演及其应用试验[J]. 地球物理学报, 55(12): 4044–4050. DOI:10.6038/j.issn.0001-5733.2012.12.016
[] 李建平. 2008. 带地形的三维复电阻率电磁场正反演研究[博士论文]. 长春: 吉林大学. http: //cdmd. cnki. com. cn/Article/CDMD-10183-2008126894. htm
[] 梁盛军. 2011. 复电阻率法三维正反演问题研究[博士论文]. 北京: 中国地质大学 (北京). http: //cdmd. cnki. com. cn/Article/CDMD-11415-1011077512. htm
[] 刘云鹤. 2011. 三维可控源电磁法非线性共轭梯度反演研究[博士论文]. 长春: 吉林大学. http: //cdmd. cnki. com. cn/Article/CDMD-10183-1011107232. htm
[] 欧鸥. 2015. 起伏地形条件下2. 5维复电阻率法数值模拟与反演成像研究[博士论文]. 成都: 成都理工大学. http: //cdmd. cnki. com. cn/Article/CDMD-10616-1015312796. htm
[] 尚晓, 李勇, 林品荣, 等. 2013. 起伏地表复电阻率法有限元数值模拟及异常特征[J]. 物探化探计算技术, 35(1): 1–9.
[] 沈金松. 2003. 用交错网格有限差分法计算三维频率域电磁响应[J]. 地球物理学报, 46(2): 281–288. DOI:10.3321/j.issn:0001-5733.2003.02.024
[] 徐凯军. 2007. 2. 5维复电阻率电磁场正反演研究[博士论文]. 长春: 吉林大学. http: //cdmd. cnki. com. cn/Article/CDMD-10183-2006109662. htm
[] 徐凯军, 李桐林. 2006. 垂直有限线源三维地电场有限差分正演研究[J]. 吉林大学学报 (地球科学版), 36(1): 137–141, 147.
[] 杨毅. 2010. 磁偶源复电阻率三维电磁模拟研究[硕士论文]. 长春: 吉林大学. http: //cdmd. cnki. com. cn/Article/CDMD-10183-2010110730. htm
[] 张辉, 李桐林, 董瑞霞. 2006. 体积分方程法模拟复电阻率三维体电磁响应[J]. 煤田地质与勘探, 34(1): 70–74.
[] 张辉. 2006. 复电阻率三维电磁场正反演研究[博士论文]. 长春: 吉林大学. http: //cdmd. cnki. com. cn/Article/CDMD-10183-2006109662. htm
[] 赵广茂. 2009. 带地形的复电阻率2. 5维电磁场正反演研究[博士论文]. 长春: 吉林大学. http: //cdmd. cnki. com. cn/Article/CDMD-10183-2008126894. htm
[] 朱成. 2016. 带地形频率域可控源电磁法三维正反演研究[硕士论文]. 长春: 吉林大学. http: //cdmd. cnki. com. cn/Article/CDMD-10491-1016312063. htm
[] 朱成, 李桐林, 范翠松. 2014. 可控源电磁场三维数值模拟中的散度校正方法研究[J]. 世界地质, 33(4): 916–924.