CN115993665A - 基于多尺度岩石物理模型的井中水合物饱和度计算方法 - Google Patents

基于多尺度岩石物理模型的井中水合物饱和度计算方法 Download PDF

Info

Publication number
CN115993665A
CN115993665A CN202310207349.1A CN202310207349A CN115993665A CN 115993665 A CN115993665 A CN 115993665A CN 202310207349 A CN202310207349 A CN 202310207349A CN 115993665 A CN115993665 A CN 115993665A
Authority
CN
China
Prior art keywords
hydrate
modulus
bulk modulus
shear modulus
model
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.)
Pending
Application number
CN202310207349.1A
Other languages
English (en)
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.)
China University of Petroleum East China
Original Assignee
China University of Petroleum East China
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 China University of Petroleum East China filed Critical China University of Petroleum East China
Priority to CN202310207349.1A priority Critical patent/CN115993665A/zh
Publication of CN115993665A publication Critical patent/CN115993665A/zh
Pending legal-status Critical Current

Links

Images

Abstract

本发明属于应用地球物理测井领域,涉及一种基于多尺度岩石物理模型的井中水合物饱和度计算方法,考虑了多个尺度的衰减机制,构建多尺度岩石物理模型;根据工区测井资料,计算固体相的体积模量和剪切模量,将固体相的体积模量和剪切模量代入多尺度岩石物理模型中,得到工区储层的纵波速度和衰减系数,与声波测井计算得出的纵波速度和衰减系数比较,确定工区储层的水合物饱和度。本发明构建的多尺度岩石物理模型能够精细刻画含水合物沉积物的纵波速度和衰减特性,计算的速度和衰减精确,并基于声学参数和多尺度岩石物理模型计算饱和度,为含水合物沉积物的饱和度定量解释提供了理论依据。

Description

