CN113189674A - 一种天然气水合物的饱和度估算方法及系统 - Google Patents

一种天然气水合物的饱和度估算方法及系统 Download PDF

Info

Publication number
CN113189674A
CN113189674A CN202110493721.0A CN202110493721A CN113189674A CN 113189674 A CN113189674 A CN 113189674A CN 202110493721 A CN202110493721 A CN 202110493721A CN 113189674 A CN113189674 A CN 113189674A
Authority
CN
China
Prior art keywords
hydrate
saturation
measured
depth point
predicted
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
CN202110493721.0A
Other languages
English (en)
Other versions
CN113189674B (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.)
China University of Petroleum Beijing
Original Assignee
China University of Petroleum Beijing
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 Beijing filed Critical China University of Petroleum Beijing
Priority to CN202110493721.0A priority Critical patent/CN113189674B/zh
Publication of CN113189674A publication Critical patent/CN113189674A/zh
Application granted granted Critical
Publication of CN113189674B publication Critical patent/CN113189674B/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
    • G01V11/00—Prospecting or detecting by methods combining techniques covered by two or more of main groups G01V1/00 - G01V9/00
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/11—Complex mathematical operations for solving equations, e.g. nonlinear equations, general mathematical optimization problems

Landscapes

  • Physics & Mathematics (AREA)
  • Engineering & Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Mathematical Physics (AREA)
  • Pure & Applied Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Data Mining & Analysis (AREA)
  • Computational Mathematics (AREA)
  • Mathematical Analysis (AREA)
  • Mathematical Optimization (AREA)
  • Geophysics (AREA)
  • Operations Research (AREA)
  • Algebra (AREA)
  • General Life Sciences & Earth Sciences (AREA)
  • Databases & Information Systems (AREA)
  • Software Systems (AREA)
  • General Engineering & Computer Science (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Geophysics And Detection Of Objects (AREA)

Abstract

本发明提供了一种天然气水合物的饱和度估算方法及系统。方法包括：获取待测地层的测井数据以及弹性参数；分别计算待测地层的预测纵波速度及预测横波速度；当饱水层段每个深度点的预测纵波速度与所述饱水层段实测纵波速度的误差、预测横波速度与饱水层段实测横波速度的误差均最小时确定临界孔隙度及在所述临界孔隙度处的配位数，重新计算待测地层的预测纵波速度及预测横波速度，并根据重新计算的预测纵波速度估算每个深度点的第一水合物饱和度，以及根据重新计算的预测横波速度估算每个深度点的第二水合物饱和度并获得该水合物饱和度。本方法估算精确度高，对水合物储层储量评估及后续勘探开发有积极的指导作用。

Description

一种天然气水合物的饱和度估算方法及系统
技术领域
本发明涉及天然气水合物储层饱和度估算领域，具体涉及一种天然气水合物的饱和度估算方法及系统。
背景技术
天然气水合物，简称水合物，又称“可燃冰”，是天然气和水在高压低温环境条件下形成的冰态、结晶状笼形化合物。作为一种非常规能源，因其巨大的资源量成为目前能源研究的热点。天然气水合物储层的定量评价是水合物开发利用和环境影响的基础，准确的预测天然气水合物的饱和度对于水合物资源定量评价至关重要。
目前在钻井资料存在的地区，应用于天然气水合物储层估算评价饱和度的方法中，电阻率测井和声波测井是评价水合物储层的有效手段。从电学角度出发，利用阿奇公式(Archie，1942)、双水模型、印度尼西亚公式由电阻率估算水合物饱和度；从声学角度出发，利用声波测井资料估算水合物饱和度。由于含水合物沉积地层的速度与水合物饱和度(含量)有关，利用波速评价水合物饱和度的关键是建立水合物含量与沉积层速度的关系，即合适的岩石物理模型。一旦建立了两者之间的数学物理方程，便可由地层速度反演预测水合物饱和度。第一种是比较常用的经验公式，如三相Wood方程、三相时间平均方程、三相加权方程；第二种是依赖于岩石物理理论，并将水合物不同的赋存状态考虑进去，基本是将天然气水合物作为骨架、流体或胶结物加入其中，即对应的三种模型：承载模型、孔隙填充模型和胶结模型，如Lee修正的Biot-Gassmann理论(BGTL)(Lee，2002)、三相Biot方程(TPBE)(Lee，2006)、Ecker胶结模型等，另外还有在孔隙尺度下基于水合物颗粒与多孔隙介质颗粒之间的相互作用而构建岩石物理模型，如等效介质理论(EMT)(Helgerud，1999)，多孔颗粒的未固结砂岩模型(Ruiz和Dvorkin，2009)，K-T方程(Zimmerman和King，1986)等。
等效介质理论(EMT)作为针对水合物储层最为常用的一种岩石物理模型，可以模拟水合物作为骨架(颗粒支撑)一部分或孔隙流体(孔隙填充)一部分这两种赋存状态下的等效弹性模量。然而，颗粒支撑模式和孔隙填充模式是物理意义上的两种极端微观分布状态。实际上，水合物在沉积物中的分布状态多种多样，且多数情况下各种分布状态可以共存，此时EMT方程适用型较差，使用此方程进行水合物储层纵、横波速度模拟与饱和度计算工作时，结果往往不准确。
发明内容
鉴于上述传统等效介质理论所采用的EMT方程进行水合物储层纵、横波速度模拟与饱和度计算结果不准确的问题，提出了本发明以便提供一种克服上述问题或者至少部分地解决上述问题的一种天然气水合物的饱和度估算方法及系统。
依据本发明的一个方面，提供一种天然气水合物的饱和度估算方法，包括：
获取待测地层的测井数据以及弹性参数；所述待测地层包括水合物层段以及饱水层段；所述测井数据包括：深度点的深度、各个深度点的孔隙度、实测纵波速度以及实测横波速度；所述弹性参数包括：体积模量、剪切模量和密度；
根据所述弹性参数分别计算待测地层的预测纵波速度及预测横波速度；
当饱水层段每个深度点的预测纵波速度与所述饱水层段实测纵波速度的误差、预测横波速度与饱水层段实测横波速度的误差均最小时确定临界孔隙度及在所述临界孔隙度处的配位数；
以实测纵波速度以及实测横波速度作为约束，根据确定临界孔隙度及在所述临界孔隙度处的配位数重新计算待测地层的预测纵波速度及预测横波速度，并根据重新计算的预测纵波速度估算每个深度点的第一水合物饱和度，以及根据重新计算的预测横波速度估算每个深度点的第二水合物饱和度；
根据每个深度点的第一水合物饱和度以及第二水合物饱和度获得该深度点的水合物饱和度。
优选的，计算地层的预测纵波速度及预测横波速度前，所述方法还包括：
模拟孔隙填充和颗粒支撑两种分布状态同时存在时的水合物等效体积参数ε；所述水合物等效体积参数ε表征水合物中孔隙填充部分的体积百分含量；其中，
ε∈[0，1]，当ε＝0时，表示水合物完全作为骨架的一部分；当ε＝1时，表示水合物完全作为孔隙流体的一部分；
基于公式一获取视孔隙度
Figure BDA0003053436920000021
所述视孔隙度
Figure BDA0003053436920000022
表征水合物对孔隙度
Figure BDA0003053436920000023
的减小作用，
Figure BDA0003053436920000031
其中，
Figure BDA0003053436920000032
为孔隙水的孔隙度，
Figure BDA0003053436920000033
为纯水合物的孔隙度；
基于公式二、公式三以及公式四分别获得待测地层的体积模量K、剪切模量μ以及密度ρ，
Figure BDA0003053436920000034
μ＝μdry 公式三，
Figure BDA0003053436920000035
其中，Kdry、μdry分别为海洋干骨架的体积模量、剪切模量，Km为含水合物颗粒时的水合物基质的体积模量，Kfl为含水合物颗粒时的孔隙流体的体积模量，ρw为孔隙水的密度，ρma为不含水合物颗粒时的岩石基质的密度，Sh为水合物饱和度，ρh为纯水合物的密度。
优选的，当饱水层段每个深度点的预测纵波速度与所述实测纵波速度的误差、预测横波速度与实测横波速度的误差均最小时确定临界孔隙度及在所述临界孔隙度处的配位数，具体包括：
在所述水合物饱和度Sh＝0时，调整临界孔隙度
Figure BDA0003053436920000036
及在所述临界孔隙度处的配位数n，其中
Figure BDA0003053436920000037
8≤n≤9.5；
基于调整后的
Figure BDA0003053436920000038
及n的值计算第i个深度点的预测纵波速度
Figure BDA0003053436920000039
及预测横波速度
Figure BDA00030534369200000310
当
Figure BDA00030534369200000311
以及
Figure BDA00030534369200000312
满足时，确定第i个深度点此时的
Figure BDA00030534369200000313
及n的值；其中，Vpi water-sat、Vsi water-sat分别为饱水层段第i个深度点的实测纵波速度、实测横波速度。
优选的，所述方法还包括：
调整水合物饱和度Sh在(0，1)之间变化；
基于调整后的水合物饱和度Sh分别重新计算待测地层第i个深度点的预测纵波速度
Figure BDA00030534369200000314
及预测横波速度
Figure BDA00030534369200000315
当
Figure BDA00030534369200000316
时，令第i个深度点的第一水合物饱和度Sh1i等于当前调整的水合物饱和度Sh；
以及当
Figure BDA0003053436920000041
时，令第i个深度点的第二水合物饱和度Sh2i等于当前调整的水合物饱和度Sh；其中，Vpi、Vsi分别为含水合物层段第i个深度点的实测纵波速度、实测横波速度。
优选的，所述方法还包括：
根据公式五获取含水合物层段第i个深度点的水合物饱和度Shi：
Figure BDA0003053436920000042
依据本发明的一个方面，提供一种天然气水合物的饱和度估算系统，包括：
数据获取模块，用于获取待测地层的测井数据以及弹性参数；所述待测地层包括水合物层段以及饱水层段；所述测井数据包括：深度点的深度、各个深度点的孔隙度、实测纵波速度以及实测横波速度；所述弹性参数包括：体积模量、剪切模量和密度；
第一计算模块，用于根据公式一及公式二分别计算待测地层的预测纵波速度及预测横波速度，
参数确定模块，用于当饱水层段每个深度点的预测纵波速度与所述实测纵波速度的误差、预测横波速度与实测横波速度的误差均最小时确定临界孔隙度及在所述临界孔隙度处的配位数；
第二计算模块，用于以实测纵波速度以及实测横波速度作为约束，根据确定的临界孔隙度及在所述临界孔隙度处的配位数重新计算待测地层的预测纵波速度及预测横波速度，并根据重新计算的预测纵波速度估算每个深度点的第一水合物饱和度，以及根据重新计算的预测横波速度估算每个深度点的第二水合物饱和度；
第三计算模块，用于根据每个深度点的第一水合物饱和度以及第二水合物饱和度获得该深度点的水合物饱和度。
优选的，所述系统包括：
模拟模块，用于模拟孔隙填充和颗粒支撑两种分布状态同时存在时的水合物等效体积参数ε；所述水合物等效体积参数ε表征水合物中孔隙填充部分的体积百分含量；其中，ε∈[0，1]，当ε＝0时，表示水合物完全作为骨架的一部分；当ε＝1时，表示水合物完全作为孔隙流体的一部分；
第四计算模块，用于基于公式一获取视孔隙度
Figure BDA0003053436920000043
所述视孔隙度
Figure BDA0003053436920000044
表征水合物对孔隙度
Figure BDA0003053436920000045
的减小作用，
Figure BDA0003053436920000051
其中，
Figure BDA0003053436920000052
为孔隙水的孔隙度，
Figure BDA0003053436920000053
为纯水合物的孔隙度；
第五计算模块，用于基于公式二、公式三以及公式四分别获得待测地层的体积模量K、剪切模量μ以及密度ρ，
Figure BDA0003053436920000054
μ＝μdry 公式三，
Figure BDA0003053436920000055
其中，Kdry、μdry分别为海洋干骨架的体积模量、剪切模量，Km为含水合物颗粒时的水合物基质的体积模量，Kfl为含水合物颗粒时的孔隙流体的体积模量，ρw为孔隙水的密度，ρma为不含水合物颗粒时的岩石基质的密度，Sh为水合物饱和度，ρh为纯水合物的密度。
优选的，所述参数确定模块还包括：
第一调整单元，用于在所述水合物饱和度Sh＝0时，调整临界孔隙度
Figure BDA0003053436920000056
及在所述临界孔隙度处的配位数n，其中
Figure BDA0003053436920000057
8≤n≤9.5；
第一预测单元，用于基于调整后的
Figure BDA0003053436920000058
及n的值计算第i个深度点的预测纵波速度
Figure BDA0003053436920000059
及预测横波速度
Figure BDA00030534369200000510
参数确定单元，用于当
Figure BDA00030534369200000511
以及
Figure BDA00030534369200000512
满足时，确定第i个深度点此时的
Figure BDA00030534369200000513
及n的值；其中，Vpi water-sat、Vsi water-sat分别为饱水层段第i个深度点的实测纵波速度、实测横波速度。
优选的，所述第二计算模块还包括：
第二调整单元，用于调整水合物饱和度Sh在(0，1)之间变化；
第二预测单元，用于基于调整后的水合物饱和度Sh分别重新计算待测地层第i个深度点的预测纵波速度
Figure BDA00030534369200000514
及预测横波速度
Figure BDA00030534369200000515
饱和度确定单元，用于当
Figure BDA00030534369200000516
时，令第i个深度点的第一水合物饱和度Sh1i等于当前调整的水合物饱和度Sh；以及当
Figure BDA00030534369200000517
时，令第i个深度点的第二水合物饱和度Sh2i等于当前调整的水合物饱和度Sh；其中，Vpi、Vsi分别为含水合物层段第i个深度点的实测纵波速度、实测横波速度。
依据本发明的一个方面，提供一种计算设备，包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序，所述处理器执行所述计算机程序时实现如上任一项所述的天然气水合物的饱和度估算方法。
本发明实施例所述的一种天然气水合物的饱和度估算方法，同时以纵、横波速度为约束进行饱和度估算，估算精度更高，为井内水合物饱和度的精确估算提供了一种比较可行的方法，对水合物储层储量评估及后续勘探开发有积极的指导作用。
上述说明仅是本发明技术方案的概述，为了能够更清楚了解本发明的技术手段，而可依照说明书的内容予以实施，并且为了让本发明的上述和其它目的、特征和优点能够更明显易懂，以下特举本发明的具体实施方式。
附图说明
为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。
图1为本发明实施例中一种天然气水合物的饱和度估算方法流程图；
图2为本发明实施例中一种天然气水合物的饱和度估算系统结构示意图；
图3为南海某含水合物工区待测地层的实测纵波速度Vp、横波速度Vs、孔隙度
Figure BDA0003053436920000061
的示意图；
图4(a)为通过本发明实施例所述的方法展现的纵波速度随水合物饱和度变化的关系示意图；
图4(b)为通过本发明实施例所述的方法展现的横波速度随水合物饱和度变化的关系示意图；
图5(a)为本发明实施例所述方法计算的饱水层的预测纵波速度与实测纵波速度对比图；
图5(b)为本发明实施例所述方法计算的饱水层的预测横波速度与实测横波速度对比图；
图6为本发明实施例中第一水合物饱和度、第二水合物饱和度以及最终的水合物饱和度值示意图；
图7为本发发明实施例中计算设备结构示意图。
具体实施方式
下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。
本发明实施例提供一种天然气水合物的饱和度估算方法，如图1所示，包括：
步骤101，获取待测地层的测井数据以及弹性参数；所述待测地层包括水合物层段以及饱水层段；所述测井数据包括：深度点的深度D、各个深度点的孔隙度
Figure BDA0003053436920000071
实测纵波速度Vp以及实测横波速度Vs；所述弹性参数包括：体积模量K、剪切模量μ和密度ρ。具体的，通过对工区测井以及岩心分析可以测得待测地层的测井数据以及弹性参数。其中，待测地层既包括含水合物层段以及不含水合物层段，不含水合物层段即为饱水层段。
步骤102，根据所述弹性参数分别计算待测地层的预测纵波速度Vp est及预测横波速度Vs est。具体而言，待测地层的预测纵波速度Vp est及预测横波速度Vs est与弹性参数有关，可以根据如下公式计算待测地层的预测纵波速度Vp est及预测横波速度Vs est：
Figure BDA0003053436920000072
步骤103，当饱水层段每个深度点的预测纵波速度Vp est与所述饱水层段实测纵波速度的误差、预测横波速度Vs est与饱水层段实测横波速度的误差均最小时，确定
Figure BDA0003053436920000073
及n；其中，
Figure BDA0003053436920000074
为临界孔隙度，n为在
Figure BDA0003053436920000075
处的配位数。其中，饱水层段意味着在此层段内水合物的饱和度为零。在本步骤中，分别判断饱水层段每个深度点的预测纵波速度Vp est与所述饱水层段实测纵波速度的误差，以及预测横波速度Vs est与饱水层段实测横波速度的误差，当在任意一个深度点的预测纵波速度Vp est与所述饱水层段实测纵波速度的误差最小，且预测横波速度Vs est与饱水层段实测横波速度的误差也为最小时，即可确定此时所对应的
Figure BDA0003053436920000076
及n的值。其中，预测纵波速度Vp est与所述饱水层段实测纵波速度的误差最小表示二者在最小二乘意义下最小；同样的，预测横波速度Vs est与饱水层段实测横波速度Vs的差最小也表示二者在最小二乘意义下最小。
步骤104，以实测纵波速度Vp以及实测横波速度Vs作为约束，根据确定的
Figure BDA0003053436920000081
及n的值重新计算待测地层的预测纵波速度Vp est及预测横波速度Vs est，并根据重新计算的预测纵波速度Vp est估算每个深度点的第一水合物饱和度Sh1，以及根据重新计算的预测横波速度Vs est估算每个深度点的第二水合物饱和度Sh2。在步骤103中，在饱水层段的纵波速度及横波速度预测较为准确的情况下确定所要使用的参数
Figure BDA0003053436920000082
和n后，根据确定好的
Figure BDA0003053436920000083
及n的值重新获得预测纵波速度Vp est及预测横波速度Vs est。
步骤105，根据每个深度点的第一水合物饱和度Sh1以及第二水合物饱和度Sh2获得该深度点的水合物饱和度。具体的，可以将每个深度点的第一水合物饱和度Sh1以及第二水合物饱和度Sh2的算术平均值作为该深度点的水合物饱和度。
本发明实施例所述的一种天然气水合物的饱和度估算方法，同时以纵、横波速度为约束进行饱和度估算，估算精度更高，为井内水合物饱和度的精确估算提供了一种比较可行的方法，对水合物储层储量评估及后续勘探开发有积极的指导作用。
本发明实施例所述的一种天然气水合物的饱和度估算方法，较佳的，计算地层的预测纵波速度Vp est及预测横波速度Vs est前，所述方法还包括：
模拟孔隙填充和颗粒支撑两种分布状态同时存在时的水合物等效体积参数ε；所述水合物等效体积参数ε表征水合物中孔隙填充部分的体积百分含量；即ε表征相对于固结作用，水合物的生成对于岩石固结的刚度的影响。其中，ε∈[0，1]，当ε＝0时，表示水合物完全作为骨架的一部分；当ε＝1时，表示水合物完全作为孔隙流体的一部分。具体的，由于现有技术的等效介质理论(EMT)可以模拟水合物作为骨架(颗粒支撑)一部分或孔隙流体(孔隙填充)一部分这两种赋存状态下的等效弹性模量，但是不包含两种分布状态同时存在的模型，因此本发明具体实施例提供一种改进的等效介质理论，即孔隙填充和颗粒支撑两种分布状态同时存在。此时，需要模拟孔隙填充和颗粒支撑两种分布状态同时存在时的水合物等效体积参数ε。根据ε的取值不同，可以模拟孔隙填充和颗粒支撑两种分布状态的不同含量。
基于公式一获取视孔隙度
Figure BDA0003053436920000084
所述视孔隙度
Figure BDA0003053436920000085
表征水合物对孔隙度
Figure BDA0003053436920000086
的减小作用，
Figure BDA0003053436920000087
其中，
Figure BDA0003053436920000088
为孔隙水的孔隙度，
Figure BDA0003053436920000089
为纯水合物的孔隙度；
较佳的，其中，
Figure BDA00030534369200000810
且
Figure BDA00030534369200000811
Sh为水合物饱和度。
较佳的实施例中，基于公式二、公式三以及公式四分别获得待测地层的体积模量K、剪切模量μ以及密度ρ，
Figure BDA0003053436920000091
μ＝μdry 公式三，
Figure BDA0003053436920000092
其中，Kdry、μdry分别为海洋干骨架的体积模量、剪切模量，Km为含水合物颗粒时的水合物基质的体积模量，Kfl为含水合物颗粒时的孔隙流体的体积模量，ρw为孔隙水的密度，ρma为不含水合物颗粒时的岩石基质的密度，Sh为水合物饱和度，ρh为纯水合物的密度。
较佳的实施例中，所述方法还包括：可以通过如下公式对公式一至公式四进行推算；
Figure BDA0003053436920000093
Figure BDA0003053436920000094
Figure BDA0003053436920000095
Figure BDA0003053436920000096
Figure BDA0003053436920000097
Figure BDA0003053436920000098
P＝(ρm-ρfl)gD
其中，KHM、μHM分别为Hertz-Mindlin模型的有效体积模量、剪切模量；ν为基质的泊松比；μm、ρm分别是含水合物颗粒时的基质的剪切模量、密度；ρfl为含水合物颗粒时的孔隙流体的密度；P是有效压力；g是重力加速度；D是深度点的深度；Z及P为中间变量。
较佳的实施例中，Km、μm、ρm分别通过以下公式获得：
Figure BDA0003053436920000101
Figure BDA0003053436920000102
ρm＝ρma(1-fh)+fhρh
Figure BDA0003053436920000103
其中，Kma、μma分别是不含水合物颗粒时的岩石基质的体积模量、剪切模量，Kh、μh、ρh、fh为纯水合物的体积模量、剪切模量、密度、体积百分比。
较佳的，不含水合物颗粒时的岩石基质为岩石组分的颗粒。其中，不含水合物颗粒时的岩石基质的体积模量Kma、剪切模量μma、密度ρma分别通过下式获得：
Figure BDA0003053436920000104
Figure BDA0003053436920000105
ρma＝∑fiρi
其中，Ki、μi、ρi和fi分别为待测地层中各种矿物颗粒的体积模量、剪切模量、密度、体积百分比，可以直接通过测井数据获得。由于待测地层包括各矿物质、孔隙水以及纯水合物，各矿物质可以包括石英、粘土、方解石等。
较佳的实施例中，含水合物颗粒时的孔隙流体的体积模量Kfl以及密度ρfl可以通过以下公式获得：
Figure BDA0003053436920000106
ρfl＝fwKw+(1-fw)Kh
Figure BDA0003053436920000107
其中，Kw为孔隙水的体积模量，Kh为纯水合物的体积模量，
Figure BDA0003053436920000111
为孔隙水的孔隙度，fw为孔隙水的体积百分比。
本发明实施例所述的一种天然气水合物的饱和度估算方法，较佳的，当饱水层段每个深度点的预测纵波速度Vp est与所述实测纵波速度Vp、预测横波速度Vs est与实测横波速度Vs的差均最小时确定
Figure BDA0003053436920000112
及n的值，具体包括：
在所述水合物饱和度Sh＝0时，调整
Figure BDA0003053436920000113
及n，其中
Figure BDA0003053436920000114
8≤n≤9.5。具体的，在本发明实施例中，饱水层段意味着水合物饱和度Sh＝0；当满足Sh＝0的前提下，在合理的范围内，即
Figure BDA0003053436920000115
8≤n≤9.5不断改变
Figure BDA0003053436920000116
及n的值，并在每次改变
Figure BDA0003053436920000117
及n的值后，重新计算第i个深度点的预测纵波速度
Figure BDA0003053436920000118
及预测横波速度
Figure BDA0003053436920000119
基于调整后的
Figure BDA00030534369200001110
及n的值计算第i个深度点的预测纵波速度
Figure BDA00030534369200001111
及预测横波速度
Figure BDA00030534369200001112
后，判断
Figure BDA00030534369200001113
以及
Figure BDA00030534369200001114
是否满足，即在饱水层段各个深度点的预测的纵横波速度与测井数据中实测的纵横波速度的误差在最小二乘意义下最小。
若
Figure BDA00030534369200001115
以及
Figure BDA00030534369200001116
同时满足，则确定该第i个深度点此时的
Figure BDA00030534369200001117
及n的值，表明此时的
Figure BDA00030534369200001118
及n的值为最终要使用的值；其中，Vpi water-sat、Vsi water-sat分别为饱水层段第i个深度点的实测纵波速度、实测横波速度。
本发明实施例所述的一种天然气水合物的饱和度估算方法，较佳的，所述方法还包括：在确定好
Figure BDA00030534369200001119
及n后，需要重新计算待测地层第i个深度点的预测纵波速度
Figure BDA00030534369200001120
及预测横波速度
Figure BDA00030534369200001121
具体的步骤为：
调整水合物饱和度Sh在(0，1)之间变化；在本步骤中，基于孔隙填充和颗粒支撑两种分布状态同时存在时的各种情况模拟下，不断改变水合物饱和度Sh在(0，1)之间变化。
基于调整后的水合物饱和度Sh分别重新计算待测地层第i个深度点的预测纵波速度
Figure BDA00030534369200001122
及预测横波速度
Figure BDA00030534369200001123
具体的计算方式同上，因此此处不再赘述。
当
Figure BDA00030534369200001124
时，令第i个深度点的第一水合物饱和度Sh1i等于当前调整的水合物饱和度Sh；
以及当
Figure BDA00030534369200001125
时，令第i个深度点的第二水合物饱和度Sh2i等于当前调整的水合物饱和度Sh；其中，Vpi、Vsi分别为含水合物层段第i个深度点的实测纵波速度、实测横波速度。
具体的，当
Figure BDA0003053436920000121
时，意味着在第i深度点的预测纵波速度
Figure BDA0003053436920000122
的误差最小，即实测纵波速度
Figure BDA0003053436920000123
与预测纵波速度
Figure BDA0003053436920000124
的差趋近于零，此时的水合物饱和度Sh记为第一水合物饱和度Sh1i；当
Figure BDA0003053436920000125
时，意味着在第i深度点的预测横波速度
Figure BDA0003053436920000126
的误差最小，即实测横波速度Vsi与预测横波速度
Figure BDA0003053436920000127
的差趋近于零，此时的水合物饱和度Sh记为第二水合物饱和度Sh2i。需要注意的是，预测纵波速度
Figure BDA0003053436920000128
的误差最小的深度点与预测横波速度
Figure BDA0003053436920000129
误差最小的深度点可以为同一深度点，也可以为不同深度点。即在确认第一水合物饱和度Sh1i时，只要是预测纵波速度的误差达到最小时，无需确认此时的预测横波速度的误差是否也为最小；同样的，在确认第二水合物饱和度Sh2i时，只需要判断第i个深度点的预测横波速度
Figure BDA00030534369200001210
是否误差最小，无需确认此时的预测纵波速度是否也为最小。
本发明实施例所述的一种天然气水合物的饱和度估算方法，较佳的，所述方法还包括：
根据公式五获取含水合物层段第i个深度点的水合物饱和度Shi：
Figure BDA00030534369200001211
通过公式五可知，将第i深度点的饱和度Sh1i值和Sh2i值取算术平均，作为该深度点最终的水合物饱和度值。
与现有技术相比，本发明的有益效果是：提出一种基于改进等效介质理论的天然气水合物饱和度估算方法，通过引入模拟孔隙填充和颗粒支撑两种分布状态同时存在时的水合物等效体积参数ε对传统的等效介质理论进行改进，可模拟多种水合物赋存状态同时存在时的水合物层等效弹性模量，因而可以实现更精确的纵横波速度预测；并结合实际测井资料，同时以纵、横波速度为约束进行饱和度估算，估算精度更高，为井内水合物饱和度精确估算提供了一种比较可行的方法，对水合物储层储量评估及后续勘探开发有积极的指导作用。
本发明实施例提供一种天然气水合物的饱和度估算系统，如图2所示，包括：
数据获取模块201，用于获取待测地层的测井数据以及弹性参数；所述待测地层包括水合物层段以及饱水层段；所述测井数据包括：深度点的深度D、各个深度点的孔隙度
Figure BDA00030534369200001212
实测纵波速度Vp以及实测横波速度Vs；所述弹性参数包括：体积模量K、剪切模量μ和密度ρ；
第一计算模块202，用于根据所述弹性参数分别计算待测地层的预测纵波速度Vp est及预测横波速度Vs est；
参数确定模块203，用于当饱水层段每个深度点的预测纵波速度Vp est与所述实测纵波速度Vp、预测横波速度Vs est与实测横波速度Vs的差均最小时确定
Figure BDA0003053436920000131
及n的值；其中，
Figure BDA0003053436920000132
为临界孔隙度，n为在
Figure BDA0003053436920000133
处的配位数；
第二计算模块204，用于根据确定的
Figure BDA0003053436920000134
及n的值重新计算待测地层的预测纵波速度Vp est及预测横波速度Vs est，并根据重新计算的预测纵波速度Vp est估算每个深度点的第一水合物饱和度Sh1，以及根据重新计算的预测横波速度Vs est估算每个深度点的第二水合物饱和度Sh2；
第三计算模块205，用于根据每个深度点的第一水合物饱和度Sh1以及第二水合物饱和度Sh2获得该深度点的水合物饱和度。
本发明实施例所述的一种天然气水合物的饱和度估算系统，较佳的，包括：
模拟模块，用于模拟孔隙填充和颗粒支撑两种分布状态同时存在时的水合物等效体积参数ε；所述水合物等效体积参数ε表征水合物中孔隙填充部分的体积百分含量；其中，ε∈[0，1]，当ε＝0时，表示水合物完全作为骨架的一部分；当ε＝1时，表示水合物完全作为孔隙流体的一部分；
第四计算模块，用于基于公式一获取视孔隙度
Figure BDA0003053436920000135
所述视孔隙度
Figure BDA0003053436920000136
表征水合物对孔隙度
Figure BDA0003053436920000137
的减小作用，
Figure BDA0003053436920000138
其中，
Figure BDA0003053436920000139
为孔隙水的孔隙度，
Figure BDA00030534369200001310
为纯水合物的孔隙度；
第五计算模块，用于基于公式二、公式三以及公式四分别获得待测地层的体积模量K、剪切模量μ以及密度ρ，
Figure BDA00030534369200001311
μ＝μdry 公式三，
Figure BDA00030534369200001312
其中，Kdry、μdry分别为海洋干骨架的体积模量、剪切模量，Km为含水合物颗粒时的水合物基质的体积模量，Kfl为含水合物颗粒时的孔隙流体的体积模量，ρw为孔隙水的密度，ρma为不含水合物颗粒时的岩石基质的密度，Sh为水合物饱和度，ρh为纯水合物的密度。
本发明实施例所述的一种天然气水合物的饱和度估算系统，较佳的，所述参数确定模块还包括：
第一调整单元，用于在所述水合物饱和度Sh＝0时，调整
Figure BDA0003053436920000141
及n的值，其中
Figure BDA0003053436920000142
8≤n≤9.5；
第一预测单元，用于基于调整后的
Figure BDA0003053436920000143
及n的值计算第i个深度点的预测纵波速度
Figure BDA0003053436920000144
及预测横波速度
Figure BDA0003053436920000145
参数确定单元，用于当
Figure BDA0003053436920000146
以及
Figure BDA0003053436920000147
满足时，确定第i个深度点此时的
Figure BDA0003053436920000148
及n的值；其中，Vpi water-sat、Vsi water-sat分别为饱水层段第i个深度点的实测纵波速度、实测横波速度。
本发明实施例所述的一种天然气水合物的饱和度估算系统，较佳的，所述第二计算模块还包括：
第二调整单元，用于调整水合物饱和度Sh在(0，1)之间变化；
第二预测单元，用于基于调整后的水合物饱和度Sh分别重新计算待测地层第i个深度点的预测纵波速度
Figure BDA0003053436920000149
及预测横波速度
Figure BDA00030534369200001410
饱和度确定单元，用于当
Figure BDA00030534369200001411
时，令第i个深度点的第一水合物饱和度Sh1i等于当前调整的水合物饱和度Sh；以及当
Figure BDA00030534369200001412
时，令第i个深度点的第二水合物饱和度Sh2i等于当前调整的水合物饱和度Sh；其中，Vpi、Vsi分别为含水合物层段第i个深度点的实测纵波速度、实测横波速度。
本发明实施例提供一种计算机设备，包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序，所述处理器执行所述计算机程序时实现如上任一项实施例所述的天然气水合物的饱和度估算方法。
下面以南海某含水合物工区为例，对本发明其中之一的实施例所述的方法进行具体的阐述，其他实施例所述的方法可以类比该具体的实施方式。
1)数据收集：
收集该工区待测地层的测井数据，并通过岩心分析，获得此研究区地层矿物组分为粘土、方解石、石英、长石，其各自的体积参数以及弹性参数中的体积百分比、体积模量、剪切模量、密度取表1给定值。孔隙水及纯水合物的体积模量、剪切模量、密度取表2给定值。如图3所示，为该工区待测地层的实测纵波速度Vp、横波速度Vs、孔隙度
Figure BDA0003053436920000151
的示意图，其中矩形框内(207-240m)为认定的水合物层段。
表1
矿物组分 体积百分比(％) 体积模量(Pa) 剪切模量(Pa) 密度(Kg/m<sup>3</sup>)
粘土 35 20.9×10<sup>9</sup> 6.85×10<sup>9</sup> 2580
方解石 26 76.8×10<sup>9</sup> 32×10<sup>9</sup> 2710
石英 25 36×10<sup>9</sup> 45×10<sup>9</sup> 2650
长石 14 75.6×10<sup>9</sup> 25.6×10<sup>9</sup> 2630
表2
组分 密度(Kg/m<sup>3</sup>) 体积模量(Pa) 剪切模量(Pa)
纯水合物 900 5.6×10<sup>9</sup> 2.4×10<sup>9</sup>
孔隙水 1032 2.5×10<sup>9</sup> 0
2)不含水合物颗粒时的岩石基质的弹性参数计算：
根据上述各矿物组分的体积模量、剪切模量、密度获得不含水合物颗粒时的岩石基质的体积模量Kma、剪切模量μma、密度ρma，以备后续计算使用，其中，Ki、μi、ρi和fi分别为各种矿物颗粒的体积模量、剪切模量、密度、体积百分比，见表1。该阶段也可以在步骤3)执行时同时执行，具体的步骤顺序不做限制。
3)预测纵波速度Vp est及预测横波速度Vs est：
根据本发明说明书具体实施例中的公式推导获得，其中，
Figure BDA0003053436920000152
为各个深度点的孔隙度，由步骤1)得到；
Figure BDA0003053436920000153
是临界孔隙度；n是在
Figure BDA0003053436920000154
处的配位数；Kma、μma、ρma由步骤2)得到；Kh、μh、ρh分别是纯水合物的体积模量、剪切模量、密度，由步骤1)得到；Kw、ρw分别是孔隙水的体积模量、密度，由步骤1)得到；g＝9.8m/s2是重力加速度；D是深度点的深度，由步骤1)得到；Sh为水合物饱和度；ε此处取ε＝0.12，为工区经验值。
如图4(a)以及图4(b)所示，分别为通过本发明具体实施例中所述的方法展现的纵波速度及横波速度随水合物饱和度变化的关系。可以发现，随着ε值的增大，颗粒支撑型水合物逐渐向孔隙填充型水合物过渡，对纵横波速度的影响也逐渐降低。该方法对模拟多赋存状态的水合物地层有着良好的优势。
4)在饱水层段的纵波速度及横波速度预测较为准确的情况下确定所要使用的参数
Figure BDA0003053436920000161
和n的值：
在合理的范围内(
Figure BDA0003053436920000162
8≤n≤9.5)，改变步骤4)中参数
Figure BDA0003053436920000163
和n的值，使得在水合物饱和度Sh＝0(即饱水层)时，重复利用步骤3)预测的各个深度点(饱水层)的纵横波速度与实测各个深度点(饱水层)的纵横波速度的误差在最小二乘意义下最小，即：
Figure BDA0003053436920000164
Figure BDA0003053436920000165
式中Vpi water-sat、Vsi water-sat表示第i个深度点(饱水层)的实测纵波速度及实测横波速度，由步骤1)得到；
Figure BDA0003053436920000166
Vsi est表示第i个深度点(饱水层)的预测纵波速度及预测横波速度，由步骤3)得到。
通过参数测试，发现
Figure BDA0003053436920000167
n＝8时，可满足上式要求，设置为最终要使用的值。具体参见图5(a)及图5(b)，分别为在
Figure BDA0003053436920000168
n＝8时，由本发明实施例所述方法计算的饱水层(不含水合物时)预测纵、横波速度与实测纵、横波速度对比图。可以发现在100-207m的饱水层段，预测的速度较为准确，说明参数
Figure BDA0003053436920000169
n＝8选择合适。
5)估算井中每个深度点的第一水合物饱和度：
将参数
Figure BDA00030534369200001610
和n＝8代入步骤3)的公式，通过调整水合物饱和度Sh的值(从0-1变化)，重新计算第i个深度点的预测纵波速度
Figure BDA00030534369200001611
并以步骤1)得到的这一深度点的实测纵波速度Vpi为约束，估算该深度点的水合物饱和度。当满足
Figure BDA00030534369200001612
时的水合物饱和度Sh的值即为第i个深度点的第一水合物饱和度，记为Sh1i。重复步骤5)得到井中各个深度点的第一水合物饱和度Sh1。
6)估算井中每个深度点的第二水合物饱和度：
将参数
Figure BDA00030534369200001613
和n＝8代入步骤3)的公式，通过调整水合物饱和度Sh的值(从0-1变化)，重新计算第i个深度点的预测横波速度
Figure BDA00030534369200001614
并以步骤1)得到的这一深度点的实测横波速度Vsi为约束，估算该深度点的水合物饱和度。当满足
Figure BDA00030534369200001615
时的水合物饱和度Sh的值即为第i个深度点的第二水合物饱和度，记为Sh2i。重复步骤6)得到井中各个深度点的第二水合物饱和度Sh2。
7)利用下式，将某个深度点的饱和度Sh1i值和Sh2i值取算术平均，作为该深度点最终的水合物饱和度值。
Figure BDA0003053436920000171
重复步骤7)得到井中各个深度点的水合物饱和度Sh。具体参见图6所示，为本实施例中第一水合物饱和度、第二水合物饱和度以及最终的水合物饱和度值。可见以实测纵波速度为约束的估算结果要大于以实测横波速度为约束的估算结果，通过取二者平均，可提供一个更为合理和准确的水合物饱和度值，可行性好。
如图7所示为本文实施例中计算设备的结构示意图，计算设备702可以包括一个或多个处理设备704，诸如一个或多个中央处理单元(CPU)，每个处理单元可以实现一个或多个硬件线程。计算设备702还可以包括任何存储资源706，其用于存储诸如代码、设置、数据等之类的任何种类的信息。非限制性的，比如，存储资源706可以包括以下任一项或多种组合：任何类型的RAM，任何类型的ROM，闪存设备，硬盘，光盘等。更一般地，任何存储资源都可以使用任何技术来存储信息。进一步地，任何存储资源可以提供信息的易失性或非易失性保留。进一步地，任何存储资源可以表示计算设备702的固定或可移除部件。在一种情况下，当处理设备704执行被存储在任何存储资源或存储资源的组合中的相关联的指令时，计算设备702可以执行相关联指令的任一操作。计算设备702还包括用于与任何存储资源交互的一个或多个驱动机构708，诸如硬盘驱动机构、光盘驱动机构等。
计算设备702还可以包括输入/输出模块710(I/O)，其用于接收各种输入(经由输入设备712)和用于提供各种输出(经由输出设备714))。一个具体输出设备可以包括呈现设备716和相关联的图形用户接口(GUI)718。在其他实施例中，还可以不包括输入/输出模块710(I/O)、输入设备712以及输出设备714，仅作为网络中的一台计算设备。计算设备702还可以包括一个或多个网络接口720，其用于经由一个或多个通信链路722与其他设备交换数据。一个或多个通信总线724将上文所描述的部件耦合在一起。
通信链路722可以以任何方式实现，例如，通过局域网、广域网(例如，因特网)、点对点连接等、或其任何组合。通信链路722可以包括由任何协议或协议组合支配的硬连线链路、无线链路、路由器、网关功能、名称服务器等的任何组合。
应理解，在本发明的各种实施例中，上述各过程的序号的大小并不意味着执行顺序的先后，各过程的执行顺序应以其功能和内在逻辑确定，而不应对本发明实施例的实施过程构成任何限定。
还应理解，在本发明实施例中，术语“和/或”仅仅是一种描述关联对象的关联关系，表示可以存在三种关系。例如，A和/或B，可以表示：单独存在A，同时存在A和B，单独存在B这三种情况。另外，本文中字符“/”，一般表示前后关联对象是一种“或”的关系。
本领域普通技术人员可以意识到，结合本文中所公开的实施例描述的各示例的模块及算法步骤，能够以电子硬件、计算机软件或者二者的结合来实现，为了清楚地说明硬件和软件的可互换性，在上述说明中已经按照功能一般性地描述了各示例的组成及步骤。这些功能究竟以硬件还是软件方式来执行，取决于技术方案的特定应用和设计约束条件。专业技术人员可以对每个特定的应用来使用不同方法来实现所描述的功能，但是这种实现不应认为超出本发明的范围。
所属领域的技术人员可以清楚地了解到，为了描述的方便和简洁，上述描述的系统、装置和模块的具体工作过程，可以参考前述方法实施例中的对应过程，在此不再赘述。
在本申请所提供的几个实施例中，应该理解到，所揭露的系统、装置和方法，可以通过其它的方式实现。例如，以上所描述的装置实施例仅仅是示意性的，例如，所述模块的划分，仅仅为一种逻辑功能划分，实际实现时可以有另外的划分方式，例如多个模块或组件可以结合或者可以集成到另一个系统，或一些特征可以忽略，或不执行。另外，所显示或讨论的相互之间的耦合或直接耦合或通信连接可以是通过一些接口、装置或模块的间接耦合或通信连接，也可以是电的，机械的或其它的形式连接。
所述作为分离部件说明的模块可以是或者也可以不是物理上分开的，作为模块显示的部件可以是或者也可以不是物理模块，即可以位于一个地方，或者也可以分布到多个网络模块上。可以根据实际的需要选择其中的部分或者全部模块来实现本发明实施例方案的目的。
另外，在本发明各个实施例中的各功能模块可以集成在一个处理模块中，也可以是各个模块单独物理存在，也可以是两个或两个以上模块集成在一个模块中。上述集成的模块既可以采用硬件的形式实现，也可以采用软件功能模块的形式实现。
所述集成的模块如果以软件功能模块的形式实现并作为独立的产品销售或使用时，可以存储在一个计算机可读取存储介质中。基于这样的理解，本发明的技术方案本质上或者说对现有技术做出贡献的部分，或者该技术方案的全部或部分可以以软件产品的形式体现出来，该计算机软件产品存储在一个存储介质中，包括若干指令用以使得一台计算机设备(可以是个人计算机，服务器，或者网络设备等)执行本发明各个实施例所述方法的全部或部分步骤。而前述的存储介质包括：U盘、移动硬盘、只读存储器(ROM，Read-OnlyMemory)、随机存取存储器(RAM，Random Access Memory)、磁碟或者光盘等各种可以存储程序代码的介质。
本发明中应用了具体实施例对本发明的原理及实施方式进行了阐述，以上实施例的说明只是用于帮助理解本发明的方法及其核心思想；同时，对于本领域的一般技术人员，依据本发明的思想，在具体实施方式及应用范围上均会有改变之处，综上所述，本说明书内容不应理解为对本发明的限制。

