CN116502557B - 基于网格分析法的储层压裂后可动流体饱和度计算方法 - Google Patents

基于网格分析法的储层压裂后可动流体饱和度计算方法 Download PDF

Info

Publication number
CN116502557B
CN116502557B CN202310485930.XA CN202310485930A CN116502557B CN 116502557 B CN116502557 B CN 116502557B CN 202310485930 A CN202310485930 A CN 202310485930A CN 116502557 B CN116502557 B CN 116502557B
Authority
CN
China
Prior art keywords
fracturing
scanning
movable fluid
core
porosity
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
Application number
CN202310485930.XA
Other languages
English (en)
Other versions
CN116502557A (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.)
Southwest Petroleum University
Original Assignee
Southwest Petroleum University
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 Southwest Petroleum University filed Critical Southwest Petroleum University
Priority to CN202310485930.XA priority Critical patent/CN116502557B/zh
Publication of CN116502557A publication Critical patent/CN116502557A/zh
Application granted granted Critical
Publication of CN116502557B publication Critical patent/CN116502557B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00Computer-aided design [CAD]
    • G06F30/20Design optimisation, verification or simulation
    • G06F30/28Design optimisation, verification or simulation using fluid dynamics, e.g. using Navier-Stokes equations or computational fluid dynamics [CFD]
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00Image analysis
    • G06T7/0002Inspection of images, e.g. flaw detection
    • G06T7/0004Industrial image inspection
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F2113/00Details relating to the application field
    • G06F2113/08Fluids
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F2119/00Details relating to the type or aim of the analysis or the optimisation
    • G06F2119/14Force analysis or force optimisation, e.g. static or dynamic forces
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/10Image acquisition modality
    • G06T2207/10072Tomographic images
    • G06T2207/10081Computed x-ray tomography [CT]
    • YGENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02ATECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
    • Y02A90/00Technologies having an indirect contribution to adaptation to climate change
    • Y02A90/30Assessment of water resources

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Mathematical Physics (AREA)
  • Evolutionary Computation (AREA)
  • Mathematical Analysis (AREA)
  • Mathematical Optimization (AREA)
  • Computing Systems (AREA)
  • Pure & Applied Mathematics (AREA)
  • Computer Hardware Design (AREA)
  • Fluid Mechanics (AREA)
  • Geometry (AREA)
  • General Engineering & Computer Science (AREA)
  • Algebra (AREA)
  • Quality & Reliability (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Analysing Materials By The Use Of Radiation (AREA)
  • Geophysics And Detection Of Objects (AREA)

Abstract

本发明公开了一种基于网格分析法的储层压裂后可动流体饱和度计算方法,包括:压裂实验与压裂前后微米CT实验设计;划分孔隙发育区与致密区;对同一块储层岩心样品压裂前后的微米CT扫描实验图像对比,划分出压裂后形成的裂缝网络;通过网格分析法分析压裂缝扩展情况;基于压裂缝经过的网格数量与孔隙发育区所包含的网格数量,确定压裂改造孔隙度;基于压裂改造孔隙度进一步确定压裂后的可动流体饱和度;与测井曲线结合,实现压裂后可动流体饱和度连续计算。本发明基于压裂前后微米CT扫描数据,通过网格分析法定量计算压裂后可动流体饱和度,并通过测井曲线计算压裂后可动流体饱和度,使得压裂后可动流体饱和度的评价更加准确高效。

Description

基于网格分析法的储层压裂后可动流体饱和度计算方法
技术领域
本发明涉及石油及天然气勘探及开发技术领域,特别是涉及一种基于网格分析法的储层压裂后可动流体饱和度计算方法。
背景技术
国家对油气资源需求日益剧增的背景下,油气勘探开发的主战场由常规转向非常规领域,致密油气成为当今及未来非常规油气勘探开发的热点和重点。由于致密储层的孔隙结构复杂,不仅影响了油气赋存,还严重制约着油气渗流和高效开采。对此,致密储层的开发一般都需要通过水力压裂实现增产,而可动流体饱和度可以有效的评价油气藏储层的储采特征及其最终的开发。
对于压裂前可动流体饱和度主要通过高压压汞实验、恒速压汞实验、核磁共振测试、岩心相渗实验等实验进行评价。但是目前国内存在的大量的非常规储层,需要进行水力压裂施工实现增产。由于岩心在压裂之后会在内部产生多条微裂缝,进行驱替或者离心实验会进一步破坏岩石的结构,导致开展常规实验得到的压裂后可动流体饱和度的实验结果与实际情况不符。所以评价压裂前可动流体饱和度的一些常规方法并不适用于对压裂后可动流体饱和度进行评价。
整体而言,目前对于压裂后流体可动性的评价还没有一套较为成熟的方法,因此开发一种储层压裂后可动流体饱和度的计算方法,在致密油气藏的勘探与开发领域具有重要的意义及广阔的前景。
发明内容
本发明的目的是针对致密油气储层勘探开发过程中,目前还没有有效的方法计算致密油气藏压裂后可动流体饱和度的现状,提供一种基于网格分析法的储层压裂后可动流体饱和度计算方法,以压裂前后微米CT扫描数据为基础,通过网格分析法定量计算压裂后可动流体饱和度,再与测井曲线特征结合开展分析,实现通过测井曲线计算压裂后可动流体饱和度,使得对致密油气藏压裂后可动流体饱和度的评价更加准确高效,更好的指导致密油气藏的勘探与开发。
为实现上述目的,本发明提供了如下方案:
一种基于网格分析法的储层压裂后可动流体饱和度计算方法,该方法包括以下步骤:
针对储层岩心样品,进行压裂实验与压裂前后微米CT实验设计;
针对储层岩心样品,划分孔隙发育区与致密区;
对同一块储层岩心样品压裂前后的微米CT扫描实验图像对比,划分出压裂后形成的裂缝网络;
基于裂缝网络,通过网格分析法分析压裂缝扩展情况;
基于压裂缝经过的网格数量与孔隙发育区所包含的网格数量,确定压裂改造孔隙度;
基于压裂改造孔隙度进一步确定压裂后的可动流体饱和度;
与测井曲线结合,实现压裂后可动流体饱和度连续计算。
进一步地,所述针对储层岩心样品,进行压裂实验与压裂前后微米CT实验设计,包括:
样品选取与制备,选择致密油气藏储层段岩心作为储层岩心样品,将致密油气藏储层段岩心钻切成设定大小的圆柱,圆柱的两端面要平整,垂直于轴线;
压裂前微米CT扫描实验,通过微米CT扫描实验获取压裂前岩心内部形态和内部微观结构;
压裂实验,通过室内三轴岩石力学实验模拟地层水力压裂实验,在应力下降到设定范围时,停止施压;
压裂后微米CT扫描实验,通过微米CT扫描实验获取压裂后岩心内部形态和内部微观结构。
进一步地,所述压裂后微米CT扫描实验满足以下条件:
a)使用与进行压裂前微米CT扫描实验的同一台仪器;
b)实验过程中岩心摆放方向、位置与压裂前扫描一致;
c)选择扫描分辨率与压裂前扫描分辨率一致。
进一步地,所述压裂实验,通过室内三轴岩石力学实验模拟地层水力压裂实验,在应力下降到设定范围时,停止施压,具体包括:
在开展三轴岩石力学实验时,首先用热塑膜包裹储层岩心样品,将储层岩心样品放入三轴力学实验测试仪,如果测试岩心所在深度段地层压力为Pf,在开展三轴岩石力学测试时,设置实验围压σ3=Pf
根据曲线特征将轴向应变曲线分为多个段,在OA段属于弹性变形;AB段应力应变基本属线性关系,卸载后可完全恢复;BC段,曲线偏离线性,出现塑性变形;CD段,岩石内部裂纹形成速度增快,裂纹密度加大,此时缓慢施加轴向压力,同时观察软件输出的应力应变曲线,轴向应变曲线在D点达到应力最大值,为岩石的最大承载能力,将此时的应力记为σmax,在应力下降1/6σmax~1/5σmax时,停止施压。
进一步地,所述针对储层岩心样品,划分孔隙发育区与致密区,包括:
对储层岩心样品的压裂前微米CT扫描实验图像进行优化处理;
将压裂前微米CT扫描实验图像转化为灰度图像,并确定分割阈值点,来提取孔隙空间与岩石骨架信息;
确定扫描区域内所有孔隙所占像素值与整个扫描区域像素值的比值为该区域的面孔率;
选择扫描区域面孔率Φap<3%的区域,划分为致密区,扫描区域面孔率Φap≥3%区域划分为孔隙发育区。
进一步地,所述将压裂前微米CT扫描实验图像转化为灰度图像,并确定分割阈值点,来提取孔隙空间与岩石骨架信息,包括:
将压裂前微米CT扫描实验图像转化为灰度图像,分为0~255共256个灰度值,将压裂前微米CT扫描实验图像中孔隙与岩石骨架和胶结物进行分割,使得分割后岩心图像中孔隙度像素占总像素的比例等于岩心孔隙度,提取孔隙空间时选择的分割阈值点的计算公式如公式(1)所示:
式中:f(n)为分割阈值判别公式,若f(n)=0,则表示n值为灰度阈值分割点;Φ为岩心测试孔隙度(%);i为灰度值;p(i)为灰度值为i的像素数。
进一步地,所述基于裂缝网络,通过网格分析法分析压裂缝扩展情况,包括:
选取压裂后的微米CT扫描的某一切面平面图;
通过边长为1mm的正方形网格将切面平面图划分为多个大小相等的微单元;
根据划分的网格,以1mm为步长,逐行对每一个网格进行扫描;扫描每个网格时,以待扫描的网格为中心,扫描周围3mm×3mm区域共9个网格,如果扫描区域内面孔率Φap>3%,则判定扫描的网格属于孔隙发育区,否则为致密区,直至完成所有区域扫描;
观察压裂后微米CT扫描图像,划出压裂缝在岩心的扩展情况。
进一步地,所述基于压裂缝经过的网格数量与孔隙发育区所包含的网格数量,确定压裂改造孔隙度,包括:
网格划分之后,边缘存在许多不完整的网格,在后续计算中,统一将面积<0.5mm2的网格记为0个,将面积>0.5mm2的网格记为1个;
压裂改造孔隙度即为压裂缝通过的网格数量与孔隙发育区所包含的网格数量之比乘以岩心测试孔隙度,表示为:
式中:Φfrac为压裂改造孔隙度(%);Nfrac为压裂缝通过的网格数量;Nnfrac为孔隙发育区所包含的网格数量;Φ为岩心测试孔隙度(%)。
进一步地,所述基于压裂改造孔隙度进一步确定压裂后的可动流体饱和度,包括:
压裂后可动流体饱和度计算公式如公式(4)、公式(5)所示:
SWMfrac=SWM+ΔSWM (5)
式中:ΔSWM为网格分析法计算得到的压裂后可动流体饱和度增量(%);SWM为核磁共振实验获取压裂前可动流体饱和度(%);SWMfrac为网格分析法计算得到的压裂后可动流体饱和度(%);由于孔隙空间内存在束缚水膜,裂缝并不能将孔隙空间内所有不可动流体转化为可动流体,所以用SWMF表示裂缝穿过时,可改造的孔隙空间内流体的比例,取SWMF=90%;Φfrac为压裂改造孔隙度(%);Φ为岩心测试孔隙度(%)。
进一步地,所述与测井曲线结合,实现压裂后可动流体饱和度连续计算,包括:
压裂前可动流体饱和度通过孔隙度、渗透率与储层品质因子进行多元拟合,如公式(6)所示:
压裂后可动流体饱和度增量与脆性指数直接相关,可通过公式(7)表示:
ΔSWM1=mBI+n (7)
压裂后可动流体饱和度可通过压裂前的可动流体饱和度与压裂后可动流体饱和度增量相加计算得到,如公式(8)所示:
SWMfrac1=SWM1+ΔSWM1 (8)
式中:SWM1为测井计算压裂前可动流体饱和度,%;Φ为岩心测试孔隙度(%);K为岩心测试渗透率(mD);BI为脆性指数;E为动态杨氏模量(GPa);ν为动态泊松比;SWMfrac1为测井计算压裂后可动流体饱和度(%);ΔSWM1为测井计算压裂后可动流体饱和度增量(%);a、b、c、d、m、n为拟合公式系数。
根据本发明提供的具体实施例,本发明公开了以下技术效果:本发明提供的基于网格分析法的储层压裂后可动流体饱和度计算方法,方法简单,包括压裂实验与压裂前后微米CT实验设计;划分孔隙发育区与致密区;对同一块岩心压裂前后的微米CT扫描实验图像对比,划分出压裂后形成的裂缝网络;通过网格分析法分析压裂缝扩展情况;通过裂缝经过的网格数量与孔隙发育区所包含的网格数量计算压裂改造孔隙度;通过压裂改造孔隙度进一步计算压裂后的可动流体饱和度;与测井曲线结合,实现压裂后可动流体饱和度连续计算。可见,本发明以压裂前后微米CT扫描数据为基础,通过网格分析法定量计算压裂后可动流体饱和度,再与测井曲线特征结合开展分析,实现通过测井曲线计算压裂后可动流体饱和度,使得对致密油气藏压裂后可动流体饱和度的评价更加准确高效,更好的指导致密油气藏的勘探与开发。
附图说明
为了更清楚地说明本发明实施例或现有技术中的技术方案,下面将对实施例中所需要使用的附图作简单地介绍,显而易见地,下面描述中的附图仅仅是本发明的一些实施例,对于本领域普通技术人员来讲,在不付出创造性劳动性的前提下,还可以根据这些附图获得其他的附图。
图1为本发明基于网格分析法的储层压裂后可动流体饱和度计算方法的流程图;
图2为本发明压裂实验与压裂前后微米CT扫描实验设计示意图;
图3为本发明图像孔隙空间提取示意图;
图4为本发明裂缝网络扩展趋势图;
图5为本发明网格分析法分析压裂缝扩展情况流程图;
图6为本发明实施例1号岩心实验流程设计图;
图7为本发明实施例1号岩心压裂缝网络扩展趋势图;
图8为本发明实施例1号岩心图像孔隙空间提取示意图;
图9为本发明实施例1号岩心网格分析法分析压裂缝扩展情况流程图;
图10为本发明实施例压裂后可动流体饱和度测井计算输出曲线示意图。
具体实施方式
下面将结合本发明实施例中的附图,对本发明实施例中的技术方案进行清楚、完整地描述,显然,所描述的实施例仅仅是本发明一部分实施例,而不是全部的实施例。基于本发明中的实施例,本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例,都属于本发明保护的范围。
本发明针对致密油气储层勘探开发过程中,目前还没有有效的方法计算致密油气藏压裂后可动流体饱和度的现状,提出一套以压裂前后微米CT扫描数据为基础,通过网格分析法定量计算压裂后可动流体饱和度,再与测井曲线特征结合开展分析,实现通过测井曲线计算压裂后可动流体饱和度,使得对致密油气藏压裂后可动流体饱和度的评价更加准确高效,更好的指导致密油气藏的勘探与开发。
为使本发明的上述目的、特征和优点能够更加明显易懂,下面结合附图和具体实施方式对本发明作进一步详细的说明。
如图1所示,本发明提供的基于网格分析法的储层压裂后可动流体饱和度计算方法,包括以下步骤:
S1,针对储层岩心样品,进行压裂实验与压裂前后微米CT实验设计;
S2,针对储层岩心样品,划分孔隙发育区与致密区;
S3,对同一块储层岩心样品压裂前后的微米CT扫描实验图像对比,划分出压裂后形成的裂缝网络;
S4,基于裂缝网络,通过网格分析法分析压裂缝扩展情况;
S5,基于压裂缝经过的网格数量与孔隙发育区所包含的网格数量,确定压裂改造孔隙度;
S6,基于压裂改造孔隙度进一步确定压裂后的可动流体饱和度;
S7,与测井曲线结合,实现压裂后可动流体饱和度连续计算。
其中,如图2所示,所述步骤S1,针对储层岩心样品,进行压裂实验与压裂前后微米CT实验设计,具体包括以下步骤:
(1)样品选取与制备
实验样品选择致密油气藏储层段岩心,将岩心钻切成直径2.5cm、长4~6cm的圆柱,圆柱的两端面要平整,垂直于轴线。
(2)压裂前微米CT扫描实验
压裂前微米CT扫描实验主要获取岩心内部形态和内部微观结构,可以观察岩心在地层原始状态下孔隙结构与天然裂缝发育情况,便于与压裂后实验结果进行对比分析。
在进行微米CT扫描实验时,需要选择扫描分辨率在8~12μm/pixel范围之内,在此范围内扫描的图像可以实现对整个样品全面精细扫描,同时可以清楚观察到压裂产生的裂缝的扩展规律,便于准确提取压裂缝,对压裂缝进行定量评价。
(3)压裂实验
通过室内三轴岩石力学实验模拟地层水力压裂实验,在开展三轴岩石力学实验时,首先用热塑膜包裹实验样品,将样品放入三轴力学实验测试仪,如果测试岩心所在深度段地层压力为Pf,在开展三轴岩石力学测试时,设置实验围压σ3=Pf,可以根据曲线特征将轴向应变曲线分为多个段,在OA段属于弹性变形;AB段应力应变基本属线性关系,卸载后可完全恢复;BC段,曲线偏离线性,出现塑性变形;CD段,岩石内部裂纹形成速度增快,裂纹密度加大,此时缓慢施加轴向压力,同时观察软件输出的应力应变曲线,轴向应变曲线在D点达到应力最大值,为岩石的最大承载能力,将此时的应力记为σmax,在应力下降1/6σmax~1/5σmax时,停止施压,在此状态下可以保证岩心内部出现压裂缝但是岩心未破损,方便进行压裂后微米CT扫描实验。
(4)压裂后微米CT扫描实验
岩心在进行压裂实验过后,即可以开展压裂后微米CT扫描实验。开展压裂后微米CT扫描需要满足以下条件:a)使用与进行压裂前微米CT扫描实验的同一台仪器;b)实验过程中岩心摆放方向、位置与压裂前扫描一致;c)选择扫描分辨率与压裂前扫描分辨率一致。在此条件下,能准确对比同一个位置压裂前后微米CT扫描图像,有利于对压裂缝的提取。
致密砂岩储层内部均质性较差,部分储层段含有大量砾石或以胶结物填充孔隙,为无效储集空间,统称为致密区。通过大量CT扫描数据观察分析发现,压裂后裂缝一般都会穿过孔隙发育区,绕过致密区扩展,而裂缝也会沟通孔隙发育区内的部分不可动流体,从而实现致密储层的有效开发。为了准确计算压裂缝改造孔隙空间体积,就需要较为准确地划分致密区域与孔隙发育区域。
因此,如图3所示,所述步骤S2,针对储层岩心样品,划分孔隙发育区与致密区,具体包括:
首先通过图像数值化处理提取孔隙,计算面孔率。具体处理步骤如下:
①CT扫描图像优化:岩心微米CT扫描获取灰度图像后,图像可由256个灰度值表示,选取数字图像中像素点及其周围临近8个像素点的像素值,将这些像素值按照升序进行排序,然后将位于中间位置的像素值作为当前像素点的像素值,让周围的像素值接近真实值,以此对图像进行优化处理,从而消除孤立的噪声点,提高图片对比度与清晰度;
②孔隙空间提取:通过阈值的确定来提取孔隙空间与岩石骨架信息,从而分割出图像中的孔隙部分。将图像转化为灰度图像,分为0~255共256个灰度值,将CT扫描灰度图像中孔隙与岩石骨架和胶结物等进行分割,使得分割后岩心图像中孔隙度像素占总像素的比例等于岩心孔隙度,保证分割结果的可靠性。提取孔隙空间时选择的分割阈值点的计算公式如公式(1)所示:
式中:f(n)为分割阈值判别公式,若f(n)=0,则表示n值为灰度阈值分割点;Φ为岩心孔隙度(%);i为灰度值(无量纲);p(i)为灰度值为i的像素数(无量纲)。
③扫描区域面孔率Φap计算:扫描区域内所有孔隙所占像素值与整个扫描区域像素值的比值为该区域的面孔率,计算公式如公式(2)所示:
式中:Φap为区域面孔率(%),N为所选区域内像素之和(无量纲),Nφ为所选区域内孔隙空间所占的像素之和(无量纲)。
④致密区与孔隙发育区划分:选择扫描区域面孔率Φap<3%的区域,划分为致密区,扫描区域面孔率Φap≥3%区域划分为孔隙发育区。
其中,如图4所示,所述步骤S3,对同一块储层岩心样品压裂前后的微米CT扫描实验图像对比,划分出压裂后形成的裂缝网络,具体包括:
微米CT扫描可以对岩心进行三维立体成像,并展示不同位置任意切面。在岩心压裂前进行微米CT扫描观察岩心在地层原始状态下的孔隙结构与微裂缝发育情况,通过三轴岩石力学实验模拟地层水力压裂,在实验结束后取下岩样,再次对岩样进行微米CT扫描,确定压裂之后的裂缝扩展情况。最后,对比岩样压裂前后的微米CT图像,判断裂缝为天然裂缝还是压裂缝,便于准确提取通过水力压裂改造形成的缝网。提取压裂缝缝网时,选择纵向五等分点 四处岩心扫描照片,选择主裂缝经过圆形切片直径附近的图像,通过网格分析法分别计算每个切片中压裂缝缝网的扩展面积,进一步计算对应点的可动流体饱和度,最后进行统计平均。
如图5所示,所述步骤S4,基于裂缝网络,通过网格分析法分析压裂缝扩展情况,具体包括:
①选取CT扫描某一切面平面图;
②通过边长为1mm的正方形网格将岩心切面划分为多个大小相等的微单元;
③根据划分的网格,以1mm为步长,逐行对每一个网格进行扫描;扫描每个网格时,以待扫描的网格为中心,扫描周围3mm×3mm区域共9个网格,如果扫描区域内面孔率Φap>3%,则判定扫描的网格属于孔隙发育区,否则为致密区,直至完成所有区域扫描;
④观察压裂后微米CT扫描图像,划出压裂缝在岩心的扩展情况。
具体地,所述步骤S5,基于压裂缝经过的网格数量与孔隙发育区所包含的网格数量,确定压裂改造孔隙度,具体包括:
网格划分之后,边缘存在许多不完整的网格,在后续计算中,统一将面积<0.5mm2的网格记为0个,将面积>0.5mm2的网格记为1个。压裂改造孔隙度即为压裂缝通过的网格数量与孔隙发育区的网格数量之比乘以岩心测试孔隙度。计算公式如公式(3)所示:
式中:Φfrac为压裂改造孔隙度(%);Nfrac为压裂缝经过的网格个数(无量纲);Nnfrac为孔隙发育区经过的网格个数(无量纲);Φ为岩心测试孔隙度(%)。
具体地,所述步骤S6,基于压裂改造孔隙度进一步确定压裂后的可动流体饱和度,包括:
压裂后可动流体饱和度计算公式如公式(4)、公式(5)所示:
SWMfrac=SWM+ΔSWM (5)
式中:ΔSWM为网格分析法计算得到的压裂后可动流体饱和度增量(%);SWM为核磁共振实验获取压裂前可动流体饱和度(%);SWMfrac为网格分析法计算得到的压裂后可动流体饱和度(%);由于孔隙空间内存在束缚水膜,裂缝并不能将孔隙空间内所有不可动流体转化为可动流体,所以用SWMF表示裂缝穿过时,可改造的孔隙空间内流体的比例,取SWMF=90%;Φfrac为压裂改造孔隙度(%);Φ为岩心测试孔隙度(%)。
具体地,所述步骤S7,与测井曲线结合,实现压裂后可动流体饱和度连续计算,包括:
通过网格分析法可以获取取心点地层压裂后的可动流体饱和度,如果要获取压裂后地层连续的可动流体饱和度,需要将实验数据与测井曲线特征联系,从而实现压裂后可动流体饱和度测井计算。
压裂前可动流体饱和度通过孔隙度、渗透率与储层品质因子进行多元拟合,如公式(6)所示:
压裂后可动流体饱和度增量与脆性指数直接相关,可通过公式(7)表示:
ΔSWM1=mBI+n (7)
压裂后可动流体饱和度可通过压裂前的可动流体饱和度与压裂后可动流体饱和度增量相加计算得到,如公式(8)所示:
SWMfrac1=SWM1+ΔSWM1 (8)
式中:SWM1为测井计算压裂前可动流体饱和度(%);Φ为岩心测试孔隙度(%);K为岩心测试渗透率(mD);BI为脆性指数(无量纲);E为动态杨氏模量(GPa);ν为动态泊松比(无量纲);SWMfrac1为测井计算压裂后可动流体饱和度(%);ΔSWM1为测井计算压裂后可动流体饱和度增量(%);a、b、c、d、m、n为拟合公式系数。
通过本发明所述方法计算得到压裂后可动流体饱和度,选择压裂后可动流体饱和度较高的井段开展压裂施工,对该井段进行压裂前后日产油气测试,测试结果显示:该层段开展压裂前基本无产能,压裂后日产油12.7方、气20680方;压裂效果显著。
在具体实施例中,以我国海上某致密油气藏为例,选择A井4592m开展相关实验与分析。岩心编号为1号,岩心基础数据如表1所示:
表1实验岩心基础数据表
(1)实验设计
①实验样品选择A井4592m储层段岩心,设置编号为1号。将岩心钻切成直径2.50cm、长4.80cm的圆柱,圆柱的两端面垂直于轴线;
②将制好样的岩心开展放入微米CT扫描仪,扫描分辨率设置为10μm/pixel,开展压裂前微米CT扫描,扫描结束后观察横切面孔隙结构特征与微裂缝发育情况;
③首先用热塑膜包裹实验样品,将样品放入三轴岩石力学实验仪,根据测压资料确定地层压力Pf=25MPa,选取实验围压σ3=Pf=25MPa,缓慢施加轴向压力,同时观察软件输出应力应变曲线,记录应变曲线最高点的差应力值σmax=181.8MPa,在差应力σ下降到约5/6σmax时(σ=151.5MPa),及时撤销轴向压力以及围压,停止轴向位移,取出岩心。此时岩心内部已出现压裂缝;
④将经过三轴岩石力学实验压裂过的样品重新放入微米CT扫描实验仪器,保证样品摆放方向、位置与压裂前扫描一致;设置扫描分辨率为10μm/pixel,开展岩心压裂后微米CT扫描。扫描完成后,观察微米CT扫描切片特征,选择与压裂前微米CT扫描切片特征相同的XY方向切片,对比观察压裂后裂缝扩展情况。1号岩心实验设计流程如图6所示。
(2)对比压裂前后微米CT扫描图像,提取压裂后形成的裂缝,准确画出压裂后产生的裂缝网络;由于压裂产生裂缝为剪切缝,压裂缝在纵向分布不均匀,此时选取压裂缝经过切片直径附近的图像(图)进行数据处理。1号岩心压裂缝网络扩展趋势图如图7所示。
(3)划分致密区与孔隙发育区
通过图像数值化处理提取孔隙,计算面孔率。具体处理步骤如下:
①CT扫描图像优化:将岩心微米CT扫描获取灰度图像用0~255的灰度值表示,选取数字图像中像素点及其周围临近8个像素点的像素值,将这些像素值按照升序进行排序,然后将位于排序中间位置的像素值作为当前像素点的像素值,以此对图像进行优化处理;
②孔隙空间提取:通过灰度阈值的确定来提取孔隙空间与岩石骨架和胶结物的信息,从而分割出图像中的孔隙部分。1号岩心图像孔隙提取过程如图8所示。
③区域面孔率计算:所选区域内所有孔隙的像素值与整个选择区域像素值的比值为该区域的面孔率计算公式如公式;
④致密区与孔隙发育区划分:选择区域面孔率<3%的区域,划分为致密区,区域面孔率≥3%区域划分为孔隙发育区。通过计算得到:1号岩心致密区比例为9.8%,孔隙发育区比例为90.2%。
(4)通过网格划分平面为多个微单元
通过网格将岩心切面划分为多个大小相等的微单元,划分致密区与孔隙发育区,观察裂缝发育形态。分析过程中将岩心切面面积通过微单元的数量进行表示,压裂后裂缝改造空间的面积即可以通过裂缝网络经过的微单元个数表示,从而实现平面内裂缝改造面积的定量计算。1号岩心网格分析法分析压裂缝扩展情况流程图如图9所示。
(5)计算压裂改造孔隙度
岩心实验测试孔隙度为9.38%,利用式(2),计算压裂缝经过网格个数与孔隙发育区所占网格总数,得到压裂改造孔隙度为3.21%。
(6)计算压裂后可动流体饱和度
压裂后未经裂缝影响区域可动流体饱和度依旧保持与压裂前一致,经过裂缝影响的孔隙空间可动流体饱和度达到理论最大值90%,利用式(3)、式(4)计算得到1号样品压裂前可动流体饱和度为56.27%。压裂后可动流体饱和度为68.06%。
(7)测井计算压裂后可动流体饱和度
通过网格分析法可以获取取心点地层压裂后的可动流体饱和度,如果要获取压裂后地层连续的可动流体饱和度,需要将实验数据与测井曲线特征联系,从而实现压裂后可动流体饱和度测井计算。
通过实验数据与测井数据分析计算得到A井:a=10.13,b=5.42,c=8.01,d=31.04,m=7.16,n=3.94。利用式(6),计算得到1号样品压裂前可动流体饱和度为57.28%;利用式(7),计算得到压裂后可动流体饱和度增量为9.60%;利用式(8),计算得到压裂后可动流体饱和度为66.88%,与网格分析法计算结果绝对误差为1.18%。
压裂后可动流体饱和度测井计算输出曲线如图10所示:在4592.0m处,网格分析法计算压裂后可动流体饱和度SWMfrac=68.06%,测井计算压裂后可动流体饱和度SWMfrac1=66.88%,绝对误差为1.18%;在4602.8m处,网格分析法计算压裂后可动流体饱和度SWMfrac=73.51%,测井计算压裂后可动流体饱和度SWMfrac1=72.36%,绝对误差为1.15%。
本文中应用了具体个例对本发明的原理及实施方式进行了阐述,以上实施例的说明只是用于帮助理解本发明的方法及其核心思想;同时,对于本领域的一般技术人员,依据本发明的思想,在具体实施方式及应用范围上均会有改变之处。综上所述,本说明书内容不应理解为对本发明的限制。

