CN109635513A - 三维TTI介质qP波波场模拟方法及装置 - Google Patents
三维TTI介质qP波波场模拟方法及装置 Download PDFInfo
- Publication number
- CN109635513A CN109635513A CN201910049386.8A CN201910049386A CN109635513A CN 109635513 A CN109635513 A CN 109635513A CN 201910049386 A CN201910049386 A CN 201910049386A CN 109635513 A CN109635513 A CN 109635513A
- Authority
- CN
- China
- Prior art keywords
- wave field
- anisotropic
- operator
- wave
- differential operator
- 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
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
- G06F30/23—Design optimisation, verification or simulation using finite element methods [FEM] or finite difference methods [FDM]
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Computer Hardware Design (AREA)
- Evolutionary Computation (AREA)
- Geometry (AREA)
- General Engineering & Computer Science (AREA)
- General Physics & Mathematics (AREA)
- Complex Calculations (AREA)
- Management, Administration, Business Operations System, And Electronic Commerce (AREA)
Abstract
本发明公开了一种三维TTI介质qP波波场模拟方法及装置，该方法包括：生成三维TTI介质qP波波场方程，其中，所述三维TTI介质qP波波场方程包括一椭圆微分算子和一标量算子；将所述三维TTI介质qP波波场方程分解成一各项异性背景波场项和一各项异性波场校正项；通过低秩分解近似计算所述椭圆微分算子；根据近似计算出的椭圆微分算子确定出所述各项异性背景波场项；根据近似计算出的椭圆微分算子以及标量算子计算出所述各项异性波场校正项；根据所述各项异性波场校正项对所述各项异性背景波场项进行校正，得到各项异性总波场。本发明解决了现有技术的三维TTI介质qP波波场模拟方法模拟数值不稳定、伪横波噪音和计算量巨大的问题。
Description
技术领域
本发明涉及波场模拟技术领域，具体而言，涉及一种三维TTI介质qP波波场模拟方法及装置。
背景技术
现有技术中，探索高效率、高精度且稳定的拟声波方程和成像算法是近几年国内和国际地球物理学家研究的热点和难点问题。对于TTI介质(倾斜对称轴各向异性介质)，Zhan et al.(2012)提出利用伪谱法单独计算波数域算子，而不是直接计算波数-空间混合域算子，并采用快速展开方法(Rapid Expansion Method)(Kosloff et al.,1989；Pestanaand Stoffa,2010；Pestana et al.,2011)延拓解耦的纯P波波场。Alkhalifah et al.,(2013)，Ibanez-Jacome et al.(2014)和Waheed和Alkhalifah(2014)提出了更为高效的算法，其主要思想是利用等效各向同性模型近似计算TTI介质纯拟声波波场。
Xu和Zhou(2014)提出了一个新的TTI介质纯qP波方程，这一方程和各向同性声波方程有相同的形式，并可采用伪谱法进行数值模拟。Xu et al.(2015)又针对Xu和Zhou(2014)提出的TTI介质纯qP波方程做了改进，他们将qP波方程中的伪微分算子解耦成两个算子：椭圆微分算子和标量算子，以此确保相位更为准确振幅更均一。然而，这种方法对于三维TTI介质计算量仍然太大。
类似的，Waheed和Alkhalifah(2015)提出利用具有垂直轴的椭圆微分算在近似TTI介质qP波频散曲线，虽然计算量比Xu et al.(2015)减小，但波场计算精度减小了。
最近，Yan和Liu(2016)利用基于低秩近似的优化拟微分算子对TTI介质纯qP波波场进行了模拟，该方法相比Fomel et al.(2012)方法更高效但模拟的波场准确度有所下降。Zhang et al.(2017)提出利用各项同性低秩近似算子高效率模拟VTI介质qP波波场，但该方法还未拓展到TTI介质。
由以上描述可以看出，现有技术的三维TTI介质qP波波场模拟方法存在数值不稳定、伪横波噪音和计算量巨大三方面问题，无法实现高效高精度稳定的三维复杂TTI介质qP波波场模拟。
发明内容
本发明的主要目的在于提供一种三维TTI介质qP波波场模拟方法及装置，以解决现有技术的三维TTI介质qP波波场模拟方法存在的模拟数值不稳定、伪横波噪音和计算量巨大的问题。
为了实现上述目的，根据本发明的一个方面，提供了一种三维TTI介质qP波波场模拟方法，该方法包括：
根据VTI介质qP波波场方程以及三维旋转因子生成三维TTI介质qP波波场方程，其中，所述三维TTI介质qP波波场方程包括一椭圆微分算子和一标量算子；
将所述三维TTI介质qP波波场方程分解成一各项异性背景波场项和一各项异性波场校正项；
通过低秩分解近似计算所述椭圆微分算子；
根据近似计算出的椭圆微分算子确定出所述各项异性背景波场项；
根据近似计算出的椭圆微分算子以及标量算子计算出所述各项异性波场校正项；
根据所述各项异性波场校正项对所述各项异性背景波场项进行校正，得到各项异性总波场。
进一步的，所述通过低秩分解近似计算所述椭圆微分算子，包括：
将所述三维TTI介质qP波波场方程转换为时间域离散表达式，其中，所述时间域离散表达式中包括拟微分算子；
通过低秩分解近似所述拟微分算子的象征，得到所述拟微分算子的近似表达式；
根据所述近似表达式计算椭圆微分算子。
进一步的，所述将所述三维TTI介质qP波波场方程转换为时间域离散表达式，包括：
将所述三维TTI介质qP波波场方程中的椭圆微分算子变换至频率-波数域得到频散关系；
根据所述频散关系以及拟解析法将所述三维TTI介质qP波波场方程转换为频率域离散表达式；
通过傅里叶反变换将所述频率域离散表达式转换为时间域离散表达式。
进一步的，所述通过低秩分解近似计算所述椭圆微分算子，包括：
通过CPU/GPU并行集群计算所述椭圆微分算子。
进一步的，所述三维TTI介质qP波波场方程的表达式为：
公式中：
Azz＝1+2εsin2θ
其中：为椭圆微分算子，Se为标量算子，x、y、z为空间坐标，u为qP波波场，v0为沿对称轴方向qP波速度，ε、θ和为各向异性参数。
进一步的，所述时间域离散表达式为：
u(x,t+Δt)-2u(x,t)+u(x,t-Δt)＝v0 2(x)Δt2FT-1{F(x,k)FFT[u(x,t)]}Se
公式中：
其中，x为空间矢量，k为空间波数矢量，Δt为时间步长，F(x,k)为优化拟拉普拉斯算子，FT-1为傅里叶反变换，FFT为快速傅里叶正变换。
进一步的，所述拟微分算子的象征的近似表达式为：
其中，aO(x,k)为拟微分算子的象征，W(x,km)是aO(x,k)的子矩阵，其由与km有关的特定的列元素组成，W(xn,k)是aO(x,k)的另一个子矩阵，其由与xn有关的特定行元素组成，amn代表系数矩阵。
为了实现上述目的，根据本发明的另一方面，提供了一种三维TTI介质qP波波场模拟装置，该装置包括：
三维TTI介质qP波波场方程计算单元，用于根据VTI介质qP波波场方程以及三维旋转因子计算出三维TTI介质qP波波场方程，其中，所述三维TTI介质qP波波场方程包括一椭圆微分算子和一标量算子；
改写计算单元，用于将所述三维TTI介质qP波波场方程分解成一各项异性背景波场项和一各项异性波场校正项；
低秩分解近似单元，用于通过低秩分解近似计算所述椭圆微分算子；
各项异性背景波场项确定单元，用于根据近似计算出的椭圆微分算子确定出所述各项异性背景波场项；
各项异性波场校正项计算单元，用于根据近似计算出的椭圆微分算子以及标量算子计算出所述各项异性波场校正项；
各项异性总波场计算单元，用于根据所述各项异性波场校正项对所述各项异性背景波场项进行校正，得到各项异性总波场。
本发明针对目前三维TTI介质qP波正演模拟和逆时偏移方法主要存在复杂介质情况下数值不稳定、伪横波成像假象噪音和计算量巨大三方面的难题，提出利用低秩分解近似椭圆微分算子进行三维TTI介质qP波波场模拟和逆时偏移成像研究，实现了如下有益效果：
(1)本发明利用椭圆微分算子近似TTI介质qP波频散曲线模拟纯qP波波场，并有效地保证qP波相位的准确性和振幅的均匀性(Xu et al.,2015)，解决目前工业界普遍采用的耦合方程(Zhou et al.,2006a,2006b；Fletcher et al.,2009；Fowler et al.,2010；Duveneck and Bakker,2011)逆时偏移成像方法在复杂TTI介质存在数值不稳定和伪横波成像假象噪音问题。
(2)本发明在(1)基础之上，利用低秩分解近似椭圆微分算子及相应的Lowrank有限差分方法，有效地减少Xu et al.(2015)所提出的基于椭圆微分算子三维TTI介质qP波伪谱法计算量，增强数值稳定性(且较大时间步长仍稳定)，大幅度减少计算时间。
(3)本发明在(1)和(2)基础之上，采用CPU/GPU并行集群加速计算低秩分解近似椭圆微分算子，为最终实现适用于工业生产的三维复杂TTI介质qP波高效率、高精度、稳定正演模拟和逆时偏移成像提供技术保障。
综上，本发明解决了目前三维复杂TTI介质qP波波场模拟数值不稳定、伪横波噪音和计算量巨大三方面问题，为实现高效高精度稳定的三维复杂TTI介质qP波波场模拟和逆时偏移成像提供了技术保障。
附图说明
为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。在附图中：
图1是本发明实施例提出的一种三维TTI介质qP波波场模拟方法的流程图；
图2是本发明实施例低秩分解近似椭圆微分算子的方法的流程图；
图3是本发明实施例将所述三维TTI介质qP波波场方程转换为时间域离散表达式的方法流程图；
图4是本发明实施例提出的一种三维TTI介质qP波波场模拟装置的结构框图；
图5是本发明实施例低秩分解近似单元的组成结构图；
图6是本发明实施例时间域离散表达式转换模块的组成结构图；
图7是本发明实施例三维TTI介质qP波波场模拟方法的流程框图；
图8a是实施例第一参数下弱各项异性参数三维TTI理论与近似慢度曲线比较图；
图8b是实施例第二参数下弱各项异性参数三维TTI理论与近似慢度曲线比较图；
图8c是实施例第三参数下弱各项异性参数三维TTI理论与近似慢度曲线比较图；
图9a是实施例第一参数下弱各项异性参数三维TTI近似慢度曲线与理论曲线误差图；
图9b是实施例第二参数下弱各项异性参数三维TTI近似慢度曲线与理论曲线误差图；
图9c是实施例第三参数下弱各项异性参数三维TTI近似慢度曲线与理论曲线误差图；
图10a是实施例第一参数下强各项异性参数三维TTI理论与近似慢度曲线比较图；
图10b是实施例第二参数下强各项异性参数三维TTI理论与近似慢度曲线比较图；
图10c是实施例第三参数下强各项异性参数三维TTI理论与近似慢度曲线比较图；
图11a是实施例第一参数下强各项异性参数三维TTI近似慢度曲线与理论曲线误差图；
图11b是实施例第二参数下强各项异性参数三维TTI近似慢度曲线与理论曲线误差图；
图11c是实施例第三参数下强各项异性参数三维TTI近似慢度曲线与理论曲线误差图；
图12a是实施例采用本发明方法正演得到的波场快照；
图12b是实施例采用Xu和Zhou(2014)方法得到的波场快照。
具体实施方式
为了使本技术领域的人员更好地理解本发明方案，下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分的实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都应当属于本发明保护的范围。
本领域内的技术人员应明白，本发明的实施例可提供为方法、系统、或计算机程序产品。因此，本发明可采用完全硬件实施例、完全软件实施例、或结合软件和硬件方面的实施例的形式。而且，本发明可采用在一个或多个其中包含有计算机可用程序代码的计算机可用存储介质(包括但不限于磁盘存储器、CD-ROM、光学存储器等)上实施的计算机程序产品的形式。
需要说明的是，本发明的说明书和权利要求书及上述附图中的术语“第一”、“第二”等是用于区别类似的对象，而不必用于描述特定的顺序或先后次序。应该理解这样使用的数据在适当情况下可以互换，以便这里描述的本发明的实施例。此外，术语“包括”和“具有”以及他们的任何变形，意图在于覆盖不排他的包含，例如，包含了一系列步骤或单元的过程、方法、系统、产品或设备不必限于清楚地列出的那些步骤或单元，而是可包括没有清楚地列出的或对于这些过程、方法、产品或设备固有的其它步骤或单元。需要说明的是，在不冲突的情况下，本发明中的实施例及实施例中的特征可以相互组合。下面将参考附图并结合实施例来详细说明本发明。
本发明针对目前三维TTI介质qP波正演模拟和逆时偏移方法主要存在复杂介质情况下数值不稳定、伪横波成像假象噪音和计算量巨大三方面的难题，提出利用低秩分解近似椭圆微分算子进行三维TTI介质qP波波场模拟和逆时偏移成像研究，本发明方案的具体推导过程如下。
在现有技术中，纯qP波在VTI介质传播可用下面解耦声波方程(Alkhalifah,2000)定义如下：
其中，u为qP波波场，v0为沿对称轴方向qP波速度，ε和δ为Thomsen各向异性参数(Thomsen,1986)。
水平方向空间微分算子定义为：
微分算子Q定义为：
在现有技术中，Xu和Zhou(2014)通过引入一个标量算子S，将方程(1)转化为：
其中，u为qP波波场，v0为沿对称轴方向qP波速度，表示拉普拉斯算子，
标量算子S写作：
其中，定义相方向矢量
由式(5)可知，在波数域计算S时，如果各项异性参数不为常数，S为拟微分算子，其计算非常昂贵。Xu和Zhou(2014)提出在空间域通过波场梯度近似计算相方向矢量。这样，求解上述VTI/TTI各项异性声波方程计算量比各项同性声波方程仅多一半，并可适应于复杂TTI介质。然而，qP波波场的准确求解取决于相方向矢量即波场梯度准确度。此外，该方法求解的qP波波场振幅不均一，虽然弹性波方程近似下的拟声波方程qP波振幅没有明确物理意义，但不均衡的振幅会造成逆时偏移成像假象。
为分析上述近似变量算子对波场相位和振幅的影响，本发明对标量算子S进行了改写，将其写为如下形式：
S＝1+ΔS(6)
并重新定于ΔS为：
这样公式(4)的波动方程可写作：
从方程(8)可知，方程右边第一项为各项同性背景波场，其为Laplacian微分算子，而方程右边第二项为各项异性波场校正项。若方程右边各项同性声波频散曲线与实际各项异性频散曲线相差越大，则第二项各项异性波场校正量越大。这样，当第二项校正量越大时，波场准确度就取决于求解校正项的相方向矢量即波场梯度准确度。然而，在实际波场模拟和逆时偏移成像过程中，若空间采样间隔大于Nyquist采样频率，那么高波数部分的波场将产生假频并影响求解低波数部分波场梯度的准确度。
为了解决上述第二项各项异性波场校正项过大问题，本发明将方程(1)中的拟微分算子解耦成两个微分算子(Xu et al.,2015)：一个椭圆微分算子和一个标量算子Se。这样，方程(1)可写为：
其中，
在公式(9)中，为椭圆微分算子，Se为标量算子。在本发明中，为了便于分析方程(9)，本发明采用公式(6)的思想。对标量算子Se进行了改写，将公式(9)其改写为：
其中，
类似于方程(8)，将方程(11)右边展开，右边第一项为各项异性背景波场项，其为椭圆微分算子，右边第二项为各项异性波场校正项。
对于计算量而言，基于方程(2)和(9)的VTI介质qP波波场模拟所需计算量大致一样。然而，对于TTI介质，特别是三维情况，求解方程(9)所需的计算量约是求解方程(2)的3倍。
本发明为了方便对比基于方程(2)和(9)的三维TTI介质qP波波场模拟所需计算量，引入三维旋转因子：
并计算出相应的三维旋转波数，结合方程(9)即可推导得到三维TTI介质时间-空间域纯qP波波场方程：
在公式(13)中，Se为标量算子，如方程(10)所定义，公式右边除Se的部分为椭圆微分算子，u为qP波波场，v0为沿对称轴方向qP波速度，A与TTI介质的各向异性参数ε、θ和有关，其表达式如下：
如前所述，qP波波场准确度就取决于求解校正项的相方向矢量即波场梯度准确度，因此，本发明采用空间上可以达到Nyquist采样极限精度的伪谱法求解方程(13)。每一次时间步长，求解方程(13)椭圆微分算子需要1次FFT(快速傅里叶变换)正变换和6次FFT反变换，而求解方程(2)Laplacian微分算子需要1次FFT正变换和3次FFT反变换；此外，求解方程(2)和(13)各项异性校正项均需要一次FFT正反变换。由此可见，采用常规的伪谱法求解基于椭圆微分算子的三维TTI介质qP波波场传播巨大。另外，时间方向二阶导数采样二阶有限差分算子近似时，伪谱法仍然存在数值误差和不稳定现象，特别是当时间步长取较大时。
然而，不难发现方程(13)系数A与Thomsen各向异性参数δ无关，且对于速度模型光滑模型波数-空间混合矩阵(椭圆微分算子)的秩将不会很大(参见实施例一)，受此启发，本发明提出利用利用低秩分解近似椭圆微分算子求解三维TTI介质qP波波场方程(13)。此方法具体步骤可分为两步，首先利用低秩分解近似椭圆微分算子求解背景波场项，然后利用标量算子Se校正更新背景波场项得到各项异性总波场。
在本发明中，首先需要对三维TTI介质qP波波场方程(13)进行改写，将方程(13)改写为一各项异性背景波场项和一各项异性波场校正项。具体的改写方法可以参照公式(11)和(12)，将公式(13)中的标量算子Se改写为Se＝(1+ΔSe)的形式，其中ΔSe为标量校正算子，在公式(12)中已给出，改写后的公式(13)为：
将改写后的公式的右侧展开，其右侧第一项：
为各项异性背景波场项，其即为椭圆微分算子。
公式展开后右边第二项：
其为各项异性波场校正项，为椭圆微分算子与标量校正算子ΔSe的乘积。由此将三维TTI介质qP波波场方程(13)改写为一各项异性背景波场项和一各项异性波场校正项。
在本发明中，低秩分解近似椭圆微分算子求解各项异性背景波场项，可先将三维TTI纯qP波波场方程(13)的椭圆微分算子，变换至频率-波数域得到频散关系：
其中：ω为圆频率，v0为沿对称轴方向qP波速度，ε、θ和为各向异性参数，kx、ky、kz分别为空间波数矢量k的三个分量。
然后，根据该频散关系，并利用拟解析法(Etgen and Brandsberg-Dahl 2009)，可以将方程(13)转换为频率域离散形式：
FT[u(x,t+Δt)-2u(x,t)+u(x,t-Δt)]＝v0 2(x)Δt2F(x,k)FFT[u(x,t)]Se (16)
其中，x为空间矢量，k为空间波数矢量，v0为沿对称轴方向qP波速度，Δt为时间步长，F(x,k)为优化拟拉普拉斯算子，FFT为快速傅里叶正变换，Se为标量算子。
进而，通过将方程(16)进行傅里叶反变换，得到三维TTI介质纯qP波波场方程的拟解析法时间域离散表达式：
u(x,t+Δt)-2u(x,t)+u(x,t-Δt)＝v0 2(x)Δt2FT-1{F(x,k)FFT[u(x,t)]}Se (17)
其中：
在公式中，F(x,k)为优化拟拉普拉斯算子(Chu and Stoffa,2010；Yan and Liu,2016)，ω＝v0|k|，x为空间矢量，k为空间波数矢量，FT-1为傅里叶反变换，FFT为快速傅里叶正变换。值得注意的是，F(x,k)与波数和空间位置有关，因而傅里叶反变换不能直接采用FFT，且其计算量为
在本发明中，为了便于计算，将拟微分算子AO记为：
AO[u(x,t)]＝FT-1{F(x,k)FFT(u(x,t))} (18)
其中，优化拟拉普拉斯算子F(x,k)也称为拟微分算子的象征。此优化拟微分算子AO不仅包括原始拟微分算子的谱估计而且还包含一个时间补偿项，此时间补偿项可在波数-空间域有效地补偿波动方程在时间方向上采用二阶差分引起的误差，并增强数值稳定性。为了减少计算量，本发明采用低秩分解近似技术加速计算拟微分算子AO。
由于拟拉普拉斯算子(拟微分算子的象征)中包含随空间坐标变化的补偿速度，若直接计算，其计算量为其中Nx为3维网格总数。为了减少计算量，本发明采用低秩分解近似优化拟微分算子的象征aO(x,k)。
根据Fomel et al.(2010,2012)，本发明将拟微分算子的象征aO(x,k)近似分解为如下形式：
其中，W(x,km)是优化拟微分算子的象征aO(x,k)的子矩阵(submatirx)，其由与km有关的特定的列元素组成；W(xn,k)是象征aO(x,k)的另一个子矩阵，其由与xn有关的特定行元素组成；amn代表中间的系数矩阵。
本发明通过用公式(19)的拟微分算子的象征aO(x,k)的近似分解形式，替换公式(18)中的优化拟拉普拉斯算子F(x,k)，重新组合公式(18)可以得到拟微分算子的近似表达式：
在TTI介质拟声波波场模拟中，由于拟微分算子的象征主要与速度有关和各向异性参数ε有关，N通常很小(特别对于光滑速度模型)。因此通过重新组合的公式(20)可以加速计算AO[u(x,t)]。
在三维TTI介质纯qP波波场方程的拟解析法时间域离散表达式(17)中，椭圆微分算子为公式右侧的第一项：
v0 2(x)Δt2FT-1{F(x,k)FFT[u(x,t)]}
通过公式(18)该椭圆微分算子可以改写为：
v0 2(x)Δt2AO[u(x,t)]
其中，Δt为时间步长，v0为沿对称轴方向qP波速度，因此通过公式(20)加速计算出AO[u(x,t)]，就可以加速计算出椭圆微分算子。
在本发明中，通过对三维TTI介质qP波波场方程(13)进行改写，将方程(13)改写为一各项异性背景波场项和一各项异性波场校正项，其中，项异性背景波场即为椭圆微分算子；各项异性波场校正项，为椭圆微分算子与标量校正算子ΔSe的乘积。因此在计算出椭圆微分算子后，可以根据椭圆微分算子分别计算出各项异性背景波场项和各项异性波场校正项，进而通过各项异性波场校正项对各项异性背景波场项进行校正，得到各项异性总波场。
在本发明实施例中，公式(13)中二阶微分算子 与优化拟微分算子的象征αO(x,k)之间的关系可表示为：
基于上述原理，本发明实施例提出了一种三维TTI介质qP波波场模拟方法。图1是本发明实施例提出的一种三维TTI介质qP波波场模拟方法的流程图，如图1所示，本实施例的三维TTI介质qP波波场模拟方法包括步骤S101至步骤S105。
步骤S101，根据VTI介质qP波波场方程以及三维旋转因子计算出三维TTI介质qP波波场方程，其中，所述三维TTI介质qP波波场方程包括一椭圆微分算子和一标量算子。
步骤S102，将所述三维TTI介质qP波波场方程分解成一各项异性背景波场项和一各项异性波场校正项；
步骤S103，通过低秩分解近似计算所述椭圆微分算子。
步骤S104，根据近似计算出的椭圆微分算子确定出所述各项异性背景波场项。
步骤S105，根据近似计算出的椭圆微分算子以及标量算子计算出所述各项异性波场校正项。
步骤S106，根据所述各项异性波场校正项对所述各项异性背景波场项进行校正，得到各项异性总波场。
图2是本发明实施例低秩分解近似椭圆微分算子的方法的流程图，如图2所示，本发明实施例低秩分解近似椭圆微分算子的方法包括步骤S201至步骤S203。
步骤S201，将所述三维TTI介质qP波波场方程转换为时间域离散表达式，其中，所述时间域离散表达式中包括拟微分算子。
步骤S202，通过低秩分解近似所述拟微分算子的象征，得到所述拟微分算子的近似表达式。
步骤S203，根据所述近似表达式计算椭圆微分算子。
图3是本发明实施例将所述三维TTI介质qP波波场方程转换为时间域离散表达式的方法流程图，如图3所示，本发明实施例将所述三维TTI介质qP波波场方程转换为时间域离散表达式的方法包括步骤S301至步骤S303。
步骤S301，将所述三维TTI介质qP波波场方程中的椭圆微分算子变换至频率-波数域得到频散关系。
步骤S302，根据所述频散关系以及拟解析法将所述三维TTI介质qP波波场方程转换为频率域离散表达式。
步骤S303，通过傅里叶反变换将所述频率域离散表达式转换为时间域离散表达式。
图7是本发明实施例三维TTI介质qP波波场模拟方法的流程框图，如图7所示，本发明一可选实施例的三维TTI介质qP波波场模拟方法的流程可以为：首先计算出三维TTI介质qP波波场方程，其中三维TTI介质qP波波场方程可以分解为一椭圆微分算子和一标量算子。进而对三维TTI介质qP波波场方程进行改写，将其分解为一各项异性背景波场项和一各项异性波场校正项。通过CPU/GPU并行集群加速计算低秩分解近似椭圆微分算子，进而根据椭圆微分算子求解各项异性背景波场项。根据椭圆微分算子通过CPU/GPU并行集群加速计算出各项异性波场校正项，以实现校正更新各项异性背景波场项得到各项异性总波场。
实施例一
本发明分别以弱Thomsen各向异性参数ε＝0.25，δ＝0.1和强Thomsen各向异性参数ε＝0.8，δ＝1.0，沿对称轴方向qP波速度Vp＝3500m/s，方位角φ＝60°以及不同TTI对称轴倾角θ＝0°，θ＝30°，θ＝75°为例，比较拉普拉斯微分算子，垂直椭圆微分算子和倾斜椭圆微分算子近似三维TTI理论慢度曲线(或频散曲线)精度。
图8a至8c以及图9a至9c比较了弱各项异性参数三维TTI理论与近似慢度曲线，图8a、图8b、图8c分别是实施例第一参数、第二参数、第三参数下弱各项异性参数三维TTI理论与近似慢度曲线比较图。在图8a、图8b、图8c的弱各项异性参数三维TTI理论与近似慢度曲线比较图中，左上为理论频散曲线，右上为拉普拉斯微分算子近似曲线，左下为垂直椭圆微分算子近似曲线，右下为椭圆微分算子近似曲线。在本实施中，图8a的弱各项异性参数三维TTI理论与近似慢度曲线比较图中的各向异性参数为：ε＝0.6，δ＝0.4，θ＝0°；图8b的弱各项异性参数三维TTI理论与近似慢度曲线比较图中的各向异性参数为：ε＝0.6，δ＝0.4，θ＝30°；图8b的弱各项异性参数三维TTI理论与近似慢度曲线比较图中的各向异性参数为：ε＝0.6，δ＝0.4，θ＝75°。
图9a、图9b、图9c分别是实施例第一参数、第二参数、第三参数下弱各项异性参数三维TTI近似慢度曲线与理论曲线误差图。在图9a、图9b、图9c的弱各项异性参数三维TTI近似慢度曲线与理论曲线误差图中，上面为拉普拉斯微分算子近似曲线，中间为垂直椭圆微分算子近似曲线，下面椭圆微分算子近似曲线。图9a的弱各项异性参数三维TTI近似慢度曲线与理论曲线误差图的各向异性参数为：ε＝0.6，δ＝0.4，θ＝0°；图9b的弱各项异性参数三维TTI近似慢度曲线与理论曲线误差图的各向异性参数为：ε＝0.6，δ＝0.4，θ＝30°；图9c的弱各项异性参数三维TTI近似慢度曲线与理论曲线误差图的各向异性参数为：ε＝0.6，δ＝0.4，θ＝75°。
由图8a至8c以及图9a至9c可以看出，本发明方法的椭圆微分算子近似精度高于拉普拉斯微分算子和垂直椭圆微分算子，特别在TTI对称轴倾角较大的情况下更为明显。
图10a至10c以及图11a至11c比较了强各项异性参数TTI理论与近似慢度曲线，图10a、图10b、图10c分别是实施例第一参数、第二参数、第三参数下强各项异性参数三维TTI理论与近似慢度曲线比较图。在图10a、图10b、图10c的强各项异性参数三维TTI理论与近似慢度曲线比较图中，左上为理论频散曲线，右上为拉普拉斯微分算子近似曲线，左下为垂直椭圆微分算子近似曲线，右下为椭圆微分算子近似曲线。在本实施中，图10a的强各项异性参数三维TTI理论与近似慢度曲线比较图中的各向异性参数为：ε＝0.8，δ＝1.0，θ＝0°；图10b的强各项异性参数三维TTI理论与近似慢度曲线比较图中的各向异性参数为：ε＝0.8，δ＝1.0，θ＝30°；图10b的强各项异性参数三维TTI理论与近似慢度曲线比较图中的各向异性参数为：ε＝0.8，δ＝1.0，θ＝75°。
图11a、图11b、图11c分别是实施例第一参数、第二参数、第三参数下强各项异性参数三维TTI近似慢度曲线与理论曲线误差图。在图11a、图11b、图11c的强各项异性参数三维TTI近似慢度曲线与理论曲线误差图中，上面为拉普拉斯微分算子近似曲线，中间为垂直椭圆微分算子近似曲线，下面椭圆微分算子近似曲线。图11a的强各项异性参数三维TTI近似慢度曲线与理论曲线误差图的各向异性参数为：ε＝0.8，δ＝1.0，θ＝0°；图11b的强各项异性参数三维TTI近似慢度曲线与理论曲线误差图的各向异性参数为：ε＝0.8，δ＝1.0，θ＝30°；图11c的强各项异性参数三维TTI近似慢度曲线与理论曲线误差图的各向异性参数为：ε＝0.8，δ＝1.0，θ＝75°。
从图10a至10c以及图11a至11c可以看出，椭圆微分算子近似精度高于拉普拉斯微分算子和垂直椭圆微分算子。而且无论是弱Thomsen各向异性参数和强Thomsen各向异性参数，倾斜椭圆微分算子近似误差都很小；但是拉普拉斯微分算子和垂直椭圆微分算子由于在不同传播角度近似理论慢度曲线误差不同，其可造成波场振幅不均一。因此，本发明可以利用椭圆微分算子进行高精度的三维TTI介质qP波波场模拟和成像。
实施例二
本发明采用均匀模型进行点源响应测试，分别采用Xu和Zhou(2014)和本发明低秩近似椭圆拟微分算子对TTI介质解耦方程纯P波波动方程进行数值模拟。
模型和正演参数为：X(水平方向)和Z(垂直方向)方向网格间距均为15m，速度Vp＝1500m/s，Thomsen各向异性参数ε＝0.25，δ＝0.1，TTI对称轴倾角θ＝60°，时间步长Δt＝2ms，雷克子波震源置于模型正中，震源子波主频为20Hz。此时精确模拟每个波长所需要的网格点数约为1500/15/50＝2个，达到了采样定理极限精度的要求。此外，伪谱法稳定性条件CFL数α＝0.3。图12a为采用本发明方法正演得到的波场快照，图12b为采用Xu和Zhou(2014)方法得到的波场快照。从图12a和图12b中可观察到，由于Xu和Zhou(2014)方法采用Laplacian微分算子，造成波场振幅严重的不均衡，而本发明方法可有效消除此问题。
由于本发明椭圆微分算子不依赖于各项异性参数δ，因此其低秩分解的秩比Fomelet al.(2010)方法要小，相应地波场正演模拟所需的FFTI/IFFT次数小。
本发明方法和Xu et al.(2015)方法对于不同时间采样间隔所需计算量(FFT/IFFT次数)和稳定性，如下表可见，本发明方法相对于Xu et al.(2015)方法计算量更小，稳定性更强。
本发明方法与Xu et al.(2015)方法对比：FFT/IFFT次数
本发明技术方案带来的有益效果：
目前三维TTI介质qP波正演模拟和逆时偏移方法主要存在复杂介质情况下数值不稳定、伪横波成像假象噪音和计算量巨大三方面的难题，探索高效率、高精度且稳定的拟声波方程和算法是近几年国内和国际地球物理学家研究重点和热点问题。如背景技术所述，目前地球物理学家已提出了各种不同的各向异性声波近似方程，这些方程可以分为两类：非解耦方程和解耦方程。
TTI介质耦合方程波场模拟主要存在三个问题：其一是存在较强的伪SV波假象，若将震源激发处介质设置为椭圆各向异性(ε＝δ)或声波各向同性介质(ε＝δ＝0)虽可消除伪SV波，但是这种人为假设会影响波场传播的运动学特征；其二是ε和δ需满足物理稳定性条件(ε≥δ)，当各向异性介质参数ε＜δ则波场传播不稳定；其三是TTI角度剧变时，即使满足ε≥δ物理稳定性条件，角度剧变区域由于SV波的影响，耦合方程也会产生不稳定的波场。
与此相反，解耦方程能解决耦合方程的上述三个问题。具体而言，解耦方程qP波完全不受伪SV波的干扰；qP波不受ε和δ关系的限制，当ε＜δ时波场也能稳定传播；对TTI角度剧变的情况qP波也有很好适应性。因此，解耦方程能在不受伪SV波干扰的同时保持稳定的传播特性。
但是，解耦VTI或TTI拟声波方程因涉及复杂的伪微分算子，采用伪谱法(Reshef,1988)虽计算比较方便，但伪谱法需要做大量的正反FFT变换，导致其计算效率极其低下。此外，尽管伪谱法在空间上可以达到Nyquist采样极限精度，但时间方向二阶导数采样二阶有限差分算子近似时，其仍然存在数值误差和不稳定现象，特别是当时间步长取较大时。
Xu et al.(2015)又针对Xu和Zhou(2014)提出的TTI介质纯qP波方程做了改进，他们将qP波方程中的伪微分算子解耦成两个算子：椭圆微分算子和标量算子，以此确保相位更为准确振幅更均一。然而，这种椭圆微分算子对于三维TTI介质计算仍然很大。
在Xu et al.(2015)工作基础上，本发明提出利用地秩分解技术近似椭圆微分算子，解决目前三维复杂TTI介质qP波波场模拟数值不稳定、伪横波噪音和计算量巨大三方面问题，实现高效高精度稳定的三维复杂TTI介质qP波波场模拟和逆时偏移成像方法。具体包括以下几点：
(1)利用椭圆微分算子近似TTI介质qP波频散曲线模拟纯qP波波场，并有效地保证qP波相位的准确性和振幅的均匀性(Xu et al.,2015)，解决目前工业界普遍采用的耦合方程(Zhou et al.,2006a,2006b；Fletcher et al.,2009；Fowler et al.,2010；Duveneckand Bakker,2011)逆时偏移成像方法在复杂TTI介质存在数值不稳定和伪横波成像假象噪音问题。
(2)利用低秩分解近似椭圆微分算子及相应的Lowrank有限差分方法，有效地减少Xu et al.(2015)所提出的基于椭圆微分算子三维TTI介质qP波伪谱法计算量，增强数值稳定性(且较大时间步长仍稳定)，大幅度减少计算时间，为达到目前工业界三维TTI介质qP波正演模拟和逆时偏移成像提供技术保障。
(3)在(1)和(2)基础之上，采用CPU/GPU并行集群加速计算低秩分解近似椭圆微分算子，为实现适用于工业生产的三维复杂TTI介质qP波高效率、高精度、稳定正演模拟和逆时偏移成像提供技术保障。
需要说明的是，在附图的流程图示出的步骤可以在诸如一组计算机可执行指令的计算机系统中执行，并且，虽然在流程图中示出了逻辑顺序，但是在某些情况下，可以以不同于此处的顺序执行所示出或描述的步骤。
基于同一发明构思，本发明实施例还提供了一种三维TTI介质qP波波场模拟装置，可以用于实现上述实施例所描述的三维TTI介质qP波波场模拟方法，如下面的实施例所述。由于三维TTI介质qP波波场模拟装置解决问题的原理与三维TTI介质qP波波场模拟方法相似，因此三维TTI介质qP波波场模拟装置的实施例可以参见三维TTI介质qP波波场模拟方法的实施例，重复之处不再赘述。以下所使用的，术语“单元”或者“模块”可以实现预定功能的软件和/或硬件的组合。尽管以下实施例所描述的装置较佳地以软件来实现，但是硬件，或者软件和硬件的组合的实现也是可能并被构想的。
图4是本发明实施例提出的一种三维TTI介质qP波波场模拟装置的结构框图，如图4所示，本发明实施例三维TTI介质qP波波场模拟装置包括：三维TTI介质qP波波场方程计算单元1、改写计算单元2、低秩分解近似单元3、各项异性背景波场项确定单元4以及各项异性波场校正项计算单元5。
三维TTI介质qP波波场方程计算单元1，用于根据VTI介质qP波波场方程以及三维旋转因子计算出三维TTI介质qP波波场方程，其中，所述三维TTI介质qP波波场方程包括一椭圆微分算子和一标量算子；
改写计算单元2，用于通过将所述标量算子改写为标量校正算子，将所述三维TTI介质qP波波场方程分解成一各项异性背景波场项和一各项异性波场校正项；
低秩分解近似单元3，用于通过低秩分解近似计算所述椭圆微分算子；
各项异性背景波场项确定单元4，用于根据近似计算出的椭圆微分算子确定出所述各项异性背景波场项；
各项异性波场校正项计算单元5，用于根据近似计算出的椭圆微分算子以及所述标量校正算子计算出所述各项异性波场校正项。
各项异性总波场计算单元6，用于根据所述各项异性波场校正项对所述各项异性背景波场项进行校正，得到各项异性总波场。
图5是本发明实施例低秩分解近似单元的组成结构图，如图5所示，本发明实施例低秩分解近似单元3包括：时间域离散表达式转换模块7、低秩分解计算模块8以及椭圆微分算子计算模块9。
时间域离散表达式转换模块7，用于将所述三维TTI介质qP波波场方程转换为时间域离散表达式，其中，所述时间域离散表达式中包括拟微分算子；
低秩分解计算模块8，用于通过低秩分解近似所述拟微分算子的象征，得到所述拟微分算子的近似表达式；
椭圆微分算子计算模块9，用于根据所述近似表达式计算椭圆微分算子。
图6是本发明实施例时间域离散表达式转换模块的组成结构图，如图6所示，本发明实施例时间域离散表达式转换模块7包括：频散关系确定子模块701、频率域离散表达式转换子模块702和傅里叶反变换子模块703。
频散关系确定子模块701，用于将所述三维TTI介质qP波波场方程中的椭圆微分算子变换至频率-波数域得到频散关系；
频率域离散表达式转换子模块702，用于根据所述频散关系以及拟解析法将所述三维TTI介质qP波波场方程转换为频率域离散表达式；
傅里叶反变换子模块703，用于通过傅里叶反变换将所述频率域离散表达式转换为时间域离散表达式。
显然，本领域的技术人员应该明白，上述的本发明的各模块或各步骤可以用通用的计算装置来实现，它们可以集中在单个的计算装置上，或者分布在多个计算装置所组成的网络上，可选地，它们可以用计算装置可执行的程序代码来实现，从而，可以将它们存储在存储装置中由计算装置来执行，或者将它们分别制作成各个集成电路模块，或者将它们中的多个模块或步骤制作成单个集成电路模块来实现。这样，本发明不限制于任何特定的硬件和软件结合。
以上所述仅为本发明的优选实施例而已，并不用于限制本发明，对于本领域的技术人员来说，本发明可以有各种更改和变化。凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。
Claims (12)
1.一种三维TTI介质qP波波场模拟方法，其特征在于，包括：
根据VTI介质qP波波场方程以及三维旋转因子生成三维TTI介质qP波波场方程，其中，所述三维TTI介质qP波波场方程包括一椭圆微分算子和一标量算子；
将所述三维TTI介质qP波波场方程分解成一各项异性背景波场项和一各项异性波场校正项；
通过低秩分解近似计算所述椭圆微分算子；
根据近似计算出的椭圆微分算子确定出所述各项异性背景波场项；
根据近似计算出的椭圆微分算子以及标量算子计算出所述各项异性波场校正项；
根据所述各项异性波场校正项对所述各项异性背景波场项进行校正，得到各项异性总波场。
2.根据权利要求1所述的三维TTI介质qP波波场模拟方法，其特征在于，所述通过低秩分解近似计算所述椭圆微分算子，包括：
将所述三维TTI介质qP波波场方程转换为时间域离散表达式，其中，所述时间域离散表达式中包括拟微分算子；
通过低秩分解近似所述拟微分算子的象征，得到所述拟微分算子的近似表达式；
根据所述近似表达式计算椭圆微分算子。
3.根据权利要求2所述的三维TTI介质qP波波场模拟方法，其特征在于，所述将所述三维TTI介质qP波波场方程转换为时间域离散表达式，包括：
将所述三维TTI介质qP波波场方程中的椭圆微分算子变换至频率-波数域得到频散关系；
根据所述频散关系以及拟解析法将所述三维TTI介质qP波波场方程转换为频率域离散表达式；
通过傅里叶反变换将所述频率域离散表达式转换为时间域离散表达式。
4.根据权利要求1所述的三维TTI介质qP波波场模拟方法，其特征在于，所述通过低秩分解近似计算所述椭圆微分算子，包括：
通过CPU/GPU并行集群计算所述椭圆微分算子。
5.根据权利要求1所述的三维TTI介质qP波波场模拟方法，其特征在于，所述三维TTI介质qP波波场方程的表达式为：
公式中：
Azz＝1+2εsin2θ
其中：为椭圆微分算子，Se为标量算子，u为qP波波场，v0为沿对称轴方向qP波速度，ε、θ和为各向异性参数。
6.根据权利要求2所述的三维TTI介质qP波波场模拟方法，其特征在于，所述时间域离散表达式为：
u(x,t+Δt)-2u(x,t)+u(x,t-Δt)＝v0 2(x)Δt2FT-1{F(x,k)FFT[u(x,t)]}Se
公式中：
其中，x为空间矢量，k为空间波数矢量，Δt为时间步长，F(x,k)为优化拟拉普拉斯算子，FT-1为傅里叶反变换，FFT为快速傅里叶正变换。
7.根据权利要求6所述的三维TTI介质qP波波场模拟方法，其特征在于，所述拟微分算子的象征的近似表达式为：
其中，aO(x,k)为拟微分算子的象征，W(x,km)是aO(x,k)的子矩阵，其由与km有关的特定的列元素组成，W(xn,k)是aO(x,k)的另一个子矩阵，其由与xn有关的特定行元素组成，amn代表系数矩阵。
8.根据权利要求3所述的三维TTI介质qP波波场模拟方法，其特征在于，所述频散关系为：
公式中：
Azz＝1+2εsin2θ
其中：ω为圆频率，v0为沿对称轴方向qP波速度，ε、θ和为各向异性参数，kx、ky、kz分别为空间波数矢量k的三个分量。
9.根据权利要求3所述的三维TTI介质qP波波场模拟方法，其特征在于，所述频率域离散表达式为：
FT[u(x,t+Δt)-2u(x,t)+u(x,t-Δt)]＝v0 2(x)Δt2F(x,k)FFT[u(x,t)]Se
其中，x为空间矢量，k为空间波数矢量，v0为沿对称轴方向qP波速度，Δt为时间步长，F(x,k)为优化拟拉普拉斯算子，FFT为快速傅里叶正变换，Se为标量算子。
10.三维TTI介质qP波波场模拟装置，其特征在于，包括：
三维TTI介质qP波波场方程计算单元，用于根据VTI介质qP波波场方程以及三维旋转因子计算出三维TTI介质qP波波场方程，其中，所述三维TTI介质qP波波场方程包括一椭圆微分算子和一标量算子；
改写计算单元，用于将所述三维TTI介质qP波波场方程分解成一各项异性背景波场项和一各项异性波场校正项；
低秩分解近似单元，用于通过低秩分解近似计算所述椭圆微分算子；
各项异性背景波场项确定单元，用于根据近似计算出的椭圆微分算子确定出所述各项异性背景波场项；
各项异性波场校正项计算单元，用于根据近似计算出的椭圆微分算子以及标量算子计算出所述各项异性波场校正项；
各项异性总波场计算单元，用于根据所述各项异性波场校正项对所述各项异性背景波场项进行校正，得到各项异性总波场。
11.根据权利要求10所述的三维TTI介质qP波波场模拟装置，其特征在于，所述低秩分解近似单元包括：
时间域离散表达式转换模块，用于将所述三维TTI介质qP波波场方程转换为时间域离散表达式，其中，所述时间域离散表达式中包括拟微分算子；
低秩分解计算模块，用于通过低秩分解近似所述拟微分算子的象征，得到所述拟微分算子的近似表达式；
椭圆微分算子计算模块，用于根据所述近似表达式计算椭圆微分算子。
12.根据权利要求11所述的三维TTI介质qP波波场模拟装置，其特征在于，所述时间域离散表达式转换模块包括：
频散关系确定子模块，用于将所述三维TTI介质qP波波场方程中的椭圆微分算子变换至频率-波数域得到频散关系；
频率域离散表达式转换子模块，用于根据所述频散关系以及拟解析法将所述三维TTI介质qP波波场方程转换为频率域离散表达式；
傅里叶反变换子模块，用于通过傅里叶反变换将所述频率域离散表达式转换为时间域离散表达式。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201910049386.8A CN109635513B (zh) | 2019-01-18 | 2019-01-18 | 三维TTI介质qP波波场模拟方法及装置 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201910049386.8A CN109635513B (zh) | 2019-01-18 | 2019-01-18 | 三维TTI介质qP波波场模拟方法及装置 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN109635513A true CN109635513A (zh) | 2019-04-16 |
| CN109635513B CN109635513B (zh) | 2020-12-15 |
Family
ID=66062170
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201910049386.8A Expired - Fee Related CN109635513B (zh) | 2019-01-18 | 2019-01-18 | 三维TTI介质qP波波场模拟方法及装置 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN109635513B (zh) |
Cited By (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN110596754A (zh) * | 2019-09-24 | 2019-12-20 | 中国矿业大学(北京) | 一种三维TTI介质qP波与qSV波波场模拟方法 |
| CN112163187A (zh) * | 2020-11-18 | 2021-01-01 | 无锡江南计算技术研究所 | 一种超长点数高性能fft计算装置 |
| CN113485385A (zh) * | 2021-07-13 | 2021-10-08 | 中国人民解放军战略支援部队信息工程大学 | 基于误差椭圆的uuv集群编队构型设计方法 |
-
2019
- 2019-01-18 CN CN201910049386.8A patent/CN109635513B/zh not_active Expired - Fee Related
Cited By (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN110596754A (zh) * | 2019-09-24 | 2019-12-20 | 中国矿业大学(北京) | 一种三维TTI介质qP波与qSV波波场模拟方法 |
| CN112163187A (zh) * | 2020-11-18 | 2021-01-01 | 无锡江南计算技术研究所 | 一种超长点数高性能fft计算装置 |
| CN112163187B (zh) * | 2020-11-18 | 2023-07-07 | 无锡江南计算技术研究所 | 一种超长点数高性能fft计算装置 |
| CN113485385A (zh) * | 2021-07-13 | 2021-10-08 | 中国人民解放军战略支援部队信息工程大学 | 基于误差椭圆的uuv集群编队构型设计方法 |
| CN113485385B (zh) * | 2021-07-13 | 2023-11-07 | 中国人民解放军战略支援部队信息工程大学 | 基于误差椭圆的uuv集群编队构型设计方法 |
Also Published As
| Publication number | Publication date |
|---|---|
| CN109635513B (zh) | 2020-12-15 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| Song et al. | Modeling of pseudoacoustic P-waves in orthorhombic media with a low-rank approximation | |
| CN109635513A (zh) | 三维TTI介质qP波波场模拟方法及装置 | |
| CN106842306B (zh) | 一种全局优化的交错网格有限差分正演模拟方法和装置 | |
| Cheng et al. | Simulating propagation of separated wave modes in general anisotropic media, Part I: qP-wave propagators | |
| Bilbao et al. | Incorporating source directivity in wave-based virtual acoustics: Time-domain models and fitting to measured data | |
| Cannon et al. | A GPU-accelerated finite-difference time-domain scheme for electromagnetic wave interaction with plasma | |
| Liu et al. | A new kind of optimal second-order symplectic scheme for seismic wave simulations | |
| CN109490956A (zh) | 一种基于交错网格的声波波动方程正演模拟方法及装置 | |
| CN110058303A (zh) | 声波各向异性逆时偏移混合方法 | |
| Yan et al. | Optimal staggered-grid finite-difference schemes by combining Taylor-series expansion and sampling approximation for wave equation modeling | |
| Jing et al. | Dispersion-relation preserving stereo-modeling method beyond Nyquist frequency for acoustic wave equation | |
| Koene et al. | Optimal finite-difference operators for arbitrarily sampled data | |
| CN113341455A (zh) | 一种粘滞各向异性介质地震波数值模拟方法、装置及设备 | |
| Zhang et al. | Efficient anisotropic quasi-P wavefield extrapolation using an isotropic low-rank approximation | |
| Liang et al. | Decoupled wave equation and forward modeling of qP wave in VTI media with the new acoustic approximation | |
| Whiting | Energy quantity estimation in radiated acoustic fields | |
| Li et al. | A generalized average-derivative optimal finite-difference scheme for 2D frequency-domain acoustic-wave modeling on continuous nonuniform grids | |
| Xu et al. | Effective modeling and reverse-time migration for novel pure acoustic wave in arbitrary orthorhombic anisotropic media | |
| Chen et al. | A time-space domain stereo finite difference method for 3D scalar wave propagation | |
| Zhang et al. | Anisotropic viscoacoustic wave modelling in VTI media using frequency-dependent complex velocity | |
| Xu et al. | Time-space-domain temporal high-order staggered-grid finite-difference schemes by combining orthogonality and pyramid stencils for 3D elastic-wave propagation | |
| Duan et al. | An adaptive node-distribution method for radial-basis-function finite-difference modeling with optimal shape parameter | |
| Jing et al. | A highly efficient time-space-domain optimized method with Lax-Wendroff type time discretization for the scalar wave equation | |
| Xu et al. | 3D acoustic wave modeling with a time-space-domain temporal high-order finite-difference scheme | |
| Wu et al. | Mesh-free radial-basis-function-generated finite differences and their application in reverse time migration |
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 | ||
| CF01 | Termination of patent right due to non-payment of annual fee | ||
| CF01 | Termination of patent right due to non-payment of annual fee |
Granted publication date: 20201215 Termination date: 20220118 |