Claims (10)

1.一种天然气水合物的饱和度估算方法，其特征在于，包括：
获取待测地层的测井数据以及弹性参数；所述待测地层包括水合物层段以及饱水层段；所述测井数据包括：深度点的深度、各个深度点的孔隙度、实测纵波速度以及实测横波速度；所述弹性参数包括：体积模量、剪切模量和密度；
根据所述弹性参数分别计算待测地层的预测纵波速度及预测横波速度；
当饱水层段每个深度点的预测纵波速度与所述饱水层段实测纵波速度的误差、预测横波速度与饱水层段实测横波速度的误差均最小时确定临界孔隙度及在所述临界孔隙度处的配位数；
以实测纵波速度以及实测横波速度作为约束，根据确定临界孔隙度及在所述临界孔隙度处的配位数重新计算待测地层的预测纵波速度及预测横波速度，并根据重新计算的预测纵波速度估算每个深度点的第一水合物饱和度，以及根据重新计算的预测横波速度估算每个深度点的第二水合物饱和度；
根据每个深度点的第一水合物饱和度以及第二水合物饱和度获得该深度点的水合物饱和度。
2.根据权利要求1所述的一种天然气水合物的饱和度估算方法，其特征在于，计算地层的预测纵波速度及预测横波速度前，所述方法还包括：
模拟孔隙填充和颗粒支撑两种分布状态同时存在时的水合物等效体积参数ε；所述水合物等效体积参数ε表征水合物中孔隙填充部分的体积百分含量；其中，
ε∈[0，1]，当ε＝0时，表示水合物完全作为骨架的一部分；当ε＝1时，表示水合物完全作为孔隙流体的一部分；
基于公式一获取视孔隙度
Figure FDA0003053436910000011
所述视孔隙度
Figure FDA0003053436910000012
表征水合物对孔隙度
Figure FDA0003053436910000013
的减小作用，
Figure FDA0003053436910000014
其中，
Figure FDA0003053436910000015
为孔隙水的孔隙度，
Figure FDA0003053436910000016
为纯水合物的孔隙度；
基于公式二、公式三以及公式四分别获得待测地层的体积模量K、剪切模量μ以及密度ρ，
Figure FDA0003053436910000017
μ＝μdry 公式三，
Figure FDA0003053436910000021
其中，Kdry、μdry分别为海洋干骨架的体积模量、剪切模量，Km为含水合物颗粒时的水合物基质的体积模量，Kfl为含水合物颗粒时的孔隙流体的体积模量，ρw为孔隙水的密度，ρma为不含水合物颗粒时的岩石基质的密度，Sh为水合物饱和度，ρh为纯水合物的密度。
3.根据权利要求2所述的一种天然气水合物的饱和度估算方法，其特征在于，当饱水层段每个深度点的预测纵波速度与所述实测纵波速度的误差、预测横波速度与实测横波速度的误差均最小时确定临界孔隙度及在所述临界孔隙度处的配位数，具体包括：
在所述水合物饱和度Sh＝0时，调整临界孔隙度
Figure FDA0003053436910000022
及在所述临界孔隙度处的配位数n，其中
Figure FDA0003053436910000023
8≤n≤9.5；
基于调整后的
Figure FDA0003053436910000024
及n的值计算第i个深度点的预测纵波速度
Figure FDA00030534369100000213
及预测横波速度
Figure FDA0003053436910000025
当
Figure FDA0003053436910000026
以及
Figure FDA0003053436910000027
满足时，确定第i个深度点此时的
Figure FDA0003053436910000028
及n的值；其中，Vpi water-sat、Vsi water-sat分别为饱水层段第i个深度点的实测纵波速度、实测横波速度。
4.根据权利要求3所述的一种天然气水合物的饱和度估算方法，其特征在于，所述方法还包括：
调整水合物饱和度Sh在(0，1)之间变化；
基于调整后的水合物饱和度Sh分别重新计算待测地层第i个深度点的预测纵波速度
Figure FDA00030534369100000211
及预测横波速度
Figure FDA0003053436910000029
当
Figure FDA00030534369100000212
时，令第i个深度点的第一水合物饱和度Sh1i等于当前调整的水合物饱和度Sh；
以及当
Figure FDA00030534369100000210
时，令第i个深度点的第二水合物饱和度Sh2i等于当前调整的水合物饱和度Sh；其中，Vpi、Vsi分别为含水合物层段第i个深度点的实测纵波速度、实测横波速度。
5.根据权利要求4所述的一种天然气水合物的饱和度估算方法，其特征在于，所述方法还包括：
根据公式五获取含水合物层段第i个深度点的水合物饱和度Shi：
Figure FDA0003053436910000031
6.一种天然气水合物的饱和度估算系统，其特征在于，包括：
数据获取模块，用于获取待测地层的测井数据以及弹性参数；所述待测地层包括水合物层段以及饱水层段；所述测井数据包括：深度点的深度、各个深度点的孔隙度、实测纵波速度以及实测横波速度；所述弹性参数包括：体积模量、剪切模量和密度；
第一计算模块，用于根据公式一及公式二分别计算待测地层的预测纵波速度及预测横波速度，
参数确定模块，用于当饱水层段每个深度点的预测纵波速度与所述实测纵波速度的误差、预测横波速度与实测横波速度的误差均最小时确定临界孔隙度及在所述临界孔隙度处的配位数；
第二计算模块，用于以实测纵波速度以及实测横波速度作为约束，根据确定的临界孔隙度及在所述临界孔隙度处的配位数重新计算待测地层的预测纵波速度及预测横波速度，并根据重新计算的预测纵波速度估算每个深度点的第一水合物饱和度，以及根据重新计算的预测横波速度估算每个深度点的第二水合物饱和度；
第三计算模块，用于根据每个深度点的第一水合物饱和度以及第二水合物饱和度获得该深度点的水合物饱和度。
7.根据权利要求6所述的一种天然气水合物的饱和度估算系统，其特征在于，包括：
模拟模块，用于模拟孔隙填充和颗粒支撑两种分布状态同时存在时的水合物等效体积参数ε；所述水合物等效体积参数ε表征水合物中孔隙填充部分的体积百分含量；其中，ε∈[0，1]，当ε＝0时，表示水合物完全作为骨架的一部分；当ε＝1时，表示水合物完全作为孔隙流体的一部分；
第四计算模块，用于基于公式一获取视孔隙度
Figure FDA0003053436910000032
所述视孔隙度
Figure FDA0003053436910000033
表征水合物对孔隙度
Figure FDA0003053436910000034
的减小作用，
Figure FDA0003053436910000035
其中，
Figure FDA0003053436910000036
为孔隙水的孔隙度，
Figure FDA0003053436910000037
为纯水合物的孔隙度；
第五计算模块，用于基于公式二、公式三以及公式四分别获得待测地层的体积模量K、剪切模量μ以及密度ρ，
Figure FDA0003053436910000041
μ＝μdry 公式三，
Figure FDA0003053436910000042
其中，Kdry、μdry分别为海洋干骨架的体积模量、剪切模量，Km为含水合物颗粒时的水合物基质的体积模量，Kfl为含水合物颗粒时的孔隙流体的体积模量，ρw为孔隙水的密度，ρma为不含水合物颗粒时的岩石基质的密度，Sh为水合物饱和度，ρh为纯水合物的密度。
8.根据权利要求7所述的一种天然气水合物的饱和度估算系统，其特征在于，所述参数确定模块还包括：
第一调整单元，用于在所述水合物饱和度Sh＝0时，调整临界孔隙度
Figure FDA0003053436910000043
及在所述临界孔隙度处的配位数n，其中
Figure FDA0003053436910000044
8≤n≤9.5；
第一预测单元，用于基于调整后的
Figure FDA0003053436910000045
及n的值计算第i个深度点的预测纵波速度
Figure FDA00030534369100000414
及预测横波速度
Figure FDA0003053436910000046
参数确定单元，用于当
Figure FDA0003053436910000047
以及
Figure FDA0003053436910000048
满足时，确定第i个深度点此时的
Figure FDA0003053436910000049
及n的值；其中，Vpi water-sat、Vsi water-sat分别为饱水层段第i个深度点的实测纵波速度、实测横波速度。
9.根据权利要求8所述的一种天然气水合物的饱和度估算系统，其特征在于，所述第二计算模块还包括：
第二调整单元，用于调整水合物饱和度Sh在(0，1)之间变化；
第二预测单元，用于基于调整后的水合物饱和度Sh分别重新计算待测地层第i个深度点的预测纵波速度
Figure FDA00030534369100000410
及预测横波速度
Figure FDA00030534369100000411
饱和度确定单元，用于当
Figure FDA00030534369100000412
时，令第i个深度点的第一水合物饱和度Sh1i等于当前调整的水合物饱和度Sh；以及当
Figure FDA00030534369100000413
时，令第i个深度点的第二水合物饱和度Sh2i等于当前调整的水合物饱和度Sh；其中，Vpi、Vsi分别为含水合物层段第i个深度点的实测纵波速度、实测横波速度。
10.一种计算设备，包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序，其特征在于，所述处理器执行所述计算机程序时实现权利要求1至5中任一项所述的天然气水合物的饱和度估算方法。
CN202110493721.0A 2021-05-07 2021-05-07 一种天然气水合物的饱和度估算方法及系统 Active CN113189674B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202110493721.0A CN113189674B (zh) 2021-05-07 2021-05-07 一种天然气水合物的饱和度估算方法及系统

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202110493721.0A CN113189674B (zh) 2021-05-07 2021-05-07 一种天然气水合物的饱和度估算方法及系统

