复电阻率法是根据岩、矿石的复电阻率频谱参数的差异来评价电性异常体,也被称为频谱激电法.岩矿石和油气藏所表现出的复电阻率性质很久以来已被人们所认识,研究结果指出可用Cole-Cole模型定量描述岩、矿石激电效应引起的复电阻率频谱特征,奠定了复电阻率法的理论基础(Pelton et al., 1978;张赛珍等,1984;罗延钟和张桂青,1988;Ulrich and Slater, 2004;何展翔和王绪本,2007;Veeken et al., 2009).随着复电阻率法的广泛应用,国内外复电阻率法正演研究得到了较大发展,一些学者开展只考虑激电效应的复电阻率法正演研究(蔡军涛等,2007;Farias et al., 2010;李勇等,2011a;Ye,2014).一些学者开展了同时考虑电磁效应和激电效应的复电阻率法正演研究(张辉等,2006;王大勇等,2010;李勇等,2011b;范翠松,2013).目前电偶源复电阻率法正演多采用结构化网格的有限元法,而基于非结构化网格的自适应有限元法可以更好地模拟复杂地形和构造(Key and Ovall, 2006;Li and Key, 2007;Ren and Tang, 2014;殷长春等,2015).本文实现了2.5维基于电偶源发射的复电阻率法自适应有限元正演计算,采用非结构化网格,可以模拟任意复杂地形和地电结构,更加符合实际应用情况,通过正演模拟分析了地形、激电参数、复杂构造对电偶源频域电磁场的影响规律.
1 复电阻率法2.5维自适应有限元正演在有限元数值模拟中,网格细化程度越高,解的精度越高.但是网格越精细,求解过程越耗时.自适应有限元法,既能自动形成具有最少网格节点数的优化网格又能够提供高精度有限元解的数值方法.
1.1 控制方程在二维地电模型中,设走向方向为x,水平电偶源位于地表,如图 1所示.假定时间因子为e-iωt,在准静态情形下,电场(E)和磁场(H)满足的控制方程为
(1) |
(2) |
这里μ0为真空中的磁导率,ω为角频率,σ为电导率.Js为电源电流分布.
由于场源是三维变量,针对二维电导率分布,利用傅里叶变换获得二维电磁方程.根据傅里叶变换
(3) |
这里kx为走向方向的波数,
(4) |
(5) |
其中
采用加权余量法,微分方程乘以电场的任意变分
(6) |
同理,微分方程乘以磁场的任意变分
(7) |
(8) |
(9) |
将求解区域Ω分解成ne个三角单元,单元编号记为e=1, 2, …, ne.积分方程则分解为各个单元积分之和.在每一个三角单元中,电磁场可表示为由三角形顶点函数值表示的关于y和z的线性函数,对每一个三角单元积分,可以得到6×6的系数矩阵.将求解区域中各三角单元系数矩阵计算出来进行叠加,可以得到总的系数矩阵,从而得到求解电磁场的线性方程组:
(10) |
K为总的系数矩阵,
对于复杂模型的有限元正演要得到高精度的数值解,需要将网格剖分的足够细,剖分的网格数越大越影响正演的计算效率.采用面向目标误差估计法(Key and Ovall, 2011)计算单元误差,从而对有限元精度进行评价,实现自适应网格剖分,保证最小网格剖分下获得高精度的数值解.
定义一个新的算子B,使得式(6)和式(7)满足:
(11) |
其中
(12) |
设定εh≈u-uh,其中u和uh分别是电磁场的精确解和有限元解,则
(13) |
定义函数G,式(11)相应的对偶方程为
(14) |
其中,B*是伴随算子,满足G(v)=B*(w, v)=B(v, w)=F(w),则有
(15) |
设定w-wh≈δh,其中w和wh分别是对偶问题的精确解和有限元解.
则
(16) |
每一个小三角单元τ的局部误差为
(17) |
自适应有限元网格细化的过程是先给定一个初始的粗糙网格,并利用有限元法计算有限元解,然后计算每一个小三角单元局部误差μτ,细化误差较大的单元,得到一个新的网格,然后重新计算小三角单元局部误差μτ,重复上述过程,直到误差满足精度为止.
1.4 Cole-Cole模型的引入为了考虑激电效应,需将地电断面中各种介质的电阻率分别换成用Cole-Cole模型表示的复电阻率,如图 2所示.
根据Pelton等(1978)的结论,对均匀岩、矿石,其复电阻率频谱可用一个Cole-Cole模型表示,即
(18) |
其中,ρ(iω)为复电阻率,ρ0为零频电阻率,m为极化率,τ为时间常数,c为频率相关系数.复电阻率各参数分布范围大致为:m=0~0.98,τ=10-3~5×103 s,c=0.1~0.6,ρ0=10-4~105 Ωm.
剖分的每个三角单元电阻率用Cole-Cole模型表示的复电阻率表示,实现了复电阻率法2.5维自适应有限元正演.
正演中发射源方位角θ为发射源天线相对x轴沿顺时针方向旋转的角度,如图 3所示.方位角为0°时,发射源沿x方向布置,为赤道装置.方位角为90°时,发射源沿y方向布置,为轴向装置.
在电偶源发射情况下,均匀大地表面电场水平分量具有解析解(底青云和王若,2008).本文计算结果与均匀半空间模型的解析解进行比较.水平电偶极子沿y轴放置,采用轴向装置,电偶极距为1 Am,发射频率为1 Hz,地下电阻率为100 Ωm.如图 4所示,正演模拟结果与解析解符合很好,除了场源附近及边界附近,相对误差不超过3%.
设计了水平层状模型,与美国Scripps海洋研究所的计算结果进行对比分析.水平电偶极子沿y轴放置在地表,中心点在0 m处,采用轴向装置,电偶极距为1 Am,发射频率为5 Hz,接收点间隔为100 m,在0~6000 m范围内均匀分布于地表.计算结果如图 5所示,可以看出两者吻合较好,除了场源及边界附近,相对误差不超过4%.
设计了一个与Mitsuhata(2000)相同的山丘起伏地表模型,模型如图 6a所示.水平电偶极子沿y轴放置在地表,中心点在-5500 m处,发射频率为50 Hz,接收点间隔为100 m沿y轴布设.卡尼亚电阻率对比结果如图 6b所示,可以看出两者计算结果一致.
通过与均匀半空间解析解、水平层状模型及起伏山丘模型计算结果比较,说明了本文的正演模拟方法是正确的,具有较高的计算精度.
复电阻率法一般采用有限长度电偶源,只观测沿场源y方向的两个测量电极之间的交流电位差
(19) |
式中,ρs为复视电阻率,K为装置系统,n为观测偶极等均剖分的个数,Ey, i为第i个等均剖分的中心电场值,ΔL为等均剖分步长,为发射电流强度.
复电阻率法观测数据同时受激电和电磁两种效应的影响,常规的数据处理需要分离激电效应和电磁效应,将电磁效应作为噪声去除.本文是基于Maxwell方程的复电阻率法正演,同时考虑了激电效应和电磁效应,对复电阻率法2.5维复杂构造电磁场进行计算,重点分析电场和磁场的变化规律.
2 地形对复电阻率法电磁场的影响特征为了分析不同地形起伏对复电阻率法电磁场的影响特征,设计了一个纯地形起伏的均匀半空间模型,均匀半空间为极化介质,其复电阻率参数为m,τ,c,ρ,模型参数及结构如图 7所示.激发源为水平电偶极子,水平电偶源方位角90°,发射频率为1 Hz,中心点放置在5000 m处,发射电流为1 A;接收点间隔为100 m, 在0~10000 m范围内均匀分布于地表.
分别计算了起伏地形和水平地形时电磁场响应曲线,通过对比分析起伏地形对电磁场的影响规律.图 8为水平地形和起伏地形时电磁场不同分量的响应曲线,从图中可以看出,起伏地形对电磁场有着非常明显的影响,整体上山脊地形处观测振幅相比水平地形时变小,而山谷地形处观测振幅相比水平地形时变大.
为了分析起伏地形对不同场值的影响程度,对振幅进行归一化处理,用起伏地形时电磁场振幅和水平地表时的振幅做比值,如电场归一化振幅可表示为
(20) |
图 9为不同电磁场分量振幅归一化曲线,从图中可以看出,山脊地形对电场影响程度明显比磁场大,山谷地形对磁场垂直分量影响稍大一些.在磁场三分量中,地形对垂直磁场影响程度较大,也可以看出Hx平行于地形构造走向方向,受地形影响小.
图 10为不同发射源方位角时电场水平分量振幅归一化曲线,从图中可以看出,采用不同发射源方位角,地形对电场水平分量都有明显的影响,整体影响程度差不多,但山谷地形对偶极装置影响稍大一点.
为了分析激电效应对复电阻率法电磁场的影响特征,设计了一个断层和背斜复合的地质模型,地层为不极化介质,地层中含有一个极化体(矿藏、油气藏等),其极化参数为m,τ,c,ρ,模型电阻率及结构分布如图 11所示.
改变发射源方位角,分别采用赤道装置(θ=0°)和轴向装置(θ=90°)进行正演计算,分析发射源方位角对激电效应影响特征.激发源为水平电偶极子,发射频率为1 Hz,中心点放置在1000 m处,发射电流为1 A;接收点间隔为100 m, 在3500~11500 m范围内均匀分布于地表.极化体的极化参数为m=0.6,τ=1 s,c=0.3, ρ=500 Ωm.为了便于对比,对振幅进行归一化处理,用极化体极化时电磁场振幅和不极化时的振幅做比值,如电场归一化振幅可表示为
图 12为发射源取不同方位角时电场水平分量振幅的归一化曲线,从图中可以看出,在相同激发电流和收发距时,相同的地下激电模型,采用轴向装置(θ=90°)时激电效应稍强于赤道装置(θ=0°).
为了分析不同激电参数对复电阻率法电磁场影响规律,模型如图 11所示,地层不极化,分别改变极化体的激电参数进行正演计算.激发源为水平电偶极子,水平电偶源方位角90°,发射频率为1 Hz,中心点放置在1000 m处,发射电流为1 A;接收点间隔为100 m, 在4500~10000 m范围内均匀分布于地表,极化体电阻率分为低阻(ρ=30 Ωm)和高阻(ρ=500 Ωm)两种情况.
3.2.1 极化率对复电阻率法电磁场影响规律为了分析极化率对复电阻率法电磁场影响规律,设置了不同的极化率值进行正演计算.极化体的极化参数为τ=1 s,c=0.3,m分别取0,0.6,0.9.
图 13为低阻极化体和高阻极化体取不同极化率时电场水平分量振幅响应曲线,从图中可以看出,极化率对电场振幅水平分量有明显的影响,而且低阻极化体的激电效应明显比高阻极化体强.相比极化体不极化时(m=0),极化率越高,激电效应越强,对电场振幅影响越大.
为了分析时间常数对复电阻率法电磁场影响规律,设置了不同的时间常数值进行正演计算.极化体的极化参数为m=0.6,c=0.3,τ分别取0.1 s,1.0 s,10 s.
图 14为取不同时间常数时电场水平分量振幅响应曲线,从图中可以看出,时间常数对电场振幅水平分量的影响没有极化率明显,整体上时间常数越高,激电效应越强,而且低阻极化体的激电效应明显比高阻极化体强.
为了分析频率相关系数对复电阻率法电磁场影响规律,设置了不同的频率相关系数值进行正演计算.极化体的极化参数为m=0.6,τ=1 s,c分别取0.1,0.3,0.5.
图 15为取不同频率相关系数时电场水平分量振幅响应曲线,从图中可以看出,频率相关系数对电磁场振幅水平分量的影响和时间常数差不多,都不如极化率明显,整体上频率相关系数越高,激电效应越强,而且低阻极化体的激电效应明显比高阻极化体强.
为了分析激电效应对复电阻率法电场和磁场影响的特征,计算了激电效应对电场和磁场影响的归一化振幅值,如图 16所示.极化体的极化参数为m=0.6,τ=1 s,c=0.3,ρ=30 Ωm.从图中可以看出,激电效应对电场的影响比磁场大,在磁场各分量中,激电效应对垂直磁场影响较大.
采用有限元自适应网格,可以模拟任意复杂地下结构,通过正演计算其电磁场分量,分析电磁场对地下构造的识别能力.对于断层和背斜复合的复杂地质模型(图 11),正演得到其电磁场分量振幅曲线,激发源为水平电偶极子,水平电偶源方位角90°,发射频率为1 Hz,中心点放置在1000 m处,发射电流为1 A;接收点间隔为100 m, 在3500~10000 m范围内均匀分布于地表.极化体参数为ρ=30 Ωm,τ=1 s,c=0.3 m分别取0,0.6,0.9.正演结果如图 17所示.
从图中可以看出电场水平振幅对构造变化反应最为敏感,在典型构造上方其振幅出现明显变化.磁场垂直分量对构造变化有一定的反应,但不如电场敏感,而磁场水平分量对构造变化基本上没有什么反应.另外从图中不同极化率振幅曲线的变化特征也可以看出,激电效应对电场的影响比磁场大,在磁场各分量中,激电效应对垂直磁场影响较大.
3.4 收发距对激电效应影响规律激发源为水平电偶极子,水平电偶源方位角90°,发射频率为1 Hz,发射电流为1 A;接收点间隔为100 m,在3500~11200 m范围内均匀分布于地表.发射源中心点位置分别位于1000 m,2000 m和2800 m处,计算了收发距R分别为2500 m,1500 m和700 m时的电场水平分量归一化振幅值,如图 18所示.极化体的极化参数为m=0.6,τ=1 s,c=0.3,ρ=30 Ωm.
从图中可以看到,整体上收发距越大时,激电效应越强,主要原因是收发距的增大,一定程度上探测的深度会变大,这时极化体极化带来的激电效应会逐步突显,并越来越大.
4 结论自适应有限元法能够提供高精度正演结果,由于采用非结构化网格,可以模拟任意地形和复杂地下结构.通过正演模拟,得到以下结论:
(1) 起伏地形对偶极源复电阻率法电磁场有着非常明显的影响,整体上山脊地形处观测振幅相比水平地形时变小,而山谷地形处观测振幅相比水平地形时变大.地形对电场影响程度比磁场大.在磁场三分量中,地形对垂直磁场影响程度较大.
(2) 激电参数对偶极源复电阻率法电磁场振幅有明显的影响.整体上极化率、时间常数和频率相关系数越大,激电效应越强,对电磁场振幅影响越大,其中极化率影响最为明显,而且低阻极化体的激电效应明显比高阻极化体强.
(3) 偶极源采用轴向装置(θ=90°)时极化体产生的激电效应稍强于赤道装置(θ=0°).激电效应对电场的影响比磁场大,在磁场各分量中,激电效应对垂直磁场影响较大.
(4) 电场水平振幅对构造变化反应最为敏感,在典型构造上方其振幅出现明显变化.磁场垂直分量对构造变化有一定的反应,但不如电场敏感.
Cai J T, Ruan B, Zhao G Z. 2007. Two-dimensional modeling of complex resistivity using finite element method. Chinese Journal of Geophysics (in Chinese), 50(6): 1869-1876. |
Di Q Y, Wang R. 2008. Forward and Inversion Method of Controlled Source Audio-Frequency MagnetoTelluric Data and its Application (in Chinese). Beijing: Science Press.
|
Farias V, Maranhao C, Rocha B, et al. 2010. Induced polarization forward modelling using finite element method and the fractal model. Applied Mathematical Modelling, 34(7): 1849-1860. DOI:10.1016/j.apm.2009.10.002 |
Fan C S. 2013. Research on complex resistivity forward and inversion with finite element method and itsapplication (in Chinese). Changchun: Jilin University.
|
He Z X, Wang X B. 2007. Geo-electrical anomalous pattern of reservoir and oil/gas detection by electromagnetic survey. Oil Geophysical Prospecting (in Chinese), 42(1): 102-106. |
Key K, Weiss C. 2006. Adaptive finite-element modeling using unstructured grids:The 2D magnetotelluric example. Geophysics, 71(6): G291-G299. DOI:10.1190/1.2348091 |
Key K, Ovall J. 2011. A parallel goal-oriented adaptive finite element method for 2.5D electromagnetic modelling. Geophysical Journal International, 186(1): 137-154. DOI:10.1111/gji.2011.186.issue-1 |
Li Y, Key K. 2007. 2D marine controlled-source electromagnetic modeling:Part 1-An adaptive finite-element algorithm. Geophysics, 72(2): WA51-WA62. DOI:10.1190/1.2432262 |
Li Y, Lin P R, Li T L. 2011a. Finite element method for solving anomalous complex potential of 2.5-D complex resistivity. Journal of Jilin University(Earth Science Edition) (in Chinese), 41(5): 1596-1604. |
Li Y, Lin P R, Xiao Y, et al. 2011b. Induced Polarization effect on frequency-domain electromagnetic sounding with electric dipole source. Chinese Journal of Geophysics (in Chinese), 54(7): 1935-1944. DOI:10.396/j.issn.0001-5733.2011.07.028 |
Luo Y Z, Zhang G Q. 1988. The Theory of Frequency-Domain Induced Polarization Method (in Chinese). Beijing: Geological Publishing House.
|
Mitsuhata Y. 2000. 2-D electromagnetic modeling by finite-element method with a dipole source and topography. Geophysics, 65(2): 465475. |
Pelton W, Ward S, Hallof P. 1978. Mineral discrimination and removal of inductive coupling with multifrequency IP. Geophysics, 43(3): 588-609. DOI:10.1190/1.1440839 |
Ren Z Y, Tang J T. 2014. A goal-oriented adaptive finite-element approach for multi-electrode resistivity system. Geophys. J. Int., 199(1): 136-145. DOI:10.1093/gji/ggu245 |
Ulrich C, Slater L. 2004. Induced polarization measurements on unsaturated, unconsolidated sands. Geophysics, 69(3): 762-771. DOI:10.1190/1.1759462 |
Veeken P, Legeydo P, Davidenko Y. 2009. Benefits of the induced polarization geoelectric method to hydrocarbon exploration. Geophysics, 74(2). |
Wang D Y, Li T L, Li J P, et al. 2010. Forward simulation of a 3D complex resistivity model. Progress in Geophysics (in Chinese), 25(1): 266-271. DOI:10.3969/j.issn.1004-2903.2010.01.035 |
Ye Y X, Li Y G, Deng J Z, et al. 2014. 5D induced polarization forward modeling using the adaptive finite-element method. Applied Geophysics, 11(4): 500-507. DOI:10.1007/s11770-014-0455-z |
Yin C C, Zhang B, Liu Y H, et al. 2015. 2.5D forward modeling of the time-domain airborne EM system in areas with topographic relief. Chinese Journal of Geophysics (in Chinese), 58(4): 1411-1424. DOI:10.6038/cjg20150427 |
Zhang H, Li T L, Dong R X. 2006. Modeling electromagnetic responses of complex resistivity 3-D body using volume integral equation method. Coal Geology & Exploration (in Chinese), 34(1): 70-73. |
Zhang S Z, Li Y X, Zhang S C, et al. 1984. The low frequency electrical phase spectra of mineralizedrocks(ores) and some factors which influence them in some sulfideore deposits in China. Chinese Journal of Geophysics (in Chinese), 27(2): 176-189. |
蔡军涛, 阮百尧, 赵国泽, 等. 2007. 复电阻率法二维有限元数值模拟. 地球物理学报, 50(6): 1869-1876. |
底青云, 王若. 2008. 可控源音频大地电磁数据正反演及方法应用. 北京: 科学出版社.
|
范翠松. 2013. 基于有限元法的复电阻率正反演研究及应用. 长春: 吉林大学. http://cdmd.cnki.com.cn/Article/CDMD-10183-1013193112.htm
|
何展翔, 王绪本. 2007. 油气藏电性异常模式与电磁法油气检测. 石油地球物理勘探, 42(1): 102-106. |
李勇, 林品荣, 李桐林, 等. 2011a. 基于异常复电位2.5维CR有限元数值模拟. 吉林大学学报(地球科学版), 41(5): 1596-1604. |
李勇, 林品荣, 肖原, 等. 2011b. 电偶源频率电磁测深激发极化效应研究. 地球物理学报, 54(7): 1935-1944. DOI:10.396/j.issn.0001-5733.2011.07.028 |
罗延钟, 张桂青. 1988. 频率域激电法原理. 北京: 地质出版社.
|
王大勇, 李桐林, 李建平, 等. 2010. 三维复电阻率模型电磁场正演模拟研究. 地球物理学进展, 25(1): 266-271. DOI:10.3969/j.issn.1004-2903.2010.01.035 |
殷长春, 张博, 刘云鹤, 等. 2015. 2.5维起伏地表条件下时间域航空电磁正演模拟. 地球物理学报, 58(4): 1411-1424. DOI:10.6038/cjg20150427 |
张辉, 李桐林, 董瑞霞. 2006. 体积分方程法模拟复电阻率三维体电磁响应. 煤田地质与勘探, 34(1): 70-73. |
张赛珍, 李英贤, 张树椿, 等. 1984. 我国几个金属矿区岩(矿)石的低频电相位频率特性及其影响因素. 地球物理学报, 27(2): 176-189. |