基于多尺度岩石物理模型的井中水合物饱和度计算方法
技术领域
本发明属于应用地球物理测井领域,具体地,涉及一种基于多尺度岩石物理模型的井中水合物饱和度计算方法。
背景技术
自然界中的天然气水合物广泛存在于永久冻土带和陆缘外围的海底沉积物中。因为天然气水合物具有分布广,资源量巨大,密度高等特点,自20世纪60年底至今一直受到许多国家广泛关注,被誉清洁、高效、储量丰富的新型潜在能源。在目前的水合物储层测井解释评价中,弹性波岩石物理模型具有重要地位。但是,当前现有的弹性波岩石物理模型无法精确描述实际水合物储层的衰减机制,影响了水合物饱和度的计算,给含水合物沉积物的地球物理定量分析带来了巨大的挑战。
波动诱导流体流动(以下简称:WIFF)是含水合物沉积物中弹性波频散和衰减的主要原因,介观和微观不均匀性是引起WIFF的主要机制。在含水合物沉积物中,介观和微观非均匀性同时存在,并可导致纵波速度的明显转变,这就意味着有必要同时考虑这两种机制对频散和衰减的影响。尽管现有模型(例如:BISQ模型)在模拟含水合物沉积物的速度和衰减特征方面取得了一定的应用效果,但这些模型还不能同时描述微观、介观和宏观三个尺度的衰减机制,导致无法精确刻画含水合物沉积物的纵波速度和衰减特性,严重影响了井中水合物饱和度的计算。
发明内容
本发明针对现有技术不能精确描述含水合物沉积物的纵波速度和衰减问题,提供了一种基于多尺度岩石物理模型的井中水合物饱和度计算方法,该方法构建了描述实际天然气水合物储层的纵波速度和衰减多尺度岩石物理模型,形成了对应的建模方法,并在此基础上基于声学参数计算饱和度。该方法构建的多尺度岩石物理模型能够精确地描述不同赋存形态的含水合物沉积物的纵波速度和衰减特性,为含水合物沉积物的饱和度定量解释提供了理论依据。
为了达到上述目的,本发明提供了一种基于多尺度岩石物理模型的井中水合物饱和度计算方法,其步骤为:
S1、建立多尺度岩石物理模型,其具体步骤为:
S11、将喷射流模型加入到Kuster-
Figure BDA0004111389350000024
模型中,计算含气/水包裹体的水合物等效体积模量和水合物等效剪切模量;
S12、将孔隙水和甲烷气体的体积模量以及步骤S11得到的水合物的等效体积模量代入wood公式,得到第一流体相的体积模量表示为:
Figure BDA0004111389350000021
式中,Kf1为第一流体相的体积模量,Kw为孔隙水的体积模块,Kg为甲烷气体的体积模量,
Figure BDA0004111389350000022
为水合物的等效体积模量,Sw为孔隙水的饱和度,Sg为甲烷气体的饱和度,SH为水合物的饱和度;
将孔隙水和甲烷气体的体积模量代入wood公式,得到第二流体相的体积模量表示为:
Figure BDA0004111389350000023
将石英颗粒的体积模量和剪切模量以及步骤S11得到的水合物的等效体积模量和等效剪切模量代入Hill平均方程,得到固体相的体积模量和剪切模量表示为:
Figure BDA0004111389350000031
Figure BDA0004111389350000032
式中,K为固体相的体积模量,G为固体相的剪切模量;Kq为石英颗粒的体积模量,Gq为石英颗粒的剪切模量,fH为在固体相中水合物的体积分数,fq为在固体相中水合物的体积分数;
S13、将步骤S11得到的水合物的等效体积模量和等效剪切模量代入胶结砂岩模型中,得到接触胶结水合物的干岩石骨架体积模量和干岩石骨架剪切模量;将石英颗粒的体积模量和剪切模量代入胶结砂岩模型中,得到颗粒涂层水合物的干岩石骨架体积模量和干岩石骨架剪切模量;将步骤S12得到的固体相的体积模量和剪切模量代入Hashin-Shtrikman模型,得到颗粒支撑水合物的干岩石骨架体积模量和干岩石骨架剪切模量;将石英颗粒的体积模量和剪切模量代入Hashin-Shtrikman模型中,得到孔隙填充水合物的干岩石骨架体积模量和干岩石骨架剪切模量;
S14、联合步骤S12中得到的流体相的体积模量和剪切模量以及步骤S13中得到的干岩石骨架的体积模量和剪切模量,根据Biot-Rayleigh理论模型,得到接触胶结水合物、颗粒涂层水合物、颗粒支撑水合物、孔隙填充水合物四种不同赋存形态的含水合物沉积物的纵波波速和衰减系数,即多尺度岩石物理模型;
S2、根据工区目标层位的自然伽马测井、孔隙度测井和岩性密度测井数据资料得到岩石成分的体积分数,根据Hill平均方程计算得到固体相的体积模量和剪切模量,根据工区的电阻率测井数据资料,得到水合物的饱和度,然后将固体相的模量与水合物的饱和度代入到S1步骤中构建的多尺度岩石物理模型,最后得到工区储层的纵波速度和衰减系数,与实测声波测井得出的纵波速度和衰减系数比较,更新水合物饱和度,直到预测和测量的纵波速度和衰减在设定的误差范围之内,得到预测的工区储层的水合物饱和度。
优选的,步骤S11中,所述喷射流模型表示为:
K′i=K1i+iωKHτi
G′i=iωμi
式中,K′i为包裹体的体积模量,G′i为包裹体的剪切模量,K1i为水合物包裹体和孔隙中的流体等效体积模量,下标i=1,2,3,4,...,N,N为不同种类包裹体;KH为水合物的体积模量,下标H为水合物,τi为弛豫时间,ω为角频率,μi为水合物中包裹体的粘度;
将喷射流模型加入到Kuster-
Figure BDA0004111389350000047
模型中,计算得到水合物的等效体积模量和等效剪切模量表示为:
Figure BDA0004111389350000041
Figure BDA0004111389350000042
Figure BDA0004111389350000043
式中,
Figure BDA0004111389350000044
为水合物的等效体积模量,
Figure BDA0004111389350000045
为水合物的等效剪切模量,GH为水合物的剪切模量,xi为第i中包裹体的体积分量;
Figure BDA0004111389350000046
为几何因子,
Figure BDA0004111389350000051
α为水和甲烷气包裹体的半径。
优选的,步骤S13中,所述胶结砂岩模型表示为:
Figure BDA0004111389350000052
Figure BDA0004111389350000053
式中,Kdry为干岩石的体积模量,Gdry为干岩石的剪切模量,Mc为胶结物或水合物的体积模量,Gc为石英颗粒或水合物的剪切模量;参数
Figure BDA0004111389350000054
为与胶结的两颗粒组合体的正向成比例,
Figure BDA0004111389350000055
与胶结的两颗粒组合体的剪切刚度成比例,参数
Figure BDA0004111389350000056
Figure BDA0004111389350000057
取决于胶结物的含量及胶结物和骨架颗粒的特性;n为颗粒配位数,在含水合物沉积物中取n=8.5;φ0为沉积物的初始孔隙度;
当孔隙度小于临界孔隙度时,所述Hashin-Shtrikman模型表示为:
Figure BDA0004111389350000058
Figure BDA0004111389350000059
当孔隙度大于等于临界孔隙度时,所述Hashin-Shtrikman模型表示为:
Figure BDA00041113893500000510
Figure BDA00041113893500000511
式中,φ是孔隙度;φc是临界孔隙度,
Figure BDA00041113893500000512
是临界孔隙度下的体积模量,
Figure BDA0004111389350000061
是临界孔隙度下的剪切模量,
Figure BDA0004111389350000062
是从K′和G'计算的矿物相的泊松比,K'是固体相或石英颗粒的体积模量,G'是固体相或石英颗粒的剪切模量,P是等效压力,
Figure BDA0004111389350000063
优选的,步骤S14中:
将步骤S12得到的第二流体相的体积模量和剪切模量以及步骤S13得到的接触胶结水合物的干岩石骨架体积模量和干岩石骨架剪切模量代入Biot-Rayleigh双孔理论模型,得到接触胶结水合物沉积物的纵波波速和衰减系数;
将步骤S12得到的第二流体相的体积模量和剪切模量以及步骤S13得到的颗粒涂层水合物的干岩石骨架体积模量和干岩石骨架剪切模量代入Biot-Rayleigh双孔理论模型,得到颗粒涂层水合物沉积物的纵波波速和衰减系数;
将步骤S12得到的第二流体相的体积模量和剪切模量以及步骤S13得到的颗粒支撑水合物的干岩石骨架体积模量和干岩石骨架剪切模量代入Biot-Rayleigh双孔理论模型,得到颗粒支撑水合物沉积物的纵波波速和衰减系数;
将步骤S12得到的第一流体相的体积模量和剪切模量以及步骤S13得到的孔隙填充水合物的干岩石骨架体积模量和干岩石骨架剪切模量代入Biot-Rayleigh双孔理论模型,得到孔隙填充水合物沉积物的纵波波速和衰减系数。
优选的,所述Biot-Rayleigh双孔理论模型的波动方程为:
Figure BDA0004111389350000064
Figure BDA0004111389350000071
Figure BDA0004111389350000072
Figure BDA0004111389350000073
式中,N、A、Qi,、Ri,i=1,2是刚性系数,u是固体平均粒子位移,
Figure BDA0004111389350000074
分别是u的一阶导数和二阶导数,ε、ζ(1)、ζ(2)表示固体、孔1流体与孔2流体的位移场散度;ρij,i=1,2,3,j=1,2,3是密度参数,ρf是主相流体的密度,φ1、φ2分别是水饱和孔隙和水合物饱和孔隙的孔隙度,φ10、φ20是两个区域的局部孔隙度,bi,i=1,2是耗散参数,
Figure BDA0004111389350000075
是局部流动过程引起的流体应变增量,
Figure BDA0004111389350000076
分别是
Figure BDA0004111389350000077
的一阶导数和二阶导数,η1是主相流体的粘度,κ1是主相流体的渗透率;若水合物作为孔隙流体,
Figure BDA0004111389350000078
是第一流相体为地层水时的位移,
Figure BDA0004111389350000079
是第二流体相为天然气水合物/游离气的位移;若水合物作为岩石骨架,
Figure BDA00041113893500000710
是第一流相体为在宿主骨架中的流体时的位移,
Figure BDA00041113893500000711
是第二流体相为在水合物包裹体中的流体时的位移;
根据平面波分析法求解上式,带入平面纵波的解析解,得到:
Figure BDA00041113893500000712
其中,k表示接触胶结水合物或颗粒涂层水合物或颗粒支撑水合物或孔隙填充水合物的波数;aij、bij表示方程的系数矩阵,i=1,2,3,j=1,2,3;
则速度和衰减的计算公式表示为:
Figure BDA0004111389350000081
Figure BDA0004111389350000082
其中,Q-1为衰减系数;v为含水合物沉积物的纵波波速;Re(k)和Im(k)分别取k的实部和虚部。
与现有技术相比,本发明的有益效果在于:
(1)本发明基于多尺度岩石物理模型的井中水合物饱和度计算方法,基于喷射流模型根据Kuster-
Figure BDA0004111389350000083
模型计算水合物等效弹性模量;根据wood公式和Hill平均方程分别计算流体相和固体相的弹性模量;利用胶结砂岩模型计算接触胶结水合物和颗粒涂层水合物的干岩石骨架弹性模量,利用Hashin-Shtrikman模型计算孔隙填充水合物和颗粒支撑水合物的干岩石骨架弹性模量;联合流体相和干岩石骨架的弹性模量,根据Biot-Rayleigh理论模型,得到4种赋存形态的含水合物沉积物的纵波速度和衰减系数,即多尺度演岩石物理模型。本发明构建的多尺度演岩石物理模型考虑了多个尺度的衰减机制,相比现有的理论模型,更加符合含水合物沉积物的实际情况,可以刻画含水合物沉积物储层的声学响应规律,进而利用声波测井进行含水合物沉积物的探测和识别。
(2)本发明提供的基于多尺度岩石物理模型的井中水合物饱和度计算方法考虑了四种不同的水合物赋存形态,与实际含水合物沉积物储层赋存形态一致,能够精确描述不同赋存形态的含水合物沉积物的纵波速度和衰减系数,为含水合物沉积物的定量解释评价提供了理论依据。
(3)本发明提供的基于多尺度岩石物理模型的井中水合物饱和度计算方法构建的多尺度岩石物理模型,由于能够描述不同赋存形态的含水合物沉积物的纵波速度和衰减系数,将构建的多尺度演示物理模型应用于水合物探测,进而提供一种利用声学参数(即纵波速度和衰减系数)进行水合物探测的新方法,可以为声波测井的改进和完善提供测量方法。
(4)本发明提供的基于多尺度岩石物理模型的井中水合物饱和度计算方法,利用了本发明构建的多尺度岩石物理模型,利用纵波速度和衰减联合确定水合物饱和度,相较于利用阿尔奇公式或变形的电阻率饱和度计算方法,本发明提供的这种新的基于声学参数的饱和度计算方法,避免了阿尔奇公式在含泥质地层或弱成岩地层的不适用问题,极大的扩展了含水合物储层的饱和度计算方法,一定程度上解决了弱成岩地层含水合物饱和度的评价解释。
附图说明
图1为本发明实施例所述基于多尺度岩石物理模型的井中水合物饱和度计算方法的流程图;
图2a为孔隙填充水合物沉积物的纵波速度的频散曲线示意图;
图2b为孔隙填充水合物沉积物的衰减系数的频散曲线示意图;
图3a为不同水合物赋存形态的纵波速度随水合物饱和度的变化曲线图;
图3b为不同水合物赋存形态的衰减系数随水合物饱和度的变化曲线图。
具体实施方式
下面,通过示例性的实施方式对本发明进行具体描述。然而应当理解,在没有进一步叙述的情况下,一个实施方式中的元件、结构和特征也可以有益地结合到其他实施方式中。
参见图1,本发明实施例提供了一种基于多尺度岩石物理模型的井中水合物饱和度计算方法,其具体步骤为:
S1、建立多尺度岩石物理模型,其具体步骤为:
S11、将喷射流模型加入到Kuster-
Figure BDA0004111389350000108
模型中,计算含气/水包裹体的水合物等效体积模量和水合物等效剪切模量。
具体地,所述喷射流模型表示为:
K′i=K1i+iωKHτi
G′i=iωμi
式中,K′i为包裹体的体积模量,G′i为包裹体的剪切模量,K1i为水合物包裹体和孔隙中的流体等效体积模量,下标i=1,2,3,4,...,N,N为不同种类包裹体;KH为水合物的体积模量,下标H为水合物,τi为弛豫时间,ω为角频率,μi为水合物中包裹体的粘度;
将喷射流模型加入到Kuster-
Figure BDA0004111389350000109
模型中,计算得到水合物的等效体积模量和等效剪切模量表示为:
Figure BDA0004111389350000101
Figure BDA0004111389350000102
Figure BDA0004111389350000103
式中,
Figure BDA0004111389350000104
为水合物的等效体积模量,
Figure BDA0004111389350000105
为水合物的等效剪切模量,GH为水合物的剪切模量,xi为第i中包裹体的体积分量;
Figure BDA0004111389350000106
为几何因子,
Figure BDA0004111389350000107
α为水和甲烷气包裹体的半径。
S12、将孔隙水和甲烷气体的体积模量以及步骤S11得到的水合物的等效体积模量代入wood公式,得到第一流体相的体积模量表示为:
Figure BDA0004111389350000111
式中,Kf1为第一流体相的体积模量,Kw为孔隙水的体积模块,Kg为甲烷气体的体积模量,
Figure BDA0004111389350000112
为水合物的等效体积模量,Sw为孔隙水的饱和度,Sg为甲烷气体的饱和度,SH为水合物的饱和度;
将孔隙水和甲烷气体的体积模量代入wood公式,得到第二流体相的体积模量表示为:
Figure BDA0004111389350000113
将石英颗粒的体积模量和剪切模量以及步骤S11得到的水合物的等效体积模量和等效剪切模量代入Hill平均方程,得到固体相的体积模量和剪切模量表示为:
Figure BDA0004111389350000114
Figure BDA0004111389350000115
式中,K为固体相的体积模量,G为固体相的剪切模量;Kq为石英颗粒的体积模量,Gq为石英颗粒的剪切模量,fH为在固体相中水合物的体积分数,fq为在固体相中水合物的体积分数。
S13、将步骤S11得到的水合物的等效体积模量和等效剪切模量代入胶结砂岩模型中,得到接触胶结水合物的干岩石骨架体积模量和干岩石骨架剪切模量;将石英颗粒的体积模量和剪切模量代入胶结砂岩模型中,得到颗粒涂层水合物的干岩石骨架体积模量和干岩石骨架剪切模量;将步骤S12得到的固体相的体积模量和剪切模量代入Hashin-Shtrikman模型,得到颗粒支撑水合物的干岩石骨架体积模量和干岩石骨架剪切模量;将石英颗粒的体积模量和剪切模量代入Hashin-Shtrikman模型中,得到孔隙填充水合物的干岩石骨架体积模量和干岩石骨架剪切模量。
需要说明的是,当水合物是孔隙流体的一部分时,将其添加到流体相中形成孔隙填充水合物。当水合物为固体颗粒时,将其添加到固体基质中形成颗粒支撑水合物。以水合物为胶结物,当胶结方式为全部水合物沉淀在颗粒接触处时形成接触胶结水合物,当胶结方式为水合物均匀地沉淀在颗粒表面时形成颗粒涂层水合物。
两种胶结方式与参数α′有关。当胶结方式为全部水合物沉淀在颗粒接触处时,
Figure BDA0004111389350000121
当胶结方式为水合物均匀地沉淀在颗粒表面时,
Figure BDA0004111389350000122
其中,S为胶结物(水合物)在孔隙空间的饱和度,n为颗粒配位数,在含水合物沉积物中取n=8.5;φ0为沉积物的初始孔隙度。
具体地,所述胶结砂岩模型表示为:
Figure BDA0004111389350000123
Figure BDA0004111389350000124
式中,Kdry为干岩石的体积模量,Gdry为干岩石的剪切模量,Mc为胶结物或水合物的体积模量,Gc为石英颗粒或水合物的剪切模量;参数
Figure BDA0004111389350000125
为与胶结的两颗粒组合体的正向成比例,
Figure BDA0004111389350000126
与胶结的两颗粒组合体的剪切刚度成比例,参数
Figure BDA0004111389350000127
Figure BDA0004111389350000128
取决于胶结物的含量及胶结物和骨架颗粒的特性;n为颗粒配位数,在含水合物沉积物中取n=8.5;φ0为沉积物的初始孔隙度。
需要说明的是,将步骤S11得到的水合物的等效体积模量和等效剪切模量代入胶结砂岩模型中,得到接触胶结水合物的干岩石骨架体积模量和干岩石骨架剪切模量时,Mc为水合物的等效体积模量
Figure BDA0004111389350000139
Gc为水合物的等效剪切模量
Figure BDA0004111389350000132
将石英颗粒的体积模量和剪切模量代入胶结砂岩模型中,得到颗粒涂层水合物的干岩石骨架体积模量和干岩石骨架剪切模量,Mc为石英颗粒的体积模量Kq,Gc为石英颗粒的剪切模量Gq
具体地,当孔隙度小于临界孔隙度时,所述Hashin-Shtrikman模型表示为:
Figure BDA0004111389350000133
Figure BDA0004111389350000134
当孔隙度大于等于临界孔隙度时,所述Hashin-Shtrikman模型表示为:
Figure BDA0004111389350000135
Figure BDA0004111389350000136
式中,φ是孔隙度;φc是临界孔隙度,
Figure BDA0004111389350000137
是临界孔隙度下的体积模量,
Figure BDA0004111389350000138
是临界孔隙度下的剪切模量,
Figure BDA0004111389350000141
是从K'和G'计算的矿物相的泊松比,K'是固体相或石英颗粒的体积模量,G'是固体相或石英颗粒的剪切模量,P是等效压力,
Figure BDA0004111389350000142
需要说明的是,将步骤S12得到的固体相的体积模量和剪切模量代入Hashin-Shtrikman模型,得到颗粒支撑水合物的干岩石骨架体积模量和干岩石骨架剪切模量时,Hashin-Shtrikman模型中的参数K'是固体相的体积模量K,G'是固体相的剪切模量G。将石英颗粒的体积模量和剪切模量代入Hashin-Shtrikman模型中,得到孔隙填充水合物的干岩石骨架体积模量和干岩石骨架剪切模量时,Hashin-Shtrikman模型中的参数K'是石英颗粒的体积模量Kq,G'是石英颗粒的剪切模量Gq
S14、联合步骤S12中得到的流体相的体积模量和剪切模量以及步骤S13中得到的干岩石骨架的体积模量和剪切模量,根据Biot-Rayleigh理论模型,得到接触胶结水合物、颗粒涂层水合物、颗粒支撑水合物、孔隙填充水合物四种不同赋存形态的含水合物沉积物的纵波波速和衰减系数,即多尺度岩石物理模型。
具体地,将步骤S12得到的第二流体相的体积模量和剪切模量以及步骤S13得到的接触胶结水合物的干岩石骨架体积模量和干岩石骨架剪切模量代入Biot-Rayleigh双孔理论模型,得到接触胶结水合物沉积物的纵波波速和衰减系数。
具体地,将步骤S12得到的第二流体相的体积模量和剪切模量以及步骤S13得到的颗粒涂层水合物的干岩石骨架体积模量和干岩石骨架剪切模量代入Biot-Rayleigh双孔理论模型,得到颗粒涂层水合物沉积物的纵波波速和衰减系数。
具体地,将步骤S12得到的第二流体相的体积模量和剪切模量以及步骤S13得到的颗粒支撑水合物的干岩石骨架体积模量和干岩石骨架剪切模量代入Biot-Rayleigh双孔理论模型,得到颗粒支撑水合物沉积物的纵波波速和衰减系数。
具体地,将步骤S12得到的第一流体相的体积模量和剪切模量以及步骤S13得到的孔隙填充水合物的干岩石骨架体积模量和干岩石骨架剪切模量代入Biot-Rayleigh双孔理论模型,得到孔隙填充水合物沉积物的纵波波速和衰减系数。
具体地,所述Biot-Rayleigh双孔理论模型的波动方程为:
Figure BDA0004111389350000151
Figure BDA0004111389350000152
Figure BDA0004111389350000153
Figure BDA0004111389350000154
式中,N、A、Qi,、Ri,i=1,2是刚性系数,u是固体平均粒子位移,
Figure BDA0004111389350000155
分别是u的一阶导数和二阶导数,ε、ζ(1)、ζ(2)表示固体、孔1流体与孔2流体的位移场散度;ρij,i=1,2,3,j=1,2,3是密度参数,ρf是主相流体的密度,φ1、φ2分别是水饱和孔隙和水合物饱和孔隙的孔隙度,φ10、φ20是两个区域的局部孔隙度,bi,i=1,2是耗散参数,
Figure BDA0004111389350000156
是局部流动过程引起的流体应变增量,
Figure BDA0004111389350000157
分别是
Figure BDA0004111389350000158
的一阶导数和二阶导数,η1是主相流体的粘度,κ1是主相流体的渗透率;若水合物作为孔隙流体,
Figure BDA0004111389350000159
是第一流相体为地层水时的位移,
Figure BDA00041113893500001510
是第二流体相为天然气水合物/游离气的位移;若水合物作为岩石骨架,
Figure BDA00041113893500001511
是第一流相体为在宿主骨架中的流体时的位移,
Figure BDA00041113893500001512
是第二流体相为在水合物包裹体中的流体时的位移;
根据平面波分析法求解上式,带入平面纵波的解析解,得到:
Figure BDA0004111389350000161
其中,k表示接触胶结水合物或颗粒涂层水合物或颗粒支撑水合物或孔隙填充水合物的波数;aij、bij表示方程的系数矩阵,i=1,2,3,j=1,2,3;
当k表示接触胶结水合物的波数k1,则速度和衰减的计算公式表示为:
Figure BDA0004111389350000162
Figure BDA0004111389350000163
其中,Q1 -1为衰减系数;v1为接触胶结水合物沉积物的纵波波速;Re(k1)和Im(k1)分别取k1的实部和虚部。
当k表示颗粒涂层水合物的波数k2,则速度和衰减的计算公式表示为:
Figure BDA0004111389350000164
Figure BDA0004111389350000165
其中,Q2 -1为衰减系数;v2为颗粒涂层水合物沉积物的纵波波速;Re(k2)和Im(k2)分别取k2的实部和虚部。
当k表示颗粒支撑水合物的波数k3,则速度和衰减的计算公式表示为:
Figure BDA0004111389350000171
Figure BDA0004111389350000172
其中,Q3 -1为衰减系数;v3为颗粒支撑水合物沉积物的纵波波速;Re(k3)和Im(k3)分别取k3的实部和虚部。
当k表示孔隙填充水合物的波数k4,则速度和衰减的计算公式表示为:
Figure BDA0004111389350000173
Figure BDA0004111389350000174
其中,Q4 -1为衰减系数;v4为孔隙填充水合物沉积物的纵波波速;Re(k4)和Im(k4)分别取k4的实部和虚部。
S2、根据工区的伽马射线测井、孔隙度测井和密度测井数据资料得到岩石成分(方解石、伊利石和石英)的体积分数,根据Hill平均方程计算得到固体相的体积模量和剪切模量,根据工区的电阻率测井数据资料,得到水合物的饱和度,然后将固体相的模量与水合物的饱和度代入到S1步骤中构建的多尺度岩石物理模型,最后得到工区储层的纵波速度和衰减系数,与实测声波测井资料得出的纵波速度和衰减系数比较,更新水合物饱和度,直到预测和测量的纵波速度和衰减在设定的误差范围之内,得到预测的工区储层的水合物饱和度。
具体地,更新水合物饱和度的具体方法为:通过多尺度岩石物理模型得到速度和衰减与实测声波测井资料计算得到的速度和衰减系数构建最小二乘目标函数,通过不断调整多尺度岩石物理模型中的各个参数,使得构建的目标函数的最小二乘达到最小的时候对应的含水饱和度即目标层位的水合物饱和度。需要说明的是,有些参数是通过区域地层的实验获得并代入多尺度岩石物理模型,有些是不同变化,如含水合物饱和度。
上述方法中,构建的含水合物沉积物的多尺度岩石物理模型,考虑了多个尺度的衰减机制,相比现有的理论模型,更加符合含水合物沉积物的实际情况,能够更加精确地描述不同赋存形态的含水合物沉积物的纵波速度和衰减系数,计算的速度和衰减更加精确,对于含水合物沉积物的声学探测和识别具有重要的意义和价值,为含水合物沉积物的饱和度定量解释提供了理论依据。将该多尺度岩石物理模型应用于水合物探测,进而提供一种利用声学参数(即纵波速度和衰减系数)进行水合物探测的新方法,可以为声波测井的改进和完善提供测量方法。
为了说明本发明上述方法构建的多尺度岩石物理模型的效果。图2a、图2b、图3a、图3b给出了含水合物沉积物中纵波速度和衰减系数的频散曲线以及4中不同水合物赋存形态的纵波速度和衰减系数随水合物饱和度的变化趋势,并与实验数据进行比较。
图2a、图2b分别显示了纵波速度和纵波衰减系数作为不同包裹体纵横比的频率函数。本发明上述方法构建的多尺度岩石物理模型预测了三个弛豫峰,即局部流(介观尺度)、全局流(宏观尺度)和喷射流(微观尺度)。通过分析图2a、图2b,可以得出以下结论:随着水(甲烷气体)包裹体纵横比的增加,第三个峰向更高的频率,纵横比控制喷射流弛豫时间,进而影响衰减峰的位置。
图3a、图3b分别展示了4中不同水合物赋存形态的纵波速度和衰减系数随水合物饱和度的变化形式。由图3a、图3b可以看出,“过量水”方法主要生成的水合物为孔隙填充水合物,“过量气”方法生成的主要是胶结水合物,并且与颗粒涂层水合物模拟结果吻合度更高,在水合物饱和度为0%-40%以下,本发明上述方法构建的多尺度岩石物理模型能够捕捉测量到纵波速度和衰减系数随水合物饱和度的变化趋势。
上述实施例用来解释本发明,而不是对本发明进行限制,在本发明的精神和权利要求的保护范围内,对本发明做出的任何修改和改变,都落入本发明的保护范围。

