基于地质统计学反演的透明化矿山岩性建模参数研究及应用

陈柏平崔凡刘波杜云飞王子昌

1.中国矿业大学(北京) 地球科学与测绘工程学院,北京 100083;2.河北工程大学 地球科学与工程学院,河北邯郸 056038

煤炭是我国目前乃至未来很长时期内能源供给的重要组成。 近年来,国家对煤炭安全愈发重视,为保障煤炭开采的安全性、高效性,煤炭开采智能化、无人化等已成为未来发展趋势[1-2]。 袁亮[3-4]于2017年提出煤炭精准开采体系,并指出精准开采和物联网架构将是实现未来无人矿山的重要组成部分。 毛善君等[5]在2018年提出透明化矿山的概念,认为透明化矿山是利用钻探、物探、煤岩层识别等技术手段,基于统一的数据管控平台和生产实时揭露信息,实现矿山井上下工程、设备、开采环境等实时信息的高精度可视化,为矿山的安全生产、培训和智能开采提供高技术手段和统一的可视化管控平台。 因此,建立精准开采和物联网架构体系,以矿井地质、钻井和测井、综合地球物理勘探数据为基础的地质保障技术至关重要,这些数据又是建立在地下空间地质体架构上的[6-7]。 精确的空间地质模型,不仅可以为透明化矿山的建设提供前期基础地质框架设计依据,还能降低后续通过综合物探技术获取的动态、实时监测等数据进行地质模型修正的工作量。

“十三五”期间,随着煤炭工业智能、绿色发展,构建精确的空间地质模型是煤炭透明、智能、精准开采方案中实现可视化、指导智能开采的重要基础[8-10]。 近年来,国内众多学者对透明化矿山建设和煤炭开采做了众多有意义的研究,并多次提到三维地质模型的重要性。 毛善君等[5]提出的透明化矿山构建原则中指出,透明化三维模型构建是关键技术之一。 李梅等[11]提到三维地质建模是矿山三维可视化系统的重要组成部分。 目前,获取井下三维地质体信息的常用方法是基于钻井、测井、三维地震等数据进行三维地质建模。 陈博等[12]基于钻井和测井资料建立煤体构造识别方法,并利用随机反演方法实现煤体结构的三维可视化。 吕堂杰等[13]结合钻井、测井地质分析,并基于层序旋回识别提高了单煤层储层厚度的解释精确度。 基于钻井、测井、地质资料的建模方法具有较高的纵向分辨率,对井点和周围区域的预测准确性较佳,但存在施工人力、物力等成本高的问题。 同时,面对较大区域时,由于井数量有限、地形和构造的影响,井间区域的预测效果已不能满足透明化矿山建设的高分辨率地质模型需求。 为解决高精度地质模型构建问题,近年来众多学者将钻井、测井、地震等数据进行结合,利用三维地震横向连续性强的特点,以钻井、测井信息为约束,形成了基于三维地震反演的储层建模方法。 胡勇等[14]基于地震反演成果进行建模,提高了建模精度;韩冰等[15]利用测井资料建立的高频地质模型和叠后地震阻抗反演建立的低频模型,实现了地层模型的三维可视化,提高了地质模型的空间精度;孙月成等[16]基于地质统计学反演方法综合钻井、测井、地震等多尺度信息建立高精度三维地质模型,进一步提高了气田三维地质模型精度。

近年来,随着储层预测精度要求的提高和三维地震反演技术的发展,常规的地震反演技术的精度已不能满足透明化矿山建设高精度数据的需求[17]。 地质统计学反演方法综合了钻井、测井、地震和地质等多方面资料,横、纵方向上都具有高分辨率的优势,已被广泛应用到油、气田、海洋等领域[18-20]。 但在煤田领域上,目前应用较少,特别是缺乏针对煤系地层环境特点的反演参数和井约束条件的研究。 为此,本文通过建立三维地震正演模型,设置多组地质统计学反演参数和井约束条件进行实验,分析地质统计学反演在煤系地层中的煤层及其厚度反演效果,并探讨反演参数和约束井选取方案,最后将其应用到实际煤田的三维岩性地质模型构建的研究中,以期获得更加精准的效果,更好地服务于透明化矿山基础模型的构建。

1 地质统计学反演和关键反演参数