Publications (2)

Publication Number Publication Date
CN113189674A true CN113189674A (zh) 2021-07-30
CN113189674B CN113189674B (zh) 2024-04-23

Family

ID=76983909

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202110493721.0A Active CN113189674B (zh) 2021-05-07 2021-05-07 一种天然气水合物的饱和度估算方法及系统

Country Status (1)

Country Link
CN (1) CN113189674B (zh)

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113466963A (zh) * 2021-08-13 2021-10-01 广州海洋地质调查局 裂缝充填ii型水合物饱和度估算方法及处理终端
CN117316329A (zh) * 2023-11-23 2023-12-29 中国石油大学(华东) 天然气水合物饱和度声电测井联合智能反演方法及系统

Citations (12)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2006031872A1 (en) * 2004-09-14 2006-03-23 The Trustees Of Columbia University In The City Of New York Systems and methods for determining in-situ gas hydrate saturation
WO2008157737A2 (en) * 2007-06-21 2008-12-24 Schlumberger Canada Limited Multi-attribute seismic characterization of gas hydrates
US20090306895A1 (en) * 2008-06-09 2009-12-10 Ohm Limited Geophysical surveying
US20110137568A1 (en) * 2009-12-07 2011-06-09 Sky Hunter Corporation Method and apparatus for detection of hydrocarbon deposits
CN103954999A (zh) * 2014-05-06 2014-07-30 成都理工大学 一种适用于低孔隙度砂泥岩地层的横波速度预测新方法
CN109212162A (zh) * 2018-08-17 2019-01-15 西南石油大学 一种估算成岩天然气水合物储层力学参数的方法
CN109298463A (zh) * 2018-09-27 2019-02-01 北京珠玛阳光科技有限公司 一种联合纵、横波速度和密度预测天然气水合物方法
CN111220522A (zh) * 2019-10-28 2020-06-02 中国石油大学(华东) 一种高泥质细粒沉积物水合物饱和度的岩心刻度测井计算方法
CN111856560A (zh) * 2020-07-08 2020-10-30 中国海洋大学 一种天然气水合物储层信息评估方法及其应用
CN112133377A (zh) * 2020-08-28 2020-12-25 中国石油天然气集团有限公司 一种海域天然气水合物赋存状态的判别方法及系统
CN112149282A (zh) * 2020-08-28 2020-12-29 中国石油天然气集团有限公司 一种井中天然气水合物饱和度岩石物理计算方法及系统
US10969514B1 (en) * 2020-09-09 2021-04-06 China University Of Petroleum (East China) Method for correcting permeability model of porous medium in dissociation process of gas hydrate, and method and system for determining permeability of hydrate-bearing porous medium

