2017, Vol. 32
2. 中交第二公路勘察设计研究院有限公司, 武汉 430056
2. CCCC Second Highway Consultants Co., Ltd, Wuhan 430056, China
可控源音频大地电磁法(CSAMT)于20世纪70年代提出并被应用到地球物理领域, 在矿产普查, 油气勘探, 水文环境等各个方面发挥了巨大的作用.该方法最大的特点是采用人工场源, 大大增加了电磁信号强度, 弥补了天然场源信号微弱、不易观测等缺点.
2.5维数值模拟指的是源为三维、模型为二维的地球物理数值模拟问题(汤井田和何继善, 2005). 2.5维数值模拟问题有着重要的实用价值, 其将三维问题简化为二维问题, 只需要对模型进行剖面离散处理而非全空间, 大大减小了矩阵的尺寸, 提高运算速率.CSAMT的有限元正演模拟.国内外都进行了研究:Unsworth等将一次场与二次场分离, 利用有限元对二维地电模型下水平电偶源的电磁场响应进行了模拟计算(Unsworth等, 1993); Mitsuhata借鉴地震模拟中的伪delta函数模拟电偶源, 只需一次性计算总场, 并利用等参有限元实现了带地形的2.5维频率域电磁正演模拟(Mitsuhata, 2000); 底青云等对2.5维频率域电磁进行了有限元模拟并给出了一组高精度波数, 在层状模型计算方面取得较好的成果(底青云等，2004；底青云和王若, 2008); 薛东川等对2.5维频率域电磁的吸收边界条件进行了讨论, 发现吸收边界条件能够有效压制边界处的反射(薛东川和戴世坤, 2008).本文作者在前人研究基础上, 重点研究基于双二次插值的CSAMT有限元正演模拟.
现阶段煤层气地球物理勘探方面, 主要运用地震勘探技术解释煤层的构造, 运用叠前AVO反演技术、叠前弹性阻抗反演技术和叠后波阻抗反演技术对煤储层含气性进行预测(彭苏萍等, 2005); 采用P波方位AVO技术识别煤层裂缝(彭晓波等, 2005); 用频谱分解技术(常锁亮等, 2009)、数字滤波法等(杨双安等, 2006)预测煤层气的富集区等.
就目前研究进展来看, 地震勘探主要对煤层储层属性进行研究, 且各方法都需要一定的条件, 如采用AVO技术对煤层气储层物性预测需要偏移距均匀分布, 且要有较高信噪比的地震资料, 纵波方位AVO技术仅适用于定向排列的角度较大的裂隙带的检测(常锁亮等, 2009)等, 各项技术都还处于探索试验阶段.就勘探成本来看, 地震勘探特别是3D地震勘探成本极高, 且在地质条件较差的地区, 很难获得高品质的地震资料.电磁勘探方法主要通过探测地下电性结构差异, 从而推断和研究地下地质构造(王绪本等, 2009).在矿产普查, 油气勘探, 水文环境等各个方面发挥了巨大的作用(底青云等, 2006; 朱金华等, 2011), 但在煤层气的勘探方面还少有研究.本文以沁水盆地北部煤层气富集与地下水、地质构造之间的关系为例, 建立相应的地球物理模型进行数值模拟, 探讨了CSAMT对薄层的识别能力, 验证了CSAMT在煤层气勘探具有一定的可行性.
1 有限元方程的推导及求解 1.1 方程的推导假定地下电性结构是二维的, 如图 1所示, 将原点取在地面上, 平行于走向方向为Y轴, 位于地面与Y轴垂直的方向为X轴, 垂直向下的方向为Z轴, 外加源的方向为Y方向, 此时只存在不随走向而变的X方向.
|
图 1 CSAMT 2.5D模拟 Figure 1 Simulation schematic diagram of CSAMT 2.5D |
CSAMT满足麦克斯韦方程.假设时谐因子为eiωt, 则Maxwell方程表示为
|
(1) |
其中:ε是介质的介电常数, σ是介质的电导率, μ是介质的磁导率, ω为圆频率, i为虚数单位, E和H分别为电场和磁场矢量, Js为外加电流源的电流密度.
对上述方程沿构造走向(Y)方向作傅立叶变换(Ky为波数)，公式为
|
(2) |
可得到关于走向电场和磁场的耦合微分方程最后通过傅氏逆变换可求得真实电磁场.这里直接给出由伽辽金有限元方法推导得到的有限元方程，公式为
|
(3) |
其中, 
由于同等计算量下双二次插值一般较双线性插值精度高, 因此, 本文采用双二次插值以提高计算精度; 先对各个单元进行分析、计算, 得到各个单元刚度矩阵; 再对所有单元总体集成得到总体刚度矩阵(刘长生等, 2010; 张继锋等, 2010).图 2所示为坐标变换关系及局部节点编号, 其中图 2a所示为母单元, 图 2b所示为子单元; a和b分别表示子单元的宽度及高度; (x0, z0)为矩形单元中心点的坐标, 采用的插值函数表示为
|
图 2 矩形单元剖分 (a)母单元；(b)子单元. Figure 2 Schematic diagram of rectangle meshing element (a) Source unit; (b)Subelement. |
|
(4) |
代入有限元方程, 则在任一单元内, 可以得到单个单元的电磁场刚度矩阵, 再对单个单元方程进行总体合成, 将方程中的各元素按节点的网格排列顺序放到求解域的整体系数矩阵中, 与对应的节点元素叠加, 最终的矩阵方程包括所有单元内所有节点的场方程.在此不再详述, 其具体过程见参考文献(徐世浙, 1994).
1.3 波数选取2.5维数值模拟是先在波数域里计算电磁场, 再通过傅里叶逆变换求得空间域的电磁场分布.由于在数值模拟时不可能取无穷多个波数来还原真实的电磁场, 因此, 需要选取若干个合适的波数以尽可能、快速地求得真实的电磁场.本文计算中, 波数范围取为10-6~1/Δ, 波数取其中的21个对数等间隔数, Δ为网格最小间距.
2 数值算例 2.1 算法验证为了验证本文算法的正确性, 首先利用本文算法计算电偶极子源A型断面响应, 并与解析解进行对比, 如图 3所示.水平电偶极子沿y向放置在坐标原点, 电极距AB=1 km, 供电电流1 A, A型断面第一层电阻率100 Ω·m, 层厚200 m, 第二层电阻率500 Ω·m, 层厚500 m, 第三层电阻率1000 Ω·m.经计算, 最大相对误差为1.326%, 平均相对误差为0.514%, 完全满足正演计算精度要求.
|
图 3 A型模型解析解与数值解对比结果 Figure 3 Comparison between analytic solution and numerical solution of A model |
图 4所示为高、阻体模型的示意图, 图 6所示为正演模拟结果.
|
图 4 高、低阻体模型 Figure 4 Schematic diagram of high/low resistivity body |
|
图 6 TE/TM模式下模型模拟结果 (a)TE模式视电阻率等值线；(b)TM模式视电阻率等值线；(c)TE模式相位等值线；(d)TM模式相位等值线. Figure 6 chematic simulation result of model of TE/TM mode (a) Apparent resistivity contour map of TE mode; (b) Apparent resistivity contour map of TM mode; (c) Phase contour map of TE mode; (d) Phase contour map of TM mode. |
图 6a和图 6c所示分别为TE模式下的视电阻率及相位等值线图, 图中都反映出了异常, 且异常的中心位置与理论模型一致, 异常都有往右下方拉长的趋势, 这是场源阴影效应的影响; 图 6b和图 6d所示分别为TM模式下的视电阻率及相位等值线图, 两个图中也都反映出了异常, 且相比TE模式异常垂向拉长明显.这是由于TM模式下测定的是电场的法向分量, 界面电荷累积造成静态效应引起的, 相位受静态效应影响较弱, 这提高了方法的灵敏度, 有利于发现小规模异常阻体, 故TM模式发现异常的能力比TE模式的强, 但是对于异常体真实形态反应, TE模式要优于TM模式.视电阻率等值线图异常最低为11.5 Ω·m, 真实异常为10 Ω·m, 接近异常体真实值; 视电阻率等值线图中异常最高为183 Ω·m, 而真实异常为1000 Ω·m, 说明CSAMT对高阻体反映不灵敏.
2.2.2 倾斜板状体模型图 5所示为倾斜板状体模型的示意图, 图 7、图 8所示为正演模拟结果.
|
图 5 倾斜板状体模型 Figure 5 Schematic model of tilt plate-like body |
|
图 7 TE/TM模式下低阻倾斜板状体体模拟结果 (a)TE模式视电阻率等值线；(b)TM模式视电阻率等值线；(c)TE模式相位等值线；(d)TM模式相位等值线. Figure 7 Schematic simulation result of high inclined resistivity body of TE/TM mode (a) Apparent resistivity contour map of TE mode; (b) Apparent resistivity contour map of TM mode; (c) Phase contour map of TE mode; (d) Phase contour map of TM mode. |
|
图 8 TE/TM模式下高阻倾斜板状体体模拟结果 (a)TE模式视电阻率等值线；(b)TM模式视电阻率等值线；(c)TE模式相位等值线；(d)TM模式相位等值线. Figure 8 Schematic simulation result of high Inclined resistivity body of TE/TM mode (a) Apparent resistivity contour map of TE mode; (b) Apparent resistivity contour map of TM mode; (c) Phase contour map of TE mode; (d) Phase contour map of TM mode. |
低阻倾斜板状体电阻率10 Ω·m, 背景电阻率100 Ω·m, 图 7为倾斜低阻体视电阻率及相位等值线图.可以看出低阻倾斜异常反映明显, TE模式和TM模式视电阻率与相位等值断面图对异常均有良好的显示, 异常形态向异常体倾向方向延伸且异常范围较大.TM模式下视电阻率等值线图异常最低为10.5 Ω·m, 接近异常体真实值.
高阻倾斜板状体电阻率1000 Ω·m, 背景电阻率100 Ω·m, 图 8为倾斜高阻体视电阻率及相位等值线图.可以看出TE模式较TM模式视电阻率与视相位等值断面图对异常形态显示较好, 异常形态向异常体倾向方向延伸, TM模式对异常体形态反应不明显.视电阻率等值线图中异常最高为380 Ω·m, 而真实异常为1000 Ω·m, 说明CSAMT对高阻体反映不灵敏.
2.3 模拟实验分析通过数值模拟, 发现CSAMT低阻体反映灵敏, 对高阻体反映不灵敏, 这是电磁法勘探方法普遍存在的一个问题.TM模式对低阻及高阻异常体反应较TE模式明显, 且横向分辨能力更强, 但是对于异常体真实形态反应, TE模式要优于TM模式.相位受静态效应影响较小.
3 CSAMT在煤气富集区勘探可行性研究 3.1 区域地质、水文地质情况沁水盆地是石炭-叠纪的华北克拉通在奥陶系地层的风化剥蚀面上接受了广泛的含煤沉积后, 经历了印支运动, 燕山运动, 使地层抬升, 并遭受大面积剥蚀, 形成的残留盆地.盆地总体为一走向北北东的宽缓复式向斜构造.整体属于阳泉—襄垣富气带, 富气带基本沿次级向斜轴部呈北北东向展布.
根据区域水文地质资料, 整个煤层气探区处于一向斜构造带上, 向斜的两翼均为地下水封堵.沿两侧煤系地层露头至深部, 地下水分区依次补给区、弱径流-缓流区、滞流区.两侧浅部的地下水补给区煤层气含气量低; 深部滞流区, 矿化度较高, 是煤层气的有利富集区.在煤层热演化过程中, 生成煤层气的同时也有大量的水生成, 这些水和地表水的补给共同形成了煤层水, 现今煤储层的压力主要表现为地下水的静压力, 合适的构造形态可以增加煤储层的静水压力.当储层压力大于煤层气的临界解吸压力时, 煤层气以吸附态保存在煤岩中, 而形成富集(王绪本等, 2013).
3.2 煤层气的富集模型本区煤层气的富集主要受构造和地下水控制作用明显.结合研究区地下水, 地质构造, 可将研究区煤层气富集模式分为以下三类, 即向斜中心部位地下水滞流区水力封闭富集模式、斜坡带弱径流水力封堵模式、陷落柱水力运移逸散模式
沁水盆地北部煤层气富集存在水力封闭模式, 而在水力封闭区一般处于宽缓向斜的中心, 其中心部位地层相对平缓, 可认为近似层状, 含水富水地层常表现为低阻特征, 上覆地层为中高阻, 其下伏地层为高阻, 主要为奥陶系灰岩, 埋深在500 m左右.为检验CSAMT是否能有效分辨低阻薄层, 并根据研究区煤层气三种富集模式, 分别设计如下物理模型进行分析:
(1) 水平层状模型, 图 9a为水平层状介质模型, 设计水平层状五层介质, 参数如图所示.图 10为层状模型在两种不同模式下CSAMT视电阻率和视相位拟断面图.
|
图 9 煤层气富集模型 Figure 9 Schematic diagram of coalbed methane enrichment model |
|
图 10 TE/TM模式下水平层状模型模拟结果 (a)TE模式视电阻率等值线；(b)TM模式视电阻率等值线；(c)TE模式相位等值线；(d)TM模式相位等值线. Figure 10 Schematic simulation result of horizontally layered model of TE/TM mode (a) Apparent resistivity contour map of TE mode; (b) Apparent resistivity contour map of TM mode; (c) Phase contour map of TE mode; (d) Phase contour map of TM mode. |
(2) 单斜地层模型, 图 9b为倾斜层状介质模型, 设计倾斜层状五层介质, 参数如图所示.图 11为单斜地层模型在两种不同模式下视电阻率和视相位拟断面图.
|
图 11 TE/TM模式下单斜地层模型模拟结果 (a)TE模式视电阻率等值线；(b)TM模式视电阻率等值线；(c)TE模式相位等值线；(d)TM模式相位等值线. Figure 11 Schematic simulation result of tilt layered model of TE/TM mode (a) Apparent resistivity contour map of TE mode; (b) Apparent resistivity contour map of TM mode; (c) Phase contour map of TE mode; (d) Phase contour map of TM mode. |
(3) 倾斜含断层模型, 图 9c为倾斜含断层层状介质模型, 设计五层介质, 参数如图所示.图 12为倾斜含断层模型在两种不同模式下视电阻率和视相位拟断面图.
|
图 12 TE/TM模式下倾斜加断层模型模拟结果 (a)TE模式视电阻率等值线；(b)TM模式视电阻率等值线；(c)TE模式相位等值线；(d)TM模式相位等值线. Figure 12 Schematic simulation result of tilt fault model of TE/TM mode (a) Apparent resistivity contour map of TE mode; (b) Apparent resistivity contour map of TM mode; (c) Phase contour map of TE mode; (d) Phase contour map of TM mode. |
通过图 10、11可以看出各层在图中均能表现在视电阻率等值线图和相位等值线图上.且对于薄层反应, 视相位等值线图中的表现比视电阻率拟等值线图表现更明显.
通过图 12可以看出, TM较TE模式在视电阻率和视相位等值线图更易识别断层, 且断层在视相位等值线图中的表现比视电阻率等值线图表现更明显.在TE模式视电阻率等值线图中, 断层仅表现为起伏地层状异常, 较难判断是否存在断层.TM模式视电阻率等值线图中, 断层下方变现为起伏地层状假异常, TM模式视相位等值线图中, 可以清楚的判断出断层的存在.
上述模型模拟试验结果表明, 利用CSAMT对地下一定埋深的薄层进行探测时, 在一定条件下可以有效识别薄层的存在.模拟结果也为CSAMT在煤层气勘探的可行性提供了理论依据.
4 结论 4.1通过A型空间模型下的有限元数值解与解析解对比, 结果显示有限元数值解与解析解基本吻合, 从而从一定程度上验证了基于双二次插值的CSAMT 2.5维正演模拟算法的正确性.后续算法能正确反映二维模型的异常, 进一步验证了该算法的可靠性.
4.2TM模式下低阻及高阻异常体视电阻率响应均有明显的静态效应, 视相位受静态效应影响较弱.TM模式对低阻及高阻异常体反应较TE模式明显, 且横向分辨能力更强, 但是对于异常体真实形态反应, TE模式要优于TM模式.
4.3可控源电磁法对低阻体反映灵敏, 对高阻体反映不灵敏.
4.4根据沁水盆地北部煤层气富集区地质构造, 建立相应地球物理模型, 通过数值模拟, 验证了CSAMT能有效识别具有一定埋深和物性条件的薄层, CSAMT方法在煤层气富集区的勘探方面具有一定的可行性, 可以作为煤层气非地震勘探的技术手段.
致谢 感谢审稿专家提出的修改意见和编辑部的大力支持!| [] | Chang S L, Liu D M, Lin Y C, et al. 2009. Application of spectral decomposition for fine seismic structural interpretation in coalfield and gas-bearing property predication of coal seam[J]. Journal of China Coal Society (in Chinese), 34(8): 1015–1021. |
| [] | Di Q Y, Unsworth M, Wang M Y. 2004. 2.5-D CSAMT modeling with the finite element method over 2-D complex earth media[J]. Chinese Journal of Geophysics (in Chinese), 47(4): 723–730. DOI:10.3321/j.issn:0001-5733.2004.04.026 |
| [] | Di Q Y, Wang G J, An Z G, et al. 2006. Geophysical exploration of a long deep tunnel on the west route of South-to-North Water Diversion Project[J]. Chinese Journal of Geophysics (in Chinese), 49(6): 1836–1842. DOI:10.3321/j.issn:0001-5733.2006.06.032 |
| [] | Liu C S, Tang J T, Ren Z Y, et al. 2010. Three-dimension magnetotellurics modeling by adaptive edge finite-element using unstructured meshes[J]. Journal of Central South University (Science and Technology) (in Chinese), 41(5): 1855–1859. |
| [] | Mitsuhata Y. 2000. 2-D electromagnetic modeling by finite-element method with a dipole source and topography[J]. Geophysics, 65(2): 465–475. DOI:10.1190/1.1444740 |
| [] | Peng S P, Gao Y F, Yang R Z, et al. 2005. Theory and application of AVO for detection of coalbed methane-A case from the Huainan coalfield[J]. Chinese Journal of Geophysics (in Chinese), 48(6): 1475–1486. DOI:10.3321/j.issn:0001-5733.2005.06.034 |
| [] | Peng X B, Peng S P, Yan G, et al. 2005. P-wave azimuthal AVO analysis of fracture detection in coalbed and its application to engineering[J]. Chinese Journal of Rock Mechanics and Engineering (in Chinese), 24(16): 2960–2965. |
| [] | Unsworth M J, Travis B J, Chave A D. 1993. Electromagnetic induction by a finite electric dipole source over a 2-D earth[J]. Geophysics (in Chinese), 58(2): 198–214. |
| [] | Wang X B, Chen J C, Guo Q S, et al. 2013. Research of the CSAMT exploration mode and experiment for the coalbed methane enrichment region in the north Qinshui basin[J]. Chinese Journal of Geophysics (in Chinese), 56(12): 4310–4323. DOI:10.6038/cjg20131233 |
| [] | Wang X B, Zhu Y T, Zhao X K, et al. 2009. Deep conductivity characteristics of the Longmen Shan, Eastern Qinghai-Tibet Plateau[J]. Chinese Journal of Geophysics (in Chinese), 52(2): 564–571. |
| [] | Xue D C, Dai S K. 2008. Absorbing boundary condition for simulating 2.5-D electromagnetic sounding in frequency domain by finite element method[J]. Journal of China University of Petroleum (in Chinese), 32(6): 57–61. |
| [] | Yang S A, Ning S N, Zhang H X, et al. 2006. Research achievements of forecasting gas using three-dimensional seismic exploration[J]. Journal of China Coal Society (in Chinese), 31(3): 334–336. |
| [] | Zhang J F, Tang J T, Wang Y, et al. 2010. Magnetotellurics fast forward based on preconditioning conjugate gradient[J]. Journal of Central South University (Science and Technology) (in Chinese), 41(5): 1877–1882. |
| [] | Zhu J H, Mao W D, Bai J L, et al. 2011. The application of CSAMT method to the evaluation of water-bearing properties of faults[J]. Geophysical and Geochemical Exploration (in Chinese), 35(4): 569–572. |
| [] | 常锁亮, 刘大锰, 林玉成, 等. 2009. 频谱分解技术在煤田精细构造解释及煤含气性预测中的应用[J]. 煤炭学报, 34(8): 1015–1021. |
| [] | 底青云, UnsworthM, 王妙月. 2004. 复杂介质有限元法2.5维可控源音频大地电磁法数值模拟[J]. 地球物理学报, 47(4): 723–730. DOI:10.3321/j.issn:0001-5733.2004.04.026 |
| [] | 底青云, 王光杰, 安志国, 等. 2006. 南水北调西线千米深长隧洞围岩构造地球物理勘探[J]. 地球物理学报, 49(6): 1836–1842. DOI:10.3321/j.issn:0001-5733.2006.06.032 |
| [] | 底青云, 王若. 2008. 可控源音频大地电磁数据正反演及方法应用[M]. 北京: 科学出版社. |
| [] | 刘长生, 汤井田, 任政勇, 等. 2010. 基于非结构化网格的三维大地电磁自适应矢量有限元模拟[J]. 中南大学学报(自然科学版), 41(5): 1855–1859. |
| [] | 彭苏萍, 高云峰, 杨瑞召, 等. 2005. AVO探测煤层瓦斯富集的理论探讨和初步实践:以淮南煤田为例[J]. 地球物理学报, 48(6): 1475–1486. DOI:10.3321/j.issn:0001-5733.2005.06.034 |
| [] | 彭晓波, 彭苏萍, 詹阁, 等. 2005. P波方位AVO在煤层裂缝探测中的应用[J]. 岩石力学与工程学报, 24(16): 2960–2965. DOI:10.3321/j.issn:1000-6915.2005.16.025 |
| [] | 汤井田, 何继善. 2005. 可控源音频大地电磁法及其应用[M]. 长沙: 中南大学出版社. |
| [] | UnsworthM J, TravisB J, ChaveA D. 1993. 2.5维偶极源CSAMT有限元正演[J]. 地球物理学, 58(2): 198–214. |
| [] | 王绪本, 陈进超, 郭全仕, 等. 2013. 沁水盆地北部煤层气富集区CSAMT勘探试验研究[J]. 地球物理学报, 56(12): 4310–4323. DOI:10.6038/cjg20131233 |
| [] | 王绪本, 朱迎堂, 赵锡奎, 等. 2009. 青藏高原东缘龙门山逆冲构造深部电性结构特征[J]. 地球物理学报, 52(2): 564–571. |
| [] | 徐世浙. 1994. 地球物理中的有限单元法[M]. 北京: 科学出版社: 98-100. |
| [] | 薛东川, 戴世坤. 2008. 频率域2.5维电磁测深有限元模拟中的吸收边界条件[J]. 中国石油大学学报(自然科学版), 32(6): 57–61. |
| [] | 杨双安, 宁书年, 张会星, 等. 2006. 三维地震勘探技术预测瓦斯的研究成果[J]. 煤炭学报, 31(3): 334–336. |
| [] | 张继锋, 汤井田, 王烨, 等. 2010. 基于预处理共轭梯度的大地电磁快速正演[J]. 中南大学学报(自然科学版), 41(5): 1877–1882. |
| [] | 朱金华, 冒我冬, 白锦琳, 等. 2011. CSAMT法在断层含水性评价中的应用[J]. 物探与化探, 35(4): 569–572. |