基于SPAC法探测松科二井深层地热储水构造

应恒成, 李洪强, 张玉敏, 靳中原, 张家铎, 符 伟, 高 磊,安栋召, 潘宗栋, 汪 伟, 侯贺晟*

1)中国地质科学院, 北京 100037; 2)黑龙江省自然资源调查院, 黑龙江 哈尔滨 150036

中国作为世界上最大的碳排放国家与工业门类最为齐全的国家, 正处在城镇化快速发展阶段,面临着经济转型、环境保护、应对气候变化等多重挑战(张希良等, 2022)。习近平总书记于 2020年 9月22日宣布, 中国将力争于2030年前完成碳达峰,争取在2060年前实现碳中和。而地热能是一种绿色低碳、可循环利用的可再生能源。在清洁能源(太阳能、风能、生物质能、地热能等)大家族中, 地热资源作为可再生能源, 是清洁能源家族中的重要成员,具有丰富的资源潜力, 利用系数高, 生命周期长,二氧化碳排放低, 潜在的低成本优势。相比于水力发电, 地热能发电能经济地建立相对较小的发电机组。并且, 我国拥有丰富的地热资源, 初步估算在3000 m的深度范围内, 地热资源可利用的热能相当于25亿t标准煤。它在未来能源供应与节能减排方面有着巨大潜力(王贵玲等, 2000; 周总瑛等, 2015;田宝卿等, 2020)。随着国家提出碳中和目标, 日后地热能源在我国总能源的占比会不断增加。目前,钻探和地球物理勘探技术是常用的地热勘探方法,其中地球物理勘探的主要作用是圈定含水破碎带和热储水层的区域分布, 是勘查地热空白区的有效手段(李安宁和吕金波, 2001; 郑东明等, 2005)。

在城市地区, 由于环境条件复杂, 常规的地质钻探和传统的物探手段在勘查中受到很大局限和约束, 如钻探工程, 施工成本较大, 并受地下的各类管线与人防工程及现有建筑物等限制, 有些钻孔无法实施, 会出现地质信息盲点; 传统物探手段则会受到噪音、电磁等各种干扰, 导致数据失真, 影响探测准确性。同时, 含水破碎带和热储水层有时埋深较大, 也加大了工作难度, 仅依靠钻探及传统物探手段很难达到勘查目的。因此, 尝试采用新的物探手段与多种传统物探方法联合勘探, 以提高地热勘查探测精度这一观点具有重要意义(付微等,2012)。

微动勘探因为其具有无需人工源、方便快捷、对环境影响小的优势, 非常合适于城市或农村进行地热或浅层速度异常勘探。Xu(2012)通过在微动速度剖面的低速异常体上打地热井的方式, 验证了低速异常体的构造意义, 得到了微动探测得到的低速异常体很大概率为断层/破碎带、含水岩层或石灰岩溶洞这一结论。本文将在这一结论与微动勘探方法的基础上, 结合在黑龙江省安达市松科二井周边的地热资源调查实践, 展示微动探测的实际应用效果。

1 微动原理与探测方法

1.1 微动定义

“微动”是地球表面无时无刻都存在的微弱震动(Okada et al., 2004)。微动是一种在时间、空间域都极不规则的震动现象, 其振动幅值大约为10-4～10-2mm, 加速度只有 1～100 μm/s2, 频率范围主要集中在0.1～50 Hz范围。既可以由诸如气压、风速、海浪、潮汐变化等自然现象产生, 也源自车辆行驶、机器运行等过程, 以及人们的日常生活、生产活动等行为(徐佩芬等, 2013a)。微动是一种由体波(P波和S波)和面波(瑞雷波和勒夫波)组成的复杂震动, 且面波能量占信号总能量的 70%以上(Toksoz and Lacoss, 1968)。

1.2 微动探测方法

要获取瑞雷波中的信息, 首先要从微动信号中提取出面波(瑞雷波)信号, 在提取面波环节, 一般是以微动信号的垂直分量作为资料, 因为勒夫波没有垂直分量, 由观测垂直分量得到的面波仅为瑞雷波, 不会混入勒夫波, 并且对于垂直分量的观测也比较容易做到(王振东, 1990)。

目前, 常用的提取面波的方法有两种, 一是Aki(1957)提出的空间自相关法(Spatial Autocorrelation Method, 简称 SPAC 法); 二是频率-波数法(Frequency-Wavenumber Spectral Method, 简称F-K法)(Ohori et al., 2002)。本文使用的就是空间自相关法。

