CN115586572A - 一种孔隙参数与储层参数的地震岩石物理解析反演方法 - Google Patents

一种孔隙参数与储层参数的地震岩石物理解析反演方法 Download PDF

Info

Publication number
CN115586572A
CN115586572A CN202211122138.XA CN202211122138A CN115586572A CN 115586572 A CN115586572 A CN 115586572A CN 202211122138 A CN202211122138 A CN 202211122138A CN 115586572 A CN115586572 A CN 115586572A
Authority
CN
China
Prior art keywords
parameters
model
pore
elastic
seismic
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN202211122138.XA
Other languages
English (en)
Other versions
CN115586572B (zh
Inventor
郭强
巴晶
雒聪
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Hohai University HHU
Original Assignee
Hohai University HHU
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Hohai University HHU filed Critical Hohai University HHU
Priority to CN202211122138.XA priority Critical patent/CN115586572B/zh
Publication of CN115586572A publication Critical patent/CN115586572A/zh
Application granted granted Critical
Publication of CN115586572B publication Critical patent/CN115586572B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
    • G01V1/28—Processing seismic data, e.g. for interpretation or for event detection
    • G01V1/30—Analysis
    • G01V1/306—Analysis for determining physical properties of the subsurface, e.g. impedance, porosity or attenuation profiles
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
    • G01V1/28—Processing seismic data, e.g. for interpretation or for event detection
    • G01V1/282—Application of seismic models, synthetic seismograms
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V2210/00—Details of seismic processing or analysis
    • G01V2210/60—Analysis
    • G01V2210/62—Physical property of subsurface
    • G01V2210/624—Reservoir parameters

Landscapes

  • Engineering & Computer Science (AREA)
  • Remote Sensing (AREA)
  • Physics & Mathematics (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Acoustics & Sound (AREA)
  • Environmental & Geological Engineering (AREA)
  • Geology (AREA)
  • General Life Sciences & Earth Sciences (AREA)
  • General Physics & Mathematics (AREA)
  • Geophysics (AREA)
  • Geophysics And Detection Of Objects (AREA)

Abstract

本申请公开了一种孔隙参数与储层参数的地震岩石物理解析反演方法，包括：采用地震道集数据，通过叠前地震反演，获取观测弹性参数模型；设置基本岩石物理参数，基于测井数据，反演得到井旁孔隙纵横比；计算模型参数的概率分布函数，构建模型参数的初始模型；利用精确正演算子模拟弹性参数，计算模拟弹性参数与观测弹性参数的误差；基于初始模型，计算精确正演算子对模型参数的偏导矩阵，构建线性正演算子；基于线性正演算子，解析计算模型参数的后验期望值；更新初始模型为后验期望值，重复步骤四至步骤六，直到误差降至预设范围内，迭代停止，输出模型参数的后验期望值。

Description

一种孔隙参数与储层参数的地震岩石物理解析反演方法
技术领域
本发明属于非常规储层地震勘探技术领域，特别涉及孔隙参数与储层参数的地震岩石物理解析反演方法。
背景技术
地震岩石物理反演结合了地震反演技术与岩石物理模型，通过地震数据或属性定量预测储层参数，旨在为油气储层的预测、勘探及开发提供技术支撑。另一方面，随着油气勘探的不断推进，常规储层日益枯竭，致密砂岩、碳酸盐岩等非常规储层已称为增储上产的重点领域，但此类储层普遍发育复杂孔隙结构，严重制约了常规地震岩石物理反演方法的适用性。
常规地震岩石物理反演方法主要预测储层参数(孔隙度、含水饱和度等)，例如文献[1]-[4]公开的现有技术，但上述方法在岩石物理建模中通常将孔隙参数(孔隙纵横比)作为已知参数，并在反演过程中采用固定的孔隙纵横比。研究发现，岩石的孔隙纵横比对地震响应具有重要影响，尤其针对发育复杂孔隙结构的储层岩石[5]-[6]；此外，非常规储层由于地质构造及沉积环境等因素，储层岩石的孔隙结构异常复杂[7]-[8]，因此基于固定孔隙参数的地震岩石物理反演方法难以适用。目前，考虑孔隙参数或结构的地震岩石物理反演方法研究起步较晚，例如Li等(2022)提出了一种预测硬孔孔隙百分比的地震反演方法[9]，但该方法仅定性考虑了孔隙结构的影响，无法直接预测孔隙参数，难以定量评估储层的孔隙结构；李宏兵等(2013)和Teillet等(2021)提出了预测孔隙纵横比的地震反演方法[10]-[11]，虽然该方法能够定量预测孔隙参数，但无法同步预测储层参数，限制了储层地震预测的准确度；李红兵等(2021)发展了孔隙扁度与含水饱和度的同步反演方法[12]，但该反演方法依赖全局优化，求解过程效率低，并且计算结果无解析表达，存在不确定性。
综上，针对具有复杂孔隙结构的非常规储层，需定量考虑孔隙结构的影响，因此，开展同步预测孔隙参数与储层参数的解析反演技术方法研究，对提高致密砂岩与碳酸盐岩类非常规储层参数定量预测精度具有重要意义。
上述引用非专利文献如下：
[1]Yin,X.,Sun,R.,Wang,B.,Zhang,G.,2014,Simultaneous inversion ofpetrophysical parameters based on geostatistical a priori information[J].Applied Geophysics,11,311–320.
[2]de Figueiredo,L.P.,Grana,D.,Santos,M.,et al.,2017,Bayesian seismicinversion based on rock-physics prior modeling for the joint estimation ofacoustic impedance,porosity and lithofacies[J].Journal of ComputationalPhysics,336,128–142.
[3]李志勇,张家树,蔡涵鹏等,2017,基于Hampel三截尾函数的储层弹性和物性参数同步反演[J].石油物探,56(2),261–272.
[4]张佳佳,印兴耀,张广智等,2020,基于线性化岩石物理反演的物性参数预测方法[J].石油勘探与开发,47(1),57–64.
[5]Weger,R.J.,Eberli,G.P.,Baechle,G.T.,et al.,2009,Quantification ofpore structure and its effect on sonic velocity and permeability incarbonates[J].AAPG Bulletin,93,1297–1317.
[6]Lehocki,I.,Avseth,P.,2021,From cradle to grave:how burial historycontrols the rock-physics properties of quartzose sandstones[J].GeophysicalProspecting,69,629–649.
[7]李王鹏,刘忠群,胡宗全等,2021,四川盆地川西坳陷新场须家河组二段致密砂岩储层裂缝发育特征及主控因素[J].石油与天然气地质,42(4),884–897.
[8]Ren,W.,Shi,W.,Xie,X.,et al.,2020,Clay mineral content,type,andtheir effects on pore throat structure and reservoir properties:Insight fromthe Permian tight sandstones in the Hangjinqi area,north Ordos Basin,China[J].Marine and Petroleum Geology,115,104281.
[9]Li,K.,Yin,X.,Zong,Z.,et al.,2022,Estimation of porosity,fluid bulkmodulus,and stiff-pore volume fraction using a multitrace Bayesian amplitude-variation-with-offset petrophysics inversion in multiporosity reservoirs[J].Geophysics,87,M25–M41.
[10]李宏兵,张佳佳,姚逢昌,2013,岩石的等效孔隙纵横比反演及其应用[J].地球物理学报,56(2),608–615.
[11]Teillet,T.,Fournier,F.,Zhao,L.,et al.,2021,Geophysical pore typeinversion in carbonate reservoir:Integration of cores,well logs,and seismicdata(Yadana field,offshore Myanmar)[J].Geophysics,86,B149–B164.
[12]李红兵,张佳佳,潘豪杰等,2021,基于弹性阻抗的孔隙结构与物性参数非线性同步反演[J].中国科学·地球科学,51,1166–1180.
发明内容
解决的技术问题：本申请主要是提出一种孔隙参数与储层参数的地震岩石物理解析反演方法，解决现有技术中存在的致密砂岩、碳酸盐岩储层普遍发育复杂孔隙结构，严重制约了常规地震岩石物理反演方法的适用性，非常规储层由于地质构造及沉积环境等因素，储层岩石的孔隙结构异常复杂，因此基于固定孔隙参数的地震岩石物理反演方法难以适用；仅定性考虑了孔隙结构的影响，无法直接预测孔隙参数，难以定量评估储层的孔隙结构；无法同步预测储层参数，限制了储层地震预测的准确度；反演方法依赖全局优化，求解过程效率低，并且计算结果无解析表达，存在不确定性等技术问题，旨在提高致密砂岩与碳酸盐岩类非常规储层参数定量预测精度。
技术方案：
一种孔隙参数与储层参数的地震岩石物理解析反演方法，具体包括如下步骤：
步骤一，采用地震道集数据，通过叠前地震反演，获取观测弹性参数模型；
步骤二，设置基本岩石物理参数，基于测井数据，反演得到井旁孔隙纵横比；
步骤三，计算模型参数的概率分布函数，构建模型参数的初始模型；
步骤四，利用地震岩石物理正演算子模拟弹性参数，计算模拟弹性参数与观测弹性参数的误差；
步骤五，基于初始模型，计算精确正演算子对模型参数的偏导矩阵，构建线性正演算子；
步骤六，基于线性正演算子，解析计算模型参数的后验期望值；
步骤七，更新模型参数的初始模型为后验期望值，重复步骤四至步骤六，直至模拟与观测弹性参数的误差降低到预设范围内，迭代停止，输出模型参数的后验期望值，即为孔隙参数与储层参数的反演结果。
作为本发明的一种优选技术方案，步骤一具体包括：
获取工区PP波地震道集数据，通过叠前地震反演，构建目标函数
Figure BDA0003846922720000031
其中，d为地震道集数据，m为弹性参数，
Figure BDA0003846922720000032
为弹性参数的平均值，Σm为弹性参数的协方差矩阵，G为Zoeppritz地震正演模型，σd为地震数据噪声的标准差，T为转置符号；弹性参数m包括纵波速度、横波速度和密度；
采用贝叶斯线性反演方法，得到满足式(1)目标函数的最优解，获取纵波速度、横波速度和密度参数模型。
作为本发明的一种优选技术方案：所述步骤二具体包括：
设置基本岩石物理参数，包括岩石基质的体积模量、岩石基质的剪切模量、岩石基质的密度、烃类的体积模量、烃类的密度、卤水的体积模量和卤水的密度；
基于测井观测的纵波速度与横波速度，通过匹配模型预测与井旁观测数据，反演得到井旁孔隙纵横比，构建目标函数
Figure BDA0003846922720000041
其中，α为孔隙纵横比，σ1与σ2为误差项权系数，
Figure BDA0003846922720000042
和
Figure BDA0003846922720000043
分别为模型预测和井旁观测的纵波速度，
Figure BDA0003846922720000044
和
Figure BDA0003846922720000045
分别为模型预测和井旁观测的横波速度，T为转置符号；
采用粒子群优化算法，得到满足式(2)目标函数的最优解，获取井旁孔隙纵横比参数。
作为本发明的一种优选技术方案：步骤三具体包括：
通过测井观测数据获得井旁孔隙度与井旁含水饱和度，结合反演得到的井旁孔隙纵横比，构建模型参数的混合高斯模型，其概率分布函数为；
Figure BDA0003846922720000046
其中，z为模型参数，Nk表示第k个高斯分布，
Figure BDA0003846922720000047
和λk分别为Nk的期望、协方差和权系数，C为高斯分布个数；模型参数包括孔隙纵横比、孔隙度与含水饱和度；
通过期望最大算法计算混合高斯模型，得到C个高斯分布的期望、协方差和权系数，再通过测井插值，得到模型参数的初始模型。
作为本发明的一种优选技术方案：步骤四具体包括：
利用解耦微分等效介质模型计算岩石骨架的弹性模量
Kd(φ)＝Km(1-φ)p (4)
μd(φ)＝μm(1-φ)q (5)
其中，φ为孔隙度，Kd和μd分别为岩石骨架的体积模量和剪切模量，Km和μm分别为岩石基质的体积模量和剪切模量，p和q为孔隙结构系数；
基于岩石骨架的弹性模量，利用Gassmann方程计算含流体岩石的弹性模量
Figure BDA0003846922720000048
μs＝μd (7)
其中，Ks和μs分别为含流体岩石的体积模量和剪切模量，Kf为混合流体的体积模量，采用Wood模型计算
Figure BDA0003846922720000051
其中，Kw和Kh分别为卤水和烃类的体积模量，Sw为含水饱和度；
根据含流体岩石的弹性模量，得到模拟弹性参数，
Figure BDA0003846922720000052
Figure BDA0003846922720000053
ρ＝φρf+(1-φ)ρm (11)
其中，VP、VS和ρ分别为纵波速度、横波速度和密度，ρf和ρm分为流体和岩石基质的密度；计算模拟弹性参数与观测弹性参数的误差ε
ε＝||mest-mobs||2 (12)
其中，mest为模拟弹性参数，mobs为观测弹性参数即由步骤一获取的观测弹性参数模型。
作为本发明的一种优选技术方案：步骤五具体包括：
基于模型参数的初始模型，计算精确地震岩石物理正演算子对模型参数的偏导矩阵
Figure BDA0003846922720000054
其中，R为精确正演算子，通过联立公式(4)–(11)获得，z为模型参数，zo为模型参数的初始模型，G为精确正演算子的偏导矩阵，模型参数包括孔隙纵横比、孔隙度与含水饱和度，其具体计算方法如下：
纵波速度对孔隙度的偏导数为：
Figure BDA0003846922720000055
其中，
Figure BDA0003846922720000061
Figure BDA0003846922720000062
纵波速度对含水饱和度的偏导数为：
Figure BDA0003846922720000063
其中，
Figure BDA0003846922720000064
纵波速度对孔隙纵横比的偏导数为：
Figure BDA0003846922720000065
其中，
Figure BDA0003846922720000066
Figure BDA0003846922720000067
其中，
Figure BDA0003846922720000068
Figure BDA0003846922720000069
Figure BDA00038469227200000610
横波速度对孔隙度的偏导数为：
Figure BDA00038469227200000611
横波速度对含水饱和度的偏导数为：
Figure BDA0003846922720000071
横波速度对孔隙纵横比的偏导数为：
Figure BDA0003846922720000072
其中，
Figure BDA0003846922720000073
其中，
Figure BDA0003846922720000074
Figure BDA0003846922720000075
Figure BDA0003846922720000076
Figure BDA0003846922720000077
Figure BDA0003846922720000078
Figure BDA0003846922720000079
Figure BDA00038469227200000710
Figure BDA00038469227200000711
Figure BDA00038469227200000712
Figure BDA0003846922720000081
密度对孔隙度的偏导数为：
Figure BDA0003846922720000082
密度对含水饱和度的偏导数为：
Figure BDA0003846922720000083
密度对孔隙纵横比的偏导数为：
Figure BDA0003846922720000084
在公式(14)–(41)中，VP、VS和ρ分别为纵波速度、横波速度和密度，Kd和μd分别为岩石骨架的体积模量和剪切模量，Ks和μs分别为含流体岩石的体积模量和剪切模量，Km和μm分别为岩石基质的体积模量和剪切模量，Kf为混合流体的体积模量，p和q为孔隙结构系数，ρf为流体的密度，ρm为岩石基质的密度，ρw为卤水的密度，ρh为烃类的密度，φ为孔隙度，α为孔隙纵横比，Sw为含水饱和度，F1、F2、F3、F4、F5、F6、F7、F8与F9为9个不同的孔隙因子系数，a和b为孔隙扁度系数，R为基质模量系数；
基于精确正演算子的偏导矩阵，构建线性正演算子
m＝Gz+R(zo)-Gzo+e (42)
其中，m为弹性参数，e为随机误差。
作为本发明的一种优选技术方案：步骤六具体包括：
利用公式(42)线性正演算子，基于贝叶斯线性反演理论，解析计算每个高斯分布中模型参数的期望值
Figure BDA0003846922720000085
其中，G为精确正演算子的偏导矩阵，m为弹性参数，Σe为随机误差的协方差，
Figure BDA0003846922720000086
和
Figure BDA0003846922720000087
为第k个高斯分布的期望和协方差，
Figure BDA0003846922720000088
为第k个高斯分布的模型参数期望值，T为转置符号；
计算模型参数的后验期望值
Figure BDA0003846922720000091
其中，Nk表示第k个高斯分布，C为高斯分布个数，λk为Nk的权系数，μpost为模型参数的后验期望值。
作为本发明的一种优选技术方案：步骤七预设范围为观测弹性参数值的百分之一。
有益效果：本申请所述孔隙参数与储层参数的地震岩石物理解析反演方法采用的以上技术方案与现有技术相比，具有以下技术效果：
1、用于同步解析反演孔隙参数及储层参数，可以提高致密砂岩与碳酸盐岩类非常规储层的储层参数定量预测精度。
2、本发明推导了精确正演算子对孔隙纵横比、孔隙度与含水饱和度的偏导矩阵，可以同步反演孔隙参数与储层参数；选用了线性反演策略以解析计算反演结果，提高了反演时效性及稳定性，相较于传统的地震岩石物理反演方法，本发明提高了复杂孔隙储层的储层参数反演精度。
附图说明：
图1为本申请一种孔隙参数与储层参数的地震岩石物理解析反演方法的流程示意图；
图2为本申请实施例中地震道集数据与观测弹性参数模型；其中a为实施例中的地震道集数据图，b为通过叠前地震反演获取的观测弹性参数图；
图3为本申请实施例中观测储层参数与井旁孔隙纵横比；其中a为实施例中的观测储层参数图，b为反演得到的井旁孔隙纵横比图；
图4为本申请实施例中利用常规方法的储层参数反演结果；
图5为本申请实施例中利用本申请方法的孔隙参数与储层参数反演结果。
具体实施方式
下面将结合本发明的附图，对本发明的技术方案进行清楚、完整地描述，显然，所描述的实施例仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动条件下所获得的所有其它实施例，都属于本发明保护的范围。
实施例1
如图1所示，本实施例提供了一种孔隙参数与储层参数的地震岩石物理解析反演方法，具体包括以下步骤：
一、采用地震道集数据，通过叠前地震反演，获取观测弹性参数模型及构建孔隙和储层参数初始模型：
获取工区PP波地震道集数据，通过叠前地震反演，其目标函数如下
Figure BDA0003846922720000101
其中，d为地震道集数据，m为弹性参数，G为Zoeppritz地震正演模型，σd为地震数据噪声的标准差，T为转置符号；
Figure BDA0003846922720000102
和Σm分别为弹性参数的平均值与协方差矩阵，可通过井旁弹性参数估算所得；弹性参数m包括纵波速度、横波速度和密度；
采用贝叶斯线性反演方法，得到满足式(1)目标函数的最优解，获取纵波速度、横波速度和密度参数模型，即观测弹性参数模型；
设置基本岩石物理参数，包括岩石基质的体积模量Km、岩石基质的剪切模量μm、岩石基质的密度ρm，烃类的体积模量Kh、烃类的密度ρh、卤水的体积模量Kw、卤水的密度ρw；实施例中基本岩石物理参数设置如下：Km＝50GPa，μm＝23GPa，Kh＝0.012GPa，Kw＝2.5GPa，ρm＝2.62g/cm3，ρh＝0.02g/cm3，ρw＝1.10g/cm3；
基于测井观测的纵波速度与横波速度，反演得到井旁孔隙纵横比，其目标函数如下
Figure BDA0003846922720000103
其中，α为孔隙纵横比，σ1与σ2为误差项权系数，
Figure BDA0003846922720000104
和
Figure BDA0003846922720000105
分别为模型预测和井旁观测的纵波速度，
Figure BDA0003846922720000106
和
Figure BDA0003846922720000107
分别为模型预测和井旁观测的横波速度，T为转置符号；
采用粒子群优化算法，实施例中设置σ1＝0.6，σ1＝0.4，得到满足式(2)目标函数的最优解，获取井旁孔隙纵横比参数；
通过测井观测数据获得井旁孔隙度与井旁含水饱和度，结合反演所得的井旁孔隙纵横比，构建混合高斯模型，其概率分布函数为
Figure BDA0003846922720000108
其中，z为模型参数，模型参数包括孔隙纵横比、孔隙度与含水饱和度，Nk表示第k个高斯分布，
Figure BDA0003846922720000109
和λk分别为Nk的期望、协方差和权系数，C为高斯分布个数；
实施例中设置C＝2，通过期望最大算法计算混合高斯模型，得到模型参数的协方差与均值，通过测井插值，得到模型参数的初始模型。
二、计算精确正演算子对模型参数的偏导矩阵：
利用解耦微分等效介质模型计算岩石骨架的弹性模量
Kd(φ)＝Km(1-φ)p (4)
μd(φ)＝μm(1-φ)q (5)
其中，φ为孔隙度，Kd和μd分别为岩石骨架的体积模量和剪切模量，Km和μm分别为岩石基质的体积模量和剪切模量，p和q为孔隙结构系数；
基于岩石骨架的弹性模量，利用Gassmann方程计算含流体岩石的弹性模量
Figure BDA0003846922720000111
μs＝μd (7)
其中，Ks和μs分别为含流体岩石的体积模量和剪切模量，Kf为混合流体的体积模量，采用Wood模型计算
Figure BDA0003846922720000112
其中，Kw和Kh分别为卤水和烃类的体积模量，Sw为含水饱和度；
基于含流体岩石的弹性模量，得到模拟弹性参数
Figure BDA0003846922720000113
Figure BDA0003846922720000114
ρ＝φρf+(1-φ)ρm (11)
其中，VP、VS和ρ分别为纵波速度、横波速度和密度，ρf和ρm分为流体和岩石基质的密度；
计算模拟弹性参数与观测弹性参数的误差ε
ε＝||mest-mobs||2 (12)
其中，mest为模拟弹性参数，mobs为观测弹性参数，即由步骤一获取的观测弹性参数模型；
基于模型参数的初始模型，计算精确地震岩石物理正演算子对模型参数的偏导矩阵
Figure BDA0003846922720000121
其中，R为精确正演算子，通过联立公式(4)–(11)获得，z为模型参数，模型参数包括孔隙纵横比、孔隙度与含水饱和度，zo为模型参数的初始模型，G为精确正演算子的偏导矩阵，其具体计算方法如下：
计算纵波速度对孔隙度的偏导数：
Figure BDA0003846922720000122
其中，
Figure BDA0003846922720000123
Figure BDA0003846922720000124
计算纵波速度对含水饱和度的偏导数：
Figure BDA0003846922720000125
其中，
Figure BDA0003846922720000126
计算纵波速度对孔隙纵横比的偏导数：
Figure BDA0003846922720000127
其中，
Figure BDA0003846922720000131
Figure BDA0003846922720000132
其中，
Figure BDA0003846922720000133
Figure BDA0003846922720000134
Figure BDA0003846922720000135
计算横波速度对孔隙度的偏导数：
Figure BDA0003846922720000136
计算横波速度对含水饱和度的偏导数：
Figure BDA0003846922720000137
计算横波速度对孔隙纵横比的偏导数：
Figure BDA0003846922720000138
其中，
Figure BDA0003846922720000139
其中，
Figure BDA00038469227200001310
Figure BDA0003846922720000141
Figure BDA0003846922720000142
Figure BDA0003846922720000143
Figure BDA0003846922720000144
Figure BDA0003846922720000145
Figure BDA0003846922720000146
Figure BDA0003846922720000147
Figure BDA0003846922720000148
Figure BDA0003846922720000149
计算密度对孔隙度的偏导数：
Figure BDA00038469227200001410
计算密度对含水饱和度的偏导数：
Figure BDA00038469227200001411
计算密度对孔隙纵横比的偏导数：
Figure BDA00038469227200001412
在公式(14)–(41)中，VP、VS和ρ分别为纵波速度、横波速度和密度，Kd和μd分别为岩石骨架的体积模量和剪切模量，Ks和μs分别为含流体岩石的体积模量和剪切模量，Km和μm分别为岩石基质的体积模量和剪切模量，Kf为混合流体的体积模量，p和q为孔隙结构系数，ρf为流体的密度，ρm为岩石基质的密度，ρw为卤水的密度，ρh为烃类的密度，φ为孔隙度，α为孔隙纵横比，Sw为含水饱和度，F1、F2、F3、F4、F5、F6、F7、F8与F9为9个不同的孔隙因子系数，a和b为孔隙扁度系数，R为基质模量系数。
基于公式(14)-(41)的计算结果，形成精确正演算子的偏导矩阵G，进而构建线性正演算子
m＝Gz+R(zo)-Gzo+e (42)
其中，m为弹性参数，e为随机误差。
三、解析计算孔隙参数与储层参数的反演结果：
利用公式(42)线性正演算子，基于贝叶斯线性反演理论，解析计算每个高斯分布中模型参数的期望值
Figure BDA0003846922720000151
其中，G为精确正演算子的偏导矩阵，m为弹性参数，Σe为随机误差的协方差，
Figure BDA0003846922720000152
和
Figure BDA0003846922720000153
为第k个高斯分布的期望和协方差，
Figure BDA0003846922720000154
为第k个高斯分布的模型参数期望值；
计算模型参数的后验期望值
Figure BDA0003846922720000155
其中，Nk表示第k个高斯分布，C为高斯分布个数，λk为Nk的权系数，μpost为模型参数的后验期望值；
更新模型参数的初始模型为后验期望值，重复上述步骤，直至模拟与观测弹性参数的误差降低到观测弹性参数值的百分之一，迭代停止，输出模型参数的后验期望值，即为孔隙参数与储层参数的反演结果。
图2(a)为本实施例中的地震道集数据，图2(b)为通过叠前地震反演获取的观测弹性参数，包括纵波速度、横波速度和密度。图3(a)为本实施例中的观测储层参数，包括孔隙度与汗水饱和，图3(b)为反演得到的井旁孔隙纵横比。图4为本实施例中利用常规方法(即固定孔隙纵横比)的储层参数反演结果。从图4中可见，孔隙度与含水饱和度的反演结果(虚线)与观测数据(实线)均有较大偏差。图5为本实施例中利用本发明方法的孔隙参数与储层参数反演结果。从图5中可见，孔隙度与含水饱和度的反演结果(虚线)与观测数据(实线)基本一致，此外孔隙纵横比的反演结果(虚线)也与井旁数据(实线)基本吻合。经计算所得，图5中孔隙度与含水饱和度的反演结果与观测数据的平均误差分别为0.056与0.258，而图4中的平均误差分别为0.092与0.517，由此可见，本发明方法可以有效提高储层参数反演结果的准确度。
尽管已经示出和描述了本发明的实施例，对于本领域的普通技术人员而言，可以理解在不脱离本发明的原理和精神的情况下可以对这些实施例进行多种变化、修改、替换和变型，本发明的范围由所附权利要求及其等同物限定。

Claims (8)

1.一种孔隙参数与储层参数的地震岩石物理解析反演方法，其特征在于，具体包括如下步骤：
步骤一，采用地震道集数据，通过叠前地震反演，获取观测弹性参数模型；
步骤二，设置基本岩石物理参数，基于测井数据，反演得到井旁孔隙纵横比；
步骤三，计算模型参数的概率分布函数，构建模型参数的初始模型；
步骤四，利用地震岩石物理正演算子模拟弹性参数，计算模拟弹性参数与观测弹性参数的误差；
步骤五，基于初始模型，计算精确正演算子对模型参数的偏导矩阵，构建线性正演算子；
步骤六，基于线性正演算子，解析计算模型参数的后验期望值；
步骤七，更新模型参数的初始模型为后验期望值，重复步骤四至步骤六，直至模拟与观测弹性参数的误差降低到预设范围内，迭代停止，输出模型参数的后验期望值，即为孔隙参数与储层参数的反演结果。
2.根据权利要求1所述的孔隙参数与储层参数的地震岩石物理解析反演方法，其特征在于，步骤一具体包括：
获取工区PP波地震道集数据，通过叠前地震反演，构建目标函数
Figure FDA0003846922710000011
其中，d为地震道集数据，m为弹性参数，
Figure FDA0003846922710000012
为弹性参数的平均值，Σm为弹性参数的协方差矩阵，G为Zoeppritz地震正演模型，σd为地震数据噪声的标准差，T为转置符号；弹性参数m包括纵波速度、横波速度和密度；
采用贝叶斯线性反演方法，得到满足式(1)目标函数的最优解，获取纵波速度、横波速度和密度参数模型。
3.根据权利要求2所述的孔隙参数与储层参数的地震岩石物理解析反演方法，其特征在于，所述步骤二具体包括：
设置基本岩石物理参数，包括岩石基质的体积模量、岩石基质的剪切模量、岩石基质的密度、烃类的体积模量、烃类的密度、卤水的体积模量和卤水的密度；
基于测井观测的纵波速度与横波速度，通过匹配模型预测与井旁观测数据，反演得到井旁孔隙纵横比，构建目标函数
Figure FDA0003846922710000013
其中，α为孔隙纵横比，σ1与σ2为误差项权系数，
Figure FDA0003846922710000014
和
Figure FDA0003846922710000015
分别为模型预测和井旁观测的纵波速度，
Figure FDA0003846922710000021
和
Figure FDA0003846922710000022
分别为模型预测和井旁观测的横波速度，T为转置符号；
采用粒子群优化算法，得到满足式(2)目标函数的最优解，获取井旁孔隙纵横比参数。
4.根据权利要求3所述的孔隙参数与储层参数的地震岩石物理解析反演方法，其特征在于，步骤三具体包括：
通过测井观测数据获得井旁孔隙度与井旁含水饱和度，结合反演得到的井旁孔隙纵横比，构建模型参数的混合高斯模型，其概率分布函数为；
Figure FDA0003846922710000023
其中，z为模型参数，Nk表示第k个高斯分布，
Figure FDA0003846922710000024
和λk分别为Nk的期望、协方差和权系数，C为高斯分布个数；模型参数包括孔隙纵横比、孔隙度与含水饱和度；
通过期望最大算法计算混合高斯模型，得到C个高斯分布的期望、协方差和权系数，再通过测井插值，得到模型参数的初始模型。
5.根据权利要求4所述的孔隙参数与储层参数的地震岩石物理解析反演方法，其特征在于，步骤四具体包括：
利用解耦微分等效介质模型计算岩石骨架的弹性模量
Kd(φ)＝Km(1-φ)p (4)
μd(φ)＝μm(1-φ)q (5)
其中，φ为孔隙度，Kd和μd分别为岩石骨架的体积模量和剪切模量，Km和μm分别为岩石基质的体积模量和剪切模量，p和q为孔隙结构系数；
基于岩石骨架的弹性模量，利用Gassmann方程计算含流体岩石的弹性模量
Figure FDA0003846922710000025
μs＝μd (7)
其中，Ks和μs分别为含流体岩石的体积模量和剪切模量，Kf为混合流体的体积模量，采用Wood模型计算
Figure FDA0003846922710000026
其中，Kw和Kh分别为卤水和烃类的体积模量，Sw为含水饱和度；
根据含流体岩石的弹性模量，得到模拟弹性参数，
Figure FDA0003846922710000031
Figure FDA0003846922710000032
ρ＝φρf+(1-φ)ρm (11)
其中，VP、VS和ρ分别为纵波速度、横波速度和密度，ρf和ρm分为流体和岩石基质的密度；计算模拟弹性参数与观测弹性参数的误差ε
ε＝||mest-mobs||2 (12)
其中，mest为模拟弹性参数，mobs为观测弹性参数即由步骤一获取的观测弹性参数模型。
6.根据权利要求5所述的孔隙参数与储层参数的地震岩石物理解析反演方法，其特征在于，步骤五具体包括：
基于模型参数的初始模型，计算精确地震岩石物理正演算子对模型参数的偏导矩阵
Figure FDA0003846922710000033
其中，R为精确正演算子，通过联立公式(4)–(11)获得，z为模型参数，zo为模型参数的初始模型，G为精确正演算子的偏导矩阵，模型参数包括孔隙纵横比、孔隙度与含水饱和度，其具体计算方法如下：
纵波速度对孔隙度的偏导数为：
Figure FDA0003846922710000034
其中，
Figure FDA0003846922710000035
Figure FDA0003846922710000041
纵波速度对含水饱和度的偏导数为：
Figure FDA0003846922710000042
其中，
Figure FDA0003846922710000043
纵波速度对孔隙纵横比的偏导数为：
Figure FDA0003846922710000044
其中，
Figure FDA0003846922710000045
Figure FDA0003846922710000046
其中，
Figure FDA0003846922710000047
Figure FDA0003846922710000048
Figure FDA0003846922710000049
横波速度对孔隙度的偏导数为：
Figure FDA00038469227100000410
横波速度对含水饱和度的偏导数为：
Figure FDA0003846922710000051
横波速度对孔隙纵横比的偏导数为：
Figure FDA0003846922710000052
其中，
Figure FDA0003846922710000053
其中，
Figure FDA0003846922710000054
Figure FDA0003846922710000055
Figure FDA0003846922710000056
Figure FDA0003846922710000057
Figure FDA0003846922710000058
Figure FDA0003846922710000059
Figure FDA00038469227100000510
Figure FDA00038469227100000511
Figure FDA00038469227100000512
Figure FDA0003846922710000061
密度对孔隙度的偏导数为：
Figure FDA0003846922710000062
密度对含水饱和度的偏导数为：
Figure FDA0003846922710000063
密度对孔隙纵横比的偏导数为：
Figure FDA0003846922710000064
在公式(14)–(41)中，VP、VS和ρ分别为纵波速度、横波速度和密度，Kd和μd分别为岩石骨架的体积模量和剪切模量，Ks和μs分别为含流体岩石的体积模量和剪切模量，Km和μm分别为岩石基质的体积模量和剪切模量，Kf为混合流体的体积模量，p和q为孔隙结构系数，ρf为流体的密度，ρm为岩石基质的密度，ρw为卤水的密度，ρh为烃类的密度，φ为孔隙度，α为孔隙纵横比，Sw为含水饱和度，F1、F2、F3、F4、F5、F6、F7、F8与F9为9个不同的孔隙因子系数，a和b为孔隙扁度系数，R为基质模量系数；
基于精确正演算子的偏导矩阵G，构建线性正演算子
m＝Gz+R(zo)-Gzo+e (42)
其中，m为弹性参数，e为随机误差。
7.根据权利要求6所述的孔隙参数与储层参数的地震岩石物理解析反演方法，其特征在于，步骤六具体包括：
利用公式(42)线性正演算子，基于贝叶斯线性反演理论，解析计算每个高斯分布中模型参数的期望值
Figure FDA0003846922710000065
其中，G为精确正演算子的偏导矩阵，m为弹性参数，Σe为随机误差的协方差，
Figure FDA0003846922710000066
和
Figure FDA0003846922710000067
为第k个高斯分布的期望和协方差，
Figure FDA0003846922710000068
为第k个高斯分布的模型参数期望值，T为转置符号；
计算模型参数的后验期望值
Figure FDA0003846922710000071
其中，Nk表示第k个高斯分布，C为高斯分布个数，λk为Nk的权系数，μpost为模型参数的后验期望值。
8.根据权利要求1所述的孔隙参数与储层参数的地震岩石物理解析反演方法，其特征在于，步骤七预设范围为观测弹性参数值的百分之一。
CN202211122138.XA 2022-09-15 2022-09-15 一种孔隙参数与储层参数的地震岩石物理解析反演方法 Active CN115586572B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202211122138.XA CN115586572B (zh) 2022-09-15 2022-09-15 一种孔隙参数与储层参数的地震岩石物理解析反演方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202211122138.XA CN115586572B (zh) 2022-09-15 2022-09-15 一种孔隙参数与储层参数的地震岩石物理解析反演方法

Publications (2)

Publication Number Publication Date
CN115586572A true CN115586572A (zh) 2023-01-10
CN115586572B CN115586572B (zh) 2023-04-07

Family

ID=84778689

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202211122138.XA Active CN115586572B (zh) 2022-09-15 2022-09-15 一种孔隙参数与储层参数的地震岩石物理解析反演方法

Country Status (1)

Country Link
CN (1) CN115586572B (zh)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116299696A (zh) * 2023-02-08 2023-06-23 中海石油(中国)有限公司深圳分公司 一种烃源岩泥质含量、toc、孔隙度、含水饱和度同时定量预测方法

Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103645509A (zh) * 2013-11-08 2014-03-19 中国石油大学(北京) 致密储层孔隙纵横比反演及横波速度预测方法
CN104516017A (zh) * 2013-09-29 2015-04-15 中国石油化工股份有限公司 一种碳酸盐岩岩石物理参数地震反演方法
CN110275202A (zh) * 2019-05-30 2019-09-24 河海大学 一种用于致密油储层脆性预测的方法
CN112965103A (zh) * 2021-02-09 2021-06-15 中国石油大学(华东) 一种多重孔隙储层叠前地震概率化多道反演方法
CN114428317A (zh) * 2020-10-15 2022-05-03 中国石油化工股份有限公司 碳酸盐岩地震岩石物理建模方法、装置、电子设备及介质
US20220236438A1 (en) * 2021-01-23 2022-07-28 Manzar Fawad Rock physics model for fluid identification and saturation estimation in subsurface reservoirs
CN114924314A (zh) * 2022-06-24 2022-08-19 中国石油大学(北京) 地震剪切横波反演模型构建方法、孔隙度预测方法及装置

Patent Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN104516017A (zh) * 2013-09-29 2015-04-15 中国石油化工股份有限公司 一种碳酸盐岩岩石物理参数地震反演方法
CN103645509A (zh) * 2013-11-08 2014-03-19 中国石油大学(北京) 致密储层孔隙纵横比反演及横波速度预测方法
CN110275202A (zh) * 2019-05-30 2019-09-24 河海大学 一种用于致密油储层脆性预测的方法
CN114428317A (zh) * 2020-10-15 2022-05-03 中国石油化工股份有限公司 碳酸盐岩地震岩石物理建模方法、装置、电子设备及介质
US20220236438A1 (en) * 2021-01-23 2022-07-28 Manzar Fawad Rock physics model for fluid identification and saturation estimation in subsurface reservoirs
CN112965103A (zh) * 2021-02-09 2021-06-15 中国石油大学(华东) 一种多重孔隙储层叠前地震概率化多道反演方法
CN114924314A (zh) * 2022-06-24 2022-08-19 中国石油大学(北京) 地震剪切横波反演模型构建方法、孔隙度预测方法及装置

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116299696A (zh) * 2023-02-08 2023-06-23 中海石油(中国)有限公司深圳分公司 一种烃源岩泥质含量、toc、孔隙度、含水饱和度同时定量预测方法

Also Published As

Publication number Publication date
CN115586572B (zh) 2023-04-07

Similar Documents

Publication Publication Date Title
CN104950331B (zh) 一种砂泥岩储层的孔隙度与泥质含量的地震预测方法
EP1820137B1 (en) Integrated anisotropic rock physics model
Vernik et al. Modeling elastic properties of siliciclastic rocks
RU2602409C2 (ru) Инверсия анизотропии многоскважинной системы
Li et al. Facies-constrained prestack seismic probabilistic inversion driven by rock physics
Stephen et al. Improved normalization of time‐lapse seismic data using normalized root mean square repeatability data to improve automatic production and seismic history matching in the Nelson field
CN111897011A (zh) 储层孔隙特征确定方法、装置及设备
Coléou et al. Petrophysical seismic inversion
CN112363226A (zh) 一种非常规油气有利区地球物理预测方法
CN115586572B (zh) 一种孔隙参数与储层参数的地震岩石物理解析反演方法
CN107831540A (zh) 储层物性参数直接提取新方法
Khadem et al. Integration of rock physics and seismic inversion for rock typing and flow unit analysis: A case study
CN110646850B (zh) 隔夹层地震预测方法及装置
Eikrem et al. Bayesian estimation of reservoir properties—effects of uncertainty quantification of 4D seismic data
CN111077578B (zh) 岩层分布预测方法和装置
Li et al. Integrated VTI model building with seismic data, geologic information, and rock-physics modeling—Part 2: Field data test
He et al. Analysis of reservoir heterogeneity‐induced amplification effect on time‐lapse seismic responses of fluid substitution: a physical modelling study
Abid et al. Unconventional reservoir characterization and sensitive attributes determination: A case study of the Eastern Sembar Formation, Lower Indus Basin, Pakistan
Bayuk Why anisotropy is important for location of microearthquake events in shale?
Wang et al. Recent advances in theory and technology of oil and gas geophysics.
CN107764697A (zh) 基于孔隙介质渐进方程非线性反演的含气性检测方法
Chi et al. Reservoir properties inversion from AVO attributes
Zhang et al. Application of pre-stack simultaneous inversion on the reservoir identification and fluid prediction: A case history of Es2 member reservoir in Shengli Oilfield
Chacon Effect of pressure depletion on hydrocarbon recovery in naturally fractured reservoirs
Kadhim Using Well logs to Predict Petrophysical Properties in Selective Wells

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant