CN112285772A - 有限差分数值模拟方法、系统、介质、计算机设备及应用 - Google Patents
有限差分数值模拟方法、系统、介质、计算机设备及应用 Download PDFInfo
- Publication number
- CN112285772A CN112285772A CN202011067839.9A CN202011067839A CN112285772A CN 112285772 A CN112285772 A CN 112285772A CN 202011067839 A CN202011067839 A CN 202011067839A CN 112285772 A CN112285772 A CN 112285772A
- Authority
- CN
- China
- Prior art keywords
- wave
- difference
- numerical simulation
- finite difference
- time
- 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
Links
- 238000004088 simulation Methods 0.000 title claims abstract description 94
- 238000000034 method Methods 0.000 title claims abstract description 88
- 239000006185 dispersion Substances 0.000 claims abstract description 28
- 238000000926 separation method Methods 0.000 claims abstract description 21
- 238000013508 migration Methods 0.000 claims abstract description 6
- 230000005012 migration Effects 0.000 claims abstract description 6
- 230000002441 reversible effect Effects 0.000 claims abstract description 6
- 238000005070 sampling Methods 0.000 claims abstract description 5
- 238000004364 calculation method Methods 0.000 claims description 7
- 230000006870 function Effects 0.000 claims description 6
- 238000009795 derivation Methods 0.000 claims description 5
- 238000004590 computer program Methods 0.000 claims description 4
- 230000010287 polarization Effects 0.000 claims description 3
- 238000013459 approach Methods 0.000 abstract description 3
- 238000012545 processing Methods 0.000 abstract description 2
- 230000002123 temporal effect Effects 0.000 description 6
- 238000007796 conventional method Methods 0.000 description 5
- 229910003460 diamond Inorganic materials 0.000 description 5
- 239000010432 diamond Substances 0.000 description 5
- 230000007547 defect Effects 0.000 description 3
- 238000010586 diagram Methods 0.000 description 3
- 238000003491 array Methods 0.000 description 2
- 230000000694 effects Effects 0.000 description 2
- 238000002474 experimental method Methods 0.000 description 2
- 238000003384 imaging method Methods 0.000 description 2
- 238000011160 research Methods 0.000 description 2
- 238000004458 analytical method Methods 0.000 description 1
- 238000004422 calculation algorithm Methods 0.000 description 1
- 238000001514 detection method Methods 0.000 description 1
- 238000012821 model calculation Methods 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 230000003287 optical effect Effects 0.000 description 1
- 238000005457 optimization Methods 0.000 description 1
- 230000008092 positive effect Effects 0.000 description 1
- 238000012827 research and development Methods 0.000 description 1
- 239000004576 sand Substances 0.000 description 1
- 239000004065 semiconductor Substances 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
- G01V1/28—Processing seismic data, e.g. for interpretation or for event detection
- G01V1/30—Analysis
- G01V1/306—Analysis for determining physical properties of the subsurface, e.g. impedance, porosity or attenuation profiles
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
- G01V1/28—Processing seismic data, e.g. for interpretation or for event detection
- G01V1/30—Analysis
- G01V1/307—Analysis for determining seismic attributes, e.g. amplitude, instantaneous phase or frequency, reflection strength or polarity
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V2210/00—Details of seismic processing or analysis
- G01V2210/60—Analysis
- G01V2210/62—Physical property of subsurface
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V2210/00—Details of seismic processing or analysis
- G01V2210/60—Analysis
- G01V2210/63—Seismic attributes, e.g. amplitude, polarity, instant phase
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)
- Management, Administration, Business Operations System, And Electronic Commerce (AREA)
Abstract
本发明属于地球物理勘探数据处理技术领域，公开了一种有限差分数值模拟方法、系统、介质、计算机设备及应用，通过时间任意偶数阶精度有限差分逼近时间偏导数、空间隐式交错网格有限差分逼近空间偏导数、纵横波波场分离实现弹性P‑SV波高精度数值模拟。本发明的弹性P‑SV波高精度有限差分数值模拟方法，采用“十字+菱形”交错网格有限差分模板逼近弹性波方程中的时间偏导数；采样隐式交错网格有限差分法逼近弹性波方程中的空间偏导数；推导了相应的频散关系；基于Taylor级数展开求取差分系数；基于纵横波分离实现了P‑SV波高精度模拟。本发明可以大幅提高P‑SV波的模拟精度，提高弹性波逆时偏移和全波形反演的精度。
Description
技术领域
本发明属于地球物理勘探数据处理技术领域，尤其涉及一种有限差分数值模拟方法、系统、介质、计算机设备及应用。
背景技术
目前：波动方程数值模拟是了解地震波在地下的传播规律、帮助解释观测数据的有效手段，也是逆时偏移成像和全波形反演的基本单元。到目前为止，波动方程数值模拟已被广泛应用于地球物理勘探领域的各个环节。常用的偏微分方程数值解法包括有限元素法、有限差分法、伪谱法等。有限差分法容易实现、计算效率高，是最常用波动方程数值模拟算法。地球物理学家对有限差分法进行过大量研究，先后提出显式有限差分和隐式有限差分；规则网格有限差分、交错网格有限差分和旋转交错网格有限差分；基于Taylor级数展开的有限差分和基于最优化的有限差分；空间域有限差分和时空域有限差分等。上述有限差分法大多是基于常规“十字”型差分模版设计的，时间精度始终只有二阶。当模拟时间较长或时间步长较大时，这些差分方法都会产生较强的时间频散，模拟精度不够。
与常规方法相比，时间高阶有限差分法可以更好地压制时间频散，且具有更好的稳定性。Lax-Wendroff方法是最早出现的时间高阶有限差分法，该方法通过空间导数来代替高阶时间导数来提高时间精度。Lax-Wendroff法需要存储多个时刻的波场。随着时间精度的增加，存储量和计算量急剧增大，现有的计算设备很难承受，而且还会出现数值不稳定问题，但Lax-Wendroff法很难推广到任意偶数阶精度、且计算量较大。研究发现有限差分法的时间模拟精度和所采用的差分模版种类密切相关。随后出现了基于“菱形”模版的声波方程规则网格有限差分法，该方法可以同时获得任意偶数阶空间和时间精度。地球物理工作者还将常规“十字型”模板和“菱形”模版相结合来获得空间(2M)阶、时间(2N)阶精度。现有的时间高阶有限差分法多数只适用于声波方程，无法直接用于弹性波方程数值模拟，且空间模拟精度仍然不够，现有的时间高阶有限差分法主要提高时间精度，空间精度仍需进一步提高。
通过上述分析，现有技术存在的问题及缺陷为：
(1)现有时间高阶有限差分法很难推广到任意偶数阶精度、且计算量较大。
(2)现有的时间高阶有限差分法多数只适用于声波方程，无法直接用于弹性波方程数值模拟，且空间模拟精度仍然不够。
解决以上问题及缺陷的难度为：弹性介质中，涉及纵波、横波和转换波等波模式。采用基于纵波或横波频散关系的差分系数进行求解只能分别得到精确的纵波或横波波场，都无法得到高精度的P-SV转换波。提高空间精度会一定程度降低时间精度和稳定性。
解决以上问题及缺陷的意义为：研发弹性P-SV波有限差分数值模拟方法可以得到大幅提高纵波、横波和P-SV转换波的模拟精度，进而改善弹性波逆时偏移和全波形反演的成像和反演效果。
发明内容
针对现有技术存在的问题，本发明提供了一种有限差分数值模拟方法、系统、介质、计算机设备及应用。
本发明是这样实现的，一种有限差分数值模拟方法，所述有限差分数值模拟方法包括：
采用十字和菱形交错网格有限差分模板逼近时间偏导数；
采用隐式交错网格有限差分法逼近空间偏导数；
频散关系推导；
差分系数求取；
得到的时间高阶差分系数和隐式空间差分系数分别代入公式进行弹性波方程数值模拟，实现基于纵横波分离的P-SV波高精度模拟。
进一步，所述十字和菱形交错网格有限差分模板逼近时间偏导数方法包括：
二维弹性波速度应力方程表示为：
其中，(vx,vz)为偏振速度矢量(τxx,τzz,τxz)，λ(x,z)和μ(x,z)为拉梅系数，ρ(x,z)为密度；采用十字+菱形交错网格有限差分模板逼近时间偏导数，公式如下：
其中，
其中dm,0和dm,n为时间差分系数，Δt为时间步长，h为空间采样间隔，M和N为差分算子长度参数。
进一步，所述隐式交错网格有限差分法逼近空间偏导数方法包括：隐式交错网格有限差分逼近空间偏导数的公式为：
其中cm和b为空间差分系数：
隐式交错网格有限差分中，空间偏导数通过求解三对角线性方程组得到。
进一步，所述频散关系推导方法包括：
平面波理论下，波场表示公式为：
其中，
为虚数单位，k为波数，ω为角频率，θ为传播角度，kx＝kcos(θ)，kz＝ksin(θ)；将波场表示公式代入公式
并化简得到：
其中，rP＝vPΔt/h和rS＝vSΔt/h为纵波和横波的库朗数，vP和vS分别为纵波和横波的传播速度：
公式
和
为纵波和横波的时间频散关系。
进一步，所述差分系数求取方法包括：
将公式
中的三角函数进行Taylor级数展开化简可得：
C+D＝E；
其中，
比较方程C+D＝E两端(kxh)2l或(kzh)2l(1≤l≤M)的系数得：
a1＝0,
比较方程C+D＝E两端
的系数得：
其中，int(x)为取整函数，通过公式：
a1＝0,
和
得到基于Taylor级数展开的纵波时间高阶差分系数，将vP换为vS得到横波时间高阶差分系数；
隐式空间交错网格有限差分系数求取公式为：
进一步，所述基于纵横波分离的P-SV波高精度模拟方法包括：将第四步中得到的时间高阶差分系数和隐式空间差分系数分别代入公式
和
进行弹性波方程数值模拟，采用波数域分离方法得到高精度的P-SV转换波精度模拟，模拟步骤为：
(1)采用纵波时间高阶差分系数和隐式空间差分系数求解弹性波方程，纵横波分离，只提取纵波场；
(2)采用横波时间高阶差分系数和隐式空间差分系数求解弹性波方程，纵横波分离，只提取横波场；
(3)将(1)中的纵波场和(2)中的横波场求和，得到最终的弹性波场；
(4)重复(1)-(3)，直到最大时刻。
本发明的另一目的在于提供一种计算机设备，所述计算机设备包括存储器和处理器，所述存储器存储有计算机程序，所述计算机程序被所述处理器执行时，使得所述处理器执行如下步骤：
采用十字和菱形交错网格有限差分模板逼近时间偏导数；
采用隐式交错网格有限差分法逼近空间偏导数；
频散关系推导；
差分系数求取；
得到的时间高阶差分系数和隐式空间差分系数分别代入公式进行弹性波方程数值模拟，实现基于纵横波分离的P-SV波高精度模拟。
本发明的另一目的在于提供一种计算机可读存储介质，存储有计算机程序，所述计算机程序被处理器执行时，使得所述处理器执行如下步骤：
采用十字和菱形交错网格有限差分模板逼近时间偏导数；
采用隐式交错网格有限差分法逼近空间偏导数；
频散关系推导；
差分系数求取；
得到的时间高阶差分系数和隐式空间差分系数分别代入公式进行弹性波方程数值模拟，实现基于纵横波分离的P-SV波高精度模拟。
本发明的另一目的在于提供一种实施所述有限差分数值模拟方法的有限差分数值模拟系统，所述有限差分数值模拟系统包括：
时间偏导数获取模块，用于采用十字+菱形交错网格有限差分模板逼近时间偏导数；
空间偏导数获取模块，用于采用隐式交错网格有限差分法逼近空间偏导数；
频散关系获取模块，用于实现频散关系推导；
差分系数求取模块，用于差分系数求取；
数值模拟模块，用于基于纵横波分离实现P-SV波高精度数值模拟。
本发明的另一目的在于提供一种所述有限差分数值模拟方法在弹性介质地震波正演、逆时偏移和全波形反演中的应用。
结合上述的所有技术方案，本发明所具备的优点及积极效果为：本发明通过时间任意偶数阶精度有限差分逼近时间偏导数、空间隐式交错网格有限差分逼近空间偏导数、纵横波波场分离实现弹性P-SV波高精度数值模拟。本发明的弹性P-SV波高精度有限差分数值模拟方法，采用“十字+菱形”交错网格有限差分模板逼近弹性波方程中的时间偏导数；采样隐式交错网格有限差分法逼近弹性波方程中的空间偏导数；推导了相应的频散关系；基于Taylor级数展开求取差分系数；基于纵横波分离实现了P-SV波高精度模拟。
本发明通过采用时间任意偶数阶精度-空间隐式交错网格有限差分法实现弹性P-SV波高精度数值模拟；提高弹性P-SV波的模拟精度，进而提高弹性波逆时偏移和全波形反演的精度。
附图说明
为了更清楚地说明本申请实施例的技术方案，下面将对本申请实施例中所需要使用的附图做简单的介绍，显而易见地，下面所描述的附图仅仅是本申请的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下还可以根据这些附图获得其他的附图。
图1是本发明实施例提供的有限差分数值模拟方法流程图。
图2是本发明实施例提供的有限差分数值模拟系统的结构示意图；
图2中：1、时间偏导数获取模块；2、空间偏导数获取模块；3、频散关系获取模块；4、差分系数求取模块；5、数值模拟模块。
图3是本发明实施例提供的有限差分数值模拟方法实现流程图。
图4是本发明实施例提供的“十字+菱形”交错网格有限差分模板示意图。
图5是本发明实施例提供的双层弹性介质模型示意图。
图6(a)是本发明实施例提供的不同方法得到的水平分量波场快照；I和II分别表示常规方法和本发明中提出的方法。
图6(b)是本发明实施例提供的不同方法得到的垂直分量波场快照；I和II分别表示常规方法和本发明中提出的方法。
具体实施方式
为了使本发明的目的、技术方案及优点更加清楚明白，以下结合实施例，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅仅用以解释本发明，并不用于限定本发明。
针对现有技术存在的问题，本发明提供了一种有限差分数值模拟方法、系统、介质、计算机设备及应用，下面结合附图对本发明作详细的描述。
如图1所示，本发明提供的有限差分数值模拟方法包括以下步骤：
S101：采用十字+菱形”交错网格有限差分模板逼近时间偏导数；
S102：采用隐式交错网格有限差分法逼近空间偏导数；
S103：频散关系推导；
S104：差分系数求取；
S105：得到的时间高阶差分系数和隐式空间差分系数分别代入公式进行弹性波方程数值模拟，实现基于纵横波分离的P-SV波高精度模拟。
本发明提供的有限差分数值模拟方法业内的普通技术人员还可以采用其他的步骤实施，图1的本发明提供的有限差分数值模拟方法仅仅是一个具体实施例而已。
如图2所示，本发明提供的有限差分数值模拟系统包括：
时间偏导数获取模块1，用于采用十字+菱形交错网格有限差分模板逼近时间偏导数；
空间偏导数获取模块2，用于采用隐式交错网格有限差分法逼近空间偏导数；
频散关系获取模块3，用于实现频散关系推导；
差分系数求取模块4，用于差分系数求取；
数值模拟模块5，用于基于纵横波分离实现P-SV波高精度数值模拟。
下面结合附图对本发明的技术方案作进一步的描述。
如图3所示，本发明提供的有限差分数值模拟具体包括以下步骤：
第一步，“十字+菱形”交错网格有限差分模板逼近时间偏导数；
第二步，隐式交错网格有限差分法逼近空间偏导数；
第三步，频散关系推导；
第四步，差分系数求取；
第五步，基于纵横波分离的P-SV波高精度数值模拟。
在本发明中，第一步中“十字+菱形”交错网格有限差分模板逼近时间偏导数方法如下：
二维弹性波速度应力方程可以表示为：
其中，(vx,vz)为偏振速度矢量(τxx,τzz,τxz)，λ(x,z)和μ(x,z)为拉梅系数，ρ(x,z)为密度。采用“十字+菱形”交错网格有限差分模板(如图4所示)逼近时间偏导数(以τxx为例)，公式如下：
其中，
其中dm,0和dm,n为时间差分系数，Δt为时间步长，h为空间采样间隔，M和N为差分算子长度参数。
在本发明中，第二步中隐式交错网格有限差分法逼近空间偏导数方法如下：隐式交错网格有限差分逼近空间偏导数的公式为(以vx为例)：
其中cm和b为空间差分系数：
隐式交错网格有限差分中，空间偏导数通过求解上述三对角线性方程组得到。
在本发明中，第三步中频散关系推导方法如下：
平面波理论假设下，波场可以表示为：
其中，
为虚数单位，k为波数，ω为角频率，θ为传播角度，kx＝kcos(θ)，kz＝ksin(θ)。将公式(6)代入公式(2)并化简得到：
其中，rP＝vPΔt/h和rS＝vSΔt/h为纵波和横波的库朗数，vP和vS分别为纵波和横波的传播速度：
公式(7)和(8)即为纵波和横波的时间频散关系。
在本发明中，第四步中差分系数求取方法如下：
将公式(7)中的三角函数进行Taylor级数展开化简可得：
C+D＝E (10)
其中，
比较方程(10)两端(kxh)2l或(kzh)2l(1≤l≤M)的系数可得：
a1＝0,
比较方程(10)两端
的系数可得：
其中，int(x)为取整函数。通过公式(14)和(15)可得到基于Taylor级数展开的纵波时间高阶差分系数，将vP换为vS得到横波时间高阶差分系数。
隐式空间交错网格有限差分系数求取公式为：
在本发明中，第五步中基于纵横波分离的P-SV波高精度模拟方法如下：
将第四步中得到的时间高阶差分系数和隐式空间差分系数分别代入公式(2)和(4)可进行弹性波方程数值模拟。但弹性介质中，涉及纵波、横波和转换波等多种波模式。采用纵波时间高阶差分系数进行求解只能得到精确的纵波场；采用横波时间高阶差分系数进行求解只能得到精确的横波场；采用纵波或横波差分系数求解都得不到高精度的P-SV转换波。为实现P-SV波高精度模拟，引入纵横波分离(本发明采用波数域分离方法)。模拟步骤变为：
(1)采用纵波时间高阶差分系数和隐式空间差分系数求解弹性波方程，纵横波分离，只提取纵波场；
(2)采用横波时间高阶差分系数和隐式空间差分系数求解弹性波方程，纵横波分离，只提取横波场；
(3)将(1)中的纵波场和(2)中的横波场求和，得到最终的弹性波场。
(4)重复上面的步骤，直到最大时刻。
下面结合实验对本发明的技术效果作详细的描述。
本实验采用如图5所示的双层弹性介质模型验证本发明中提出的方法。模型计算区域有201*201个网格点，空间间隔为10.0m，时间步长为1.5ms，最大记录时间为1.5s。震源为40Hz的雷克子波，位于(1000m，800m)，201个检波点均匀分布于地表。差分算子长度M＝8，N＝4。图6(a)和图6(b)给出不同有限差分方法得到的波场快照。为了便于比较，将常规方法和本发明的新方法的结果拼接显示(分别位于左右两侧)。由图可知，常规方法结果中出现了明显的时间和空间频散，本发明中提出的方法有效压制了这些数值频散。因此，本发明提出的弹性P-SV波高精度有限差分数值模拟方法方法可以得到精确的纵波、横波和P-SV转换波。
应当注意，本发明的实施方式可以通过硬件、软件或者软件和硬件的结合来实现。硬件部分可以利用专用逻辑来实现；软件部分可以存储在存储器中，由适当的指令执行系统，例如微处理器或者专用设计硬件来执行。本领域的普通技术人员可以理解上述的设备和方法可以使用计算机可执行指令和/或包含在处理器控制代码中来实现，例如在诸如磁盘、CD或DVD-ROM的载体介质、诸如只读存储器(固件)的可编程的存储器或者诸如光学或电子信号载体的数据载体上提供了这样的代码。本发明的设备及其模块可以由诸如超大规模集成电路或门阵列、诸如逻辑芯片、晶体管等的半导体、或者诸如现场可编程门阵列、可编程逻辑设备等的可编程硬件设备的硬件电路实现，也可以用由各种类型的处理器执行的软件实现，也可以由上述硬件电路和软件的结合例如固件来实现。
以上所述，仅为本发明的具体实施方式，但本发明的保护范围并不局限于此，任何熟悉本技术领域的技术人员在本发明揭露的技术范围内，凡在本发明的精神和原则之内所作的任何修改、等同替换和改进等，都应涵盖在本发明的保护范围之内。
Claims (10)
1.一种有限差分数值模拟方法，其特征在于，所述有限差分数值模拟方法包括：
采用十字和菱形交错网格有限差分模板逼近时间偏导数；
采用隐式交错网格有限差分法逼近空间偏导数；
频散关系推导；
差分系数求取；
得到的时间高阶差分系数和隐式空间差分系数分别代入公式进行弹性波方程数值模拟，实现基于纵横波分离的P-SV波高精度模拟。
2.如权利要求1所述的有限差分数值模拟方法，其特征在于，所述十字和菱形交错网格有限差分模板逼近时间偏导数方法包括：
二维弹性波速度应力方程表示为：
其中，(vx,vz)为偏振速度矢量(τxx,τzz,τxz)，λ(x,z)和μ(x,z)为拉梅系数，ρ(x,z)为密度；采用十字+菱形交错网格有限差分模板逼近时间偏导数，公式如下：
其中，
其中dm,0和dm,n为时间差分系数，Δt为时间步长，h为空间采样间隔，M和N为差分算子长度参数。
3.如权利要求1所述的有限差分数值模拟方法，其特征在于，所述隐式交错网格有限差分法逼近空间偏导数方法包括：隐式交错网格有限差分逼近空间偏导数的公式为：
其中cm和b为空间差分系数：
隐式交错网格有限差分中，空间偏导数通过求解三对角线性方程组得到。
4.如权利要求1所述的有限差分数值模拟方法，其特征在于，所述频散关系推导方法包括：
平面波理论下，波场表示公式为：
其中，
为虚数单位，k为波数，ω为角频率，θ为传播角度，kx＝k cos(θ)，kz＝ksin(θ)；将波场表示公式代入公式
并化简得到：
其中，rP＝vPΔt/h和rS＝vSΔt/h为纵波和横波的库朗数，vP和vS分别为纵波和横波的传播速度：
公式
和
为纵波和横波的时间频散关系。
5.如权利要求1所述的有限差分数值模拟方法，其特征在于，所述差分系数求取方法包括：
将公式
中的三角函数进行Taylor级数展开化简可得：
C+D＝E；
其中，
比较方程C+D＝E两端(kxh)2l或(kzh)2l(1≤l≤M)的系数得：
比较方程C+D＝E两端
的系数得：
其中，int(x)为取整函数，通过公式：
a1＝0，
和
得到基于Taylor级数展开的纵波时间高阶差分系数，将vP换为vS得到横波时间高阶差分系数；
隐式空间交错网格有限差分系数求取公式为：
6.如权利要求1所述的有限差分数值模拟方法，其特征在于，所述基于纵横波分离的P-SV波高精度模拟方法包括：将第四步中得到的时间高阶差分系数和隐式空间差分系数分别代入公式
和
进行弹性波方程数值模拟，采用波数域分离方法得到高精度的P-SV转换波精度模拟，模拟步骤为：
(1)采用纵波时间高阶差分系数和隐式空间差分系数求解弹性波方程，纵横波分离，只提取纵波场；
(2)采用横波时间高阶差分系数和隐式空间差分系数求解弹性波方程，纵横波分离，只提取横波场；
(3)将(1)中的纵波场和(2)中的横波场求和，得到最终的弹性波场；
(4)重复(1)-(3)，直到最大时刻。
7.一种计算机设备，其特征在于，所述计算机设备包括存储器和处理器，所述存储器存储有计算机程序，所述计算机程序被所述处理器执行时，使得所述处理器执行如下步骤：
采用十字和菱形交错网格有限差分模板逼近时间偏导数；
采用隐式交错网格有限差分法逼近空间偏导数；
频散关系推导；
差分系数求取；
得到的时间高阶差分系数和隐式空间差分系数分别代入公式进行弹性波方程数值模拟，实现基于纵横波分离的P-SV波高精度模拟。
8.一种计算机可读存储介质，存储有计算机程序，所述计算机程序被处理器执行时，使得所述处理器执行如下步骤：
采用十字和菱形交错网格有限差分模板逼近时间偏导数；
采用隐式交错网格有限差分法逼近空间偏导数；
频散关系推导；
差分系数求取；
得到的时间高阶差分系数和隐式空间差分系数分别代入公式进行弹性波方程数值模拟，实现基于纵横波分离的P-SV波高精度模拟。
9.一种实施权利要求1～6任意一项所述有限差分数值模拟方法的有限差分数值模拟系统，其特征在于，所述有限差分数值模拟系统包括：
时间偏导数获取模块，用于采用十字+菱形交错网格有限差分模板逼近时间偏导数；
空间偏导数获取模块，用于采用隐式交错网格有限差分法逼近空间偏导数；
频散关系获取模块，用于实现频散关系推导；
差分系数求取模块，用于差分系数求取；
数值模拟模块，用于基于纵横波分离实现P-SV波高精度数值模拟。
10.一种如权利要求1～6任意一项所述有限差分数值模拟方法在弹性介质地震波正演、逆时偏移和全波形反演中的应用。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202011067839.9A CN112285772B (zh) | 2020-10-07 | 2020-10-07 | 有限差分数值模拟方法、系统、介质、计算机设备及应用 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202011067839.9A CN112285772B (zh) | 2020-10-07 | 2020-10-07 | 有限差分数值模拟方法、系统、介质、计算机设备及应用 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN112285772A true CN112285772A (zh) | 2021-01-29 |
| CN112285772B CN112285772B (zh) | 2023-05-19 |
Family
ID=74422316
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202011067839.9A Active CN112285772B (zh) | 2020-10-07 | 2020-10-07 | 有限差分数值模拟方法、系统、介质、计算机设备及应用 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN112285772B (zh) |
Citations (9)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20030130796A1 (en) * | 2002-01-04 | 2003-07-10 | Westerngeco L.L.C. | Method for computing finite-frequency seismic migration traveltimes from monochromatic wavefields |
| EP2184620A1 (en) * | 2008-08-26 | 2010-05-12 | PGS Geophysical AS | Fourier finite-difference migration for three dimensional tilted transverse isotropic media |
| CN102323613A (zh) * | 2011-06-01 | 2012-01-18 | 西南石油大学 | 一种基于有理切比雪夫逼近优化系数有限差分偏移方法 |
| CN106814390A (zh) * | 2015-11-27 | 2017-06-09 | 中国石油化工股份有限公司 | 基于时空域优化的交错网格正演模拟方法 |
| CN107526105A (zh) * | 2017-08-09 | 2017-12-29 | 西安交通大学 | 一种波场模拟交错网格有限差分方法 |
| WO2018047111A1 (en) * | 2016-09-08 | 2018-03-15 | King Abdullah University Of Science And Technology | Methods for efficient wavefield solutions |
| CN107976710A (zh) * | 2017-11-17 | 2018-05-01 | 河海大学 | 一种基于声波方程的线性优化隐式时空域有限差分数值模拟方法 |
| CN108279437A (zh) * | 2018-01-17 | 2018-07-13 | 中国石油大学(华东) | 变密度声波方程时间高阶精度交错网格有限差分方法 |
| US20180356547A1 (en) * | 2017-06-12 | 2018-12-13 | Saudi Arabian Oil Company | Modeling angle domain common image gathers from reverse time migration |
-
2020
- 2020-10-07 CN CN202011067839.9A patent/CN112285772B/zh active Active
Patent Citations (9)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20030130796A1 (en) * | 2002-01-04 | 2003-07-10 | Westerngeco L.L.C. | Method for computing finite-frequency seismic migration traveltimes from monochromatic wavefields |
| EP2184620A1 (en) * | 2008-08-26 | 2010-05-12 | PGS Geophysical AS | Fourier finite-difference migration for three dimensional tilted transverse isotropic media |
| CN102323613A (zh) * | 2011-06-01 | 2012-01-18 | 西南石油大学 | 一种基于有理切比雪夫逼近优化系数有限差分偏移方法 |
| CN106814390A (zh) * | 2015-11-27 | 2017-06-09 | 中国石油化工股份有限公司 | 基于时空域优化的交错网格正演模拟方法 |
| WO2018047111A1 (en) * | 2016-09-08 | 2018-03-15 | King Abdullah University Of Science And Technology | Methods for efficient wavefield solutions |
| US20180356547A1 (en) * | 2017-06-12 | 2018-12-13 | Saudi Arabian Oil Company | Modeling angle domain common image gathers from reverse time migration |
| CN107526105A (zh) * | 2017-08-09 | 2017-12-29 | 西安交通大学 | 一种波场模拟交错网格有限差分方法 |
| CN107976710A (zh) * | 2017-11-17 | 2018-05-01 | 河海大学 | 一种基于声波方程的线性优化隐式时空域有限差分数值模拟方法 |
| CN108279437A (zh) * | 2018-01-17 | 2018-07-13 | 中国石油大学(华东) | 变密度声波方程时间高阶精度交错网格有限差分方法 |
Non-Patent Citations (3)
| Title |
|---|
| ENJIANG WANG ET AL.: "Effective finite-difference modelling methods with 2D acoustic wave equation using a combination of cross and rhombus stencils", 《GEOPHYSICAL JOURNAL INTERNATIONAL》 * |
| ZHIMING REN ET AL.: "Modeling of the Acoustic Wave Equation by Staggered-Grid Finite-Difference Schemes with High-Order Temporal and Spatial Accuracy", 《BULLETIN OF THE SEISMOLOGICAL SOCIETY OF AMERICA》 * |
| 白亚东等: "交错网格高阶差分算法的改进", 《物探化探计算技术》 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN112285772B (zh) | 2023-05-19 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN112327358B (zh) | 一种粘滞性介质中声波地震数据正演模拟方法 | |
| CN101614826A (zh) | 三维地震数据处理中实现面元均化的方法和装置 | |
| CN105467443A (zh) | 一种三维各向异性弹性波数值模拟方法及系统 | |
| CN103699798A (zh) | 一种实现地震波场数值模拟方法 | |
| Cao et al. | Joint deblending and data reconstruction with focal transformation | |
| Symes et al. | Kirchhoff simulation migration and inversion using finite-difference travel-times and amplitudes | |
| US11686870B2 (en) | Interpretive-guided velocity modeling seismic imaging method and system, medium and device | |
| Petrov et al. | Higher-order grid-characteristic schemes for the acoustic system | |
| CN112285772A (zh) | 有限差分数值模拟方法、系统、介质、计算机设备及应用 | |
| Treister et al. | A multigrid solver to the Helmholtz equation with a point source based on travel time and amplitude | |
| Zhou et al. | Optimizing orthogonal-octahedron finite-difference scheme for 3D acoustic wave modeling by combination of Taylor-series expansion and Remez exchange method | |
| Wang et al. | Memory optimization in RNN-based full waveform inversion using boundary saving wavefield reconstruction | |
| Fomel | Traveltime computation with the linearized eikonal equation | |
| Song et al. | Implicit finite difference wave equation forward modeling based on partial fraction expansion | |
| CN107462925B (zh) | 一种三维孔隙弹性介质中快速波场模拟方法 | |
| Shao‐Lin et al. | Symplectic RKN schemes for seismic scalar wave simulations | |
| CN114924313B (zh) | 基于行波分离的声波最小二乘逆时偏移梯度求取方法 | |
| CN114624766B (zh) | 基于行波分离的弹性波最小二乘逆时偏移梯度求取方法 | |
| CN113341460B (zh) | 一种基于连续算子分裂的循环式极小化地震数据重建方法 | |
| Liu et al. | Using pseudo-spectral method on curved grids for SH-wave modeling of irregular free-surface | |
| Zhou et al. | Local cosine similarity attribute based on fast-streaming computation | |
| Ross | Multiple suppression: beyond 2-D. Part I: theory | |
| Calandra et al. | Recent advances in numerical methods for solving the wave equation in the context of seismic depth imaging | |
| He et al. | A practical implementation of 3D full-waveform inversion on heterogeneous HPC systems | |
| Cai et al. | A mesh-free finite-difference scheme for frequency-domain acoustic wave simulation with topography |
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 |