1.3 基于SPAC法的S波速度结构的估算

SPAC法基于两点基本假设: (1)微动在时空上符合平稳随机过程; (2)微动所包含的各种成分的波中, 基阶面波占优。其基本原理如下:

圆心台阵圆心O(0, 0)和半径为r的圆周上A(r,θ)点的微动信号垂直分量表达式为(Okada, 2006):

该式表明, 波的相速度 c(ω)可以通过第一类零阶贝塞尔函数和空间自相关系数求得。

下一步将圆心和圆周上两点的微动信号通过傅里叶变换转化到频率域的表示, 得到空间自相关系数在频率域下的表示, 将其代入式(6)可得相速度c(ω), 获得相速度频散曲线(徐佩芬等, 2013b)。

获得频散曲线之后, 利用全局寻优算法——分歧型遗传算法(fGA)(Cho et al., 1999)反演地下S波速度结构。

2 应用实例

2.1 测区概况

本次测区位于黑龙江安达市南部松科二井及周边地区(如图1)。松科二井完钻井深7018 m, 是亚洲国家实施的最深的大陆科学钻井, 也是国际大陆科学钻探计划(ICDP)成立 22年来实施的最深钻井(侯贺晟等, 2018)。松科二井位于中国东北的松辽盆地, 松辽盆地是一个整体形状呈现为菱形的中—新生代沉积盆地, 是沉积盆地型地热资源在我国的主要分布区之一(王贵玲等, 2017)。松辽盆地沉积了7000 m厚的中—新生代陆源碎屑岩, 是良好的热储层。具体到地层上, 姚家组、青二段、青三段、泉四段及泉三段都是可能的储层。而青山口组的底部青一段、嫩江组的底部嫩一、二段为盆地分别处于最大沉降期时所形成的半深湖-深湖相沉积, 沉积环境稳定, 泥岩纯且厚度大, 分布面积广, 封盖性能好, 是非常好的区域性隔水层(朱焕来, 2011)。在热源方面, 根据中国大陆地区大地热流数据汇编(第四版)显示, 松辽盆地在统计上属于我国的高热流带。地表热流数据反映了岩石圈浅部热状态和能量均衡信息(罗凡等, 2022)。在具体的热流数据上,中国大陆地区(含渤海海域)平均热流值为(61.5±13.9) mW·m-2, 而松辽盆地的平均热流值高达(70.9±14.4) mW·m-2(姜光政等, 2016)。同时, 松辽盆地中部位于我国地温梯度最高值分布区域(张薇等, 2019), 表明该区域具有很高的地温背景。松科二井的多次井温测井结果显示, 盆地具有较大的地温梯度, 大地热流值达到84 mW·m-2, 2000 m以深地温梯度达 3.7 ℃/100 m, 具有良好的深部热传导和浅部热聚集构造机制, 地热资源赋存条件较好(侯贺晟等, 2018)。

图1 微动测区以及探测点位置图Fig. 1 Survey location and layout of the microtremor station point site

2.2 仪器一致性测试

正式进行数据采集之前, 为保证实验结果的可靠性, 需要保证每台仪器间具有良好的一致性。因此, 本次实验使用的10台SmartSolo IGU-16HR在正式采集数据前, 由深圳面元智能科技有限公司统一进行一致性检测, 使仪器一致性满足微动探测的需求。深圳面元智能科技有限公司提供的仪器的一致性表现图在文章最后的附录中展示。

2.3 微动台阵布设

空间自相关法和频率-波数法对台站布设的要求不同, 本次实验采用空间自相关法, 组成台阵一般需要十台拾震器, 至少为四台。将一台拾震器布置在圆心, 剩下的仪器等间距布置在半径分别为 r,2r, 4r的圆周上(半径r由测深需求和现场环境决定),在本次实验中三个同心圆的半径分别为 200 m、400 m和800 m, 台阵布设如图2。仪器布设完成后对每个观测点进行6小时的连续观测, 共10个观测点, 每个观测点上使用10台节点式地震仪, 最终获得100份三分量地震数据。微动台阵观测对布局的要求很严格, 需要对各台站进行精确定位, 且台阵中的各台站要同时观测, 才能保证数据的准确。

台阵的半径对探测深度起主要决定作用。一般来说, 在使用图 2中的圆形台阵观测时, 勘探深度是圆形台阵半径r的3～5倍。