地质统计学反演是基于贝叶斯框架,通过地震、测井等验证性地质信息拟合得到概率密度函数(Probability Density Function,PDF)来表示储层的岩性或物性特征的概率分布。 反演得到的计算公式是一个多维函数,较为复杂,需要采用马尔科夫链蒙特卡洛算法(Markov Chain Monte Carlo,MCMC)构建马尔科夫链进行随机采样和蒙特卡洛模拟,从而获得需要的高分辨率数据体[21-22]。 贝叶斯判别是一种利用概率统计原理对多种信息进行整合并判别的方法,能够公平地将多种不确定信息源进行统一分析判断,常被用来在已知先验信息概率(E)和假设条件概率(H)的前提下,计算某一事件(X)发生的后验概率[23]。 数学表达式如下:

式中,P(X|H,E)为待求解的后验函数(后验概率分布);P(X|H)为在已知假设条件(H)下事件X的先验函数(先验概率分布);P(E|X)为似然函数;P(H|E)为近似全概率模型。

本文涉及的先验信息概率是通过测井的岩性和纵波阻抗信息拟合得到,主要包括概率密度函数和变程函数两个主要的反演参数。

概率密度函数表示储层参数或属性在空间上的分布规律[24],可表示为

式中,m为岩性;f为纵波阻抗;i为采样点;P(f)为岩性的概率;P(m1|f1)为顶部样本点岩性对应下阻抗的概率的分布;P(mi|m(i-1),f(i-1),fi)为第i处样本的岩性对应下阻抗的概率的分布。

可以看出,用来拟合纵波阻抗的高斯分布函数的均值和方差直接影响到先验函数的准确性[25]。

变程函数表示储层属性在空间位置上的延展性,可以通过式(3)进行变差分析后拟合得到[18]。

式中,σf为样本点对应的纵波阻抗拟合的方差;Δt为采样间隔时间;R为通过变差分析得到的变程。

变程具有一定的主观性,变程越大相关性越强。 对于纵波阻抗,变程一般采用高斯函数进行拟合,变程的大小会影响反演储层参数的横、纵向的延伸程度;对于岩性,横向变程越大,同一岩性的在横向方向的延展程度越长,因此横向变程的取值至关重要[26]。

2 楔状煤层模型正演和反演分析

基于MCMC 算法的地质统计学反演过程中,先验信息的概率密度函数和变程参数以及约束井条件对最终的反演结果至关重要[27]。 目前地质统计学反演方法已被应用到油田、气田领域中,对薄储层的预测效果良好。 但在煤田方面,尤其是透明化矿山建设中煤系地层地质建模的需求,反演参数的设置和约束井的选择,目前尚缺少具有针对性的研究。

本文主要基于GeoSoftware 地学软件平台中的Jsaon 软件进行相关试验和研究。 煤系地层地质统计学反演参数和约束井条件选择方法如下:

(1) 建立楔状煤层纵波阻抗数值模型,通过波阻抗模型正演技术[28]获得地震数据。

(2) 从阻抗模型中提取伪测井数据来模拟实际测井数据,包括纵波阻抗和岩性数据。

(3) 根据伪井的纵波阻抗和岩性分布信息,设置不同的概率密度函数均值和横向变程参数进行地质统计学反演,并分析最佳的反演参数组合。

(4) 在最佳反演参数条件下,进一步分析不同井约束条件下反演的煤层厚度和分布情况。

(5) 获得反演参数和井约束,选取最佳方案。

2.1 楔状煤层纵波阻抗数值模型介绍

楔状煤层纵波阻抗数值模型如图1 所示,模型尺寸设置为101 条线和101 个数据道,长度为1 000 m×1 000 m(线道距为10 m×10 m),总时深设置为1 s,采样时间间隔为1 ms,由中部楔状煤层和上部砂岩、下部泥岩构成,煤层的厚度从1 ～25 m均匀变化。 砂岩、泥岩和煤层的密度、纵波传播速度设置见表1。 图1 中不同的颜色代表对应岩性的纵波阻抗值。 利用纵波阻抗模型正演方法进行正演,正演过程中通过纵波阻抗与子波进行褶积,采用的子波类型为雷克子波,频率为50 Hz,波长为80 ms,得到的地震数据结果如图2 所示。

图2 正演数据剖面Fig.2 Forward data section

表1 模型各岩性的密度、纵波速度参数Table 1 Parameter table of density and p-wave velocity of each rock in the model

图1 楔状煤层纵波阻抗模型示意图Fig.1 Model diagram of p-wave impedance of wedge-shaped coal seam

2.2 地质统计学反演参数和井约束条件

地质统计学反演方法流程如图3 所示,主要包括数据准备、反演参数设置、反演3 个阶段。 数据准备阶段中,需要准备优化好的三维地震、测井、岩性、子波等数据。 反演参数设置是影响最终反演结果的重要环节。 本文重点对岩性和纵波阻抗的概率密度、横向变程以及井约束条件进行了分析。