Claims (8)

1.一种基于网格分析法的储层压裂后可动流体饱和度计算方法,其特征在于,包括以下步骤:
针对储层岩心样品,进行压裂实验与压裂前后微米CT实验设计;
针对储层岩心样品,划分孔隙发育区与致密区;
对同一块储层岩心样品压裂前后的微米CT扫描实验图像对比,划分出压裂后形成的裂缝网络;
基于裂缝网络,通过网格分析法分析压裂缝扩展情况;
基于压裂缝经过的网格数量与孔隙发育区所包含的网格数量,确定压裂改造孔隙度;压裂改造孔隙度即为压裂缝通过的网格数量与孔隙发育区所包含的网格数量之比乘以岩心测试孔隙度,表示为:
式中:Φfrac为压裂改造孔隙度,%;Nfrac为压裂缝通过的网格数量;Nnfrac为孔隙发育区所包含的网格数量;Φ为岩心测试孔隙度,%;
基于压裂改造孔隙度进一步确定压裂后的可动流体饱和度;包括:压裂后可动流体饱和度计算公式如公式(4)、公式(5)所示:
SWMfrac=SWM+ΔSWM (5)
式中:ΔSWM为网格分析法计算得到的压裂后可动流体饱和度增量,%;SWM为核磁共振实验获取压裂前可动流体饱和度,%;SWMfrac为网格分析法计算得到的压裂后可动流体饱和度,%;由于孔隙空间内存在束缚水膜,裂缝并不能将孔隙空间内所有不可动流体转化为可动流体,所以用SWMF表示裂缝穿过时,可改造的孔隙空间内流体的比例,取SWMF=90%;Φfrac为压裂改造孔隙度,%;Φ为岩心测试孔隙度,%;
与测井曲线结合,实现压裂后可动流体饱和度连续计算,包括:
压裂前可动流体饱和度通过孔隙度、渗透率与储层品质因子进行多元拟合,如公式(6)所示:
压裂后可动流体饱和度增量与脆性指数直接相关,可通过公式(7)表示:
ΔSWM1=mBI+n (7)
压裂后可动流体饱和度可通过压裂前的可动流体饱和度与压裂后可动流体饱和度增量相加计算得到,如公式(8)所示:
SWMfrac1=SWM1+ΔSWM1 (8)
式中:SWM1为测井计算压裂前可动流体饱和度,%;Φ为岩心测试孔隙度,%;K为岩心测试渗透率,mD;BI为脆性指数;E为动态杨氏模量,GPa;ν为动态泊松比;SWMfrac1为测井计算压裂后可动流体饱和度,%;ΔSWM1为测井计算压裂后可动流体饱和度增量,%;a、b、c、d、m、n为拟合公式系数。
2.根据权利要求1所述的基于网格分析法的储层压裂后可动流体饱和度计算方法,其特征在于,所述针对储层岩心样品,进行压裂实验与压裂前后微米CT实验设计,包括:
样品选取与制备,选择致密油气藏储层段岩心作为储层岩心样品,将致密油气藏储层段岩心钻切成设定大小的圆柱,圆柱的两端面要平整,垂直于轴线;
压裂前微米CT扫描实验,通过微米CT扫描实验获取压裂前岩心内部形态和内部微观结构;
压裂实验,通过室内三轴岩石力学实验模拟地层水力压裂实验,在应力下降到设定范围时,停止施压;
压裂后微米CT扫描实验,通过微米CT扫描实验获取压裂后岩心内部形态和内部微观结构。
3.根据权利要求2所述的基于网格分析法的储层压裂后可动流体饱和度计算方法,其特征在于,所述压裂后微米CT扫描实验满足以下条件:
a)使用与进行压裂前微米CT扫描实验的同一台仪器;
b)实验过程中岩心摆放方向、位置与压裂前扫描一致;
c)选择扫描分辨率与压裂前扫描分辨率一致。
4.根据权利要求2所述的基于网格分析法的储层压裂后可动流体饱和度计算方法,其特征在于,所述压裂实验,通过室内三轴岩石力学实验模拟地层水力压裂实验,在应力下降到设定范围时,停止施压,具体包括:
在开展三轴岩石力学实验时,首先用热塑膜包裹储层岩心样品,将储层岩心样品放入三轴力学实验测试仪,如果测试岩心所在深度段地层压力为Pf,在开展三轴岩石力学测试时,设置实验围压σ3=Pf
根据曲线特征将轴向应变曲线分为多个段,在OA段属于弹性变形;AB段应力应变基本属线性关系,卸载后可完全恢复;BC段,曲线偏离线性,出现塑性变形;CD段,岩石内部裂纹形成速度增快,裂纹密度加大,此时缓慢施加轴向压力,同时观察软件输出的应力应变曲线,轴向应变曲线在D点达到应力最大值,为岩石的最大承载能力,将此时的应力记为σmax,在应力下降1/6σmax~1/5σmax时,停止施压。
5.根据权利要求1所述的基于网格分析法的储层压裂后可动流体饱和度计算方法,其特征在于,所述针对储层岩心样品,划分孔隙发育区与致密区,包括:
对储层岩心样品的压裂前微米CT扫描实验图像进行优化处理;
将压裂前微米CT扫描实验图像转化为灰度图像,并确定分割阈值点,来提取孔隙空间与岩石骨架信息;
确定扫描区域内所有孔隙所占像素值与整个扫描区域像素值的比值为该区域的面孔率;
选择扫描区域面孔率Φap<3%的区域,划分为致密区,扫描区域面孔率Φap≥3%区域划分为孔隙发育区。
6.根据权利要求5所述的基于网格分析法的储层压裂后可动流体饱和度计算方法,其特征在于,所述将压裂前微米CT扫描实验图像转化为灰度图像,并确定分割阈值点,来提取孔隙空间与岩石骨架信息,包括:
将压裂前微米CT扫描实验图像转化为灰度图像,分为0~255共256个灰度值,将压裂前微米CT扫描实验图像中孔隙与岩石骨架和胶结物进行分割,使得分割后岩心图像中孔隙度像素占总像素的比例等于岩心孔隙度,提取孔隙空间时选择的分割阈值点的计算公式如公式(1)所示:
式中:f(n)为分割阈值判别公式,若f(n)=0,则表示n值为灰度阈值分割点;Φ为岩心测试孔隙度,%;i为灰度值;p(i)为灰度值为i的像素数。
7.根据权利要求1所述的基于网格分析法的储层压裂后可动流体饱和度计算方法,其特征在于,所述基于裂缝网络,通过网格分析法分析压裂缝扩展情况,包括:
选取压裂后的微米CT扫描的某一切面平面图;
通过边长为1mm的正方形网格将切面平面图划分为多个大小相等的微单元;
根据划分的网格,以1mm为步长,逐行对每一个网格进行扫描;扫描每个网格时,以待扫描的网格为中心,扫描周围3mm×3mm区域共9个网格,如果扫描区域内面孔率Φap>3%,则判定扫描的网格属于孔隙发育区,否则为致密区,直至完成所有区域扫描;
观察压裂后微米CT扫描图像,划出压裂缝在岩心的扩展情况。
8.根据权利要求7所述的基于网格分析法的储层压裂后可动流体饱和度计算方法,其特征在于,所述基于压裂缝经过的网格数量与孔隙发育区所包含的网格数量,确定压裂改造孔隙度,包括:
网格划分之后,边缘存在许多不完整的网格,在后续计算中,统一将面积<0.5mm2的网格记为0个,将面积>0.5mm2的网格记为1个。
CN202310485930.XA 2023-05-04 2023-05-04 基于网格分析法的储层压裂后可动流体饱和度计算方法 Active CN116502557B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202310485930.XA CN116502557B (zh) 2023-05-04 2023-05-04 基于网格分析法的储层压裂后可动流体饱和度计算方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202310485930.XA CN116502557B (zh) 2023-05-04 2023-05-04 基于网格分析法的储层压裂后可动流体饱和度计算方法