图2 微动台阵布设以及本次实验采用的布设参数Fig. 2 Layout of the microtremor station points and the designed parameters

2.4 数据处理流程

一致性测试达标后, 就要布置台阵并采集数据。本次实验中, 台阵半径分别设置为200 m、400 m、800 m, 观测系统中共10个探测点, 探测点间距为1000 m。利用RTK(Real-time kinematic)在实地找到预设的点位, 由于本次实验预设的点位有少数在田间小道旁, 甚至有个别接近乡道, 因此实验中在满足接近预设点位的同时, 对部分点位的仪器布设进行点位偏移, 以确保尽量能够避免不必要的干扰,比如道路上的行车, 村民的家畜擅自触碰仪器等情况, 导致数据采集不准确。采集时间设定为 6个小时。使用SmartSolo IGU-16HR进行采集。获得原始数据后的数据处理流程(何正勤等, 2007)如图 3所示。①对原始数据进行预处理, 并分析其频谱, 检查是否存在明显的干扰或数据缺失; ②对预处理后的数据进行数字滤波, 计算其空间自相关系数, 并进行方位平均; ③经过以上处理后可以获得瑞雷波的相速度频散, 手动或自动选取频散谱的峰值以获得频散曲线; ④对频散曲线进行反演得到S波速度结构并绘制成图。

图3 微动台阵探测数据处理流程图Fig. 3 Flow chart of data processing for microtremor array exploration

图4展示了2021年7月10日探测点记录的野外记录。可以看出十台拾震仪接受到的垂直分量信号在连续性上表现较好, 表明仪器在观测期间一直处于正常记录状态。但是原始记录上存在较多毛刺,即地震记录中的突兀的高峰, 这是地震仪器偶尔出现问题导致的, 可以在数据预处理部分改正。此外,从图中还能看到 st.04记录的数据存在一个较小的倾斜趋势, 这也可以在预处理步骤中进行改正。

图4 2021年7月10日探测点原始微动数据Fig. 4 Microtremor recording acquired on July 10, 2021

获得原始数据后, 首先需要对数据进行预处理。利用SAC(Seismic Analysis Code)读入原始数据,成图并检查每天获得的地震波波形, 排除受人为因素而严重干扰的数据或存在仪器没能正常记录的数据。对质量较好的原始数据进行预处理, 进行如下操作: 一是去均值、去线性趋势和波形尖灭, 因为波形数据总会存在一个非零的均值或者存在一个长周期的线性趋势, 这会影响到数据的分析, 必须在数据分析前去除。另一方面, 在对数据进行谱域操作(如FFT、滤波等)时, 若数据的两端不为零, 则会出现谱域假象, 因而实际数据经常需要做尖灭处理,使得数据两端在短时间窗内逐渐变成零值。二是去毛刺, 地震仪器偶尔会出现问题, 导致连续地震数据流中出现尖锋或者数据丢失。这些所谓的毛刺,肉眼很容易识别, 但是在使用程序自动处理数据时却很容易被误认为是地震信号, 因而需要在数据分析之前将毛刺去除。三是对原始数据进行滤波, 本文使用的是10 Hz的低通滤波器。最后要对数据进行分段, 目的是在提取频散曲线时将不同时间段的频散数据进行叠加, 获得更为准确的结果。本次的数据原长6 h, 按45 min为标准分成若干段。

完成数据预处理后, 需要①定义 10个检波器的排列参数, 排列参数按照仪器的实际排列情况给出; ②进行频散曲线的提取, 以2021年7月 10日采集的数据为例, 在获得瑞雷波频散谱后, 手动选取频散谱的峰值, 此时提取的频散点不宜太过密集,因为密集的频散点对最终结果并无太大意义, 且会极大的降低频散曲线的反演效率。下面是经过叠加后得到的瑞利波相速度频散曲线(图5)。③使用分歧型遗传算法(fGA)反演获得观测点地下 S波速度结构, 如图6。图中的绿线为初始模型, 初始模型可以自动创建或者手动给定, 这里选择利用松科二井的地层分层资料给定。④将10个观测点的地下S波速度结构整合, 获得测区内一条横跨 9 km,深度到地下4000 m的地下二维S波速度结构剖面,结果如图7b。

图5 瑞利波频散曲线Fig. 5 Rayleigh wave dispersion curve

图6 横波速度模型Fig. 6 Shear wave velocity mode

3 结果与讨论

