CN116559941A - 一种基于Norris-KG模型的地震纵波模拟分析方法 - Google Patents

一种基于Norris-KG模型的地震纵波模拟分析方法 Download PDF

Info

Publication number
CN116559941A
CN116559941A CN202310364518.2A CN202310364518A CN116559941A CN 116559941 A CN116559941 A CN 116559941A CN 202310364518 A CN202310364518 A CN 202310364518A CN 116559941 A CN116559941 A CN 116559941A
Authority
CN
China
Prior art keywords
wave
model
longitudinal
norris
seismic
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN202310364518.2A
Other languages
English (en)
Other versions
CN116559941B (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.)
Oil & Gas Survey Cgs
Original Assignee
Oil & Gas Survey Cgs
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 Oil & Gas Survey Cgs filed Critical Oil & Gas Survey Cgs
Priority to CN202310364518.2A priority Critical patent/CN116559941B/zh
Publication of CN116559941A publication Critical patent/CN116559941A/zh
Application granted granted Critical
Publication of CN116559941B publication Critical patent/CN116559941B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • GPHYSICS
    • G01MEASURING; TESTING
    • G01VGEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V1/00Seismology; Seismic or acoustic prospecting or detecting
    • G01V1/28Processing seismic data, e.g. for interpretation or for event detection
    • 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
    • Y02A10/00TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE at coastal zones; at river basins
    • Y02A10/40Controlling or monitoring, e.g. of flood or hurricane; Forecasting, e.g. risk assessment or mapping

Landscapes

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

Abstract

本发明提供了一种基于Norris‑KG模型的地震纵波模拟分析方法,其基于Norris‑KG裂缝‑孔隙等效介质模型,通过qP波相速度公式及纵、横波分离的一阶速度应力方程对介质的波场控制方程进行推导,并利用基于PML边界条件的交错网格高阶有限差分算法对地震波场进行正演数值模拟,最终,在给定参数下,对目的层的qP波相速度方位各向异性及其顶、底界面的反射波振幅方位各向异性分析,对目的层的各向异性及其影响参数进行优选分析,得到裂缝‑孔隙介质地震纵波相速度和波场响应数据,最后进行纵波数值模拟分析。本发明构思合理,构思合理,能从地震纵波响应特征中直接对储层的物性参数进行识别,分析效果好且精确。

Description