图3 地质统计学反演流程Fig.3 Geostatistical inversion process

2.2.1 概率密度函数及变程参数设置

基于正演数据的地质统计学反演中,反演参数的大小设定主要针对煤层,其他岩性参数设置为表1 纵波阻抗模型参数的对应值,且反演过程中加入井约束。 由于实际地层的密度和声波的测井数据不是定值,一般可以用高斯分布进行拟合,因此在伪井信息统计的基础上加入高斯白噪声,使得伪井的纵波阻抗值接近实际分布。 伪井的提取遵循均匀分布原则,从模型平面范围上均匀选取9 个点进行提取,如图4 所示。 在实际反演参数设置时,一般通过测井、地震数据统计拟合得到。 在拟合过程中,由于数据自身和人为经验影响,参数的选取会存在一定不确定性[29]。 为了进一步针对煤层探讨反演参数对反演结果的影响,设置了不同概率密度函数的均值ε 和方差σ2以及横向变程H进行实验分析,详细参数设置见表2。 针对概率密度函数,参照表1 中煤的纵波阻抗值设置3种不同的值。 针对横向变程,参照阻抗数值模型尺寸大小设置4 种不同的值,见表2 中3、7、8、9号实验。 其中,第3 组实验参数与实际加噪后的伪井数据一致,理论上第3 组实验参数下效果最好。

表2 地质统计学反演参数Table 2 Geostatistical inversion parameter

2.2.2 井约束条件参数设置

此外,反演过程采用井约束的数量和位置条件对反演结果也有重要影响[30]。 为分析不同数量和位置条件下井约束对地质统计学反演的效果,基于上面第3 组实验的地质统计学反演参数,并加入不同数量和不同位置的井约束条件进行了分析研究。9 组不同井约束条件的实验参数见表3,表中代表约束井位置的编号与图4(a)对应。

表3 不同井约束条件反演参数Table 3 Parameter of different well constraint conditions inversion

图4 伪井提取位置和纵波阻抗测井曲线(加噪)Fig.4 Pseudo-well extraction location and p-wave impedance log (imnoised)

2.3 地质统计学反演结果分析

表2 不同反演参数条件下煤层厚度反演及其相对误差如图5 和图6 所示。 其中实验组1、3、5和实验组2、4、6 的厚度和相对误差结果表明,在其他条件一样时,均值的增大和减小都会造成最终的反演结果与实际误差增大;对比实验组1、2,实验组3、4 和实验组5、6 的厚度和相对误差可以看出,在其他条件一样时,方差的选取不当也会造成最终的反演结果与实际误差增大;对比实验组3、7、8、9的厚度和相对误差可以看出,横向变程为200 m 的误差最大,1 200 m 的误差最小,并随着变程的增大,反演结果的相对误差逐渐降低。

图5 不同反演参数条件下煤层厚度地质统计学反演平面图Fig.5 Geostatistical inversion plane of coal seam thickness under different inversion parameters

图6 不同反演参数条件下煤层厚度地质统计学反演与实际模型厚度的相对误差平面图Fig.6 The relative error of the thickness of coal seam thickness and actual model thickness under different inversion parameters

概率密度函数的均值和方差影响较小,这反映出基于MCMC 法的地质统计学反演结果倾向于地震信息。 在本次实验条件下,地震数据在反演过程中,对概率密度函数均值和方差的取值具有一定的约束能力,反演结果趋近地震信息。

横向变程取值对反演结果影响较大。 这是由于横向变程是通过地震信息变差分析得到的,进一步反映地震数据在地质统计学反演中的重要性。 本次实验条件下,煤层横向上是连续分布的,横向变差取值越接近实际煤层横向延展大小,反演结果越好,但考虑实际煤层受沉积环境、构造等影响造成局部变薄、尖灭现象,实际反演过程中变程取值需要综合考虑。

表3 不同井约束条件下煤层厚度反演及其相对误差如图7 和图8 所示。 对比实验1 ～9 结果可以看出,约束井数量越多,反演的厚度整体误差越小,越接近实际值。 约束井的分布对反演结果存在影响,在约束井附近,反演的厚度相对误差较小,越近越准确;远离约束井的位置,准确性降低。 此外,均匀分散选取约束井,可以降低约束井的数量要求,如实验2 相比较实验1,选取的井减少1/3,但由于选取的井分散在整个区域,实验2 的反演结果基本接近实验1。

图7 不同井约束条件下煤层厚度地质统计学反演平面图Fig.7 Coal seam thickness geostatistical inversion floor plan under different well constraints