本次微动实验, 获得了一条横跨9 km, 深度到地下4000 m的地下二维S波速度结构剖面。根据剖面图可以发现地下S波速度总体呈现西侧速度高,东侧速度低的趋势; 1000 m以上速度基本成层分布,但是在松科二井(图7b中星形所在位置)西侧100 m、东侧2300 m、东侧3700 m处, 深度在1400～2800 m范围内存在三个明显的低速异常体。根据Xu(2012)的研究与实验结果表明, 隐伏断裂破碎带在微动视S波速度剖面上显示出明显的低速异常, 这成为微动剖面解释隐伏地热构造的重要标志, 为地热井的选址提供重要依据。除此之外, 此类速度异常还可能是含水岩层或石灰岩溶洞。根据这一结论, 说明图中显示的低速异常体很可能存在断层或隐伏断裂构造。结合松科二井地层分层资料(侯贺晟等, 2018),如图 7a所示, 在 1675～2 964.92 m的深度范围内,地层以泉头组和登娄库组为主, 泉头组部分(1675～2540 m), 泉三段岩性为砂岩和泥岩的互层;泉四段, 泉二段和泉一段岩性以泥岩为主, 泉一段底部存在粉砂质泥岩层。登娄库组(2540～2 964.92 m)部分, 除了登二段的底部存在部分砂砾岩层外, 其余部分均为不同粒径的砂岩与泥岩的互层。因此从整体上看, S波速度剖面存在低速异常的深度范围内, 岩性整体以砂岩和泥岩互层为主, 而砂岩孔隙度大, 渗透率高, 适合水热的储存, 说明此处具有良好的热储存条件。其次, 结合测区概况部分提及的松辽盆地的地热背景, 相关研究显示松辽盆地热流值远高于全球平均热流值, 且其下地层如姚二、三段、姚一段、青二、三段等地热水自然产能丰富。因此这些位置可能存在丰富的地热资源。

图7 松科二井地下1500～3100 m范围的岩性剖面(a)与微动横波速度剖面(b)Fig. 7Lithologic profile of the SK-2 area between1500and 3100 mbelow the surface (a) andmicrotremor Vs cross section(b)

相较于传统的地震勘探方法, 本次实验使用的微动探测方法,该方法不受城市中强电磁环境、人文工业活动的干扰; 其无需主动源, 因此对城市环境没有损害(田宝卿和丁志峰, 2021)。本次研究区在东北广阔的玉米地中, 对周围环境的影响极小, 即使测线附近均为当地村民种植的玉米, 仪器依旧能按照预期位置布置而不影响居民的日常生活及其农业生产。同时, 仪器的布置以及数据的采集效率高,每天能采集一到两个测点的数据, 并能及时对数据质量做出检查, 若是仪器及人员充足, 每天能采集的数据数量也能进一步增加。

数据处理方面, 目前相关软件已经较为成熟,在保证数据质量的前提下, 可以通过及时的数据处理结果显示, 对相关参数进行改正, 以获取更准确、直观的结果。

4 结论

本文通过微动探测方法, 获得了一条长 9 km,探测深度可达4000 m的二维S波地下速度剖面。本次实验很好地发挥了微动探测效率高, 抗干扰能力强的优点, 野外施工能够高效、且没有对测区内农民的生产生活造成影响的情况下完成。应用目前已经成熟的各种软件, 在保证数据质量的前提下,可以及时的通过显示的数据处理结果, 对相关参数进行改正, 以获取更准确、直观的结果。在获得的S波反演速度剖面上识别出了三个低速异常体, 依据前人进行相关工作得到的结论(Xu et al., 2012),以及结合测区的地热背景(姜光政等, 2016; 王贵玲等, 2017; 张薇等, 2019), 结合松科二井获得的地下分层数据(侯贺晟等, 2018), 判断这三个低速异常体很有可能对应了深部的地热储水构造。

在现今世界, 碳中和已经是世界关注的热点,地球科学可以在其中发挥巨大作用(马冰等, 2021)。微动方法凭借其应用区域广、成本低、探测效率高的优势, 相信能在地热探测方面做出贡献, 为东北经济振兴助力。

致谢: 感谢深圳面元智能科技有限公司提供的仪器与相关技术支持。

Acknowledgements:

This study was supported by China Geological Survey (Nos. DD20190010 and DD20221643-2), and Institute of Natural Resources of Heilongjiang Province (No. HX2020-23).