Claims (5)

1.一种基于多尺度岩石物理模型的井中水合物饱和度计算方法,其特征在于,其步骤为:
S1、建立多尺度岩石物理模型,其具体步骤为:
S11、将喷射流模型加入到
Figure FDA0004111389320000016
模型中,计算含气/水包裹体的水合物等效体积模量和水合物等效剪切模量;
S12、将孔隙水和甲烷气体的体积模量以及步骤S11得到的水合物的等效体积模量代入wood公式,得到第一流体相的体积模量表示为:
Figure FDA0004111389320000011
式中,Kf1为第一流体相的体积模量,Kw为孔隙水的体积模块,Kg为甲烷气体的体积模量,
Figure FDA0004111389320000012
为水合物的等效体积模量,Sw为孔隙水的饱和度,Sg为甲烷气体的饱和度,SH为水合物的饱和度;
将孔隙水和甲烷气体的体积模量代入wood公式,得到第二流体相的体积模量表示为:
Figure FDA0004111389320000013
将石英颗粒的体积模量和剪切模量以及步骤S11得到的水合物的等效体积模量和等效剪切模量代入Hill平均方程,得到固体相的体积模量和剪切模量表示为:
Figure FDA0004111389320000014
Figure FDA0004111389320000015
式中,K为固体相的体积模量,G为固体相的剪切模量;Kq为石英颗粒的体积模量,Gq为石英颗粒的剪切模量,fH为在固体相中水合物的体积分数,fq为在固体相中水合物的体积分数;
S13、将步骤S11得到的水合物的等效体积模量和等效剪切模量代入胶结砂岩模型中,得到接触胶结水合物的干岩石骨架体积模量和干岩石骨架剪切模量;将石英颗粒的体积模量和剪切模量代入胶结砂岩模型中,得到颗粒涂层水合物的干岩石骨架体积模量和干岩石骨架剪切模量;将步骤S12得到的固体相的体积模量和剪切模量代入Hashin-Shtrikman模型,得到颗粒支撑水合物的干岩石骨架体积模量和干岩石骨架剪切模量;将石英颗粒的体积模量和剪切模量代入Hashin-Shtrikman模型中,得到孔隙填充水合物的干岩石骨架体积模量和干岩石骨架剪切模量;
S14、联合步骤S12中得到的流体相的体积模量和剪切模量以及步骤S13中得到的干岩石骨架的体积模量和剪切模量,根据Biot-Rayleigh理论模型,得到接触胶结水合物、颗粒涂层水合物、颗粒支撑水合物、孔隙填充水合物四种不同赋存形态的含水合物沉积物的纵波波速和衰减系数,即多尺度岩石物理模型;
S2、根据工区目标层位的自然伽马测井、孔隙度测井和岩性密度测井数据资料得到岩石成分的体积分数,根据Hill平均方程计算得到固体相的体积模量和剪切模量,根据工区的电阻率测井数据资料,得到水合物的饱和度,然后将固体相的模量与水合物的饱和度代入到S1步骤中构建的多尺度岩石物理模型,最后得到工区储层的纵波速度和衰减系数,与声波测井得出的纵波速度和衰减系数比较,更新水合物饱和度,直到预测和测量的纵波速度和衰减在设定的误差范围之内,得到预测的工区储层的水合物饱和度。
2.如权利要求1所述的基于多尺度岩石物理模型的井中水合物饱和度计算方法,其特征在于,步骤S11中,所述喷射流模型表示为:
K′i=K1i+iωKHτi
G′i=iωμi
式中,K′i为包裹体的体积模量,G′i为包裹体的剪切模量,K1i为水合物包裹体和孔隙中的流体等效体积模量,下标i=1,2,3,4,...,N,N为不同种类包裹体;KH为水合物的体积模量,下标H为水合物,τi为弛豫时间,ω为角频率,μi为水合物中包裹体的粘度;
将喷射流模型加入到
Figure FDA0004111389320000031
模型中,计算得到水合物的等效体积模量和等效剪切模量表示为:
Figure FDA0004111389320000032
Figure FDA0004111389320000033
Figure FDA0004111389320000034
式中,
Figure FDA0004111389320000035
为水合物的等效体积模量,
Figure FDA0004111389320000036
为水合物的等效剪切模量,GH为水合物的剪切模量,xi为第i中包裹体的体积分量;
Figure FDA0004111389320000037
为几何因子,
Figure FDA0004111389320000038
α为水和甲烷气包裹体的半径。
3.如权利要求1所述的基于多尺度岩石物理模型的井中水合物饱和度计算方法,其特征在于,步骤S13中,所述胶结砂岩模型表示为:
Figure FDA0004111389320000041
Figure FDA0004111389320000042
式中,Kdry为干岩石的体积模量,Gdry为干岩石的剪切模量,Mc为胶结物或水合物的体积模量,Gc为石英颗粒或水合物的剪切模量;参数
Figure FDA0004111389320000043
为与胶结的两颗粒组合体的正向成比例,
Figure FDA0004111389320000044
与胶结的两颗粒组合体的剪切刚度成比例,参数
Figure FDA0004111389320000045
Figure FDA0004111389320000046
取决于胶结物的含量及胶结物和骨架颗粒的特性;n为颗粒配位数,在含水合物沉积物中取n=8.5;φ0为沉积物的初始孔隙度;
当孔隙度小于临界孔隙度时,所述Hashin-Shtrikman模型表示为:
Figure FDA0004111389320000047
Figure FDA0004111389320000048
当孔隙度大于等于临界孔隙度时,所述Hashin-Shtrikman模型表示为:
Figure FDA0004111389320000049
Figure FDA00041113893200000410
式中,φ是孔隙度;φc是临界孔隙度,
Figure FDA00041113893200000411
是临界孔隙度下的体积模量,
Figure FDA00041113893200000412
是临界孔隙度下的剪切模量,
Figure FDA00041113893200000413
是从K'和G'计算的矿物相的泊松比,K'是固体相或石英颗粒的体积模量,G'是固体相或石英颗粒的剪切模量,P是等效压力,
Figure FDA0004111389320000051
4.如权利要求2所述的基于多尺度岩石物理模型的井中水合物饱和度计算方法,其特征在于,步骤S14中:
将步骤S12得到的第二流体相的体积模量和剪切模量以及步骤S13得到的接触胶结水合物的干岩石骨架体积模量和干岩石骨架剪切模量代入Biot-Rayleigh双孔理论模型,得到接触胶结水合物沉积物的纵波波速和衰减系数;
将步骤S12得到的第二流体相的体积模量和剪切模量以及步骤S13得到的颗粒涂层水合物的干岩石骨架体积模量和干岩石骨架剪切模量代入Biot-Rayleigh双孔理论模型,得到颗粒涂层水合物沉积物的纵波波速和衰减系数;
将步骤S12得到的第二流体相的体积模量和剪切模量以及步骤S13得到的颗粒支撑水合物的干岩石骨架体积模量和干岩石骨架剪切模量代入Biot-Rayleigh双孔理论模型,得到颗粒支撑水合物沉积物的纵波波速和衰减系数;
将步骤S12得到的第一流体相的体积模量和剪切模量以及步骤S13得到的孔隙填充水合物的干岩石骨架体积模量和干岩石骨架剪切模量代入Biot-Rayleigh双孔理论模型,得到孔隙填充水合物沉积物的纵波波速和衰减系数。
5.如权利要求4所述的基于多尺度岩石物理模型的井中水合物饱和度计算方法,其特征在于,所述Biot-Rayleigh双孔理论模型的波动方程为:
Figure FDA0004111389320000052
Figure FDA0004111389320000061
Figure FDA0004111389320000062
Figure FDA0004111389320000063
式中,N、A、Qi,、Ri,i=1,2是刚性系数,u是固体平均粒子位移,
Figure FDA0004111389320000064
Figure FDA0004111389320000065
分别是u的一阶导数和二阶导数,ε、ζ(1)、ζ(2)表示固体、孔1流体与孔2流体的位移场散度;ρij,i=1,2,3,j=1,2,3是密度参数,ρf是主相流体的密度,φ1、φ2分别是水饱和孔隙和水合物饱和孔隙的孔隙度,φ10、φ20是两个区域的局部孔隙度,bi,i=1,2是耗散参数,
Figure FDA0004111389320000066
是局部流动过程引起的流体应变增量,
Figure FDA0004111389320000067
分别是
Figure FDA0004111389320000068
的一阶导数和二阶导数,η1是主相流体的粘度,κ1是主相流体的渗透率;若水合物作为孔隙流体,
Figure FDA0004111389320000069
是第一流相体为地层水时的位移,
Figure FDA00041113893200000610
是第二流体相为天然气水合物/游离气的位移;若水合物作为岩石骨架,
Figure FDA00041113893200000611
是第一流相体为在宿主骨架中的流体时的位移,
Figure FDA00041113893200000612
是第二流体相为在水合物包裹体中的流体时的位移;
根据平面波分析法求解上式,带入平面纵波的解析解,得到:
Figure FDA00041113893200000613
其中,k表示接触胶结水合物或颗粒涂层水合物或颗粒支撑水合物或孔隙填充水合物的波数;aij、bij表示方程的系数矩阵,i=1,2,3,j=1,2,3;
则速度和衰减的计算公式表示为:
Figure FDA0004111389320000071
Figure FDA0004111389320000072
其中,Q-1为衰减系数;v为含水合物沉积物的纵波波速;Re(k)和Im(k)分别取k的实部和虚部。
CN202310207349.1A 2023-03-06 2023-03-06 基于多尺度岩石物理模型的井中水合物饱和度计算方法 Pending CN115993665A (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202310207349.1A CN115993665A (zh) 2023-03-06 2023-03-06 基于多尺度岩石物理模型的井中水合物饱和度计算方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202310207349.1A CN115993665A (zh) 2023-03-06 2023-03-06 基于多尺度岩石物理模型的井中水合物饱和度计算方法

Publications (1)

Publication Number Publication Date
CN115993665A true CN115993665A (zh) 2023-04-21

Family

ID=85995161

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202310207349.1A Pending CN115993665A (zh) 2023-03-06 2023-03-06 基于多尺度岩石物理模型的井中水合物饱和度计算方法

Country Status (1)

Country Link
CN (1) CN115993665A (zh)

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN117055114A (zh) * 2023-10-09 2023-11-14 中国石油大学(华东) 储层沉积物游离气饱和度定量分析方法
CN117316329A (zh) * 2023-11-23 2023-12-29 中国石油大学(华东) 天然气水合物饱和度声电测井联合智能反演方法及系统

Cited By (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN117055114A (zh) * 2023-10-09 2023-11-14 中国石油大学(华东) 储层沉积物游离气饱和度定量分析方法
CN117055114B (zh) * 2023-10-09 2023-12-29 中国石油大学(华东) 储层沉积物游离气饱和度定量分析方法
CN117316329A (zh) * 2023-11-23 2023-12-29 中国石油大学(华东) 天然气水合物饱和度声电测井联合智能反演方法及系统
CN117316329B (zh) * 2023-11-23 2024-03-29 中国石油大学(华东) 天然气水合物饱和度声电测井联合智能反演方法及系统

Similar Documents

Publication Publication Date Title
CN115993665A (zh) 基于多尺度岩石物理模型的井中水合物饱和度计算方法
CN103713320B (zh) 一种富有机质泥页岩岩石物理模型的建立方法
CN102508296B (zh) 一种非饱和双重孔隙介质地震波频散衰减分析方法及装置
CN103954999B (zh) 一种适用于低孔隙度砂泥岩地层的横波速度预测方法
CN105089615B (zh) 一种基于油藏模型的测井数据历史回归处理方法
CN106154351A (zh) 一种低孔渗储层渗透率的估算方法
WO2006062612A2 (en) Integrated anisotropic rock physics model
Chen et al. Pore‐scale modeling of hydromechanical coupled mechanics in hydrofracturing process
CN103412336A (zh) 一种非均质油藏中岩石系统的纵波速度预测方法
Van Der Kolk et al. The 3D shear experiment over the Natih field in Oman: the effect of fracture‐filling fluids on shear propagation
Liu et al. Quantitative multiparameter prediction of fractured tight sandstone reservoirs: a case study of the Yanchang Formation of the Ordos Basin, Central China
CN113075728A (zh) 一种建立致密砂岩多尺度三维岩石物理图板的方法
Li et al. Nonlinear simultaneous inversion of pore structure and physical parameters based on elastic impedance
CN106869915B (zh) 一种水平井井间隔夹层预测方法及装置
CN106443772B (zh) 一种去底辟原始地层厚度恢复方法
CN105301642B (zh) 非均匀孔隙岩石及其固态有机质体积含量确定方法及装置
CN114114409A (zh) 海域天然气水合物的岩石物理建模方法、电子设备及介质
Pan et al. A unified effective medium modeling framework for quantitative characterization of hydrate reservoirs
Barros-Galvis et al. Analytical modeling and contradictions in limestone reservoirs: Breccias, vugs, and fractures
CN114675328A (zh) 基于变化孔隙纵横比Xu-white模型的富有机质页岩横波预测方法
Parra et al. Characterization of fractured low Q zones at the Buena Vista Hills reservoir, California
CN113281825A (zh) 岩石物理模型构建方法及装置
Zhang et al. Construction and validation of the upscaling model of organic-rich shale by considering water-sensitivity effects
CN110716235A (zh) 一种砂泥岩测井孔隙结构反演方法
CN117233845B (zh) 多尺度岩石物理模型估算水合物饱和度的方法及系统

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