图8 不同井约束条件下煤层厚度地质统计学反演与实际模型厚度相对误差平面图Fig.8 The thickness of coal seam thickness geostatistical inversion and actual model thickness relative error plan for different well constraints

综上分析,在针对煤系地层的地质统计学反演中,煤层的概率密度函数均值应接近实际值,方差不易过大,横向变程应考虑实际煤层分布情况进行选取,由于地震波的振幅信息可以一定程度上反映地层的延展情况,可以考虑通过煤层界面的地震数据均方根振幅值的平面分布范围选取合适的横向变程。 此外,约束井应均匀分散在整个研究区。

3 实例应用

以山西省孝义市曙光煤矿井田某工作面区域为例,基于正演实验结果对该研究区的目标煤系地层,进行地质统计学反演并获取三维岩性模型。 研究区地层主要有长城系、寒武系、奥陶系、石炭系、二叠系、三叠系及新生界的第三系、第四系,主要目标煤层为二叠系下统山西组的2 号煤层,煤层的厚度为0 ～3.25 m,平均1.30 m,煤层纵波平均速度为2 450 m/s。 研究区面积约为3.85 km2,叠后地震的采样间隔为1 ms,线道间距为5 m×5 m,频带范围在10 ～120 Hz,主频约50 Hz,区内共有5 个分布较均匀的钻井,其中05 井不参与最终的反演,作为反演结果的盲井验证,如图9 所示。

图9 研究区范围及三维地震数据Fig.9 Research area range and three-dimensional seismic data

为保障反演的准确性,前期对地震数据、测井数据进行了优化和校正处理,以确保先验信息的准确性。 5 个井的测井曲线如图10 所示,煤层厚度分别为1.51 m、2.08 m、2.2 m、0.98 m、1.24 m,煤层的平均纵波阻抗最低,砂岩最高。 由于煤层纵波阻抗相对周围岩性较低,在地震数据中煤层顶底界面表现为强反射轴,通过同相轴追踪获取煤层的地震界面,并将其作为约束条件加入到地质统计学反演中,进一步提高反演的正确性。

图10 测井岩性和纵波阻抗曲线Fig.10 Logging lithology and p-wave impedance curves

为消除子波的旁瓣影响,选取的反演范围为研究区2 号煤层顶底界面上下各扩大20 ms,反演参数取值见表4。 根据2 号煤层地震的均方根振幅信息,煤层的横向变程设置为600 m,概率密度拟合结果如图11 所示。

表4 反演参数取值Table 4 Value of inversion parameters

图11 地质统计学反演概率密度函数拟合结果Fig.11 Fitting results of geostatistical inversion probability density function

2 号煤地质统计学反演的厚度如图12 所示,厚度范围为0 ～2.6 m。 地质统计学反演厚度在测井处的煤层厚度对比如图13 所示,结合表5 的井上厚度误差分析可以看出,地质统计学反演的厚度与测井的实际厚度在井点处吻合度较好,反演的相对误差为1.82% ～12.24% ,05 井(盲井)的误差最大。 此外,煤层厚度越大,反演的误差越小。 基于反演结果对2 号煤及其顶、底板部分岩性进行了三维岩性地质体建模(图14)。 可以看出,地质统计学反演在横向和纵向上都具有很好的分辨能力,能够精确地刻画出各岩性厚度及其空间展布,尤其是针对薄煤层也有比较准确的刻画能力。

图12 2 号煤地质统计学反演的厚度平面图Fig.12 Thickness plane of no.2 coal from geostatistical inversion

图13 地质统计学反演与测井处的煤层厚度对比剖面Fig.13 Coal seam thickness contrast profile at geostatistical inversion and logging

表5 地质统计学反演煤层厚度井点处相对误差分析Table 5 Analysis of relative error at well point in geostatistical inversion of coal seam thickness

图14 基于地质统计学反演岩性结果建立的三维地质岩性模型Fig.14 A 3d geological lithology model based on geostatistical inversion

4 结论与展望

(1) 基于楔状煤层纵波阻抗数值模型正演数据的多组反演实验结果对比分析,概率密度函数、横向变程以及约束井条件对地质统计学反演结果均有影响,其中,概率密度函数影响较小,横向变程和约束井条件影响较大。 井的数量越多,对全局的约束力和范围越大,反演的准确性越高。 此外,约束井位置分布均匀,可以适当降低井数量的要求。

(2) 实际应用表明,地质统计学反演的厚度与测井数据基本一致,相对误差为1.82% ～12.24%。能够较好地刻画出煤系地层的空间展布信息,可以为透明化矿山前期的地质模型构建提供可靠的岩性体数据。