Publications (2)

Publication Number Publication Date
CN116502557A CN116502557A (zh) 2023-07-28
CN116502557B true CN116502557B (zh) 2024-01-30

Family

ID=87319836

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202310485930.XA Active CN116502557B (zh) 2023-05-04 2023-05-04 基于网格分析法的储层压裂后可动流体饱和度计算方法

Country Status (1)

Country Link
CN (1) CN116502557B (zh)

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101929973A (zh) * 2009-06-22 2010-12-29 中国石油天然气股份有限公司 裂缝储层含油气饱和度定量计算方法
CN112459754A (zh) * 2020-11-04 2021-03-09 中国石油天然气集团有限公司 一种干法压裂焖井后co2与储层流体置换规律实验方法
CN115434698A (zh) * 2022-09-01 2022-12-06 西南石油大学 基于光电吸收截面指数的地层钙质胶结物含量计算方法
CN115713049A (zh) * 2022-11-25 2023-02-24 西南石油大学 一种耦合页岩水化膨胀与致裂作用的焖井时间优化方法

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101929973A (zh) * 2009-06-22 2010-12-29 中国石油天然气股份有限公司 裂缝储层含油气饱和度定量计算方法
CN112459754A (zh) * 2020-11-04 2021-03-09 中国石油天然气集团有限公司 一种干法压裂焖井后co2与储层流体置换规律实验方法
CN115434698A (zh) * 2022-09-01 2022-12-06 西南石油大学 基于光电吸收截面指数的地层钙质胶结物含量计算方法
CN115713049A (zh) * 2022-11-25 2023-02-24 西南石油大学 一种耦合页岩水化膨胀与致裂作用的焖井时间优化方法

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
Fracture-matrix interactions during immiscible three-phase flow;Elfeel, MA等;《JOURNAL OF PETROLEUM SCIENCE AND ENGINEERING》;171-186 *
松辽盆地古龙页岩油人工油藏的动态模拟及预测;赵国忠等;《大庆石油地质与开发》;第40卷(第5期);170-180 *