一种基于Norris-KG模型的地震纵波模拟分析方法
技术领域
本发明涉及油气勘探技术领域,具体涉及一种基于Norris-KG模型的地震纵波模拟分析方法。
背景技术
裂缝性储层内油、水或气、水两相流体共存时,储层地震响应具有频变方位各向异性特征,如何准确建立地震响应特征与流体饱和度参数之间的关系,是储层预测的关键和难点。常规方法主要基于地震解释、属性提取等手段对二者的关系进行推测和估计,预测结果的多解性强,可靠性差。基于裂缝-孔隙等效介质建模的地震纵波模拟方法,能够通过裂缝-孔隙等效介质模型计算裂缝性储层内不同流体饱和度岩石的弹性参数与储层物性参数之间的联系,通过正演数值模拟给定岩石地震响应特征与岩石岩性参数的关系,从而确定纵波地震响应特征与储层物性参数之间的确定性关系,是实现多相流体饱和裂缝性储层高精度预测的基础。
周期成层理论是目前裂缝-孔隙等效介质建模的基础理论。Norris(1993)指出,周期成层两套孔隙地层的弹性属性可以通过两套地层各自的属性联合得到,并给出垂直入射时周期成层两套孔隙介质纵波模量的通用表达式。Brajanovshi等(2005)在Norris模型的基础上,将两套地层分别赋予基质孔隙地层和裂缝地层的物理含义,给出了单一流体饱和的周期成层裂缝-孔隙介质在垂直裂缝面方向纵波模型的近似表达式,但由于两套地层均为体积模量相对较大的同一种流体(比如含水时)所饱和,导致该模型在流体体积含量较小(比如含气时)时并不适用。Kong等(2013)对Brajanovski等的模型进行了改进,通过在Norris模型中引入裂缝流体指示因子F,突破了Brajanovski等(2005)模型关于流体体积模量的限制,使得Norris模型真正转化为具有实际物理意义的裂缝-孔隙等效介质模型。但Kong等(2013)模型(以下称KG模型)假设孔隙和裂缝仍然是单一流体饱和的两套地层,虽然两套地层内所包含的流体可以不同,并没有实现真正意义上建立整体模型地震波速度与流体饱和度的定量关系。孔丽云等(2023)借助于KG模型对流体无任何假设的优势,结合Norris周期层状模型的普适性,建立部分饱和裂缝-孔隙等效介质模型(Norris-KG模型),并将模型推广至任意入射角的情况。通过数值计算,模拟分析了含水饱和度对裂缝性储层频率依赖纵波相速度及其各向异性的影响,并通过实验室测量数据的拟合分析了Norris-KG模型的有效性,能够对海相碳酸盐岩裂缝性储层的各向异性和频散衰减特性进行有效描述,并在此基础上对介质的AVO响应特征进行分析。
地震资料的数值模拟在地震勘探处理方法的研究中起着非常重要的作用。一方面可以利用数值模拟地震记录检验处理结果的好坏、处理方法的有效性,另一方面正演模拟可以作为反演研究的基础(吴国忱,2006;梁锴,2009)。地震波场正演模拟方法主要可以分为三类:直接法,积分方程和射线追踪方法;其中,动力学孔隙弹性波动方程的求解最常用的是直接方法(Carcione et al.,2010),主要有三种类型:伪谱法,有限元法和有限差分法。有限元法是以变分原理和插值技术为基础的地震波方程数值模拟方法(Lysmer andDrake,1972),对计算机内存要求很高,计算量大。伪谱法在地震领域的应用则大大加快了计算速度,并节省了计算机存储空间(Gazdag,1981;Kosloff and Baysal,1982)[39-40],但该方法对边界的处理较难。有限差分法是目前声波或弹性波场数值模拟中应用最广的方法(Alterman and Karal,1968;Boore,1972;Alford et al.,1974;Carcione,1996;牟永光,1996;孙卫涛,2009),其实现简单,运行速度快,但难以克服空间离散带来的数值频散问题。高阶差分方法(Virieux,1986)和一阶弹性波交错网格方法(Dablain,1986)的提出与结合,形成了一阶交错网格高阶有限差分方法(董良国等,2000),该方法有力地减弱了常规有限差分法的数值频散问题,因而在弹性介质(裴正林,2004)、粘弹性介质(李军等,2009;韩令贺等,2010)、非均匀介质(裴正林,牟永光,2003)、各向异性介质(裴正林,2006)、孔隙介质(王秀明等,2003)、裂缝-孔隙介质(杜启振等,2009;孔丽云等,2013)的正演模拟中取得了相对成功的应用。
纵横波波场分离方面,目前常用的方法是亥姆霍兹分解(Yan和Sava,2008)。该方法所依据的基本原理是由于纵横波偏振方向的不同(纵波偏振方向平行于传播方向,横波偏振方向垂直于传播方向),可以用波场的散度和旋度来分别对纵波和横波进行求取。这在各向同性的情况下是正确的,但是在各向异性介质中,地震波传播的方向发生了偏离,其振动方向并不再是纵波与其相同、横波与其垂直的关系,所以使用亥姆霍兹分解并不能很好地将各向异性介质中的耦合纵横波分离(Yan和Sava,2009)。另外,由于在分离过程中采用散度场代表纵波,旋度场代表横波,这就意味着,纵波仅代表整体位移的输入输出,而横波则代表旋转变换的角速度,其物理意义与原波场并不相同,且其振幅和相位信息都发生了改变。qP波近似方程方面,Alkhalifah从VTI介质的精确qP-qSV波频散关系出发,假设垂向剪切波速度Vs0=0,推导了VTI介质四阶拟声波方程。为降低四阶偏微分方程计算的复杂性(Alkhalifah,2008),Zhou等(2009)和Du等(2008)从声学近似的频散关系出发,通过引入不同的辅助波场函数,分别推导了两种VTI介质二阶耦合拟声波方程。为结合交错网格技术实现精确的有限差分数值模拟,Hestholm推导了一阶偏微分形式的VTI介质拟声波方程(Hestholm,2009)。Duveneck等(2008)从声学近似的胡克定律出发,也推导了二阶形式和一阶应力-速度形式的VTI介质拟声波方程。与ALkhalifah方程一样,几种拟声波方程均存在低速度、低振幅的qSV人为干扰波问题(Grechka等,2004),而这种干扰波的存在会影响正演模拟和偏移成像的最终结果。同时,由于无法考虑到转换波的信息,利用qP波近似方程得到的各向异性介质纵波数据并不完整。为得到振幅和相位信息相对准确且波场信息相对完整的各向异性介质纵波数据,Zhang等提出“方程解耦”的分离方法(Zhang等,2007)。在各向同性介质中,Zhang等(2007)得到精确的纵横波分离的一阶速度应力方程;Zhou等(2016)将这种方法进一步拓展到三维各向同性模型的情况。孔丽云等(2022)通过将弹性刚度系数分量中的耦合项C13近似为纵波模量和横波模量的线性函数,推导得到了各向异性介质的纵横波分离的一阶速度-应力方程,从而能够实现各向异性介质的纵横波波场分离。
综上可知,目前的研究大多将等效介质建模、正演数值模拟、纵横波场分离这三方面的技术分开研究,无法直观体现裂缝性储层纵波地震响应特征与储层物性参数之间的联系。
发明内容
针对上述背景技术中存在的问题,本发明提出了一种基于Norris-KG模型的地震纵波模拟分析方法,其构思合理,以储层物性参数为输入,通过三方面理论的融合、计算,能够直接给出储层的纵波地震响应,并通过分析研究纵波地震响应特征随储层物性参数的变化规律,有助于揭示裂缝性储层地震响应特征的产生机理,为非常规储层的勘探开发提供可靠的理论依据。
为解决上述技术问题,本发明提供的一种基于Norris-KG模型的地震纵波模拟分析方法,基于Norris-KG裂缝-孔隙等效介质模型,通过qP波相速度公式及纵、横波分离的一阶速度应力方程对介质的波场控制方程进行推导,并利用基于PML边界条件的交错网格高阶有限差分算法对地震波场进行正演数值模拟,最终,在给定参数下,对目的层的qP波相速度方位各向异性及其顶、底界面的反射波振幅方位各向异性分析,对目的层的各向异性及其影响参数进行优选分析,得到裂缝-孔隙介质地震纵波相速度和波场响应数据,最后进行纵波数值模拟分析。
所述基于Norris-KG模型的地震纵波模拟分析方法,其中,所述波场控制方程的推导过程为:
1.1)构建qP波相速度方程
假设裂缝倾角为θ0,地震波的传播方向为n=(nx,ny,nz)T=(sinθcosφ,sinθsinφ,cosθ)T,其中,θ为传播方向与观测系统z轴的夹角,φ为传播方向与x轴的夹角,任意入射角的纵波相速度表达式为:
其中,
上式(2)建立了裂缝-多孔隙介质的纵波相速度Vqp与储层弹性参数C之间的数学联系,结合Norris-KG模型所给出的介质弹性模量表达式(1),通过数值计算,得到介质的纵波相速度随储层物性参数的变化规律;
1.2)构建纵横波分离的一阶速度-应力方程
基于“方程分离”的思想,将弹性刚度系数分量中的耦合项C13近似为纵波模量和横波模量的线性函数,推导得到了各向异性介质的纵、横波分离的一阶速度-应力方程如下:
其中,上述方程(3)中下角标p和s分别代表纵波和横波,G为偏微分算子矩阵,Cp为与纵波相关的弹性系数矩阵,Cs为与横波相关的弹性系数矩阵;且Cp与Cs表达式分别如下:
上述的Cp与Cs表达式中,ε和γ为Thomsen各向异性参数,Vp0为qP波沿VTI介质对称轴方向的相速度,Vs0为qSV波沿VTI介质对称轴方向的相速度;和B=-2分别为纵波和横波模量的系数;
上述方程(3)建立了裂缝-多孔隙介质的地震纵波波场v与储层物性参数C之间的数学联系,通过对波动方程的数值求解,可得到介质的纵波地震响应vp,进而对储层纵波地震响应特征随物性参数的变化规律进行分析。
所述基于Norris-KG模型的地震纵波模拟分析方法,其中,所述纵波数值模拟分析的具体过程为:
2.1)qP波相速度分析
2.1.1)基准模型结果分析
基于任意入射角的纵波相速度公式(2),根据给定井段的岩石物理测量数据和测井数据,对目的层qP波的相速度及储层物性参数的影响进行分析;
2.1.2)储层参数的影响分析
在给定参数的基础上,分别改变裂缝倾角、裂缝线密度、基质孔隙度、含气饱和度等参数,以分析不同储层物性参数对qP波相速度的影响;
2.2)波场传播特征分析
2.2.1)基准模型结果分析
设计三层理论模型对目的层顶、底界面的地震波场传播特征进行分析,然后基于所述三层模型的参数,根据一阶速度-应力方程(3),采用基于PML边界条件的交错网格高阶有限差分方法,对所述三层模型进行正演数值模拟;
2.2.2)储层参数的影响分析
在给定参数的基础上,保持其他参数不变,分别改变裂缝倾角、裂缝线密度、基质孔隙度和含气饱和度,以分析不同的参数对地震波场响应特征的影响。
所述基于Norris-KG模型的地震纵波模拟分析方法,其中,所述正演数值模拟是采用交错网格高阶有限差分法对介质的一阶速度-应力方程进行数值求解;具体包括以下步骤:
3.1)方程的离散化处理
采用交错网格高阶有限差分算法,其时间、空间高阶差分近似公式分别为:
其中,上式(4)-(5)中Δt为时间步长,2M为时间差分精度,2L为空间差分精度;Δx为差分网格间距,为不同阶的交错网格空间差分系数,可以通过如下方程组求取:
将该方程组带入到上述方程(3)中得到一阶速度应力方程的差分近似格式;
3.2)PML边界条件处理
采用PML边界条件处理方法,对上述方程(3)中的总体波场v进行变量分离,每一个波场变量可以分为纵波分量vp和横波分量vs两部分,T是总的应力场,上角标x、y、z代表三个坐标轴方向:
然后对于所给定的三维空间模型设置边界衰减系数。
所述基于Norris-KG模型的地震纵波模拟分析方法,其中,所述Norris-KG裂缝-孔隙等效介质模型的具体构建过程为:基于Norris周期成层模型和KG裂缝-孔隙等效介质模型建立Norris-KG等效介质模型;基于Norris-KG等效介质模型,假设裂缝为非常薄且松软的无限大裂缝平面,将Norris-KG等效介质模型看成是周期成层两套裂缝-孔隙地层,根据Norris理论及Heuristic假设,Norris-KG等效介质模型的弹性属性表达式为:
上式(1)中,为垂直裂缝面方向纵波模量,/>为模型高频极限弹性系数,/>为垂直裂缝面方向的高频极限弹性系数,/>为模型低频极限弹性系数;式(1)给出的是水平裂缝性储层的弹性系数矩阵,可计算任意入射角情况的频变纵波速度随流体饱和度的变化情况。
采用上述技术方案,本发明具有如下有益效果:
本发明基于Norris-KG模型的地震纵波模拟分析方法是在已知储层物性参数(颗粒类型、孔隙形态、流体属性等)的情况下,给出地震纵波响应特征(速度、振幅等)与储层物性参数之间的对应关系,让我们有能力从地震纵波响应特征中直接对储层的物性参数进行识别。
本发明将Norris-KG等效介质模型、纵横波分离算法、数值模拟计算三者相融合,实现了基于Norris-KG模型的地震纵波模拟分析。在给定参数下,通过对某研究区目的层的qP波相速度方位各向异性及其顶、底界面的反射波振幅方位各向异性分析,对目的层的各向异性及其影响参数进行优选。结果表明,在给定的岩石物理参数和测井参数条件下,该研究区目的层的观测各向异性较为明显,且主要受裂缝倾角和裂缝线密度的影响。
附图说明
为了更清楚地说明本发明具体实施方式或现有技术中的技术方案下面将对具体实施方式或现有技术描述中所需要使用的附图作简单地介绍,显而易见地,下面描述中的附图是本发明的一些实施方式,对于本领域普通技术人员来讲,在不付出创造性劳动的前提下,还可以根据这些附图获得其他的附图。
图1为本发明基于Norris-KG模型的地震纵波模拟分析方法中涉及的Norris-KG模型示意图;
图2为本发明基于Norris-KG模型的地震纵波模拟分析方法中涉及的三维交错网格示意图;
图3为本发明基于Norris-KG模型的地震纵波模拟分析方法中涉及的三维完全匹配层吸收边界示意图;
图4为本发明基于Norris-KG模型的地震纵波模拟分析方法中涉及的qP波相速度随入射角和观测方位角的变化(红、紫、蓝线代表入射角分别为17°,32°和47°的情况)示意图;
图5为本发明基于Norris-KG模型的地震纵波模拟分析方法中涉及的不同裂缝倾角、地震波入射角下的qP波相速度随方位角的变化规律((a)-(i)分别对应裂缝倾角10°-90°,间隔为10°;红、紫、蓝线代表入射角分别为17°,32°和47°的情况)示意图;
图6为本发明基于Norris-KG模型的地震纵波模拟分析方法中涉及的裂缝线密度分别为(a)1条/m、(b)15条/m和(c)30条/m时的qP波相速度随方位角的变化规律(红、紫、蓝线代表入射角分别为17°,32°和47°的情况)示意图;
图7为本发明基于Norris-KG模型的地震纵波模拟分析方法中涉及的基质孔隙度分别为(a)1.12%、(b)3.40%和(c)5.81%时的qP波相速度随方位角的变化规律(红、紫、蓝线代表入射角分别为17°,32°和47°)示意图;
图8为本发明基于Norris-KG模型的地震纵波模拟分析方法中涉及的含气饱和度分别为(a)0%、(b)25%、(c)50%、(d)75%和(e)100%时的qP波相速度(红、紫、蓝线代表入射角分别为17°,32°和47°的情况)示意图;
图9为本发明基于Norris-KG模型的地震纵波模拟分析方法中涉及的层顶界面(a)、底界面(b)的叠前全方位角道集模拟数据示意图;
图10为本发明基于Norris-KG模型的地震纵波模拟分析方法中涉及的目的层顶界面(a)、底界面(b)的振幅方位各向异性示意图;
图11为本发明基于Norris-KG模型的地震纵波模拟分析方法中涉及的裂缝倾角分别为1条/m(左)、15条/m(中)和30条/m(右)时,目的层顶界面蜗牛道集((a)-(c))、目的层底界面蜗牛道集((d)-(f))、目的层顶界面反射波振幅随方位角的变化((g)-(i))、目的层底界面反射波振幅随方位角的变化((j)-(l))示意图;
图12为本发明基于Norris-KG模型的地震纵波模拟分析方法中涉及的裂缝线密度分别为(a)1条/m、(b)15条/m和(c)30条/m时,目的层顶界面蜗牛道集((a)-(c))、目的层底界面蜗牛道集((d)-(f))、目的层顶界面反射波振幅随方位角的变化((g)-(i))、目的层底界面反射波振幅随方位角的变化((j)-(l))示意图;
图13为本发明基于Norris-KG模型的地震纵波模拟分析方法中涉及的基质孔隙度分别为(a)1.12%、(b)3.40%和(c)5.81%时,目的层顶界面蜗牛道集((a)-(c))、目的层底界面蜗牛道集((d)-(f))、目的层顶界面反射波振幅随方位角的变化((g)-(i))、目的层底界面反射波振幅随方位角的变化((j)-(l))示意图;
图14为本发明基于Norris-KG模型的地震纵波模拟分析方法中涉及的含气饱和度分别为(a)0%、(b)50%和(c)100%时,目的层顶界面蜗牛道集((a)-(c))、目的层底界面蜗牛道集((d)-(f))、目的层顶界面反射波振幅随方位角的变化((g)-(i))、目的层底界面反射波振幅随方位角的变化((j)-(l))示意图。
具体实施方式
下面将结合附图对本发明的技术方案进行清楚、完整地描述,显然,所描述的实施例是本发明一部分实施例,而不是全部的实施例。基于本发明中的实施例,本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例,都属于本发明保护的范围。
下面结合具体的实施方式对本发明做进一步的解释说明。
本实施例提供的基于Norris-KG模型的地震纵波模拟分析方法,基于Norris-KG裂缝-孔隙等效介质模型,通过qP波相速度公式及纵、横波分离的一阶速度应力方程对介质的波场控制方程进行推导,并利用基于PML边界条件的交错网格高阶有限差分算法对地震波场进行正演数值模拟,最终,在给定参数下,对目的层的qP波相速度方位各向异性及其顶、底界面的反射波振幅方位各向异性分析,对目的层的各向异性及其影响参数进行优选分析,得到裂缝-孔隙介质地震纵波相速度和波场响应数据,最后进行纵波数值模拟分析。
本实施例提供的基于Norris-KG模型的地震纵波模拟分析方法,主要包括以下步骤:
S100)Norris-KG等效介质建模
基于Norris周期成层模型和KG裂缝-孔隙等效介质模型,建立图1所示的部分饱和裂缝-孔隙等效介质模型(以下称Norris-KG模型)。图中,模型周期为H,两套地层均由干燥的KG模型所组成,蓝色代表的是水饱和KG模型(流体体积模量和粘度系数分别为Kwater、ηwater),灰色代表的是气饱和KG模型(流体体积模量和粘度系数分别为Kgas、ηgas)。这样,含水饱和度Sw和含气饱和度Sg就分别代表水饱和和气饱和地层的空间占比。
基于图1所示的模型,假设裂缝为非常薄且松软的无限大裂缝平面(符合KG模型假设条件),将Norris-KG模型看成是周期成层两套裂缝-孔隙地层(单一流体饱和),根据Norris理论及Heuristic假设(水平层状介质中,流体的能量传播方向始终垂直水平地层),模型的弹性属性表达式为:
其中,为垂直裂缝面方向纵波模量,/>分别为模型高、低频极限弹性系数。
式(1)给出的是水平裂缝性储层(即VTI介质,Vertical Transverse Isotropy)的弹性系数矩阵,能够计算任意入射角情况的频变纵波速度随流体饱和度的变化情况。
S200)波场控制方程推导
S210)qP波相速度方程
假设裂缝倾角为θ0,地震波的传播方向为n=(nx,ny,nz)T=(sinθcosφ,sinθsinφ,cosθ)T,其中,θ为传播方向与观测系统z轴的夹角,φ为传播方向与x轴的夹角,任意入射角的纵波相速度表达式:
其中,
式(2)建立了裂缝-多孔隙介质的纵波相速度Vqp与储层弹性参数C之间的数学联系,结合Norris-KG模型所给出的介质弹性模量表达式(1),通过数值计算,可以得到介质的纵波相速度随储层物性参数的变化规律。
S220)纵横波分离的一阶速度-应力方程
基于“方程分离”的思想,将弹性刚度系数分量中的耦合项C13近似为纵波模量和横波模量的线性函数,推导得到了各向异性介质的纵、横波分离的一阶速度-应力方程:
其中,下角标p和s分别代表纵波和横波,G为偏微分算子矩阵;Cp和Cs分别为与纵波和横波相关的弹性系数矩阵,Cp和Cs的表达式如下:
ε和γ为Thomsen各向异性参数;Vp0、Vs0分别为qP波和qSV波沿VTI介质对称轴方向的相速度;和B=-2分别为纵波和横波模量的系数。
表达式(3)建立了裂缝-多孔隙介质的地震纵波波场v与储层物性参数C之间的数学联系,通过对波动方程的数值求解,可以得到介质的纵波地震响应vp,进而对储层纵波地震响应特征随物性参数的变化规律进行分析。
S300)正演数值模拟
交错网格高阶有限差分法是通过将速度(应力)对时间的奇数阶高阶导数转化为应力(速度)对空间的导数,运用时间和空间差分精度均可达到任意阶的高阶差分法,采用交错网格技术,对介质的一阶速度-应力方程进行数值求解的过程。
S310)方程的离散化处理
采用交错网格高阶有限差分算法,其时间、空间高阶差分近似公式分别为:
其中,上式(4)-(5)中Δt为时间步长,2M为时间差分精度,2L为空间差分精度;Δx为差分网格间距,为不同阶的交错网格空间差分系数,可以通过如下方程组求取:
将上述表达式带入到方程(3)中,就得到一阶速度应力方程的差分近似格式。
进一步,时间上,将速度场放在主时间网格上,应力场放在次时间网格上;空间上,速度场和应力场的位置如图2所示。这样,就可以得到差分近似格式的离散化形式。
S320)PML边界条件处理
采用PML边界条件处理方法,对上述方程(3)中的总体波场v进行变量分离,每一个波场变量可以分为纵波分量vp和横波分量vs两部分,T是总的应力场,上角标x、y、z代表三个坐标轴方向:
对于三维模型,按照图3设置边界衰减系数。边界面区:A,dx=0,dy>0,dz=0;B,dx=0,dy=0,dz>0;C,dx>0,dy=0,dz=0。在边界棱区:D,dx=0,dy>0,dz>0;E,dx>0,dy=0,dz>0;F,dx>0,dy>0,dz=0。在边界角点区:G,dx>0,dy>0,dz>0。
最终,通过步骤S100-S300,得到裂缝-孔隙介质地震纵波相速度和波场响应数据。
S400)纵波数值模拟分析
S410)qP波相速度分析
S4101)基准模型结果分析
基于任意入射角的纵波相速度公式(2),根据给定井段的岩石物理测量数据和测井数据,对目的层qP波的相速度及储层物性参数的影响进行分析。
岩石固体颗粒为白云石,体积模量、剪切模量和密度分别为94.9GPa,45GPa和2870Kg/m3;储层内所含流体为气和水,其中,气的体积模量和密度分别为0.213GPa和258Kg/m3,水的体积模量和密度分别为2.5GPa和1020Kg/m3;储层中的基质孔隙度为5.81%;储层中裂缝的线密度为30条/m,裂缝平均宽度为0.0024m,裂缝倾角为90°。
将以上参数带入到qP波相速度表达式(2)中,初始含气饱和度为100%,可以得到在给定参数下的qP波相速度随入射角和观测方位角的变化如图4所示。图中展示了入射角分别为17°(红线),32°(紫线)和47°(蓝线)时qP波相速度随观测方位角从0°-360°的变化规律。从图中可以看出,当裂缝倾角为90°时,随着观测方位角从0°到360°变化,qP波相速度呈现明显的方位各向异性特征,且入射角越大,方位各向异性越明显。
S4102)储层参数的影响
在给定参数的基础上,分别改变裂缝倾角、裂缝线密度、基质孔隙度、含气饱和度等参数,以分析不同储层物性参数对qP波相速度的影响。
S41021)裂缝倾角的影响
对不同裂缝倾角、地震波入射角以及观测方位角下的qP波相速度进行计算,结果如图5所示。其中,图5(a)-(i)分别对应裂缝倾角10°-90°,间隔为10°;每幅图都展示了入射角分别为17°(红线),32°(紫线)和47°(蓝线)时qP波相速度随观测方位角从0°-360°的变化规律。
从图5中可以看出,当入射角θ=17°时(红线),qP波相速度在所有裂缝倾角的情况下,几乎都关于观测方位角φ=180°对称,这表明当入射角较小时,各向异性对qP波相速度的影响很小。当入射角θ=47°时(紫线),qP波相速度在大部分情况下关于φ=180°是不对称的,且呈现以下规律:(1)当φ<180°时,裂缝倾角较小(<40°)时的qP波相速度随方位角出现先增大后减小的趋势(图5(a)-(d)),且在φ=90°时出现最大值,表明此时地震波传播方向与各向同性面(裂缝面)的平行方向最为接近;而当裂缝倾角较大时,地震波随方位角的增大逐渐透过裂缝面而后又返回,因此会出现qP波相速度先增大后减小然后又增大的趋势(图5(e)-(h)));(2)当φ>180°时,qP波相速度随方位角的增大出现先增大后减小的趋势(图5(a)-(i)),且在φ=270°时出现最小值,表明此时地震波传播方向与各向同性面(裂缝面)的垂直方向最为接近;(3)只有当裂缝倾角为90°时,qP波相速度是关于φ=180°完全对称的(图5(i))。因此,对于各向异性较强的海相碳酸盐岩储层,其实际地震资料的采集和处理应尽可能考虑大角度(大于30°)入射时的全方位(0°-360°)地震信息。
S41022)裂缝线密度的影响
图6给出了裂缝线密度分别为1条/m(图6(a))、15条/m(图6(b))和30条/m(图6(c))时的qP波相速度计算结果。从图中可以看到,随着裂缝线密度的增大,qP波相速度整体变小,且其方位各向异性逐渐增大。这是因为,随着裂缝线密度增大,等效介质的体积模量变小,会导致速度减小;而裂缝的增多会导致介质的等效各向异性增强。
S41023)基质孔隙度的影响
图7给出了基质孔隙度分别为1.12%(图7(a))、3.40%(图7(b))和5.81%(图7(c))时的qP波相速度计算结果。从图中可以看到,随着基质孔隙度的增大,qP波相速度整体变小,但其方位各向异性变化不大。这是由于,基质孔隙度增大,等效介质的弹性模量减小,会导致速度减小;但是各向异性不受基质孔隙度的影响。
S41024)含气饱和度的影响
图8给出了含气饱和度分别为0%(图8(a))、25%(图8(b))、50%(图8(c))、75%(图8(d))和100%(图8(e))时的qP波相速度计算结果。从图中可以看到,随着含气饱和度的增大,qP波相速度整体变小,但其方位各向异性变化不大。这是由于,随着当储层含气时,介质的等效体积模量减小,会导致速度减小;当储层含气后,混合流体的体积模量近似等于气体的体积模量,而等效介质的密度变小,从而导致整体的速度随着含气饱和度的增大而增大。但由于基质孔隙度较小(5.81%),含气饱和度对速度的影响并不明显。
S420)波场传播特征分析
S4201)基准模型结果分析
为了对目的层顶、底界面的地震波场传播特征进行分析,特设计三层理论模型。中间层是目的层,其参数与3.1中所给定的初始参数相同;模型顶层和底层均为各向同性弹性介质,通过测井资料给定其纵、横波速度和密度分别为:顶层,4729m/s,2606m/s,2600kg/m3;底层,5660m/s,3364m/s,2730kg/m3
基于上述模型参数,根据一阶速度-应力方程(3),采用基于PML边界条件的交错网格高阶有限差分方法,对上述三层模型进行正演数值模拟。震源采用爆炸震源,主频25Hz,时间间隔为0.4ms;模型大小为704m×704m×528m,网格间距为8m;交错网格差分法精度为O(Δt2,Δx16,Δy16,Δz16)。计算入射角分别为:0°、9°、17°、25°、32°、37°、43°、47°、51°、54°、57°、60°时的模型叠前全方位角道集数据,对于每一个入射角,方位角从0°到360°变化。图9给出了目的层顶界面(图9(a))、底界面(图9(b))的叠前全方位角道集数据,图10给出入射角分别为17°(红色)、32°(紫色)、47°(蓝色)时中间层顶界面(图10(a))、底界面(图10(b))的反射波振幅信息随观测方位角的变化曲线。
从图9和10中可以看出,当入射角较小时(小于30°),目的层顶、底界面反射波振幅的各向异性并不明显(图9,图10中红线),随着入射角的增大,目的层顶界面反射波振幅呈现明显各向异性的特点(图9,图10紫线),而底界面反射波的走时则呈现较强的各向异性特征(图9(b))。这表明,裂缝-多孔隙介质模型中纵波速度的各向异性会影响裂缝性储层的顶、底界面反射波走时和振幅的方位各向异性,可以利用qP波地震资料的这一特点对储层地震资料的各向异性时间、振幅属性进行分析,为储层的裂缝预测提供相应的理论依据。
S4202)储层参数的影响
在给定参数的基础上,保持其他参数不变,分别改变裂缝倾角、裂缝线密度、基质孔隙度和含气饱和度,以分析不同的参数对地震波场响应特征的影响。图11-14中,每组图的12幅图中,从上到下分别为目的层顶界面蜗牛道集、目的层底界面蜗牛道集、目的层顶界面反射波振幅随方位角的变化、目的层底界面反射波振幅随方位角的变化。从左向右,图11展示的是裂缝倾角分别为30°、60°、90°时的数值模拟结果,图12展示的是裂缝线密度分别为1条/m、15条/m、30条/m的数值模拟结果,图13展示的是基质孔隙度分别为1.12%、3.40%和5.81%的数值模拟结果,图14展示的是含气饱和度分别为0%、50%、100%时的数值模拟结果。
S42021)裂缝倾角的影响
从图11中可以看出,当裂缝倾角增大时,目的层顶、底界面的反射振幅的方位各向异性呈现增大的趋势(从左到右,反射波振幅随方位角变化的椭圆程度增强),且不再关于180°方位角对称,这是由于qP波相速度的变化引起的。
S42022)裂缝线密度的影响
从图12中可以看出,当裂缝线密度变化时,目的层顶、底界面反射波振幅的方位各向异性随裂缝线密度的减小而减小(从右至左,反射波振幅随方位角变化的椭圆程度降低)。
S42023)基质孔隙度的影响
从图13中可以看出,当基质孔隙度变化时,目的层顶、底界面反射波振幅的方位各向异性基本不变,但是反射波振幅强度增强。
S42024)含气饱和度的影响
从图14中可以看出,当含气饱和度增大时,目的层顶、底界面反射波振幅的方位各向异性基本不变,但是反射波振幅强度增强。
本发明将Norris-KG等效介质模型、纵横波分离算法、数值模拟计算三者相融合,实现了基于Norris-KG模型的地震纵波模拟分析。在给定参数下,通过对某研究区目的层的qP波相速度方位各向异性及其顶、底界面的反射波振幅方位各向异性分析,对目的层的各向异性及其影响参数进行优选。结果表明,在给定的岩石物理参数和测井参数条件下,该研究区目的层的观测各向异性较为明显,且主要受裂缝倾角和裂缝线密度的影响。
最后应说明的是:以上各实施例仅用以说明本发明的技术方案,而非对其限制;尽管参照前述各实施例对本发明进行了详细的说明,本领域的普通技术人员应当理解:其依然可以对前述各实施例所记载的技术方案进行修改,或者对其中部分或者全部技术特征进行等同替换;而这些修改或者替换,并不使相应技术方案的本质脱离本发明各实施例技术方案的范围。

