CN106054242B - 三维各向异性衰减介质波场模拟方法 - Google Patents

三维各向异性衰减介质波场模拟方法 Download PDF

Info

Publication number
CN106054242B
CN106054242B CN201610289461.4A CN201610289461A CN106054242B CN 106054242 B CN106054242 B CN 106054242B CN 201610289461 A CN201610289461 A CN 201610289461A CN 106054242 B CN106054242 B CN 106054242B
Authority
CN
China
Prior art keywords
wave
dimensional
indicate
field
attenuation medium
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
CN201610289461.4A
Other languages
English (en)
Other versions
CN106054242A (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 Geosciences Beijing
Original Assignee
China University of Geosciences 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 Geosciences Beijing filed Critical China University of Geosciences Beijing
Priority to CN201610289461.4A priority Critical patent/CN106054242B/zh
Publication of CN106054242A publication Critical patent/CN106054242A/zh
Application granted granted Critical
Publication of CN106054242B publication Critical patent/CN106054242B/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
    • G01V1/282Application of seismic models, synthetic seismograms
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01VGEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V2210/00Details of seismic processing or analysis
    • G01V2210/60Analysis
    • G01V2210/67Wave propagation modeling

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

本发明公开一种三维各向异性衰减介质波场模拟方法,包括:采用三维插值法,建模三维地质模型;定义震源类型;定义三维三分量观测***;根据所述三维地质模型、所述震源类型和所述三维三分量观测***,计算各向异性衰减波场,以产生地震炮记录和波场快照;输出所述地震炮记录和波场快照。通过本发明,以解决现有技术存在的介质各向异性衰减三维模拟的问题。

Description

三维各向异性衰减介质波场模拟方法
技术领域
本发明涉及地震勘探的技术领域,尤其涉及一种三维各向异性衰减介质波场模拟方法。
背景技术
地下蕴藏的油气大部分是裂缝型油气藏,裂缝不仅提供了油气储集空间,也提供了油气运移的通道,因此裂缝检测对油气勘探来说极其重要,近年来,油气勘探的重点方向,例如,裂缝型油气藏的油水检测以及页岩气等非常规资源勘探,都离不开对裂缝的精确分析,因此现阶段对于裂缝的研究显得尤为重要和迫切。采用等效介质理论,水平或垂直裂缝可以采用VTI/HTI介质描述,该类介质仅考虑了地震波传播的速度各向异性,而当考虑流体填充的裂缝储层时,VTI/HTI各向异性衰减是一种更为近似的描述方法,目前国内外在各向异性衰减三维波场计算和分析方面的研究尚待进一步发展。因此,针对各向异性衰减的波场模拟显得尤为重要。
发明内容
本发明的主要目的在于提供一种三维各向异性衰减介质波场模拟方法,以解决现有技术存在的介质各向异性衰减三维模拟的问题。
为解决上述问题,本发明实施例提供一种三维各向异性衰减介质波场模拟方法,包括:采用三维插值法,建模三维地质模型;定义震源类型;定义三维三分量观测***;根据所述三维地质模型、所述震源类型和所述三维三分量观测***,计算各向异性衰减波场,以产生地震炮记录和波场快照;输出所述地震炮记录和波场快照。
根据本发明的技术方案,通过建模三维地质模型、定义震源类型及三维三分量观测***,再根据所述三维地质模型、所述震源类型和所述三维三分量观测***,计算各向异性衰减波场,以产生地震炮记录和波场快照。如此一来,有效地解决了介质各向异性衰减三维模拟的问题,尤其适用于油气水识别分析、非常规油气资源模型正演分析工作,还能够快速实现三维各向异性衰减波场模拟,准确记录数值模拟结果,并提高地质解释的精度。
附图说明
此处所说明的附图用来提供对本发明的进一步理解,构成本申请的一部分,本发明的示意性实施例及其说明用于解释本发明,并不构成对本发明的不当限定。在附图中:
图1是根据本发明实施例的三维各向异性衰减介质波场模拟方法的流程图;
图2是根据本发明实施例的三维地质模型的示意图;
图3是根据本发明实施例的地震子波的示意图;
图4是根据本发明实施例的三维观测***的示意图;
图5a、图5b及图5c分别是根据本发明实施例的三维三分量各向异性模型模拟记录的示意图;
图6a、图6b及图6c分别是根据本发明实施例的三个分量波场快照的示意图。
具体实施方式
本发明的主要思想在于,基于建模三维地质模型、定义震源类型及三维三分量观测***,再根据所述三维地质模型、所述震源类型和所述三维三分量观测***,计算各向异性衰减波场,以产生地震炮记录和波场快照。如此一来,有效地解决了介质各向异性衰减三维模拟的问题,尤其适用于油气水识别分析、非常规油气资源模型正演分析工作,还能够快速实现三维各向异性衰减波场模拟,准确记录数值模拟结果,并提高地质解释的精度。
为使本发明的目的、技术方案和优点更加清楚,以下结合附图及具体实施例,对本发明做进一步地详细说明。
根据本发明的实施例,提供了一种三维各向异性衰减介质波场模拟方法。
图1是根据本发明实施例的三维各向异性衰减介质波场模拟方法的流程图。
在步骤S102中,采用三维插值法,建模三维地质模型。在本实施例中,三维地质模型可包括速度模型或弹性系数矩阵。其中,所述速度模型包括模型参数,例如纵波速度、横波速度、密度、Thomson各向异性参数、纵波品质因子和横波品质因子。而弹性系数矩阵例如包括描述各向异性衰减的6*6弹性系数矩阵,且弹性系数矩阵的模型参数由HTI衰减的弹性系数给定,如图2所示。
进一步的,所述三维地质模型满足如下公式:
其中,V表示介质速度,ρ是密度,L是散度运算符,F表示体力,σ是应力矢量,ε是应变矢量,E是记忆变量,是各向异性衰减介质的弹性系数矩阵,
其中,E包括:
其中,E={e1,e11,e22,e23,e13,e12}T,上标”.”表示时间的一阶导数,上标“..”表示时间的二阶导数,表示变量的方向导数,函数φα表达式如下:
其中,τσα和τεα分别表示应力和应变松弛时间,上角标1表示纵波,2表示横波,vi,i=x,y,z表示质点速度分量,函数Θ表达式如下:
在步骤S104中,定义震源类型。在本实施例中,地震子波是采用Ricker子波,并可自定义子波的主频、传播时间和震源类型。其中,震源类型例如可以定义为***震源或者点震源,而震源的选取可以根据实际数值模拟需求来定义。震源类型满足如下公式:
其中,表Zk示体力项,Z表示体力的垂向分j量,表示垂向单位向量,φ表示一标量函数,所述标量函数是关于时间的函数, 分别为x和y方向的单位向量,x、y、z分别表示相互垂直的三维坐标轴。
在步骤S106中,定义三维三分量观测***。在本实施例中,三维三分量观测***是各向异性衰减研究的关键,根据目标体的特点,可以预先设计合理的观测***,计算在预设观测点观测到的数值地震记录后,采集地震信息,制作出与实际地震资料品质相适应的地震记录,便于对比分析。進一步的,三维三分量观测***可在主测线(Inline)和联络测线(Crossline)方向分别定义,给定接收线起始位置、线间隔、以及接收线条数后,可以确定接收线位置。并且,接收线道间距和接收点个数也可以在此时给定,以确定检波器关系。
在步骤S108中,根据所述三维地质模型、所述震源和所述三维三分量观测***,计算各向异性衰减波场,以产生地震炮记录和波场快照。在本实施例中,计算各向异性衰减波场是采用傅里叶伪谱法、时间四阶龙格-库塔法。也就是说,将地质模型、震源、观测***输入,作为模拟***的初始条件,采用傅里叶伪谱法、时间四阶龙格-库塔法计算各向异性衰减波场,以产生地震炮记录和波场快照。其中,傅里叶虚谱法的运算分两步,首先利用傅氏变换将波场函数表示成傅里叶级数的展开形式,然后在时间-波数域(或时间-频率域)中对波动方程进行数值求解。
另外,四阶龙格-库塔法可以实现各向异性衰减介质波动方程的时间导数计算,其目的是用来逼近泰勒级数解的一种算法,和泰勒展开式不同之处在于,龙格一库塔法除了一阶导数外,不用计算其高阶导数,代替的是算几个一阶导数,用这几个一阶导数的线性规律组合来逼近泰勒展开方程,去掉了计算高阶导数的繁琐。即可以通过简单的计算来实现。该方法是稳定的和递推的,也就是说,计算时只要用到前面一个点就能计算出该点的函数值。四阶龙格-库塔法满足如下公式:
其中,
H1=NVn+Dn
H2=N(Vn+dt/2H1)+Dn+1/2
H3=N(Vn+dt/2H2)+Dn+1/2
H4=N(Vn+dtH3)+Dn+1
Vn表示第n阶波场,Vn+1表示第n+1阶波场,dt表示时间采样间隔,H1表示时间段开始时的斜率,H2表示时间段中点的斜率,H3表示时间段中点的斜率,H4表示时间段终点的斜率,N表示速度在空间三个方向的导数算子,Dn表示各向异性衰减方程中常量在开始时刻(第n时刻)的数值,Dn+1/2表示常量在时间段中点(第n+1/2时刻)的数值,Dn+1表示常量在时间段终点(第n+1时刻)的数值。
进一步的,可将卷积完全匹配层(CPML)吸收边界引入到三维各向异性衰减介质数值模拟,且通过设置适当的参数值CPML。该方法引入复频域变换使离散方程与模拟的介质类型无关,是一种无需显式卷积计算的不***的迭代格。该方法通过坐标扩展函数将实数空间坐标变换为复数坐标,并将扩展函数的极点移动到复平面的虚轴上,对常规完全匹配层进行了改善。该方法引入虚拟的介质层,该层介质的波阻抗与相邻介质的波阻抗完全匹配,因而理论上人射波可以无反射的进人该层,并被完全吸收。该方法不仅吸收效果明显、易于实现,而且比传统PML明显节省存储量。
在步骤S110中,输出所述地震炮记录和波场快照。如此一来,本实施例能够快速实现三维各向异性衰减波场模拟,准确记录数值模拟结果。
上述已说明了三维各向异性衰减介质波场模拟方法,以下将提供一些实例来验证上述方法的适用性及准确性。
假设震源类型为***震源,且地震子波主频设定为50Hz,时间延迟70ms,如图3所示。进一步,根据HTI衰减各向异性特性,可以定义出三维三分量观测***,如图4所示。图4给出了0度和90度的检波线,两条检波线(即主侧线和联络测线)相互垂直,可以接收垂直裂缝和平行裂缝方向的快慢横波,检波器间隔是10m。将设置好的三维地质模型、震源类型和三维三分量观测***运行模拟,其中每次运行一个时间采样,模拟会在终端显示运行情况,给出最直观的运行情况。
可选地,根据模拟需求可以设置时间采样间隔,默认值是1ms,时间采样大小可根据情况设置。同样,可以设置波场快照的间隔,***默认值是100ms,即每隔100ms,模拟***自动存储地震波的波场快照,存储文件以波场快照时间来命名,显示快捷方便。模拟完成后,可以显示三维三分量各向异性模型模拟记录,其为Z分量、R分量和T分量的三分量地震炮记录,如图5a、图5b和图5c所示。并且,波场快照包含了时间起始时间到终止时间的全部结果,利用X分量、Y分量、Z分量,可以显示地震波场在某时刻的传播情况,如图6a、图6b和图6c所示。
综上所述,根据本发明的技术方案,通过建模三维地质模型、定义震源类型及三维三分量观测***,再根据所述三维地质模型、所述震源类型和所述三维三分量观测***,计算各向异性衰减波场,以产生地震炮记录和波场快照。如此一来,有效地解决了介质各向异性衰减三维模拟的问题,尤其适用于油气水识别分析、非常规油气资源模型正演分析工作,还能够快速实现三维各向异性衰减波场模拟,准确记录数值模拟结果,并提高地质解释的精度。
显然,本领域的技术人员应该明白,上述的本发明的三维各向异性衰减波场模拟方法和各步骤可以用通用的计算装置来实现,它们可以集中在单个的计算装置上,或者分布在多个计算装置所组成的网络上,可选地,它们可以用计算装置可执行的程序代码来实现,从而,可以将它们存储在存储装置中由计算装置来执行。这样,本发明不限制于任何特定的硬件和软件结合。存储装置为非易失性存储器,如:ROM/RAM、闪存、磁碟、光盘等。
以上所述仅为本发明的实施例而已,并不用于限制本发明,对于本领域的技术人员来说,本发明可以有各种更改和变化。凡在本发明的精神和原则之内,所作的任何修改、等同替换、改进等,均应包含在本发明的权利要求范围之内。

Claims (5)

1.一种三维各向异性衰减介质波场模拟方法,其特征在于,包括以下步骤:
采用三维插值法,建模三维地质模型,所述三维地质模型包括速度模型或弹性系数矩阵;
定义震源类型,所述震源类型满足如下公式:
其中,表示体力项,Z表示体力的垂向分量,表示垂向单位向量,φ表示一标量函数,所述标量函数是关于时间的函数,分别为x和y方向的单位向量,x、y、z分别表示相互垂直的三维坐标轴;
定义三维三分量观测***;
根据所述三维地质模型、所述震源类型和所述三维三分量观测***,计算各向异性衰减波场,以产生地震炮记录和波场快照;
输出所述地震炮记录和波场快照。
2.根据权利要求1所述的三维各向异性衰减介质波场模拟方法,其特征在于,所述计算各向异性衰减波场是采用傅里叶伪谱法、时间四阶龙格-库塔法。
3.根据权利要求1所述的三维各向异性衰减介质波场模拟方法,其特征在于,所述速度模型的模型参数包括:纵波速度、横波速度、密度、Thomson各向异性参数、纵波品质因子和横波品质因子。
4.根据权利要求1所述的三维各向异性衰减介质波场模拟方法,其特征在于,所述弹性系数矩阵的模型参数由HT I衰减的弹性系数给定。
5.根据权利要求1所述的三维各向异性衰减介质波场模拟方法,其特征在于,所述三维地质模型满足如下公式:
其中,V表示介质速度,ρ是密度,L是散度运算符,F表示体力,σ是应力矢量,ε是应变矢量,E是记忆变量,是各向异性衰减介质的弹性系数矩阵,
其中,E包括:
其中,E={e1,e11,e22,e23,e13,e12}T,上标“ .”表示时间的一阶导数,上标“..”表示时间的二阶导数,表示变量的方向导数,函数φα表达式如下:
其中,τσα和τεα分别表示应力和应变松弛时间,上角标(1)表示纵波,上角标(2)表示横波,vi,i=x,y,z表示质点速度分量,函数Θ表达式如下:
CN201610289461.4A 2016-05-04 2016-05-04 三维各向异性衰减介质波场模拟方法 Active CN106054242B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201610289461.4A CN106054242B (zh) 2016-05-04 2016-05-04 三维各向异性衰减介质波场模拟方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201610289461.4A CN106054242B (zh) 2016-05-04 2016-05-04 三维各向异性衰减介质波场模拟方法

Publications (2)

Publication Number Publication Date
CN106054242A CN106054242A (zh) 2016-10-26
CN106054242B true CN106054242B (zh) 2018-11-13

Family

ID=57176680

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201610289461.4A Active CN106054242B (zh) 2016-05-04 2016-05-04 三维各向异性衰减介质波场模拟方法

Country Status (1)

Country Link
CN (1) CN106054242B (zh)

Families Citing this family (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN106959469A (zh) * 2017-04-14 2017-07-18 中国石油天然气股份有限公司 地震波的速度及衰减模拟分析方法及装置
CN110297272A (zh) * 2019-07-05 2019-10-01 中南大学 三维速度模型生成方法、装置、设备及存储介质
CN113341455B (zh) * 2021-06-24 2024-02-09 中国石油大学(北京) 一种粘滞各向异性介质地震波数值模拟方法、装置及设备
CN113960663A (zh) * 2021-10-25 2022-01-21 中国地质大学(北京) 基于并行计算的三维各向异性衰减正演模拟方法及***

Family Cites Families (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN1230260A (zh) * 1996-09-13 1999-09-29 Pgs张量公司 时滞储层监控方法
CN101369024A (zh) * 2008-09-09 2009-02-18 中国石油天然气集团公司 一种地震波动方程生成方法及***
CN103364835A (zh) * 2013-07-01 2013-10-23 西安交通大学 一种地层结构自适应中值滤波方法
CN104280768B (zh) * 2013-07-12 2017-03-15 中国石油天然气集团公司 一种适用于逆时偏移的吸收边界条件方法
CN104166159B (zh) * 2014-07-15 2015-10-21 刘改成 四维微地震监测的裂缝形态处理方法和***
CN105158797B (zh) * 2015-10-16 2018-02-09 成都理工大学 一种基于实际地震资料的交错网格波动方程正演的方法

Also Published As

Publication number Publication date
CN106054242A (zh) 2016-10-26

Similar Documents

Publication Publication Date Title
CN108181652B (zh) 一种海底节点地震资料上下行波场数值模拟方法
Regone et al. Geologic model building in SEAM Phase II—Land seismic challenges
CN106054242B (zh) 三维各向异性衰减介质波场模拟方法
CN103149585B (zh) 一种弹性偏移地震波场构建方法及装置
CN106772578B (zh) 一种合成地震记录的方法和装置
CN108139499A (zh) Q-补偿的全波场反演
CN106772579B (zh) 一种薄煤层中地震叠前反演方法和装置
CN104570082B (zh) 一种基于格林函数表征的全波形反演梯度算子的提取方法
CN106526674A (zh) 一种三维全波形反演能量加权梯度预处理方法
CN106597545B (zh) 一种水平裂缝地震叠前反演方法和装置
CN104459770B (zh) 一种高维地震数据规则化方法
CN109725345A (zh) 一种初至波正演模拟方法及装置
CN107390270B (zh) 一种基于弹性波逆时偏移ADCIGs的AVA分析方法
CN106415321A (zh) 用于储层建模的基于瞬时等时属性的地质体识别
CN106501872B (zh) 一种裂缝储层地应力特征的计算方法及装置
CN103758511A (zh) 一种井下逆时偏移成像识别隐蔽储层的方法及装置
CN105911584A (zh) 一种隐式交错网格有限差分弹性波数值模拟方法及装置
AU2011380936B2 (en) Seismic imaging systems and methods employing correlation-based stacking
CN110850469A (zh) 一种基于克希霍夫积分解的地震槽波深度偏移的成像方法
CN105807317B (zh) 基于切比雪夫伪谱法的各向异性衰减面波模拟方法
CN104597489A (zh) 一种震源子波优化设置方法和装置
CN110687598A (zh) 一种加速微震数值模拟的方法及装置
CN109738944A (zh) 基于广角反射的地震采集参数确定方法及装置
Ribeiro et al. Spatial analysis of oil reservoirs using detrended fluctuation analysis of geophysical data
CN107807392A (zh) 一种自适应抗频散的分块时空双变逆时偏移方法

Legal Events

Date Code Title Description
C06 Publication
PB01 Publication
C10 Entry into substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant