CN112149282B - 一种井中天然气水合物饱和度岩石物理计算方法及系统 - Google Patents
一种井中天然气水合物饱和度岩石物理计算方法及系统 Download PDFInfo
- Publication number
- CN112149282B CN112149282B CN202010889675.1A CN202010889675A CN112149282B CN 112149282 B CN112149282 B CN 112149282B CN 202010889675 A CN202010889675 A CN 202010889675A CN 112149282 B CN112149282 B CN 112149282B
- Authority
- CN
- China
- Prior art keywords
- hydrate
- well
- saturation
- porosity
- density
- 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.)
- Active
Links
- 238000004364 calculation method Methods 0.000 title claims abstract description 65
- NMJORVOYSJLJGU-UHFFFAOYSA-N methane clathrate Chemical compound C.C.C.C.O.O.O.O.O.O.O.O.O.O.O.O.O.O.O.O.O.O.O.O.O.O.O NMJORVOYSJLJGU-UHFFFAOYSA-N 0.000 title claims abstract description 38
- 239000011435 rock Substances 0.000 title claims description 22
- 230000015572 biosynthetic process Effects 0.000 claims abstract description 93
- 238000007596 consolidation process Methods 0.000 claims abstract description 71
- 238000000034 method Methods 0.000 claims abstract description 71
- 238000012937 correction Methods 0.000 claims abstract description 50
- 238000013459 approach Methods 0.000 claims abstract description 8
- 230000006870 function Effects 0.000 claims description 41
- 230000004044 response Effects 0.000 claims description 35
- 239000011159 matrix material Substances 0.000 claims description 20
- 239000011148 porous material Substances 0.000 claims description 19
- 238000004590 computer program Methods 0.000 claims description 17
- 239000012530 fluid Substances 0.000 claims description 13
- 239000002245 particle Substances 0.000 claims description 13
- 229910052500 inorganic mineral Inorganic materials 0.000 claims description 12
- 239000011707 mineral Substances 0.000 claims description 12
- XLYOFNOQVPJJNP-UHFFFAOYSA-N water Substances O XLYOFNOQVPJJNP-UHFFFAOYSA-N 0.000 claims description 12
- 239000007787 solid Substances 0.000 claims description 11
- 239000004927 clay Substances 0.000 claims description 7
- 150000001875 compounds Chemical class 0.000 claims description 7
- 238000012545 processing Methods 0.000 claims description 7
- 230000001174 ascending effect Effects 0.000 claims description 4
- 238000005259 measurement Methods 0.000 claims description 4
- 239000013535 sea water Substances 0.000 claims description 4
- 239000000758 substrate Substances 0.000 claims description 4
- 238000004458 analytical method Methods 0.000 claims description 3
- VNWKTOKETHGBQD-UHFFFAOYSA-N methane Natural products C VNWKTOKETHGBQD-UHFFFAOYSA-N 0.000 claims description 3
- 238000011835 investigation Methods 0.000 claims description 2
- 239000003345 natural gas Substances 0.000 claims description 2
- -1 natural gas hydrates Chemical class 0.000 claims 1
- 230000000452 restraining effect Effects 0.000 abstract 1
- 238000010586 diagram Methods 0.000 description 13
- 238000011160 research Methods 0.000 description 7
- VEXZGXHMUGYJMC-UHFFFAOYSA-M Chloride anion Chemical compound [Cl-] VEXZGXHMUGYJMC-UHFFFAOYSA-M 0.000 description 4
- 238000011161 development Methods 0.000 description 4
- 238000004519 manufacturing process Methods 0.000 description 4
- 230000000694 effects Effects 0.000 description 3
- 238000011158 quantitative evaluation Methods 0.000 description 2
- 238000006467 substitution reaction Methods 0.000 description 2
- 229910021532 Calcite Inorganic materials 0.000 description 1
- 239000004568 cement Substances 0.000 description 1
- 230000008859 change Effects 0.000 description 1
- 150000004683 dihydrates Chemical class 0.000 description 1
- 238000005553 drilling Methods 0.000 description 1
- 230000007613 environmental effect Effects 0.000 description 1
- 238000011156 evaluation Methods 0.000 description 1
- 239000010433 feldspar Substances 0.000 description 1
- 150000004677 hydrates Chemical class 0.000 description 1
- 230000003993 interaction Effects 0.000 description 1
- 239000010445 mica Substances 0.000 description 1
- 229910052618 mica group Inorganic materials 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 239000008239 natural water Substances 0.000 description 1
- 230000003287 optical effect Effects 0.000 description 1
- 230000008569 process Effects 0.000 description 1
- 239000010453 quartz Substances 0.000 description 1
- 229920006395 saturated elastomer Polymers 0.000 description 1
- VYPSYNLAJGMNEJ-UHFFFAOYSA-N silicon dioxide Inorganic materials O=[Si]=O VYPSYNLAJGMNEJ-UHFFFAOYSA-N 0.000 description 1
- 238000000638 solvent extraction Methods 0.000 description 1
- 239000002023 wood Substances 0.000 description 1
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/20—Design optimisation, verification or simulation
-
- 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
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06Q—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES; SYSTEMS OR METHODS SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES, NOT OTHERWISE PROVIDED FOR
- G06Q10/00—Administration; Management
- G06Q10/04—Forecasting or optimisation specially adapted for administrative or management purposes, e.g. linear programming or "cutting stock problem"
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06Q—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES; SYSTEMS OR METHODS SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES, NOT OTHERWISE PROVIDED FOR
- G06Q50/00—Information and communication technology [ICT] specially adapted for implementation of business processes of specific business sectors, e.g. utilities or tourism
- G06Q50/06—Energy or water supply
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- Business, Economics & Management (AREA)
- Theoretical Computer Science (AREA)
- General Physics & Mathematics (AREA)
- Economics (AREA)
- Strategic Management (AREA)
- Mathematical Physics (AREA)
- Human Resources & Organizations (AREA)
- Operations Research (AREA)
- Computational Mathematics (AREA)
- Tourism & Hospitality (AREA)
- Marketing (AREA)
- General Business, Economics & Management (AREA)
- Data Mining & Analysis (AREA)
- Pure & Applied Mathematics (AREA)
- Health & Medical Sciences (AREA)
- Mathematical Optimization (AREA)
- Mathematical Analysis (AREA)
- General Engineering & Computer Science (AREA)
- Computer Hardware Design (AREA)
- General Health & Medical Sciences (AREA)
- Evolutionary Computation (AREA)
- Public Health (AREA)
- Primary Health Care (AREA)
- Water Supply & Treatment (AREA)
- Algebra (AREA)
- Geometry (AREA)
- Databases & Information Systems (AREA)
- Software Systems (AREA)
- Development Economics (AREA)
- Game Theory and Decision Science (AREA)
- Entrepreneurship & Innovation (AREA)
- Quality & Reliability (AREA)
- Geophysics And Detection Of Objects (AREA)
Abstract
本发明提出了一种井中天然气水合物饱和度岩石物理计算方法及系统，该方法包括：根据工区资料，计算得到井中地层孔隙度；经过反演得到当前深度点的地层最佳固结系数；对井中每个深度点进行井中地层孔隙度及固结系数的计算，得到整个井段的地层孔隙度及固结系数；利用三相Biot模型预测地层纵波速度及地层密度，并以修正系数为约束计算井中每个深度点的水合物饱和度；通过修改权重因子调整修正系数，利用新的修正系数为约束计算井中不含水合物层位的饱和度平均值，直至不含水合物层位的饱和度平均值趋近于0时，得到本井的权重因子；利用本井的权重因子，重新计算井中每个深度点的水合物饱和度，得到整个井段的水合物饱和度值。
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)等。
然而，水合物的存在不仅会改变地层声速，还会对地层密度、固结系数等参数产生影响，因此仅仅利用地层声速去估算饱和度是不够的，会影响结果的准确性；另外，不同的水合物饱和度估算方法或岩石物理模型适用于不同的地区和赋存状态，且模型中众多参数的设定都会对饱和度估算结果产生影响。
综上来看，亟需一种可以结合测井资料对模型参数进行约束，并根据不同地区、水合物赋存状态等因素综合考虑来计算水合物饱和度的技术方案。
发明内容
本发明提出了一种井中天然气水合物饱和度岩石物理计算方法及系统，该方法及系统利用改进的岩石物理模型，综合考虑水合物的存在对地层固结系数的影响，通过合理控制计算参数，实现对孔隙填充型和颗粒支撑型两种赋存状态水合物储层的饱和度估算，并引入修正系数对估算结果进行修正，进一步提高海域天然气水合物饱和度的估算精度。
具体的，在本发明实施例的第一方面，提出了一种井中天然气水合物饱和度岩石物理计算方法，该方法包括：
获取工区资料；
根据所述工区资料，分别对当前深度点的含水合物层位及不含水合物层位进行井中地层孔隙度计算，得到井中地层孔隙度；
根据所述工区资料及井中地层孔隙度，经过反演得到当前深度点的地层最佳固结系数；
对井中每个深度点进行井中地层孔隙度及固结系数的计算，得到整个井段的地层孔隙度及固结系数；
根据地层孔隙度及固结系数，利用三相Biot模型预测地层纵波速度及地层密度，并以修正系数为约束计算井中每个深度点的水合物饱和度；
通过修改权重因子调整修正系数，利用新的修正系数为约束计算井中不含水合物层位的饱和度平均值，直至所述不含水合物层位的饱和度平均值趋近于0时，得到本井的权重因子；
利用所述本井的权重因子，重新计算井中每个深度点的水合物饱和度，得到整个井段的水合物饱和度值。
在本发明实施例的第二方面，提出了一种井中天然气水合物饱和度岩石物理计算系统，该系统包括：
工区资料获取模块，用于获取工区资料；
孔隙度计算模块，用于根据所述工区资料，分别对当前深度点的含水合物层位及不含水合物层位进行井中地层孔隙度计算，得到井中地层孔隙度；
固结系数计算模块，用于根据所述工区资料及井中地层孔隙度，经过反演得到当前深度点的地层最佳固结系数；
井段数据处理模块，用于对井中每个深度点进行井中地层孔隙度及固结系数的计算，得到整个井段的地层孔隙度及固结系数；
饱和度计算模块，用于根据地层孔隙度及固结系数，利用三相Biot模型预测地层纵波速度及地层密度，并以修正系数为约束计算井中每个深度点的水合物饱和度；
权重因子调整模块，用于通过修改权重因子调整修正系数，利用新的修正系数为约束计算井中不含水合物层位的饱和度平均值，直至所述不含水合物层位的饱和度平均值趋近于0时，得到本井的权重因子；
所述饱和度计算模块，还用于利用所述本井的权重因子，重新计算井中每个深度点的水合物饱和度，得到整个井段的水合物饱和度值。
在本发明实施例的第三方面，提出了一种计算机设备，包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序，所述处理器执行所述计算机程序时实现井中天然气水合物饱和度岩石物理计算方法。
在本发明实施例的第四方面，提出了一种计算机可读存储介质，所述计算机可读存储介质存储有计算机程序，所述计算机程序被处理器执行时实现井中天然气水合物饱和度岩石物理计算方法。
本发明提出的井中天然气水合物饱和度岩石物理计算方法及系统，将基于岩石物理模型的井内水合物饱和度估算理论进行改进，并利用提出的修正系数对估算结果进行修正，估算精度更高；为井内水合物饱和度精确估算提供了一种可行性高的方案，对水合物的精确定量和加快实际生产开发有积极的推进作用。
附图说明
为了更清楚地说明本申请实施例技术方案，下面将对实施例描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图是本申请的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其它的附图。
图1是本发明一实施例的井中天然气水合物饱和度岩石物理计算方法的流程示意图。
图2是本发明一具体实施例的井内水合物饱和度估算的详细流程图。
图3是本发明一具体实施例的利用密度与中子交会法估算含水合物层孔隙度流程图。
图4是本发明一具体实施例的南海某井测井曲线示意图。
图5是本发明一具体实施例的南海某井地层固结系数散点图。
图6是本发明一具体实施例的权重因子为0时南海某井模型预测速度基线与实测速度对比、饱和度估算结果与氯离子法实测饱和度对比的示意图。
图7是本发明一具体实施例的权重因子为1.2时南海某井模型预测速度基线与实测速度对比、饱和度估算结果与氯离子法实测饱和度对比的示意图。
图8是本发明一具体实施例的是本发明一实施例的井中天然气水合物饱和度岩石物理计算方法的流程示意图。
图9是本发明一实施例的计算机设备结构示意图。
具体实施方式
下面将参考若干示例性实施方式来描述本发明的原理和精神。应当理解，给出这些实施方式仅仅是为了使本领域技术人员能够更好地理解进而实现本发明，而并非以任何方式限制本发明的范围。相反，提供这些实施方式是为了使本公开更加透彻和完整，并且能够将本公开的范围完整地传达给本领域的技术人员。
本领域技术人员知道，本发明的实施方式可以实现为一种系统、装置、设备、方法或计算机程序产品。因此，本公开可以具体实现为以下形式，即：完全的硬件、完全的软件(包括固件、驻留软件、微代码等)，或者硬件和软件结合的形式。
根据本发明的实施方式，提出了一种井中天然气水合物饱和度岩石物理计算方法及系统。
下面参考本发明的若干代表性实施方式，详细阐释本发明的原理和精神。
图1是本发明一实施例的井中天然气水合物饱和度岩石物理计算方法的流程示意图。如图1所示，该该方法包括：
步骤S1，获取工区资料；
步骤S2，根据所述工区资料，分别对当前深度点的含水合物层位及不含水合物层位进行井中地层孔隙度计算，得到井中地层孔隙度；
步骤S3，根据所述工区资料及井中地层孔隙度，经过反演得到当前深度点的地层最佳固结系数；
步骤S4，对井中每个深度点进行井中地层孔隙度及固结系数的计算，得到整个井段的地层孔隙度及固结系数；
步骤S5，根据地层孔隙度及固结系数，利用三相Biot模型预测地层纵波速度及地层密度，并以修正系数为约束计算井中每个深度点的水合物饱和度；
步骤S6，通过修改权重因子调整修正系数，利用新的修正系数为约束计算井中不含水合物层位的饱和度平均值，直至所述不含水合物层位的饱和度平均值趋近于0时，得到本井的权重因子；
步骤S7，利用所述本井的权重因子，重新计算井中每个深度点的水合物饱和度，得到整个井段的水合物饱和度值。
为了对上述井中天然气水合物饱和度岩石物理计算方法进行更为清楚的解释，下面结合一个具体的实施例对每一步进行详细说明。
参考图2，为本发明一具体实施例的井中水合物饱和度计算的详细流程示意图。如图2所示，具体流程为：
步骤S1：
步骤S11，通过工区测井和岩心分析，获得测井数据及岩石矿物数据；其中，
所述测井数据包括：研究区的深度、纵波速度、电阻率、密度、自然伽马及声波速度数据；
所述岩石矿物数据包括：井中地层矿物颗粒的体积模量、剪切模量、体积百分比、密度及中子响应值；
步骤S12，根据岩石矿物数据，得到固体成分的体积模量、剪切模量、密度及中子响应值，计算公式如下：
ρma＝∑fiρi；
其中，Kma、μma、ρma、分别为基质的体积模量、剪切模量、密度、中子响应值；Ki、μi、ρi、/>fi分别为每种矿物颗粒的体积模量、剪切模量、密度、中子响应值、体积百分比。
步骤S13，根据所述工区资料，利用含水合物层在测井曲线上表现为自然伽马值低于第一阈值、电阻率值高于第二阈值及声波速度值高于第三阈值的特征，划分得到含水合物层，并将其它层位划分为不含水合物层。
步骤S2：
参考图3，为本发明一实施例的密度与中子交会法估算含水合物层孔隙度的流程示意图。如图3所示，具体流程为：
步骤S21，利用含泥质修正的密度与中子交会法计算井中地层孔隙度；其中，
利用含泥质修正的密度法计算密度井中地层孔隙度计算公式为：
其中，ρma为基质密度；ρb为测井实测密度；ρsh为泥质密度；ρf为孔隙流体密度；ρh为纯水合物密度；Sh为预估的水合物饱和度；Vsh为泥质含量；
利用含泥质修正的中子测井方法计算中子井中地层孔隙度计算公式为：
其中，为基质的中子响应；/>为测井实测中子响应；/>为泥质中子响应；/>为孔隙流体中子响应；/>为纯水合物中子响应；Sh为预估的水合物饱和度；Vsh为泥质含量；
井中地层孔隙度取密度和中子井中地层孔隙度的平均值，计算公式为：
在本实施例中，地层泥质含量Vsh用以下公式计算：
其中，GRmin＝15API和GRmax＝90API分别是纯砂岩和泥岩层的自然伽马测井值；GR是实测自然伽马测井值；IGR为自然伽马相对值；C＝3.7，为Hilchie指数。
步骤S22，对于含水合物层位和不含合物层位分别计算得到井中地层孔隙度，其中，
不含水合物层位，令式中Sh＝0，井中地层孔隙度为
含水合物层位，采用迭代法计算：
设定水合物饱和度Sh＝Sh1＝0，计算得到水合物饱和度Sh1对应的井中地层孔隙度将/>代入下式反计算出饱和度值Sh2：
令Sh＝Sh2，计算新的孔隙度反复进行孔隙度计算，判断计算出的两次相邻的孔隙度是否满足目标函数：
在满足目标函数时，含水合物层位的预估水合物饱和度值为Sh(i+1)，孔隙度为
步骤S3：
步骤S31，对于含水合物层位和不含水合物层位分别计算对应的固体成分体积模量Km和剪切模量μm，计算公式如下：
对于不含水合物层位：
Km＝Kma；
μm＝μma；
其中，Kma、μma分别是基质的体积模量和剪切模量；
对于含水合物层位：
其中，Kma、μma是基质的体积模量和剪切模量；为井中地层孔隙度，取值为/>Sh为预估的水合物饱和度，取值为Sh(i+1)；Kh、μh分别是纯水合物体积模量和剪切模量；
步骤S32，设定α＝0，利用下式计算固结系数对应的Vp *2ρ*：
其中，α为固结系数；Km、μm是固体成分体积模量和剪切模量；为井中地层孔隙度，对于含水合物层位取值为/>Kf为孔隙流体体积模量；
判断Vp *2ρ*与实测Vp mesure2ρmesure(α)是否满足目标函数：
Vp *2ρ*(α)-Vp mesure2ρmesure(α)→0；
其中，Vp mesure、ρmesure分别是测井实测纵波速度和密度；
在满足目标函数时，对应的α值为当前深度点的地层最佳固结系数；若不满足，通过梯度上升法逐步搜索满足目标函数的α的地层最佳固结系数。
步骤S4：
重复上述步骤S2及S3，对井中每个深度点进行井中地层孔隙度及固结系数的计算，得到整个井段的地层孔隙度及固结系数。
步骤S5：
步骤S51，根据地层孔隙度及固结系数，利用三相Biot模型预测地层纵波速度Vp *及地层密度ρ*，计算公式如下：
其中，
K＝Kma(1-βp)+βp 2Kav；
G＝μma(1-βs)；
γ＝(1+2α)/(1+α)；
其中，Kma、μma、ρma分别是基质的体积模量、剪切模量和密度；是井中地层孔隙度；ρf是海水密度；ρh是纯水合物密度；Sh是预估水合物饱和度；α是最佳固结系数；ε取值范围为0-1，当ε＝1时适用于孔隙填充型水合物，当ε＝0时适用于颗粒支撑型水合物，根据工区实际水合物赋存状态调整取值。
步骤S52：
定义修正系数M＝Vpρw，其中，式中，w为权重因子，初始值为0。
以修正系数为约束，利用地层纵波速度Vp *及地层密度ρ*，通过调整Sh的值，使预测的修正系数M*＝Vp *(Sh)·ρ(Sh)w与测井实测系数Mmeaure＝Vp mesureρmesurew的误差满足目标函数：
M*-Msesure→0；
在满足目标函数时，得到当前深度点含水合物层位及不含水合物层位的水合物饱和度值；
其中，Sh取值范围为0-1；w为权重因子，初始值为0；Vp *(Sh)、ρ*(Sh)分别为在预估水合物饱和度Sh下，由三相Biot模型预测得到的纵波速度及密度；Vp mesure、ρmesure分别为测井实测纵波速度及密度。
步骤S6：
修改权重因子w的值，使(由步骤S5计算得到的)不含水合物层位的饱和度平均值满足目标函数：
Sh no-hydrate(w*)→0；
在满足目标函数时，得到本井的权重因子w*；其中，Sh no-hydrate为不含水合物层位(饱水层段)的饱和度平均值。
步骤S7：
利用所述本井的权重因子w*，(利用步骤S5)重新计算井中每个深度点的水合物饱和度，得到整个井段的水合物饱和度值。
需要说明的是，尽管在上述实施例及附图中以特定顺序描述了本发明方法的操作，但是，这并非要求或者暗示必须按照该特定顺序来执行这些操作，或是必须执行全部所示的操作才能实现期望的结果。附加地或备选地，可以省略某些步骤，将多个步骤合并为一个步骤执行，和/或将一个步骤分解为多个步骤执行。
为了对上述井中天然气水合物饱和度岩石物理计算方法进行更为清楚的解释，下面参考图1所述的步骤结合一个具体的实施例来进行说明，然而值得注意的是该实施例仅是为了更好地说明本发明，并不构成对本发明不当的限定。
步骤S1，采集研究区基本数据。
研究区为南海某含水合物井位，通过岩心分析，获得此研究区地层矿物组分为方解石、石英、长石、云母和粘土，体积分数分别为14％、28％、12％、26％、20％，其各自弹性模量、密度、中子响应值取下表经验值。进一步，利用步骤S12中的计算公式，计算得到岩石基质模量、密度、中子响应值。
固体成分的体积模量、剪切模量、密度及中子响应值，计算公式如下：
ρma＝∑fiρi；
其中，Kma、μma、ρma、分别为基质的体积模量、剪切模量、密度、中子响应值；Ki、μi、ρi、/>fi分别为每种矿物颗粒的体积模量、剪切模量、密度、中子响应值、体积百分比。
表1弹性模量及中子响应值参数表
如表1所示，为弹性模量及中子响应值参数表。其中，研究区测井曲线为海底下50m-250m，包含自然伽马、密度、中子孔隙度、纵波速度、电阻率曲线。基于各测井曲线特征划分含水合物层与不含水合物层。
在一般含水合物层在测井曲线上表现为低自然伽马值、高电阻率值、高声波速度值。参考图4，为南海某井测井曲线及预估的含水合物层的示意图。如图4所示，图中标记框401为划分的含水合物层。据此划分[190m，220m]段为含水合物储层，其它层段为不含水合物层。
步骤S2，使用密度-中子交会法进行地层孔隙度估算。
结合步骤S21，利用含泥质修正的密度法计算密度井中地层孔隙度计算公式为：
其中，ρma为基质密度，由步骤S1得到；ρb为测井实测密度；ρsh＝2.58g/cm3为泥质密度；ρf＝1.032为孔隙流体密度；ρh为纯水合物密度；Sh为预估的水合物饱和度；Vsh为泥质含量。
利用含泥质修正的中子测井方法计算中子井中地层孔隙度计算公式为：
其中，为基质的中子响应；/>为测井实测中子响应；/>为泥质中子响应；为孔隙流体中子响应；/>为纯水合物中子响应；Sh为预估的水合物饱和度；Vsh为泥质含量。
地层泥质含量Vsh用以下公式计算：
其中，GRmin＝15API和GRmax＝90API分别是纯砂岩和泥岩层的自然伽马测井值；GR是实测自然伽马测井值；IGR为自然伽马相对值；C＝3.7。
结合步骤S22，对于含水合物层位和不含合物层位分别计算得到井中地层孔隙度，其中，
不含水合物层位，令式中Sh＝0，井中地层孔隙度为
含水合物层位，采用迭代法计算：
设定水合物饱和度Sh＝Sh1＝0，计算得到水合物饱和度Sh1对应的井中地层孔隙度将/>代入下式反计算出饱和度值Sh2：
令Sh＝Sh2，计算新的孔隙度反复进行孔隙度计算，判断计算出的两次相邻的孔隙度是否满足目标函数：
在满足目标函数时，含水合物层位的预估水合物饱和度值为Sh(i+1)，孔隙度为
步骤S3，反演得到地层最佳固结系数。
结合步骤S31，对于含水合物层位和不含水合物层位分别计算对应的固体成分体积模量Km和剪切模量μm，计算公式如下：
对于不含水合物层位：
Km＝Kma；
μm＝μma；
其中，Kma、μma分别是基质的体积模量和剪切模量；
对于含水合物层位：
其中，Kma、μma是基质的体积模量和剪切模量；为井中地层孔隙度，取值为/>Sh为预估的水合物饱和度，取值为Sh(i+1)；Kh、μh分别是纯水合物体积模量和剪切模量；
结合步骤S32，设定α＝0，利用下式计算固结系数对应的Vp *2ρ*：
其中，α为固结系数；Km、μm是固体成分体积模量和剪切模量；为井中地层孔隙度，对于含水合物层位取值为/>Kf为孔隙流体体积模量；
判断Vp *2ρ*与实测Vp mesure2ρmesure(α)是否满足目标函数：
Vp *2ρ*(α)-Vp mesure2ρmesure(α)→0；
其中，Vp mesure、ρmesure分别是测井实测纵波速度和密度；
在满足目标函数时，对应的α值为当前深度点的地层最佳固结系数。若不满足，通过梯度上升法逐步搜索满足目标函数的α的地层最佳固结系数。参考图5，为本发明一具体实施例的井内地层固结系数散点图，如图5所示，可以看出固结系数的散点分布。
步骤S4，重复上述步骤S2及S3，对井中每个深度点进行井中地层孔隙度及固结系数的计算，得到整个井段的地层孔隙度及固结系数。
步骤S5：将步骤S4得到的地层孔隙度值、固结系数值，代入简化的三相Biot方程中预测地层纵波速度Vp *及密度ρ*，并以修正系数为约束计算井中每个深度点的水合物饱和度。
步骤S51，根据地层孔隙度及固结系数，利用三相Biot模型预测地层纵波速度Vp *及地层密度ρ*，计算公式如下：
其中，
K＝Kma(1-βp)+βp 2Kav；
G＝μma(1-βs)；
γ＝(1+2α)/(1+α)；
其中，Kma、μma、ρma分别是基质的体积模量、剪切模量和密度；是井中地层孔隙度；ρf是海水密度；ρh是纯水合物密度；Sh是预估水合物饱和度；α是最佳固结系数；ε取值范围为0-1，当ε＝1时适用于孔隙填充型水合物，当ε＝0时适用于颗粒支撑型水合物，根据工区实际水合物赋存状态调整取值，本研究区ε取值为0.12。
步骤S52：
定义修正系数M＝Vpρw，其中，式中，w为权重因子，初始值为0。
以修正系数为约束，利用地层纵波速度Vp *及地层密度ρ*，通过调整Sh的值，使预测的修正系数M*＝Vp *(Sh)·ρ(Sh)w与测井实测系数Mmeaure＝Vp mesureρmesurew的误差满足目标函数：
M*-Msesure→0；
在满足目标函数时，得到当前深度点含水合物层位及不含水合物层位的水合物饱和度值；
其中，Sh取值范围为0-1；w为权重因子，初始值为0；Vp *(Sh)、ρ*(Sh)分别为在预估水合物饱和度Sh下，由三相Biot模型预测得到的纵波速度及密度；Vp mesure、ρmesure分别为测井实测纵波速度及密度。
步骤S6，修改并得到适合本井的权重因子。
修改权重因子w的值，使(由步骤S5计算得到的)不含水合物层位的饱和度平均值满足目标函数：
Sh no-hydrate(w*)→0；
在满足目标函数时，得到本井的权重因子w*，w*＝1.2；其中，Sh no-hydrate为不含水合物层位的饱和度平均值。
参考图6，为权重因子为0时南海某井模型预测速度基线与实测速度对比、饱和度估算结果与氯离子法实测饱和度对比示意图。
参考图7，为权重因子为1.2时南海某井模型预测速度基线与实测速度对比、饱和度估算结果与氯离子法实测饱和度对比示意图。
经过对比图6及图7可知，适合本井的权重因子为w*＝1.2。
步骤S7，计算修正后的饱和度值。
利用所述本井的权重因子w*＝1.2，(利用步骤S5)重新计算井中每个深度点的水合物饱和度，得到整个井段的水合物饱和度值。
本发明提出的井中天然气水合物饱和度岩石物理计算方法及系统，将基于岩石物理模型的井内水合物饱和度估算理论进行改进，并利用提出的修正系数对估算结果进行修正，估算精度更高；为井内水合物饱和度精确估算提供了一种可行性高的方案，对水合物的精确定量和加快实际生产开发有积极的推进作用。
在介绍了本发明示例性实施方式的方法之后，接下来，参考图8对本发明示例性实施方式的井中天然气水合物饱和度岩石物理计算系统进行介绍。
井中天然气水合物饱和度岩石物理计算系统的实施可以参见上述方法的实施，重复之处不再赘述。以下所使用的术语“模块”或者“单元”，可以是实现预定功能的软件和/或硬件的组合。尽管以下实施例所描述的装置较佳地以软件来实现，但是硬件，或者软件和硬件的组合的实现也是可能并被构想的。
基于同一发明构思，本发明还提出了一种井中天然气水合物饱和度岩石物理计算系统，如图8所示，该系统包括：
工区资料获取模块810，用于获取工区资料；
孔隙度计算模块820，用于根据所述工区资料，分别对当前深度点的含水合物层位及不含水合物层位进行井中地层孔隙度计算，得到井中地层孔隙度；
固结系数计算模块830，用于根据所述工区资料及井中地层孔隙度，经过反演得到当前深度点的地层最佳固结系数；
井段数据处理模块840，用于对井中每个深度点进行井中地层孔隙度及固结系数的计算，得到整个井段的地层孔隙度及固结系数；
饱和度计算模块850，用于根据地层孔隙度及固结系数，利用三相Biot模型预测地层纵波速度及地层密度，并以修正系数为约束计算井中每个深度点的水合物饱和度；
权重因子调整模块860，用于通过修改权重因子调整修正系数，利用新的修正系数为约束计算井中不含水合物层位的饱和度平均值，直至所述不含水合物层位的饱和度平均值趋近于0时，得到本井的权重因子；
所述饱和度计算模块850，还用于利用所述本井的权重因子，重新计算井中每个深度点的水合物饱和度，得到整个井段的水合物饱和度值。
应当注意，尽管在上文详细描述中提及了井中天然气水合物饱和度岩石物理计算系统的若干模块，但是这种划分仅仅是示例性的并非强制性的。实际上，根据本发明的实施方式，上文描述的两个或更多模块的特征和功能可以在一个模块中具体化。反之，上文描述的一个模块的特征和功能可以进一步划分为由多个模块来具体化。
基于前述发明构思，如图9所示，本发明还提出了一种计算机设备900，包括存储器910、处理器920及存储在存储器910上并可在处理器920上运行的计算机程序930，所述处理器920执行所述计算机程序930时实现前述井中天然气水合物饱和度岩石物理计算方法。
基于前述发明构思，本发明提出了一种计算机可读存储介质，所述计算机可读存储介质存储有计算机程序，所述计算机程序被处理器执行时实现前述井中天然气水合物饱和度岩石物理计算方法。
本发明提出的井中天然气水合物饱和度岩石物理计算方法及系统，将基于岩石物理模型的井内水合物饱和度估算理论进行改进，并利用提出的修正系数对估算结果进行修正，估算精度更高；为井内水合物饱和度精确估算提供了一种可行性高的方案，对水合物的精确定量和加快实际生产开发有积极的推进作用。
本领域内的技术人员应明白，本发明的实施例可提供为方法、系统、或计算机程序产品。因此，本发明可采用完全硬件实施例、完全软件实施例、或结合软件和硬件方面的实施例的形式。而且，本发明可采用在一个或多个其中包含有计算机可用程序代码的计算机可用存储介质(包括但不限于磁盘存储器、CD-ROM、光学存储器等)上实施的计算机程序产品的形式。
本发明是参照根据本发明实施例的方法和计算机程序产品的流程图和/或方框图来描述的。应理解可由计算机程序指令实现流程图和/或方框图中的每一流程和/或方框、以及流程图和/或方框图中的流程和/或方框的结合。可提供这些计算机程序指令到通用计算机、专用计算机、嵌入式处理机或其他可编程数据处理设备的处理器以产生一个机器，使得通过计算机或其他可编程数据处理设备的处理器执行的指令产生用于实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能的装置。
这些计算机程序指令也可存储在能引导计算机或其他可编程数据处理设备以特定方式工作的计算机可读存储器中，使得存储在该计算机可读存储器中的指令产生包括指令装置的制造品，该指令装置实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能。
这些计算机程序指令也可装载到计算机或其他可编程数据处理设备上，使得在计算机或其他可编程设备上执行一系列操作步骤以产生计算机实现的处理，从而在计算机或其他可编程设备上执行的指令提供用于实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能的步骤。
最后应说明的是：以上所述实施例，仅为本发明的具体实施方式，用以说明本发明的技术方案，而非对其限制，本发明的保护范围并不局限于此，尽管参照前述实施例对本发明进行了详细的说明，本领域的普通技术人员应当理解：任何熟悉本技术领域的技术人员在本发明揭露的技术范围内，其依然可以对前述实施例所记载的技术方案进行修改或可轻易想到变化，或者对其中部分技术特征进行等同替换；而这些修改、变化或者替换，并不使相应技术方案的本质脱离本发明实施例技术方案的精神和范围，都应涵盖在本发明的保护范围之内。因此，本发明的保护范围应所述以权利要求的保护范围为准。
Claims (7)
1.一种井中天然气水合物饱和度岩石物理计算方法，其特征在于，该方法包括：
获取工区资料；
根据所述工区资料，分别对当前深度点的含水合物层位及不含水合物层位进行井中地层孔隙度计算，得到井中地层孔隙度；
根据所述工区资料及井中地层孔隙度，经过反演得到当前深度点的地层最佳固结系数；
对井中每个深度点进行井中地层孔隙度及固结系数的计算，得到整个井段的地层孔隙度及固结系数；
根据地层孔隙度及固结系数，利用三相Biot模型预测地层纵波速度及地层密度，并以修正系数为约束计算井中每个深度点的水合物饱和度；
通过修改权重因子调整修正系数，利用新的修正系数为约束计算井中不含水合物层位的饱和度平均值，直至所述不含水合物层位的饱和度平均值趋近于0时，得到本井的权重因子；
利用所述本井的权重因子，重新计算井中每个深度点的水合物饱和度，得到整个井段的水合物饱和度值；
其中，根据所述工区资料，分别对当前深度点的含水合物层位及不含水合物层位进行井中地层孔隙度计算，得到井中地层孔隙度，包括：
利用含泥质修正的密度与中子交会法计算井中地层孔隙度；其中，
利用含泥质修正的密度法计算密度井中地层孔隙度计算公式为：
其中，ρma为基质密度；ρb为测井实测密度；ρsh为泥质密度；ρf为孔隙流体密度；ρh为纯水合物密度；Sh为预估的水合物饱和度；Vsh为泥质含量；
利用含泥质修正的中子测井方法计算中子井中地层孔隙度计算公式为：
其中，为基质的中子响应；/>为测井实测中子响应；/>为泥质中子响应；/>为孔隙流体中子响应；/>为纯水合物中子响应；Sh为预估的水合物饱和度；Vsh为泥质含量；
井中地层孔隙度取密度和中子井中地层孔隙度的平均值，计算公式为：
其中，地层泥质含量Vsh用以下公式计算：
其中，GRmin＝15API和GRmax＝90API分别是纯砂岩和泥岩层的自然伽马测井值；GR是实测自然伽马测井值；IGR为自然伽马相对值；C为Hilchie指数，取值为3.7；
对于含水合物层位和不含合物层位分别计算得到井中地层孔隙度，其中，
不含水合物层位，令式中Sh＝0，井中地层孔隙度为
含水合物层位，采用迭代法计算：
设定水合物饱和度Sh＝Sh1＝0，计算得到水合物饱和度Sh1对应的井中地层孔隙度将代入下式反计算出饱和度值Sh2：
令Sh＝Sh2，计算新的孔隙度反复进行孔隙度计算，判断计算出的两次相邻的孔隙度是否满足目标函数：
在满足目标函数时，含水合物层位的预估水合物饱和度值为Sh(i+1)，孔隙度为
其中，根据所述工区资料及井中地层孔隙度，经过反演得到当前深度点的地层最佳固结系数，包括：
对于含水合物层位和不含水合物层位分别计算对应的固体成分体积模量Km和剪切模量μm，计算公式如下：
对于不含水合物层位：
Km＝Kma；
μm＝μma；
其中，Kma、μma分别是基质的体积模量和剪切模量；
对于含水合物层位：
其中，Kma、μma是基质的体积模量和剪切模量；为井中地层孔隙度，取值为/>Sh为预估的水合物饱和度，取值为Sh(i+1)；Kh、μh分别是纯水合物体积模量和剪切模量；
设定α＝0，利用下式计算固结系数对应的Vp *2ρ*：
其中，α为固结系数；Km、μm是固体成分体积模量和剪切模量；为井中地层孔隙度，对于含水合物层位取值为/>Kf为孔隙流体体积模量；
判断Vp *2ρ*与实测Vp mesure2ρmesure(α)是否满足目标函数：
Vp *2ρ*(α)-Vp mesure2ρmesure(α)→0；
其中，Vp mesure、ρmesure分别是测井实测纵波速度和密度；
在满足目标函数时，对应的α值为当前深度点的地层最佳固结系数；若不满足，通过梯度上升法逐步搜索满足目标函数的α的地层最佳固结系数；
其中，根据地层孔隙度及固结系数，利用三相Biot模型预测地层纵波速度及地层密度，并以修正系数为约束计算井中每个深度点的水合物饱和度，包括：
根据地层孔隙度及固结系数，利用三相Biot模型预测地层纵波速度Vp *及地层密度ρ*，计算公式如下：
其中，
K＝Kma(1-βp)+βp 2Kav；
G＝μma(1-βs)；
γ＝(1+2α)/(1+α)；
其中，Kma、μma、ρma分别是基质的体积模量、剪切模量和密度；是井中地层孔隙度；ρf是海水密度；ρh是纯水合物密度；Sh是预估水合物饱和度；α是最佳固结系数；ε取值范围为0-1，当ε＝1时适用于孔隙填充型水合物，当ε＝0时适用于颗粒支撑型水合物，根据工区实际水合物赋存状态调整取值；
以修正系数为约束，利用地层纵波速度Vp *及地层密度ρ*，通过调整Sh的值，使预测的修正系数M*＝Vp *(Sh)·ρ(Sh)w与测井实测系数Mmeaure＝Vp mesureρmesurew的误差满足目标函数：
M*-Msesure→0；
在满足目标函数时，得到当前深度点含水合物层位及不含水合物层位的水合物饱和度值；
其中，Sh取值范围为0-1；w为权重因子，初始值为0；Vp *(Sh)、ρ*(Sh)分别为在预估水合物饱和度Sh下，由三相Biot模型预测得到的纵波速度及密度；Vp mesure、ρmesure分别为测井实测纵波速度及密度。
2.根据权利要求1所述的井中天然气水合物饱和度岩石物理计算方法，其特征在于，获取工区资料，包括：
通过工区测井和岩心分析，获得测井数据及岩石矿物数据；其中，
所述测井数据包括：研究区的深度、纵波速度、电阻率、密度、自然伽马及声波速度数据；
所述岩石矿物数据包括：井中地层矿物颗粒的体积模量、剪切模量、体积百分比、密度及中子响应值；
根据岩石矿物数据，得到固体成分的体积模量、剪切模量、密度及中子响应值，计算公式如下：
ρma＝∑fiρi；
其中，Kma、μma、ρma、分别为基质的体积模量、剪切模量、密度、中子响应值；Ki、μi、ρi、fi分别为每种矿物颗粒的体积模量、剪切模量、密度、中子响应值、体积百分比。
3.根据权利要求2所述的井中天然气水合物饱和度岩石物理计算方法，其特征在于，该方法还包括：
根据所述工区资料，利用含水合物层在测井曲线上表现为自然伽马值低于第一阈值、电阻率值高于第二阈值及声波速度值高于第三阈值的特征，划分得到含水合物层，并将其它层位划分为不含水合物层。
4.根据权利要求3所述的井中天然气水合物饱和度岩石物理计算方法，其特征在于，通过修改权重因子调整修正系数，利用新的修正系数为约束计算井中不含水合物层位的饱和度平均值，直至所述不含水合物层位的饱和度平均值趋近于0时，得到本井的权重因子，包括：
修改权重因子w的值，使不含水合物层位的饱和度平均值满足目标函数：
Sh no-hydrate(w*)→0；
在满足目标函数时，得到本井的权重因子w*；其中，Sh no-hydrate为不含水合物层位的饱和度平均值。
5.一种井中天然气水合物饱和度岩石物理计算系统，其特征在于，该系统包括：
工区资料获取模块，用于获取工区资料；
孔隙度计算模块，用于根据所述工区资料，分别对当前深度点的含水合物层位及不含水合物层位进行井中地层孔隙度计算，得到井中地层孔隙度；
固结系数计算模块，用于根据所述工区资料及井中地层孔隙度，经过反演得到当前深度点的地层最佳固结系数；
井段数据处理模块，用于对井中每个深度点进行井中地层孔隙度及固结系数的计算，得到整个井段的地层孔隙度及固结系数；
饱和度计算模块，用于根据地层孔隙度及固结系数，利用三相Biot模型预测地层纵波速度及地层密度，并以修正系数为约束计算井中每个深度点的水合物饱和度；
权重因子调整模块，用于通过修改权重因子调整修正系数，利用新的修正系数为约束计算井中不含水合物层位的饱和度平均值，直至所述不含水合物层位的饱和度平均值趋近于0时，得到本井的权重因子；
所述饱和度计算模块，还用于利用所述本井的权重因子，重新计算井中每个深度点的水合物饱和度，得到整个井段的水合物饱和度值；
其中，孔隙度计算模块，用于根据所述工区资料，分别对当前深度点的含水合物层位及不含水合物层位进行井中地层孔隙度计算，得到井中地层孔隙度，包括：
利用含泥质修正的密度与中子交会法计算井中地层孔隙度；其中，
利用含泥质修正的密度法计算密度井中地层孔隙度计算公式为：
其中，ρma为基质密度；ρb为测井实测密度；ρsh为泥质密度；ρf为孔隙流体密度；ρh为纯水合物密度；Sh为预估的水合物饱和度；Vsh为泥质含量；
利用含泥质修正的中子测井方法计算中子井中地层孔隙度计算公式为：
其中，为基质的中子响应；/>为测井实测中子响应；/>为泥质中子响应；/>为孔隙流体中子响应；/>为纯水合物中子响应；Sh为预估的水合物饱和度；Vsh为泥质含量；
井中地层孔隙度取密度和中子井中地层孔隙度的平均值，计算公式为：
其中，地层泥质含量Vsh用以下公式计算：
其中，GRmin＝15API和GRmax＝90API分别是纯砂岩和泥岩层的自然伽马测井值；GR是实测自然伽马测井值；IGR为自然伽马相对值；C为Hilchie指数，取值为3.7；
对于含水合物层位和不含合物层位分别计算得到井中地层孔隙度，其中，
不含水合物层位，令式中Sh＝0，井中地层孔隙度为
含水合物层位，采用迭代法计算：
设定水合物饱和度Sh＝Sh1＝0，计算得到水合物饱和度Sh1对应的井中地层孔隙度将代入下式反计算出饱和度值Sh2：
令Sh＝Sh2，计算新的孔隙度反复进行孔隙度计算，判断计算出的两次相邻的孔隙度是否满足目标函数：
在满足目标函数时，含水合物层位的预估水合物饱和度值为Sh(i+1)，孔隙度为
其中，固结系数计算模块，用于根据所述工区资料及井中地层孔隙度，经过反演得到当前深度点的地层最佳固结系数，包括：
对于含水合物层位和不含水合物层位分别计算对应的固体成分体积模量Km和剪切模量μm，计算公式如下：
对于不含水合物层位：
Km＝Kma；
μm＝μma；
其中，Kma、μma分别是基质的体积模量和剪切模量；
对于含水合物层位：
其中，Kma、μma是基质的体积模量和剪切模量；为井中地层孔隙度，取值为/>Sh为预估的水合物饱和度，取值为Sh(i+1)；Kh、μh分别是纯水合物体积模量和剪切模量；
设定α＝0，利用下式计算固结系数对应的Vp *2ρ*：
其中，α为固结系数；Km、μm是固体成分体积模量和剪切模量；为井中地层孔隙度，对于含水合物层位取值为/>Kf为孔隙流体体积模量；
判断Vp *2ρ*与实测Vp mesure2ρmesure(α)是否满足目标函数：
Vp *2ρ*(α)-Vp mesure2ρmesure(α)→0；
其中，Vp mesure、ρmesure分别是测井实测纵波速度和密度；
在满足目标函数时，对应的α值为当前深度点的地层最佳固结系数；若不满足，通过梯度上升法逐步搜索满足目标函数的α的地层最佳固结系数；
其中，饱和度计算模块，用于根据地层孔隙度及固结系数，利用三相Biot模型预测地层纵波速度及地层密度，并以修正系数为约束计算井中每个深度点的水合物饱和度，包括：
根据地层孔隙度及固结系数，利用三相Biot模型预测地层纵波速度Vp *及地层密度ρ*，计算公式如下：
其中，
K＝Kma(1-βp)+βp 2Kav；
G＝μma(1-βs)；
γ＝(1+2α)/(1+α)；
其中，Kma、μma、ρma分别是基质的体积模量、剪切模量和密度；是井中地层孔隙度；ρf是海水密度；ρh是纯水合物密度；Sh是预估水合物饱和度；α是最佳固结系数；ε取值范围为0-1，当ε＝1时适用于孔隙填充型水合物，当ε＝0时适用于颗粒支撑型水合物，根据工区实际水合物赋存状态调整取值；
以修正系数为约束，利用地层纵波速度Vp *及地层密度ρ*，通过调整Sh的值，使预测的修正系数M*＝Vp *(Sh)·ρ(Sh)w与测井实测系数Mmeaure＝Vp mesureρmesurew的误差满足目标函数：
M*-Msesure→0；
在满足目标函数时，得到当前深度点含水合物层位及不含水合物层位的水合物饱和度值；
其中，Sh取值范围为0-1；w为权重因子，初始值为0；Vp *(Sh)、ρ*(Sh)分别为在预估水合物饱和度Sh下，由三相Biot模型预测得到的纵波速度及密度；Vp mesure、ρmesure分别为测井实测纵波速度及密度。
6.一种计算机设备，包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序，其特征在于，所述处理器执行所述计算机程序时实现权利要求1至4任一所述方法。
7.一种计算机可读存储介质，其特征在于，所述计算机可读存储介质存储有计算机程序，所述计算机程序被处理器执行时实现权利要求1至4任一所述方法。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202010889675.1A CN112149282B (zh) | 2020-08-28 | 2020-08-28 | 一种井中天然气水合物饱和度岩石物理计算方法及系统 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202010889675.1A CN112149282B (zh) | 2020-08-28 | 2020-08-28 | 一种井中天然气水合物饱和度岩石物理计算方法及系统 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN112149282A CN112149282A (zh) | 2020-12-29 |
| CN112149282B true CN112149282B (zh) | 2024-05-28 |
Family
ID=73890566
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202010889675.1A Active CN112149282B (zh) | 2020-08-28 | 2020-08-28 | 一种井中天然气水合物饱和度岩石物理计算方法及系统 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN112149282B (zh) |
Families Citing this family (8)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN112946737B (zh) * | 2021-01-20 | 2023-10-31 | 中国地质大学（北京） | 一种利用纵横波速度增量交汇图识别天然气水合物的方法 |
| CN112946783B (zh) * | 2021-01-29 | 2024-02-06 | 中国石油大学(北京) | 一种水合物饱和度确定方法、装置及设备 |
| CN113189674B (zh) * | 2021-05-07 | 2024-04-23 | 中国石油大学(北京) | 一种天然气水合物的饱和度估算方法及系统 |
| CN113341464B (zh) * | 2021-06-04 | 2024-01-26 | 中国石油大学(北京) | 天然气水合物储层的识别方法、装置、设备和存储介质 |
| CN113376709B (zh) * | 2021-06-21 | 2022-05-13 | 西南石油大学 | 一种利用测井数据预测储层天然气水合物饱和度的方法 |
| CN113466963B (zh) * | 2021-08-13 | 2022-02-01 | 广州海洋地质调查局 | 裂缝充填ii型水合物饱和度估算方法及处理终端 |
| CN114783531B (zh) * | 2022-04-02 | 2022-11-22 | 广州海洋地质调查局 | 计算砂岩型水合物中胶结型和孔隙充填型饱和度的方法 |
| CN114879259B (zh) * | 2022-04-15 | 2023-03-31 | 中国科学院海洋研究所 | 裂隙充填型天然气水合物层的各向异性地震饱和度预测方法 |
Citations (8)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| SU1056109A1 (ru) * | 1981-10-14 | 1983-11-23 | Саратовский Ордена Трудового Красного Знамени Государственный Университет Им.Н.Г.Чернышевского | Способ определени коэффициента нефтегазонасыщени пород |
| CN104134101A (zh) * | 2014-07-23 | 2014-11-05 | 中国石油集团川庆钻探工程有限公司 | 低渗透储层天然气产能预测方法 |
| CN106680878A (zh) * | 2016-12-29 | 2017-05-17 | 中国石油大学(华东) | 一种基于改进Biot系数分析的横波速度估算方法 |
| CN107290782A (zh) * | 2016-03-30 | 2017-10-24 | 中国石油化工股份有限公司 | 储层孔隙度、含水饱和度和泥质含量参数同时反演新方法 |
| CN109212162A (zh) * | 2018-08-17 | 2019-01-15 | 西南石油大学 | 一种估算成岩天然气水合物储层力学参数的方法 |
| CN110927256A (zh) * | 2019-04-22 | 2020-03-27 | 广东石油化工学院 | 一种基于Wood波阻抗法的天然气水合物饱和度计算方法 |
| CN111220522A (zh) * | 2019-10-28 | 2020-06-02 | 中国石油大学(华东) | 一种高泥质细粒沉积物水合物饱和度的岩心刻度测井计算方法 |
| CN111368245A (zh) * | 2020-03-12 | 2020-07-03 | 成都理工大学 | 一种考虑多影响因素的定量计算岩石Biot系数方法 |
-
2020
- 2020-08-28 CN CN202010889675.1A patent/CN112149282B/zh active Active
Patent Citations (8)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| SU1056109A1 (ru) * | 1981-10-14 | 1983-11-23 | Саратовский Ордена Трудового Красного Знамени Государственный Университет Им.Н.Г.Чернышевского | Способ определени коэффициента нефтегазонасыщени пород |
| CN104134101A (zh) * | 2014-07-23 | 2014-11-05 | 中国石油集团川庆钻探工程有限公司 | 低渗透储层天然气产能预测方法 |
| CN107290782A (zh) * | 2016-03-30 | 2017-10-24 | 中国石油化工股份有限公司 | 储层孔隙度、含水饱和度和泥质含量参数同时反演新方法 |
| CN106680878A (zh) * | 2016-12-29 | 2017-05-17 | 中国石油大学(华东) | 一种基于改进Biot系数分析的横波速度估算方法 |
| CN109212162A (zh) * | 2018-08-17 | 2019-01-15 | 西南石油大学 | 一种估算成岩天然气水合物储层力学参数的方法 |
| CN110927256A (zh) * | 2019-04-22 | 2020-03-27 | 广东石油化工学院 | 一种基于Wood波阻抗法的天然气水合物饱和度计算方法 |
| CN111220522A (zh) * | 2019-10-28 | 2020-06-02 | 中国石油大学(华东) | 一种高泥质细粒沉积物水合物饱和度的岩心刻度测井计算方法 |
| CN111368245A (zh) * | 2020-03-12 | 2020-07-03 | 成都理工大学 | 一种考虑多影响因素的定量计算岩石Biot系数方法 |
Non-Patent Citations (2)
| Title |
|---|
| 刘洁 ; 张建中 ; 孙运宝 ; 赵铁虎 ; .南海神狐海域天然气水合物储层参数测井评价.天然气地球科学.(01),168-176. * |
| 应用改进的Biot-Gassmann模型估算天然气水合物的饱和度;高红艳，等;海洋地质与第四纪地质;20120831;第32卷(第4期);83-89 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN112149282A (zh) | 2020-12-29 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN112149282B (zh) | 一种井中天然气水合物饱和度岩石物理计算方法及系统 | |
| Keys et al. | An approximation for the Xu-White velocity model | |
| US10605951B2 (en) | Method of determining CEC and other properties from multi-frequency dielectric measurements | |
| Harris et al. | A two-dimensional, time-dependent model of suspended sediment transport and bed reworking for continental shelves | |
| CN103954999B (zh) | 一种适用于低孔隙度砂泥岩地层的横波速度预测方法 | |
| CN112133377B (zh) | 一种海域天然气水合物赋存状态的判别方法及系统 | |
| CN112946783B (zh) | 一种水合物饱和度确定方法、装置及设备 | |
| CN112182840A (zh) | 一种天然气水合物地层衰减建模方法与装置 | |
| AU2012390305B2 (en) | A dual porosity methodology for mineral volume calculations in source rock | |
| CN113189674B (zh) | 一种天然气水合物的饱和度估算方法及系统 | |
| CN113466963B (zh) | 裂缝充填ii型水合物饱和度估算方法及处理终端 | |
| CN114280670A (zh) | 多物理模型测井横波速度曲线重构方法、系统及电子设备 | |
| CN109655936B (zh) | 一种碎屑岩岩性替换的弹性参数计算方法及系统 | |
| CN114114409A (zh) | 海域天然气水合物的岩石物理建模方法、电子设备及介质 | |
| CN105301642A (zh) | 非均匀孔隙岩石及其固态有机质体积含量确定方法及装置 | |
| CN110927787A (zh) | 基于虚拟岩性的横波速度计算方法及系统 | |
| CN116930023A (zh) | 致密砂岩相控分类孔隙度测井精细解释方法和装置 | |
| Spichak et al. | Electromagnetic prediction of rock thermal properties beyond boreholes: Soultz-sous-Forets (France) case study | |
| Glazunov et al. | Large-eddy simulation of interaction of ocean and atmospheric boundary layers | |
| CN113917531B (zh) | 基于扰动体积模量的压力系数预测方法及电子设备 | |
| CN113189645B (zh) | 一种基质矿物模量确定方法、装置、电子设备和存储介质 | |
| CN113791457B (zh) | 一种天然气水合物储层岩石骨架模量的计算方法及装置 | |
| Pan et al. | Joint Interpretation of elastic and electrical data for petrophysical properties of gas-hydrate-bearing sediments using inverse rock physics modeling method | |
| CN114076987B (zh) | 基于致密砂岩油气藏有效储层的测井识别方法及装置 | |
| CN114076983B (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 |