Patent Citations (12)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2006031872A1 (en) * 2004-09-14 2006-03-23 The Trustees Of Columbia University In The City Of New York Systems and methods for determining in-situ gas hydrate saturation
WO2008157737A2 (en) * 2007-06-21 2008-12-24 Schlumberger Canada Limited Multi-attribute seismic characterization of gas hydrates
US20090306895A1 (en) * 2008-06-09 2009-12-10 Ohm Limited Geophysical surveying
US20110137568A1 (en) * 2009-12-07 2011-06-09 Sky Hunter Corporation Method and apparatus for detection of hydrocarbon deposits
CN103954999A (zh) * 2014-05-06 2014-07-30 成都理工大学 一种适用于低孔隙度砂泥岩地层的横波速度预测新方法
CN109212162A (zh) * 2018-08-17 2019-01-15 西南石油大学 一种估算成岩天然气水合物储层力学参数的方法
CN109298463A (zh) * 2018-09-27 2019-02-01 北京珠玛阳光科技有限公司 一种联合纵、横波速度和密度预测天然气水合物方法
CN111220522A (zh) * 2019-10-28 2020-06-02 中国石油大学(华东) 一种高泥质细粒沉积物水合物饱和度的岩心刻度测井计算方法
CN111856560A (zh) * 2020-07-08 2020-10-30 中国海洋大学 一种天然气水合物储层信息评估方法及其应用
CN112133377A (zh) * 2020-08-28 2020-12-25 中国石油天然气集团有限公司 一种海域天然气水合物赋存状态的判别方法及系统
CN112149282A (zh) * 2020-08-28 2020-12-29 中国石油天然气集团有限公司 一种井中天然气水合物饱和度岩石物理计算方法及系统
US10969514B1 (en) * 2020-09-09 2021-04-06 China University Of Petroleum (East China) Method for correcting permeability model of porous medium in dissociation process of gas hydrate, and method and system for determining permeability of hydrate-bearing porous medium

