页岩储层游离气含气量地震预测方法及系统
技术领域
[0001] 本
发明涉及油气地球物理勘探领域，具体地，涉及一种页岩储层游离气含气量地震预测方法及系统。
背景技术
[0002] 页岩气是蕴藏于页岩层可供开采的
天然气能源，具有开采寿命长、生产周期长等优点，是一种清洁、高效的能源资源。游离状态的页岩气存在于页岩的孔隙或裂隙中，其数量的多少决定于页岩内的孔隙裂隙空间，而页岩含气量对页岩储层评价、有利区优选具有重要的意义。现有的预测含气量的方法有现场
解吸法、
测井解释法、TOC拟合法等。
[0003]
发明人发现，现有方法由于忽略了其它影响因素，具有估算量偏小、，误差较大等问题，导致获得的含气量值差别较大，给利用含气量进行资源储量的预测和有利区的优选带来困难。因此，有必要开发一种能够高
精度预测页岩储层游离气含气量的方法及系统。
[0004] 公开于本发明背景技术部分的信息仅仅旨在加深对本发明的一般背景技术的理解，而不应当被视为承认或以任何形式暗示该信息构成已为本领域技术人员所公知的
现有技术。
发明内容
[0005] 本发明提出了一种页岩储层游离气含气量地震预测方法及系统，其能够通过
地层温度、地
层压力、孔隙度及含气
饱和度等因素约束，实现页岩储层游离气含气量的高精度预测。
[0006] 根据本发明的一方面，提出了一种页岩储层游离气含气量地震预测方法，所述方法可以包括：基于工区叠前地震数据体与深度域的测井数据，分别获取工区弹性参数与时间域的测井数据；基于所述工区弹性参数与所述时间域的测井数据，建立
岩石物理模型，结合贝叶斯反演
算法，获得孔隙度与含气饱和度；基于地震资料与
地温梯度，分别获得地层压力与地层温度，进而获得体积系数；以及基于所述孔隙度、所述含气饱和度、所述体积系数与所述地层压力，建立游离气含气量的计算模型，获得所述游离气含气量，其中，所述工区弹性参数包括纵波阻抗、横波阻抗及
密度。
[0007] 根据本发明的另一方面，提出了一种页岩储层游离气含气量地震预测系统，所述系统可以包括：用来基于工区叠前地震数据体与深度域的测井数据，分别获取工区弹性参数与时间域的测井数据的单元；用来基于所述工区弹性参数与所述时间域的测井数据，建立岩石物理模型，结合贝叶斯反演算法，获得孔隙度与含气饱和度的单元；用来基于地震资料与地温梯度，分别获得地层压力与地层温度，进而获得体积系数的单元；以及用来基于所述孔隙度、所述含气饱和度、所述体积系数与所述地层压力，建立游离气含气量的计算模型，获得所述游离气含气量的单元，其中，所述工区弹性参数包括纵波阻抗、横波阻抗及密度。
[0008] 本发明的方法和装置具有其它的特性和优点，这些特性和优点从并入本文中的
附图和随后的具体实施方式中将是显而易见的，或者将在并入本文中的附图和随后的具体实施方式中进行详细陈述，这些附图和具体实施方式共同用于解释本发明的特定原理。
附图说明
[0009] 通过结合附图对本发明示例性实施方式进行更详细的描述，本发明的上述以及其它目的、特征和优势将变得更加明显，其中，在本发明示例性实施方式中，相同的参考标号通常代表相同部件。
[0010] 图1示出了根据本发明的页岩储层游离气含气量地震预测方法的步骤的
流程图。
[0011] 图2示出了根据本发明的一个实施方式的游离气含气量的示意图。
具体实施方式
[0012] 下面将参照附图更详细地描述本发明。虽然附图中显示了本发明的优选实施方式，然而应该理解，可以以各种形式实现本发明而不应被这里阐述的实施方式所限制。相反，提供这些实施方式是为了使本发明更加透彻和完整，并且能够将本发明的范围完整地传达给本领域的技术人员。
[0013] 实施方式1
[0014] 图1示出了根据本发明的页岩储层游离气含气量地震预测方法的步骤的流程图。
[0015] 在该实施方式中，根据本发明的页岩储层游离气含气量地震预测方法可以包括：步骤101，基于工区叠前地震数据体与深度域的测井数据，分别获取工区弹性参数与时间域的测井数据；步骤102，基于工区弹性参数与时间域的测井数据，建立岩石物理模型，结合贝叶斯反演算法，获得孔隙度与含气饱和度；步骤103，基于地震资料与地温梯度，分别获得地层压力与地层温度，进而获得体积系数；以及步骤104，基于孔隙度、含气饱和度、体积系数与地层压力，建立游离气含气量的计算模型，获得游离气含气量，其中，工区弹性参数包括纵波阻抗、横波阻抗及密度。
[0016] 该实施方式通过地层温度、地层压力、孔隙度及含气饱和度等因素约束，实现页岩储层游离气含气量的高精度预测。
[0017] 下面详细说明根据本发明的页岩储层游离气含气量地震预测方法的具体步骤。
[0018] 在一个示例中，基于工区叠前地震数据体与深度域的测井数据，分别获取工区弹性参数与时间域的测井数据，其中，工区弹性参数包括纵波阻抗、横波阻抗及密度。
[0019] 具体地，基于工区叠前地震数据体，通过地震反演的方法，可以利用商业
软件，比如Jason、Strata等，选择合适的子波、低频模型以及合适的反演参数，可以获得纵向
分辨率较高的纵波阻抗Ip(x,y,t)、横波阻抗Is(x,y,t)与密度ρ(x,y,t)。
[0020] 具体地，深度域的测井数据Well(z)可以包括：纵波阻抗Ip(z)、横波阻抗Is(z)、密度ρ(z)、孔隙度φ(z)、含气饱和度Sg(z)、纵
波速度vp(z)。
[0021] 具体地，对深度域的测井数据进行时深转换，获得时间域的测井数据，时深转换的关系式可以为：
[0022]
[0023] 其中，t为时间域，z为深度域，vp(z)为纵波速度。
[0024] 具体地，时间域的测井数据可以包括：时间域的纵波阻抗曲线Ip(t)、时间域的横波阻抗Is(t)、时间域的密度曲线ρ(t)、时间域的孔隙度φ(t)、时间域的含气饱和度Sg(t)。
[0025] 在一个示例中，基于工区弹性参数与时间域的测井数据，建立岩石物理模型，结合贝叶斯反演算法，获得孔隙度与含气饱和度。
[0026] 具体地，基于工区弹性参数与时间域的测井数据，建立页岩储层
各向异性岩石物理模型，进而构建孔隙度、含气饱和度与工区弹性参数之间的函数关系可以为：
[0027] {Ip(t),Is(t),ρ(t)}＝f(φ(t),Sg(t)…) (2)
[0028] 在一个示例中，基于贝叶斯反演算法，构建目标函数，并确定孔隙度与含气饱和度的先验分布。
[0029] 具体地，本发明可以采用贝叶斯分类算法，目标函数表示为在已知工区弹性参数条件下，孔隙度φ和含气饱和度Sg的最大后验概率分布：
[0030] {φ,Sg}＝arg Maxp(φj,Sgj…|Ip,Is,ρ)j＝1,2…,Nc (3)
[0031] 其中，{φ,Sg}表示孔隙度和含气饱和度的解，p(φj,Sgj…|Ip,Is,ρ)表示孔隙度与含气饱和度的后验概率，Nc表示分类的类别数，j表示类别编号，本领域技术人员可以通过具体情况设定Nc用来控制预测精度。
[0032] 即，如果：
[0033] p(φj,Sgj…|Ip,Is,ρ)>p(φm,Sgm…|Ip,Is,ρ)1≤j,m≤Nc,j≠m (4)[0034] 则：
[0035] φ＝φj，Sg＝Sgj (5)
[0036] 其中，(4)与(5)的意义为，在已知纵波阻抗Ip，横波阻抗Is，密度ρ的情况下，求取孔隙度φ和含气饱和度Sg为不同值时的后验概率，其中，最大后验概率的值φ＝φj，Sg＝Sgj为孔隙度和含气饱和度的解。
[0037] 在一个示例中，目标函数可以为：
[0038] {φ,Sg}＝arg Max{p(Ip,Is,ρ|φ,Sg…)·p(φ,Sg…)} (6)
[0039] 其中，{φ,Sg}表示孔隙度和含气饱和度的解，Ip表示纵波阻抗，Is表示横波阻抗，ρ表示密度，p(φ,Sg…)表示孔隙度与含气饱和度的先验分布，p(Ip,Is,ρ|φ,Sg…)表示似然函数。
[0040] 具体地，基于贝叶斯公式，公式(3)可以改写为公式(6)。假设孔隙度与含气饱和度服从多维高斯分布，即孔隙度与含气饱和度的先验分布可以为：
[0041]
[0042] 其中，F代表多维高斯分布， 为多维高斯分布的均值， 为多维高斯分布的方差，Nr为岩石物理模型中，孔隙度与含气饱和度的变量的个数，权重系数αk满足[0043] 在一个示例中，基于孔隙度与含气饱和度的先验分布，结合抽样算法，获得孔隙度、含气饱和度与工区弹性参数的联合随机模拟结果。
[0044] 具体地，基于孔隙度与含气饱和度的先验分布，利用抽样算法对孔隙度与含气饱和度进行随机模拟，可以得到孔隙度与含气饱和度的随机模拟结果{φi,Sgi…}i＝1…N。其中，可以采用MCMC抽样Metropolis Hastings抽样算法等现有抽样算法进行随机模拟。
[0045] 可以将{φi,Sgi…}i＝1…N带入(2)，得到工区弹性参数、孔隙度与含气饱和度的联合随机模拟结果为：
[0046] {Ipi,Isi,ρi,φi,Sgi…}i＝1…N (8)
[0047] 其中，N为随机抽样次数，i表示第i次随机抽样，i＝1…N，Ipi为第i次随机抽样的纵波阻抗，Isi为第i次随机抽样的横波阻抗，ρi为第i次随机抽样的密度，φi为第i次随机抽样的孔隙度，Sgi为第i次随机抽样的含气饱和度。
[0048] 在一个示例中，基于贝叶斯分类算法与孔隙度、含气饱和度与工区弹性参数的联合随机模拟结果，求解目标函数。
[0049] 在一个示例中，求解目标函数可以为：
[0050]
[0051] 其中，N表示随机抽样次数，n表示计数统计，n(Ip∩φi,Sgi,…)表示统计{Ipi,Isi,ρi,φi,Sgi,…}i＝1…N中φ的值等于φi，Sg的值为Sgi，且纵波阻抗为Ip的随机
采样样点的个数，n(Is∩φi,Sgi,…)表示统计{Ipi,Isi,ρi,φi,Sgi,…}i＝1…N中φ的值等于φi，Sg的值为Sgi，且横波阻抗为Is的随机采样样点的个数，n(ρ∩φi,Sgi,…)表示统计{Ipi,Isi,ρi,φi,Sgi,…}i＝1…N中φ的值等于φi，Sg的值为Sgi，且密度为ρ的随机采样样点的个数。
[0052] 具体地，基于贝叶斯分类算法的条件独立性假设，目标函数的计算形式可以改写为：
[0053] p(φj,Sgj…|Ip,Is,ρ)＝p(Ip|φj,Sgj,…)·p(Is|φj,Sgj,…)·p(ρ|φj,Sgj,…)·p(φj,Sgj,…) (10)
[0054] 将(8)带入(10)可以得到：
[0055]
[0056] 在一个示例中，将工区弹性参数带入目标函数，获得孔隙度与含气饱和度。
[0057] 具体地，可以将纵波阻抗Ip(x,y,t)、横波阻抗Is(x,y,t)、密度ρ(x,y,t)带入(11)，求得孔隙度φ(x,y,t)与含气饱和度Sg(x,y,t)。
[0058] 在一个示例中，基于地震资料与地温梯度，分别获得地层压力与地层温度，进而获得体积系数。
[0059] 在一个示例中，体积系数可以为：
[0060]
[0061] 其中，Bg表示体积系数，Z表示天然气压缩因子，T表示地层温度；P表示地层压力。
[0062] 具体地，地层温度可以为
[0063] T＝T0+d*ΔT (13)
[0064] 其中，T表示地层温度，T0表示恒温带温度，单位均为℃，d表示地层埋深，ΔT表示地温梯度，地温梯度是指在恒温带以下每增加100米深地温的变化情况，在实际应用过程以实际地温测试结果为准。
[0065] 具体地，由于求得的孔隙度φ(x,y,t)与含气饱和度Sg(x,y,t)均为时间域的数据，因此将(12)转换为时间域可以为：
[0066]
[0067] 其中，Bg(x,y,t)表示时间域的体积系数，Z表示天然气压缩因子，T(x,y,t)表示时间域的地层温度，P(x,y,t)表示时间域的地层压力。
[0068] 在一个示例中，基于孔隙度、含气饱和度、体积系数与地层压力，建立游离气含气量的计算模型，获得游离气含气量。
[0069] 在一个示例中，游离气含气量可以为：
[0070]
[0071] 其中，Gg表示游离气含气量，φ表示孔隙度，Sg表示含气饱和度，ρ表示密度。
[0072] 具体地，由于求得的孔隙度φ(x,y,t)、含气饱和度Sg(x,y,t)与体积系数Bg(x,y,t)均为时间域的数据，因此将(13)转换为时间域可以为：
[0073]
[0074] 其中，Gg(x,y,t)表示时间域的游离气含气量。
[0075] 应用示例
[0076] 为便于理解本发明实施方式的方案及其效果，以下给出一个具体应用示例。本领域技术人员应理解，该示例仅为了便于理解本发明，其任何具体细节并非意在以任何方式限制本发明。
[0077] 101)、首先，基于工区叠前地震数据体，通过地震反演的方法，利用Jason反演软件的RockTrace模
块，也可以替换成其他可以得到纵波阻抗、横波阻抗和密度的叠前反演方法或者软件，选择合适的子波、低频模型以及合适的反演参数，获得纵向分辨率较高的纵波阻抗Ip(x,y,t)、横波阻抗Is(x,y,t)、密度ρ(x,y,t)。
[0078] 对深度域的测井数据进行时深转换，获得时间域的测井数据，时深转换的关系式为：
[0079]
[0080] 其中，t为时间域，z为深度域，vp(z)为纵波速度。获得的时间域的测井数据可以包括：时间域的纵波阻抗曲线Ip(t)、时间域的横波阻抗Is(t)、时间域的密度曲线ρ(t)、时间域的孔隙度φ(t)、时间域含气饱和度Sg(t)。
[0081] 102)、基于工区弹性参数与时间域测井数据，建立页岩储层各向异性岩石物理模型，进而构建孔隙度、含气饱和度与工区弹性参数之间的函数关系为：
[0082] {Ip(t),Is(t),ρ(t)}＝f(φ(t),Sg(t)…) (2)
[0083] 采用贝叶斯分类算法，目标函数表示为在已知工区弹性参数条件下，孔隙度φ和含气饱和度Sg的最大后验概率分布：
[0084] {φ,Sg}＝arg Maxp(φj,Sgj…|Ip,Is,ρ)j＝1,2…,Nc (3)
[0085] 其中，{φ,Sg}表示孔隙度和含气饱和度的解，p(φj,Sgj…|Ip,Is,ρ)表示孔隙度与含气饱和度的后验概率，Nc表示分类的类别数，j表示类别编号，本领域技术人员可以通过具体情况设定Nc用来控制预测精度。
[0086] 即，如果：
[0087] p(φj,Sgj…|Ip,Is,ρ)>p(φm,Sgm…|Ip,Is,ρ)1≤j,m≤Nc,j≠m (4)[0088] 则：
[0089] φ＝φj，Sg＝Sgj (5)
[0090] (4)与(5)的意义为，在已知纵波阻抗Ip，横波阻抗Is，密度ρ的情况下，求取孔隙度φ和含气饱和度Sg为不同值时的后验概率，其中，最大后验概率的值φ＝φj，Sg＝Sgj为孔隙度和含气饱和度的解。
[0091] 基于贝叶斯公式，(3)可以改写为：
[0092] {φ,Sg}＝arg Max{p(Ip,Is,ρ|φ,Sg…)·p(φ,Sg…)} (6)
[0093] 其中，Ip表示纵波阻抗，Is表示横波阻抗，ρ表示密度，p(φ,Sg…)表示孔隙度与含气饱和度的先验分布，p(Ip,Is,ρ|φ,Sg…)表示似然函数。
[0094] 假设孔隙度与含气饱和度服从多维高斯分布，即孔隙度与含气饱和度的先验分布可以为：
[0095]
[0096] 其中，F代表多维高斯分布， 为多维高斯分布的均值， 为多维高斯分布的方差，Nr为岩石物理模型中，孔隙度与含气饱和度的变量的个数，权重系数αk满足[0097] 基于孔隙度与含气饱和度的先验分布，利用MCMC抽样Metropolis Hastings抽样算法对孔隙度与含气饱和度进行随机模拟，可以得到孔隙度与含气饱和度的随机模拟结果{φi,Sgi…}i＝1…N。
[0098] 可以将{φi,Sgi…}i＝1…N带入(2)，得到工区弹性参数、孔隙度与含气饱和度的联合随机模拟结果为：
[0099] {Ipi,Isi,ρi,φi,Sgi…}i＝1…N (8)
[0100] 其中，N为随机抽样次数，i表示第i次随机抽样，i＝1…N，Ipi为第i次随机抽样的纵波阻抗，Isi为第i次随机抽样的横波阻抗，ρi为第i次随机抽样的密度，φi为第i次随机抽样的孔隙度，Sgi为第i次随机抽样的含气饱和度。
[0101] 基于贝叶斯分类算法的条件独立性假设，目标函数的计算形式改写为：
[0102] p(φj,Sgj…|Ip,Is,ρ)＝p(Ip|φj,Sgj,…)·p(Is|φj,Sgj,…)·p(ρ|φj,Sgj,…)·p(φj,Sgj,…) (10)
[0103] 将(8)带入(10)可以得到：
[0104]
[0105] 其中，N表示随机抽样次数，n表示计数统计，n(Ip∩φi,Sgi,…)表示统计{Ipi,Isi,ρi,φi,Sgi,…}i＝1…N中φ的值等于φi，Sg的值为Sgi，且纵波阻抗为Ip的随机采样样点的个数，n(Is∩φi,Sgi,…)表示统计{Ipi,Isi,ρi,φi,Sgi,…}i＝1…N中φ的值等于φi，Sg的值为Sgi，且横波阻抗为Is的随机采样样点的个数，n(ρ∩φi,Sgi,…)表示统计{Ipi,Isi,ρi,φi,Sgi,…}i＝1…N中φ的值等于φi，Sg的值为Sgi，且密度为ρ的随机采样样点的个数。
[0106] 将纵波阻抗Ip(x,y,t)、横波阻抗Is(x,y,t)、密度ρ(x,y,t)带入(11)，求得孔隙度φ(x,y,t)与含气饱和度Sg(x,y,t)。
[0107] 103)、基于地震资料与地温梯度，获得地层温度为
[0108] T＝T0+d*ΔT (13)
[0109] 其中，T表示地层温度，T0表示恒温带温度，d表示地层埋深，ΔT表示地温梯度。
[0110] 基于地层压力与地层温度，获得体积系数为：
[0111]
[0112] 其中，Bg(x,y,t)表示时间域的体积系数，Z表示天然气压缩因子，T(x,y,t)表示时间域的地层温度，P(x,y,t)表示时间域的地层压力。
[0113] 104)、图2示出了根据本发明的一个实施方式的游离气含气量的示意图。基于孔隙度、含气饱和度、体积系数与地层压力，建立游离气含气量的计算模型，获得游离气含气量为：
[0114]
[0115] 其中，Gg(x,y,t)表示时间域的游离气含气量。
[0116] 本发明在贝叶斯
框架下，基于纵波阻抗、横波阻抗和密度参数多属性约束的孔隙度、含气饱和度参数地震预测方法，结合地层压力和地层温度预测结果，实现对页岩储层游离气含气量的高精度预测，能够获得更加稳定可靠的游离气含气量预测结果。
[0117] 本领域技术人员应理解，上面对本发明的实施方式的描述的目的仅为了示例性地说明本发明的实施方式的有益效果，并不意在将本发明的实施方式限制于所给出的任何示例。
[0118] 实施方式2
[0119] 根据本发明的实施方式，提供了一种页岩储层游离气含气量地震预测系统，所述系统可以包括：用来基于工区叠前地震数据体与深度域的测井数据，分别获取工区弹性参数与时间域的测井数据的单元；用来基于工区弹性参数与时间域的测井数据，建立岩石物理模型，结合贝叶斯反演算法，获得孔隙度与含气饱和度的单元；用来基于地震资料与地温梯度，分别获得地层压力与地层温度，进而获得体积系数的单元；以及用来基于孔隙度、含气饱和度、体积系数与地层压力，建立游离气含气量的计算模型，获得游离气含气量的单元，其中，工区弹性参数包括纵波阻抗、横波阻抗及密度。
[0120] 该实施方式通过地层温度、地层压力、孔隙度及含气饱和度等因素约束，实现页岩储层游离气含气量的高精度预测。
[0121] 在一个示例中，体积系数可以为：
[0122]
[0123] 其中，Bg表示体积系数，Z表示天然气压缩因子，T表示地层温度，P表示地层压力，[0124] 在一个示例中，游离气含气量可以为：
[0125]
[0126] 其中，Gg表示游离气含气量，φ表示孔隙度，Sg表示含气饱和度，ρ表示密度。
[0127] 在一个示例中，孔隙度与含气饱和度可以包括：基于贝叶斯反演算法，构建目标函数，并确定孔隙度与含气饱和度的先验分布；基于孔隙度与含气饱和度的先验分布，结合抽样算法，获得孔隙度、含气饱和度与工区弹性参数的联合随机模拟结果；基于贝叶斯分类算法与孔隙度、含气饱和度与工区弹性参数的联合随机模拟结果，求解目标函数；以及将工区弹性参数带入目标函数，获得孔隙度与含气饱和度。
[0128] 在一个示例中，目标函数可以为：
[0129] {φ,Sg}＝arg Max{p(Ip,Is,ρ|φ,Sg…)·p(φ,Sg…)} (6)
[0130] 其中，{φ,Sg}表示孔隙度和含气饱和度的解，Ip表示纵波阻抗，Is表示横波阻抗，ρ表示密度，p(φ,Sg…)表示反演目标参数先验分布，p(Ip,Is,ρ|φ,Sg…)表示似然函数。
[0131] 在一个示例中，求解目标函数可以为：
[0132]
[0133] 其中，N表示随机抽样次数，n表示计数统计，n(Ip∩φi,Sgi,…)表示统计{Ipi,Isi,ρi,φi,Sgi,…}i＝1…N中φ的值等于φi，Sg的值为Sgi，且纵波阻抗为Ip的随机采样样点的个数，n(Is∩φi,Sgi,…)表示统计{Ipi,Isi,ρi,φi,Sgi,…}i＝1…N中φ的值等于φi，Sg的值为Sgi，且横波阻抗为Is的随机采样样点的个数，n(ρ∩φi,Sgi,…)表示统计{Ipi,Isi,ρi,φi,Sgi,…}i＝1…N中φ的值等于φi，Sg的值为Sgi，且密度为ρ的随机采样样点的个数。
[0134] 本领域技术人员应理解，上面对本发明的实施方式的描述的目的仅为了示例性地说明本发明的实施方式的有益效果，并不意在将本发明的实施方式限制于所给出的任何示例。
[0135] 以上已经描述了本发明的各实施方式，上述说明是示例性的，并非穷尽性的，并且也不限于所披露的各实施方式。在不偏离所说明的各实施方式的范围和精神的情况下，对于本技术领域的普通技术人员来说许多
修改和变更都是显而易见的。本文中所用术语的选择，旨在最好地解释各实施方式的原理、实际应用或对市场中的技术的改进，或者使本技术领域的其它普通技术人员能理解本文披露的各实施方式。