Claims (5)

1.一种基于Norris-KG模型的地震纵波模拟分析方法,其特征在于,基于Norris-KG裂缝-孔隙等效介质模型,通过qP波相速度公式及纵、横波分离的一阶速度应力方程对介质的波场控制方程进行推导,并利用基于PML边界条件的交错网格高阶有限差分算法对地震波场进行正演数值模拟,最终,在给定参数下,对目的层的qP波相速度方位各向异性及其顶、底界面的反射波振幅方位各向异性分析,对目的层的各向异性及其影响参数进行优选分析,得到裂缝-孔隙介质地震纵波相速度和波场响应数据,最后进行纵波数值模拟分析。
2.如权利要求1所述的基于Norris-KG模型的地震纵波模拟分析方法,其特征在于,所述波场控制方程的推导过程为:
1.1)构建qP波相速度方程
假设裂缝倾角为θ0,地震波的传播方向为n=(nx,ny,nz)T=(sinθcosφ,sinθsinφ,cosθ)T,其中,θ为传播方向与观测系统z轴的夹角,φ为传播方向与x轴的夹角,任意入射角的纵波相速度表达式为:
其中,
上式(2)建立了裂缝-多孔隙介质的纵波相速度Vqp与储层弹性参数C之间的数学联系,结合Norris-KG模型所给出的介质弹性模量表达式(1),通过数值计算,得到介质的纵波相速度随储层物性参数的变化规律;
1.2)构建纵横波分离的一阶速度-应力方程
基于“方程分离”的思想,将弹性刚度系数分量中的耦合项C13近似为纵波模量和横波模量的线性函数,推导得到了各向异性介质的纵、横波分离的一阶速度-应力方程如下:
其中,上述方程(3)中下角标p和s分别代表纵波和横波,G为偏微分算子矩阵,Cp为与纵波相关的弹性系数矩阵,Cs为与横波相关的弹性系数矩阵;且Cp与Cs表达式分别如下:
上述的Cp与Cs表达式中,ε和γ为Thomsen各向异性参数,Vp0为qP波沿VTI介质对称轴方向的相速度,Vs0为qSV波沿VTI介质对称轴方向的相速度;和B=-2分别为纵波和横波模量的系数;
上述方程(3)建立了裂缝-多孔隙介质的地震纵波波场v与储层物性参数C之间的数学联系,通过对波动方程的数值求解,可得到介质的纵波地震响应vp,进而对储层纵波地震响应特征随物性参数的变化规律进行分析。
3.如权利要求2所述的基于Norris-KG模型的地震纵波模拟分析方法,其特征在于,所述纵波数值模拟分析的具体过程为:
2.1)qP波相速度分析
2.1.1)基准模型结果分析
基于任意入射角的纵波相速度公式(2),根据给定井段的岩石物理测量数据和测井数据,对目的层qP波的相速度及储层物性参数的影响进行分析;
2.1.2)储层参数的影响分析
在给定参数的基础上,分别改变裂缝倾角、裂缝线密度、基质孔隙度、含气饱和度等参数,以分析不同储层物性参数对qP波相速度的影响;
2.2)波场传播特征分析
2.2.1)基准模型结果分析
设计三层理论模型对目的层顶、底界面的地震波场传播特征进行分析,然后基于所述三层模型的参数,根据一阶速度-应力方程(3),采用基于PML边界条件的交错网格高阶有限差分方法,对所述三层模型进行正演数值模拟;
2.2.2)储层参数的影响分析
在给定参数的基础上,保持其他参数不变,分别改变裂缝倾角、裂缝线密度、基质孔隙度和含气饱和度,以分析不同的参数对地震波场响应特征的影响。
4.如权利要求2所述的基于Norris-KG模型的地震纵波模拟分析方法,其特征在于,所述正演数值模拟是采用交错网格高阶有限差分法对介质的一阶速度-应力方程进行数值求解;具体包括以下步骤:
3.1)方程的离散化处理
采用交错网格高阶有限差分算法,其时间、空间高阶差分近似公式分别为:
其中,上式(4)-(5)中Δt为时间步长,2M为时间差分精度,2L为空间差分精度;Δx为差分网格间距,为不同阶的交错网格空间差分系数,可以通过如下方程组求取:
将该方程组带入到上述方程(3)中得到一阶速度应力方程的差分近似格式;
3.2)PML边界条件处理
采用PML边界条件处理方法,对上述方程(3)中的总体波场v进行变量分离,每一个波场变量可以分为纵波分量vp和横波分量vs两部分,T是总的应力场,上角标x、y、z代表三个坐标轴方向:
然后对于所给定的三维空间模型设置边界衰减系数。
5.如权利要求1所述的基于Norris-KG模型的地震纵波模拟分析方法,其特征在于,所述Norris-KG裂缝-孔隙等效介质模型的具体构建过程为:基于Norris周期成层模型和KG裂缝-孔隙等效介质模型建立Norris-KG等效介质模型;基于Norris-KG等效介质模型,假设裂缝为非常薄且松软的无限大裂缝平面,将Norris-KG等效介质模型看成是周期成层两套裂缝-孔隙地层,根据Norris理论及Heuristic假设,Norris-KG等效介质模型的弹性属性表达式为:
上式(1)中,为垂直裂缝面方向纵波模量,/>为模型高频极限弹性系数,/>为垂直裂缝面方向的高频极限弹性系数,/>为模型低频极限弹性系数;式(1)给出的是水平裂缝性储层的弹性系数矩阵,可计算任意入射角情况的频变纵波速度随流体饱和度的变化情况。
CN202310364518.2A 2023-04-07 2023-04-07 一种基于Norris-KG模型的地震纵波模拟分析方法 Active CN116559941B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202310364518.2A CN116559941B (zh) 2023-04-07 2023-04-07 一种基于Norris-KG模型的地震纵波模拟分析方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202310364518.2A CN116559941B (zh) 2023-04-07 2023-04-07 一种基于Norris-KG模型的地震纵波模拟分析方法

Publications (2)

Publication Number Publication Date
CN116559941A true CN116559941A (zh) 2023-08-08
CN116559941B CN116559941B (zh) 2024-03-12

Family

ID=87495530

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202310364518.2A Active CN116559941B (zh) 2023-04-07 2023-04-07 一种基于Norris-KG模型的地震纵波模拟分析方法

Country Status (1)

Country Link
CN (1) CN116559941B (zh)

Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20110007604A1 (en) * 2009-07-10 2011-01-13 Chevron U.S.A. Inc. Method for propagating pseudo acoustic quasi-p waves in anisotropic media
US20170059742A1 (en) * 2014-02-19 2017-03-02 Repsol, S.A. Method Implemented In A Computer For The Numerical Simulation Of A Porous Medium
CN106772608A (zh) * 2017-02-16 2017-05-31 甘肃省地震局 等效孔隙裂缝介质的弹性阻抗及广义流体因子分析方法
CN109946742A (zh) * 2019-03-29 2019-06-28 中国石油大学(华东) 一种TTI介质中纯qP波地震数据模拟方法
CN110596754A (zh) * 2019-09-24 2019-12-20 中国矿业大学(北京) 一种三维TTI介质qP波与qSV波波场模拟方法
RU2732035C1 (ru) * 2020-01-28 2020-09-10 Публичное акционерное общество «Татнефть» имени В.Д. Шашина Способ определения трещинной пористости пород
CN115600373A (zh) * 2022-09-15 2023-01-13 山东科技大学(Cn) 黏滞各向异性介质qP波模拟方法、系统、设备及应用