Also Published As

Publication number Publication date
CN116502557A (zh) 2023-07-28

Similar Documents

Publication Publication Date Title
Ju et al. Laboratory in situ CT observation of the evolution of 3D fracture networks in coal subjected to confining pressures and axial compressive loads: a novel approach
AU2011258594B2 (en) Method for obtaining consistent and integrated physical properties of porous media
CN110644980B (zh) 一种超低渗透油藏储层综合分类评价方法
US8170799B2 (en) Method for determining in-situ relationships between physical properties of a porous medium from a sample thereof
AU2009316427B2 (en) Method for determining rock physics relationships using computer tomographic images thereof
Weniger et al. Characterizing coal cleats from optical measurements for CBM evaluation
AU2009316347B2 (en) Method for determining elastic-wave attenuation of rock formations using computer tomograpic images thereof
CN113609696B (zh) 基于图像融合的多尺度多组分数字岩心构建方法及系统
Jiang et al. Three-dimensional visualization of the evolution of pores and fractures in reservoir rocks under triaxial stress
US20150331145A1 (en) Method for producing a three-dimensional characteristic model of a porous material sample for analysis of permeability characteristics
CN108956424A (zh) 一种页岩中孔隙定量表征的方法
CN111425193A (zh) 一种基于聚类分析测井岩石物理相划分的储层可压性评价方法
CN110208487A (zh) 一种基于ct扫描的页岩水化损伤测试方法
Duan et al. An investigation of the evolution of the internal structures and failure modes of Longmaxi shale using novel X-ray microscopy
Wong Strain-induced anisotropy in fabric and hydraulic parameters of oil sand in triaxial compression
Yang et al. Quantifying the impact of 2D and 3D fractures on permeability in wellbore cement after uniaxial compressive loading
CN114897767A (zh) 一种致密混积岩储层储集空间多尺度表征与储层分类方法
CN116502557B (zh) 基于网格分析法的储层压裂后可动流体饱和度计算方法
CN112487620A (zh) 一种页岩油可动资源量评价模型、评价方法、应用
CN115147340A (zh) 岩石裂缝变形信息获取方法、装置、设备和存储介质
Zhang et al. Fine characterization of pore structures in coral reef limestones based on three-dimensional geometrical reconstruction
Crandall et al. Foamed cement analysis with computed tomography
CN114089421B (zh) 一种油气储层非均质性分析方法
Lukas An Experimental Investigation of the Influence of Sampling on the behavior of Intermediate Soils
Wang Ultradeep Carbonate Gas Reservoirs: Reservoir Characteristics and Percolation Mechanism

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
CB03 Change of inventor or designer information

Inventor after: Wu Feng

Inventor after: Chen Chunchao

Inventor after: Wang Ao

Inventor after: Long Toujing

Inventor after: Chen Siyuan

Inventor after: Luo Yingying

Inventor before: Chen Chunchao

Inventor before: Wu Feng

Inventor before: Wang Ao

Inventor before: Long Toujing

Inventor before: Chen Siyuan

Inventor before: Luo Yingying

CB03 Change of inventor or designer information