Non-Patent Citations (4)

* Cited by examiner, † Cited by third party
Title
刘欣欣;印兴耀;栾锡武;: "天然气水合物地层岩石物理模型构建", 中国科学:地球科学, no. 09, 20 September 2018 (2018-09-20) *
孟大江;文鹏飞;张如伟;赵斌;李延;: "天然气水合物横波速度等效介质模型预测方法", 石油地球物理勘探, no. 01, 15 February 2020 (2020-02-15) *
李红星;陶春辉;周建平;邓居智;邓显明;方根显;: "非胶结含水合物沉积物修正等效介质速度模型及其地震波场特征研究", 物理学报, no. 11, 15 November 2009 (2009-11-15) *
高红艳;钟广法;梁金强;郭依群;: "应用改进的Biot-Gassmann模型估算天然气水合物的饱和度", 海洋地质与第四纪地质, no. 04, 28 August 2012 (2012-08-28) *

Cited By (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113466963A (zh) * 2021-08-13 2021-10-01 广州海洋地质调查局 裂缝充填ii型水合物饱和度估算方法及处理终端
CN113466963B (zh) * 2021-08-13 2022-02-01 广州海洋地质调查局 裂缝充填ii型水合物饱和度估算方法及处理终端
CN117316329A (zh) * 2023-11-23 2023-12-29 中国石油大学(华东) 天然气水合物饱和度声电测井联合智能反演方法及系统
CN117316329B (zh) * 2023-11-23 2024-03-29 中国石油大学(华东) 天然气水合物饱和度声电测井联合智能反演方法及系统

Also Published As

Publication number Publication date
CN113189674B (zh) 2024-04-23

Similar Documents

Publication Publication Date Title
CN102353989B (zh) 基于自适应基质矿物等效弹性模量反演的横波速度估算方法
CN112149282B (zh) 一种井中天然气水合物饱和度岩石物理计算方法及系统
Priest Determination of discontinuity size distributions from scanline data
Piggott et al. Anisotropic mesh adaptivity for multi-scale ocean modelling
CN109655903B (zh) 页岩层横波速度预测方法及系统
CN113189674B (zh) 一种天然气水合物的饱和度估算方法及系统
CN112133377B (zh) 一种海域天然气水合物赋存状态的判别方法及系统
CN105425280A (zh) 一种矿物模量和孔隙结构预测方法
CN112946783B (zh) 一种水合物饱和度确定方法、装置及设备
CN103954999A (zh) 一种适用于低孔隙度砂泥岩地层的横波速度预测新方法
WO2013158873A2 (en) System and method for calibrating permeability for use in reservoir modeling
CN107194104A (zh) 一种采用测井曲线计算页岩储层中无机质孔隙度、总孔隙度以及组分的方法
CN105804737A (zh) 一种基于迭代算法求解地层孔隙度的方法
CN105317429A (zh) 一种确定页岩地层干酪根含量的方法
CN106285652B (zh) 确定页岩游离气体饱和度的方法
CN105275459A (zh) 一种确定页岩地层可动水体积含量的方法
CN109839676B (zh) 一种基于记忆模拟退火的基质模量估算方法及电子设备
EP2885660A1 (en) A dual porosity methodology for mineral volume calculations in source rock
CN104316962A (zh) 岩石孔隙结构参数反演方法及其系统
CN112343574B (zh) 一种溶蚀型储层的孔隙度录井计算方法
CN109655936B (zh) 一种碎屑岩岩性替换的弹性参数计算方法及系统
CN111208566B (zh) 基于sca模型的孔缝参数反演方法、装置及存储介质
CN110344822A (zh) 含钙质致密储层含水饱和度的确定方法
CN113466962B (zh) 一种基于高斯过程回归的横波预测方法及系统
CN110954942A (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
GR01 Patent grant
GR01 Patent grant