Patent Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20110007604A1 (en) * 2009-07-10 2011-01-13 Chevron U.S.A. Inc. Method for propagating pseudo acoustic quasi-p waves in anisotropic media
US20170059742A1 (en) * 2014-02-19 2017-03-02 Repsol, S.A. Method Implemented In A Computer For The Numerical Simulation Of A Porous Medium
CN106772608A (zh) * 2017-02-16 2017-05-31 甘肃省地震局 等效孔隙裂缝介质的弹性阻抗及广义流体因子分析方法
CN109946742A (zh) * 2019-03-29 2019-06-28 中国石油大学(华东) 一种TTI介质中纯qP波地震数据模拟方法
CN110596754A (zh) * 2019-09-24 2019-12-20 中国矿业大学(北京) 一种三维TTI介质qP波与qSV波波场模拟方法
RU2732035C1 (ru) * 2020-01-28 2020-09-10 Публичное акционерное общество «Татнефть» имени В.Д. Шашина Способ определения трещинной пористости пород
CN115600373A (zh) * 2022-09-15 2023-01-13 山东科技大学(Cn) 黏滞各向异性介质qP波模拟方法、系统、设备及应用

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
XIAOHUI YANG 等: "Frequency-Dependent Amplitude Versus Offset Variations in Porous Rocks with Aligned Fractures", 《PURE AND APPLIED GEOPHYSICS》, vol. 174, pages 1043 - 1059, XP036183445, DOI: 10.1007/s00024-016-1423-8 *
孔丽云 等: "基于两级尺度粗化的裂缝-多孔隙介质地震响应分析方法", 《地球物理学报》, vol. 61, no. 9, pages 3791 - 3799 *
韩令贺;何兵寿;: "VTI介质一阶准P波方程正演模拟及边界条件", 石油地球物理勘探, vol. 45, no. 06, pages 819 - 825 *

Also Published As

Publication number Publication date
CN116559941B (zh) 2024-03-12

Similar Documents

Publication Publication Date Title
Ladopoulos Petroleum and gas reserves exploration by real-time expert seismology and non-linear seismic wave motion'
Virieux P-SV wave propagation in heterogeneous media: Velocity-stress finite-difference method
Carcione et al. Computational poroelasticity—A review
Kostin et al. Local time–space mesh refinement for simulation of elastic wave propagation in multi-scale media
Shragge Reverse time migration from topography
CN110058303B (zh) 声波各向异性逆时偏移混合方法
EP1820137A2 (en) Integrated anisotropic rock physics model
CN108845351A (zh) 一种vsp地震资料转换波全波形反演方法
CN110286410B (zh) 基于绕射波能量的裂缝反演方法和装置
Wang et al. Massively parallel structured direct solver for equations describing time-harmonic qP-polarized waves in TTI media
CN111505714B (zh) 基于岩石物理约束的弹性波直接包络反演方法
CN113536638B (zh) 基于间断有限元和交错网格的高精度地震波场模拟方法
CN114895351A (zh) 模拟地震波在任意不连续界面处传播的介质模型化的方法及装置
CN111273341B (zh) 根据裂缝空间分布的含裂缝储层岩石物理建模方法
Masson Distributional finite-difference modelling of seismic waves
CN105093351A (zh) 识别储层微裂缝的方法
NO20190489A1 (en) Seismic modeling
Hanyga et al. Wave field simulation for heterogeneous transversely isotropic porous media with the JKD dynamic permeability
Ladopoulos Non-linear Seismic Wave Motion in Elastodynamics with Application to Real-Time Expert Seismology
Ladopoulos Non-linear Elastodynamics by Seismic Wave Motion in Real-Time Expert Seismology
Mesgouez et al. Transient solution for multilayered poroviscoelastic media obtained by an exact stiffness matrix formulation
CN116559941B (zh) 一种基于Norris-KG模型的地震纵波模拟分析方法
CN109521470B (zh) 分析地质构造对地震反演裂缝密度影响的方法
CN104483702B (zh) 一种适用于非均匀运动水体的地震正演模拟方法
Käser et al. Wavefield modeling in exploration seismology using the discontinuous Galerkin finite-element method on HPC infrastructure

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