CN113536638B - 基于间断有限元和交错网格的高精度地震波场模拟方法 - Google Patents
基于间断有限元和交错网格的高精度地震波场模拟方法 Download PDFInfo
- Publication number
- CN113536638B CN113536638B CN202110831354.0A CN202110831354A CN113536638B CN 113536638 B CN113536638 B CN 113536638B CN 202110831354 A CN202110831354 A CN 202110831354A CN 113536638 B CN113536638 B CN 113536638B
- Authority
- CN
- China
- Prior art keywords
- finite element
- wave field
- seismic
- intermittent
- region
- 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
Links
- 238000000034 method Methods 0.000 title claims abstract description 223
- 238000004088 simulation Methods 0.000 title claims abstract description 91
- 238000004364 calculation method Methods 0.000 claims abstract description 43
- 238000000354 decomposition reaction Methods 0.000 claims abstract description 24
- 238000005516 engineering process Methods 0.000 claims abstract description 16
- 230000001629 suppression Effects 0.000 claims abstract description 4
- 238000004422 calculation algorithm Methods 0.000 claims description 18
- 238000010521 absorption reaction Methods 0.000 claims description 14
- 239000013598 vector Substances 0.000 claims description 14
- 239000011159 matrix material Substances 0.000 claims description 13
- 238000006243 chemical reaction Methods 0.000 claims description 11
- 238000007781 pre-processing Methods 0.000 claims description 10
- 238000004458 analytical method Methods 0.000 claims description 9
- 238000013459 approach Methods 0.000 claims description 8
- 238000003825 pressing Methods 0.000 claims description 7
- 230000008569 process Effects 0.000 claims description 6
- 238000005070 sampling Methods 0.000 claims description 6
- 230000004907 flux Effects 0.000 claims description 5
- 239000002245 particle Substances 0.000 claims description 5
- 238000002156 mixing Methods 0.000 claims description 4
- 238000001514 detection method Methods 0.000 claims description 3
- 238000012545 processing Methods 0.000 claims description 3
- 230000010354 integration Effects 0.000 claims description 2
- 238000012821 model calculation Methods 0.000 claims description 2
- 230000002238 attenuated effect Effects 0.000 claims 1
- 230000008859 change Effects 0.000 abstract description 4
- 239000010410 layer Substances 0.000 description 13
- 239000006185 dispersion Substances 0.000 description 6
- 230000005540 biological transmission Effects 0.000 description 5
- 238000010276 construction Methods 0.000 description 5
- 238000010586 diagram Methods 0.000 description 5
- 230000007246 mechanism Effects 0.000 description 4
- 239000004215 Carbon black (E152) Substances 0.000 description 3
- BVKZGUZCCUSVTD-UHFFFAOYSA-L Carbonate Chemical compound [O-]C([O-])=O BVKZGUZCCUSVTD-UHFFFAOYSA-L 0.000 description 3
- 230000009286 beneficial effect Effects 0.000 description 3
- 230000007547 defect Effects 0.000 description 3
- 230000000694 effects Effects 0.000 description 3
- 229930195733 hydrocarbon Natural products 0.000 description 3
- 150000002430 hydrocarbons Chemical class 0.000 description 3
- 239000007787 solid Substances 0.000 description 3
- 238000011161 development Methods 0.000 description 2
- 238000002224 dissection Methods 0.000 description 2
- 230000005284 excitation Effects 0.000 description 2
- 230000004044 response Effects 0.000 description 2
- 230000003595 spectral effect Effects 0.000 description 2
- 238000012876 topography Methods 0.000 description 2
- 230000002159 abnormal effect Effects 0.000 description 1
- 230000015572 biosynthetic process Effects 0.000 description 1
- 230000008878 coupling Effects 0.000 description 1
- 238000010168 coupling process Methods 0.000 description 1
- 238000005859 coupling reaction Methods 0.000 description 1
- 238000011438 discrete method Methods 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 239000003208 petroleum Substances 0.000 description 1
- 239000011295 pitch Substances 0.000 description 1
- 230000000644 propagated effect Effects 0.000 description 1
- 150000003839 salts Chemical class 0.000 description 1
- 238000001228 spectrum Methods 0.000 description 1
- 239000002344 surface layer Substances 0.000 description 1
- 238000012795 verification Methods 0.000 description 1
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]
-
- 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
-
- 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
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Life Sciences & Earth Sciences (AREA)
- Remote Sensing (AREA)
- Theoretical Computer Science (AREA)
- Geometry (AREA)
- General Engineering & Computer Science (AREA)
- Evolutionary Computation (AREA)
- Computer Hardware Design (AREA)
- Acoustics & Sound (AREA)
- Environmental & Geological Engineering (AREA)
- Geology (AREA)
- General Life Sciences & Earth Sciences (AREA)
- Geophysics (AREA)
- Geophysics And Detection Of Objects (AREA)
Abstract
本发明公布了一种基于间断有限元和交错网格的高精度地震波场模拟方法，通过有限差分计算、分解边界波场交换、间断有限元计算和人工边界反射压制等技术，实现在地质模型复杂构造中地震波场的高精度传播。通过区域分解将地质模型分解为构造简单区域和构造复杂区域，在构造简单区域应用交错网格方法，在构造复杂区域应用间断有限元方法，本发明对地质模型的全波场模拟具有更大的灵活性和高精度性。本发明可以提高交错网格方法数值模拟的精度，避免规则网格在地层剧烈变化地区产生虚假地震波；又可以在保证精度的同时，进一步提高间断有限元方法数值模拟的效率。
Description
技术领域
本发明属于高精度地震波场正演数值模拟技术领域，涉及间断有限元方法数值模拟、交错网格方法数值模拟、混合算法的完全匹配层技术以及波场参数转换和高精度插值方法，尤其涉及一种基于间断有限元和交错网格的高精度地震波场模拟方法，可用于复杂地表、剧烈构造变化地区高精度的地震波波场数值模拟并获取数据。
背景技术
中国的西部地区蕴含着大量的油气储层，生产的油气支撑着国家经济的发展。但是，西部地区的地表地下环境剧烈变化，构造非常复杂。未探明的油气资源储量依旧相当惊人。可是，勘探难度也相当巨大，四川盆地的地表是丘陵，而塔里木盆地的地表被沙漠所覆盖；这些地区的地下储藏着丰富的石油资源，拥有大量的储层，尤其是碳酸盐岩储层。地表条件异常复杂，剧烈变化，地下构造也非常复杂，断层和孔洞发育，褶皱强烈，地层变化大等。这些复杂的地形环境使得西部地区的地震资料普遍信噪比低，地震资料处理和解释困难，不利于勘探开发。
四川地区、鄂尔多斯地区和塔里木地区有大量的碳酸盐岩储层。裂缝和孔洞是碳酸盐岩储层油气主要的储集体。但是，对于这些地区的地震资料和储层特征有待研究，通过地震波波场正演模拟方法理解地震波的传播规律和响应机制对该地区油气勘探变得非常重要，尤其是在复杂的构造区域。从根本上理解地震波的传播机制，地震波波场数值模拟是一种非常有效的技术手段。
其中，Alterman和Karal(1968)首先将有限差分方法引入地震学领域求解层状介质地震波传播[1]。Virieux(1986)提出有限差分的交错网格形式求解在不均匀介质中地震波传播[2]。有限差分方法有计算效率高，存储量小，编程简单等优点；但是有限差分使用统一的规则网格，对复杂的地形条件缺乏灵活性。而有限元方法是在复杂地表和复杂构造地区进行地震波场数值模拟最有效的技术，但是有限元方法需要进行质量矩阵求逆[3]，以及求解一个大型的线性代数方程，需要巨大的计算存储空间。间断有限元方法是将有限体积法与有限元方法结合的一类新型的高精度数值模拟方法[4]。利用有限体积方法的数值流和有限元方法的三角网格对地震波进行数值模拟。具有高精度、低频散的优点，同时间断有限元方法也具有许多其它独特的优势：一是允许基函数在单元界面之间可以不连续；二是只需要相邻单元的信息就可以进行波场更新；三是对大型的线性代数矩阵解耦计算[5][6]。
在实际的地质模型当中使用有限差分方法，虽然效率高，但是在孔洞、起伏地表等复杂构造地区模拟精度显然是较低的；只使用间断有限元方法，数值模拟精度高，但是需要对每个三角单元逐个计算，效率较低。尤其在油气勘探领域纯粹使用一种数值方法进行正演模拟是不够的。既要考虑模拟的精度，又要考虑效率因素的影响。为了解决这个问题，在有限差分方法、有限元法、伪谱法、谱元法等基础上发展了多种混合方法。其中有一种混合算法是将不同方法在原理上进行结合，在横向上使用有限元方法或者谱元法，纵向上使用有限差分方法[7]。这种方法提高了数值模拟的精度，但是使用的依然是规则网格，缺乏灵活性，也会产生绕射现象。还有一种可行的方法是变网格技术，通过组装不同间距的网格覆盖不同区域，对交错网格采用非统一网格的技术，但是该方法在不同的网格下需要推导新的差分算子，不利于勘探领域大规模计算。
现有技术还基于区域分解思想，在不同的区域应用不同的计算方法。例如，Moczo等(1997)首先提出将有限元与二阶交错网格方法结合[8]；之后，Ma等(2004)将其推广到四阶精度[9]。黄自萍等(2004)也提出相似的混合算法[10]；马德堂等(2004)将有限元与伪谱法结合[11]；在2016年，杨顶辉等将二阶近似解析离散方法和间断有限元方法结合[12]。
在油气勘探领域，对地质模型进行正演模拟时，由于低速地表层和储层的存在，使用二阶精度或者是四阶精度数值模拟方法得到结果的精度是不够的。现有的地震波场数值模拟技术很难有效正确地进行地震资料处理和解释，对复杂构造地区的数值模拟精度和效率不足。
参考文献：
[1]Alterman,Z.,and Karal,F.C..Propagation of elastic wave in layeredmedia by finite difference method.Bulletin of the Seismological Society ofAmerica,1968,58(1),367-398.
[2]Virieux J..P-SV wave propagation in heterogeneous media:Velocity-stress finite-difference method.Geophysics,1986,51,889–901.
[3]Padovani E.,Priolo E.and Seriani G..Low and high order finiteelement method:experience in seismic modeling.Journal of ComputationalAcoustics,1994,2,371-422.
[4]Cockburn B.and Shu C.W..The Runge-Kutta discontinuous Galerkinmethod for conservation laws.V:Multidimensional systems.Journal ofComputational Physics,1998,141,199–224.
[5]Dumbser,M.,M..An arbitrary high-order discontinuous Galerkinmethod for elastic waves on unstructured meshes—II.The three-dimensionalisotropic case.Geophysical Journal International,2006,167,319–336.
[6]He X.J.,Yang D.H.,Ma X.and Lang C..Dispersion-dissipation analysesof the triangle-based discontinuous Galerkin method for scalar waveequation.Geophysical Journal International,2019,218,1174-1198.
[7]Liu T.,Hu T.Y.Sen M.K.et al.A hybrid scheme for seismic modellingbased on Galerkin method.Geophysical Journal International,2008,186,1165-1178.
[8]Moczo P.,Bystricky E.,Kristek J.,Carcione J.M.and BouchonM..Hybrid modeling of P-SV seismic motion at inhomogeneous viscoelastictopographic structures.Bulletin of the Seismological Society of America,1997,87,1305-1323.
[9]Ma S.,Archuleta R.J.and Liu P..Hybrid modeling of elastic P-SVwave motion:A combined finite-element and staggered-grid finite-differenceapproach.Bulletin of the Seismological Society of America,2004,94,1557–1563.
[10]黄自萍,张铭,吴文青,董良国.弹性波传播数值模拟的区域分裂法.地球物理学报,2004,6,1094-1100.
[11]马德堂,朱光明.有限元法与伪谱法混合求解弹性波动方程.地球科学与环境学报,2004,1,61-64.
[12]Yang D.H.,He X.J.,Ma,X.,Zhou Y.J.and Li,J.S..An optimal nearlyanalytic discrete-weighted Runge-Kutta discontinuous Galerkin hybrid methodfor acoustic wavefield modeling.Geophysics,2016,81,T251-T263.
发明内容
针对上述现有技术存在的问题，克服地震波场数值模拟精度和效率的不足，本发明提供一种基于间断有限元方法和交错网格方法结合的高精度地震波场数值模拟方法，既具有交错网格方法的高效性又具有间断有限元方法的灵活性，还能有效提高地震波场的数值模拟精度。
现有的方法中，有限差分方法和间断有限元方法均为单独适用的方法，两种方法不可兼容使用，无法同时满足效率和精度的要求。本发明提出使用高阶插值方式，使得间断有限元方法和交错网格方法耦合，直接进行数据交换，在新构建的同一个地质模型中执行地震波场模拟。而且，本发明根据频散特性的精度要求，在相应的计算区域里使用不同的阶数，可以提高在复杂地质模型中地震波的传播精度。根据模型的复杂程度调整网格模式，本发明对地质模型的网格灵活性比间断有限元方法更好，又能够提高计算效率，有利于帮助理解地震波在复杂构造中的传播规律和响应机制。
本发明是一种全波场数值模拟混合方法。通过考虑模型的复杂程度，将地质模型(本说明书中表示速度模型和密度模型)分解为构造简单区域和构造复杂区域进行计算，构造复杂区域主要有低速的起伏近地表层，孔洞和断层，高速盐丘体等。引入3阶总变差不变的龙格-库塔格式作为时间离散格式，增强该混合算法的稳定性，本方法具有很大的灵活性和高精度性；并采用完全匹配层技术对人工边界进行吸收，以压制边界反射效应。同时，通过高阶拉格朗日插值公式对混合算法的边界进行插值处理，使得每个时刻波场变量在该边界稳定地传播，不发生边界反射现象。在提高数值模拟精度同时，又可以进一步提高全波场模拟的计算效率。将本发明方法具体应用到实际的地质模型当中，与实际的地震记录相比，本发明具有很高的数值模拟精度。
本发明的核心是：本发明充分考虑地震波场数值模拟的精度和效率问题，根据地质模型的复杂程度，将初始地质模型进行分解，更加灵活自由地进行全波场地震数值模拟；并通过采用完全匹配层技术解决混合算法人工边界难吸收问题。本发明在构造简单区域应用交错网格方法提高计算效率，在构造复杂区域应用间断有限元方法提高数值精度，对复杂的地质模型应用高精度的数值模拟算法。和常规交错网格模拟方法相比，本发明可以提高地震波场数值模拟的精度，避免交错网格在地层剧烈变化地区产生的虚假地震波；和间断有限元模拟方法相比，本发明可以在保证精度的同时，进一步提高地震波场数值模拟的效率，对复杂构造的地震波场进行更加精细的模拟。
本发明提供的技术方案如下：
一种基于间断有限元和交错网格的高精度地震波场模拟方法，构建本发明所需的地质模型，通过区域分解方法，将地质模型分解为构造简单区域和构造复杂区域，在构造简单区域应用交错网格方法，在构造复杂区域应用间断有限元方法；模型输入数据包括交错网格地质模型数据、间断有限元地质模型数据，以及模型区域分解参数。在模型计算过程中，通过交错网格计算、分解边界波场交换、间断有限元计算、数值模拟阶数设置、数值流设置和人工边界反射压制等技术，实现在地质模型复杂构造中地震波场的高精度传播模拟，并获得高精度的勘探地震波传播数据。包括以下步骤(图1)：
A.构建本发明地震波场模拟所需的地质模型，包括：
A1.根据初始地质模型(速度模型和密度模型)，获取交错网格和间断有限元地质模型数据，以及模型区域分解参数；
A2.构建震源子波模块，根据时间采样率、最大模拟时间、震源子波的类型及震源位置，获取地震子波数据；
A3.构建观测系统模块，接收地震波场，用于建立地面接收地震数据记录，确定炮点和检波点的位置关系；
B.进行数据前处理和混合算法稳定性分析；
B1.根据模型区域分解参数，对分解边界进行处理；
包括间断有限元方法和交错网格方法边界确定，间断有限元方法和交错网格方法之间参数变量的转化关系以及时间插值；
B2.对间断有限元地质模型数据进行预处理，包括计算质量矩阵、寻找每个单元的相邻单元、以及每个单元边的法向量求取等；
B3.数值流设置和时间演化格式选择；选择在间断有限元计算区域使用数值通量的种类，以及时间离散格式的选择；
B4.混合算法精度和稳定性分析；对间断有限元计算区域和交错网格计算区域分别进行精度和稳定性分析，确定两种方法对应计算区域的精度，选择各自对应的阶数，使得本发明在两个计算子域中能够保持精度一致性；确定采样时间大小，保持数值模拟的稳定性。
C.人工边界反射压制：由于有限的传播区域，在人工边界处会产生反射，在不同的计算区域应用完全匹配层吸收边界条件；
D.进行波场模拟演化计算：
D1.建立初始化模块，用于将所有的变量参数，包括波场、地震记录等初始化为0；
D2.加载震源子波模块、观测系统模块和人工边界反射压制模块；在相应的地质模型中输入人工震源，激发地震波场，对模型边界产生的人工反射波进行吸收，避免影响真实地震波场的传播，在表面对地震波进行接收，获取表面地震观测数据；
D3.计算模块：
从n时刻到n+1时刻波场演化过程如下(图2中黑色表示n+1时刻，灰色表示n时刻)：在n时刻，交错网格(构造简单)区域内的变量值表示为：(P)n+1/2、(Vx)n、(Vz)n，采用高阶交错网格方法计算n+1时刻的波场变量值表示为：(P)n+3/2、(Vx)n+1、(Vz)n+1，如图2(a)所示。其中，(P)nn+1/2和(P)n+3/2分别表示n+1/2和n+3/2时刻的正应力，(Vx)n和(Vx)n+1分别表示n和n+1时刻在x方向上的质点速度，(Vz)n和(Vz)n+1分别表示n和n+1时刻在z方向上的质点速度。
接着根据上一步的信息，利用时间转换公式和高精度插值方法得到n+1时刻构造复杂(间断有限元)区域边界处的变量值(P)n+1、(Vx)n+1、(Vz)n+1，如图2(b)所示。根据得到的间断有限元区域n+1时刻的边界信息，计算构造复杂区域中所有的变量值，得到n+1时刻的(P)n+1、(Vx)n+1、(Vz)n+1，如图2(c)所示。最后，根据图2(c)中得到的间断有限元区域n+1时刻的信息，再次通过时间转换公式和高精度插值公式反向得到交错网格区域边界处n+1时刻的变量值(P)n+3/2、(Vx)n+1、(Vz)n+1，完成时间迭代更新，如图2(d)所示。其中，(P)n+1和(P)n+3/2分别表示n+1和n+3/2时刻的正应力，(Vx)n+1和(Vz)n+1表示n+1时刻在x和z方向上的质点速度。此时整个区域上所有的变量值已经从n时刻更新到n+1时刻(计算区域全部变成黑色)；
D4.完成波场更新，得到交错网格波场和间断有限元波场系数；将间断有限元波场系数转化为相应区域对应的波场，获得该地质模型某一时刻的波场快照。
E.建立地震波场生成模块，用于生成某一时刻的地震波场切片，当达到最大模拟时间时输出地震记录结果，即实现基于间断有限元法和交错网格法的高精度地震波场的模拟。
与现有技术相比，本发明的有益效果是：
交错网格方法和间断有限元方法都有各自的优缺点。而本发明基于区域分解的思想，通过在构造简单区域使用交错网格方法，在构造复杂区域使用间断有限元方法，通过插值方式让本不可能兼容的两种数值方法进行有效融合，实现使用交错网格方法和间断有限元方法对同一个地质模型进行高精度数值模拟，且在模型分解边界直接进行波场信息交换。摒弃交错网格方法和间断有限元方法这两种方法的缺点，将两种方法的优点同时作用在地震波场数值模拟当中，提高复杂构造模型数值模拟的精度，又提高了地震波场的计算效率。本发明的技术优势在于：(一)本发明提出使用高阶插值方式，使得间断有限元方法和交错网格方法耦合，直接进行数据交换，在新构建的同一个地质模型中执行地震波场模拟；设计完全匹配层技术应用于混合方法，成功解决人工边界反射和不同算法的兼容问题；通过交错网格方法计算构造简单区域，提高计算效率；在构造复杂区域，如孔洞、断层等，通过应用间断有限元方法，提高数值模拟的精度；
(二)在交错网格或间断有限元计算区域可以分别选择各自所用的阶数，进一步提高数值模拟的精度；
(三)通过高精度插值方式，使得不同方法之间的应力和速度变量稳定地传播；(四)间断有限元使用的三角形单元边可以耦合多个交错网格网格点，进一步提高计算效率。
附图说明
图1是本发明提供的一种基于间断有限元和交错网格的高精度地震波场模拟方法的程序框图。
图2是本发明提供的波场演化示意图：
其中，黑色表示n+1时刻波场，灰色表示n时刻波场；(a)表示交错网格区域首先完成计算更新；(b)为利用交错网格(构造简单)区域波场得到分解边界的波场值；(c)表示间断有限元区域完成波场更新；(d)表示利用间断有限元区域的波场得到交错网格区域的边界值。
图3是本发明进行验证用的起伏地表均匀介质模型。
其中，黑色五角星表示震源；黑色实线表示区域分解线。
图4是起伏地表进行三角网格剖分示意图。
图5是本发明进行实施例数值模拟的波场快照图：
其中，白色实线表示区域分界线，黑色实线表示起伏地表边界线；(a)是0.3秒时刻波场快照图；(b)是0.4秒时刻波场快照图；(c)是0.5秒时刻波场快照图；(d)是0.6秒时刻波场快照图。
图6是本发明进行实施例和其它数值方法的波场快照对比图：
其中，波场快照对应时刻为0.7秒；(a)表示本发明方法得到的结果；(b)表示10阶交错网格方法结果；(c)表示2阶间断有限元方法得到结果。
图7是本发明在实施例中用到的某油田地区速度模型构造图：
其中，虚线框内表示构造复杂区域，虚线框外表示构造简单区域。
图8是利用本发明在实施例中某油田地区得到的炮集记录。
图9是利用本发明在实施例中某油田地区得到的炮集记录的局部放大图：
其中，(a)表示区域为在图8中的位置1虚线框内；(b)表示区域为在图8中的位置2的虚线框内。
图10是利用本发明在实施例中某油田地区得到的炮集记录与实际地震记录对比图：
其中，(a)表示实际地震数据的局部区域；(b)表示本发明得到的模拟地震记录局部区域。
具体实施方式
下面结合附图，通过实施例进一步描述本发明，但不以任何方式限制本发明的范围。
本发明是一种基于间断有限元和交错网格的高精度地震波场模拟方法，为了解决交错网格数值方法对复杂模型适应性差的问题，以及间断有限元方法计算效率低等问题，通过区域分解方法，将地质模型分解为构造简单区域和构造复杂区域，在构造简单区域应用交错网格方法，在构造复杂区域应用间断有限元方法，本发明对地质模型的全波场模拟具有更大的灵活性。引入3阶总变差不变的龙格-库塔时间格式，使得本发明方法具有更大的稳定性，同时提出完全匹配层技术对混合算法的人工边界反射进行压制，并使得间断有限元方法和交错网格方法的吸收边界能够兼容；本发明可以根据频散分析的精度需要，在各自计算区域应用不同阶数的数值模拟方法，使得交错网格方法计算的区域和间断有限元方法计算的区域保持相同的数值模拟精度。通过实施例验证本发明的正确性，并将本发明应用到某油田地区的实际模型当中，与实际地震数据进行了对比，对该油田地区地震波的传播机制有更好的解释。本发明需要输入地质模型，复杂构造区域三角网格数据和分解边界参数，通过对地质模型进行数据前处理，包括边界波场转化，质量矩阵，寻找每个单元的相邻单元、以及每个单元边的法向量求取，实现全波场数值模拟。
本发明采用正演模拟的计算方程为标量应力-速度方程，如公式(1)：
上式中，P表示正应力，Vx和Vx分别是在x方向和z方向上的质点速度。κ＝ρc2(c表示波速)表示体积模量，ρ表示介质密度。
利用间断有限元方法的基本原理推导标量应力-速度方程，应用梯度算子 公式(1)可以重新被写成下列的向量形式：
其中，U＝(P，Vx，Vz)T表示在声波方程中的未知数向量，H(AU，BU)是数值流向量，A和B表示雅克比矩阵，其中
公式(2)对应的积分方程可以简写为：
其中n表示单元边界的外法向量，其值为(nx，nz)；/>表示Ωm的边界。Uint表示单元Ωm内部趋近于边界/>的值，Uext表示单元Ωm外部趋近于边界/>的值。H(Uint，Uext，n)是数值流函数，可以写成下列的统一形式：
对于中心数值流：
局部Lax-Friedrichs数值流:
Godunov数值流:
在这里，是一个常数，|An|是由对雅克比矩阵An＝nxA+nzB的特征值矩阵取绝对值后获得。
利用交错网格方法对公式(1)进行空间离散，可以得到半离散形式为：
其中，Δx和Δz是交错网格离散空间步长。交错网格方法拥有2M阶精度，其中an是交错网格系数：
图1所示是本发明提供的一种基于间断有限元和交错网格的高精度地震波场模拟方法步骤，方法包含以下步骤:
1)将初始地质模型数据按照本发明的要求：交错网格区域需要网格点数据，间断有限元计算区域需要三角网格数据，将模型数据划分为交错网格地质模型数据和间断有限元地质模型数据。
2)确定模型区域分解参数，使间断有限元和交错网格两种数值方法分别计算相应的区域，不相互产生串扰。
3)根据时间采样率Δt、最大模拟时间tmax，震源子波的类型f及震源位置，定义震源子波；根据实际地震数据在地表激发和接收规则，将震源子波放置在地表激发地震波信号，同时设置地震波场的演化网格间距和接收检波器间隔，检波器的深度位置等参数，确定地表检波点数量，定义好地面接收观测系统。
4)对数据进行前处理：
本发明是基于两种方法结合的混合算法，在波场的传递过程中，波场需要相互转化。通过拉格朗日插值法将交错网格区域规则网格点上的波场值插值到间断有限元边界上。需要提前使用变量存储好它们之间的位置和转化关系。具体插值公式为：
其中，P为构造简单区域网格点上的正应力波场；Pj-1、Pj、Pj+1分别为在网格点j-1，j，j+1上的波场值。
进行前处理包括如下：
通量形式的间断有限元公式为：
在同一时间层面，对间断有限元模型数据中的每个单元执行方程(12)时，只有含变量Uext的边界积分项需要外部波场来完成与相邻元素的信息交换。此外，本发明在间断有限元和交错网格的接口中采用三点高斯积分方法计算方程(12)中的含有Uext的边界积分项。也就是说，当波场从构造简单(交错网格)区域传播到构造复杂(间断有限元)区域时，本发明只需要知道变量Uext在界面的高斯积分点上的波场值，就可以执行方程(12)。这也意味着完成了间断有限元和交错网格方法在界面中的信息交换。
本发明应用间断有限元方法，在进行波场传递时，需要相邻单元信息。输入的数据是没有相应的信息，需要设置变量提前寻找并存储相邻单元信息，计算本单元边的法向量信息，计算间断有限元方法需要的质量矩阵、刚度矩阵以及通量矩阵等。
为了压制震源噪声、网格噪声等假波问题，除采用更加严格的数值流方式，有局部Lax-Friedrichs数值流和Godunov数值流，同时也为了保持两种方法结合的稳定性，统一使用3阶龙格库塔时间格式，相应公式如下所示。
C(1)＝C(n)+ΔtL(C(n))，
上式中，C表示波场向量系数，具体表示式(1)中标量应力-速度方程的应力和速度分量，L表示离散空间计算算子，n表示时刻。
5)人工边界反射压制：利用完全匹配层技术使得两种方法的吸收边界可以兼容。
在混合计算区域外面设置一个人工区域，当地震波传播到吸收区域时对波进行吸收衰减。通过边界波场衰减达到吸收边界波场的效果，能很大程度上减少人工边界反射的干扰，边界条件的生成主要包括如下步骤：首先设置衰减系数得到衰减因子，本发明采用衰减因子如下：
其中，L为匹配层宽度，R为理想边界层反射系数(取0～10-12)，vmax为地质模型中最大的纵波速度，x表示相应位置到吸收边界的距离，d(x)和d(z)分别表示横向和纵向上的衰减因子。当d(x)≠0，d(z)≠0表示地震波衰减状态；当d(x)＝0，d(z)＝0时表示不衰减。
本发明将公式(14)中的衰减因子直接加载在交错网格与间断有限元方法的边界进行人工边界反射压制。将衰减因子加载到两种数值方法对应的方程当中，交错网格方法与间断有限元方法应用的标量应力-速度方程形式不同，需要分别进行加载处理。
间断有限元对应的吸收方程变形为：
交错网格方法对应的吸收方程为：
对本发明的稳定性分析也是需要关注的标准之一。对于本发明来说，需要同时考虑两种算法的稳定性。
对于交错网格方法：
其中，ΔtSSM为交错网格方法时间步长；h为交错网格的网格尺度；vmax为最大速度。
对于间断有限元方法：
上式中，ΔtDGM为间断有限元对应的时间步长；αSSM和ΔDGM分别对应交错网格方法和间断有限元方法的最大库朗数。因此，对于本发明提出的混合方法，时间步长Δt必须满足：
Δt≤min{ΔtSSM，ΔtDGM} 式(19)
6)进行波场演化计算：
在进行波场演化的初始时刻，应力和速度分量需要全局初始化为零，才不会影响到接下来时间迭代。在每一时刻，将地震子波加载到相应的位置，可以是间断有限元区域，也可以是交错网格区域。
不同计算区域之间波场传递是本发明最重要的核心内容，以n时刻到n+1时刻波场演化为例，对本发明的波场传递做更细致的说明，首先在n时刻，交错网格区域内的变量值为(P)n+1/2、(Vx)n、(Vz)n，采用高阶交错网格方法计算n+1时刻的变量值(P)n+3/2、(Vx)n+1、(Vz)n+1，如图2(a)所示。接着根据上一步的信息，利用时间转换公式和高精度插值方法得到n+1时刻间断有限元区域边界处(分解区域交界)的变量值(P)n+1、(Vx)n+1、(Vz)n+1，如图2(b)中黑色实线所示。根据得到的间断有限元区域tn+1时刻的边界信息，计算间断有限元区域的波场值，得到n+1时刻的(P)n+1、(Vx)n+1、(Vz)n+1，如图2(c)所示，三角形单元网格全部变成黑色。根据(c)中得到的间断有限元区域n+1时刻的信息，再次通过时间转换公式和高精度插值公式得到交错网格区域边界处n+1时刻的变量值(P)n+3/2、(Vx)n+1、(Vz)n+1，完成时间迭代更新，如图2(d)所示。此时整个区域上所有的变量值已经从n时刻更新到n+1时刻，在图2中即表示从灰色全部变成黑色。
7)输出波场时间切片结果和地面地震记录结果。
本发明中，高阶交错网格方法具有较高的精度和较好的数值稳定性；间断有限元拥有低频散、高精度的优点。将两种方法进行结合对复杂的地质模型进行高精度，高效率的数值模拟。
为了让本发明的目的、技术和优势更加清晰，下面将通过实施例做进一步的说明。本发明在一个简单的起伏地表均匀介质模型中进行数值正演模拟，并同时比较了只使用交错网格和间断有限元的波场。实施例地质模型的水平长度为2500米，深度也为2500米，其速度模型结构如图3所示。模型的顶部是一个起伏地表的地形，下面是均匀介质。将地质模型分解为两个区域(图中用黑线分开)，地表地形起伏，且是低速层，因此上面应用间断有限元方法进行数值模拟，下面应用交错网格方法进行数值模拟。顶部的地形用distmesh软件进行网格剖分，其剖分结果示意图如图4所示，三角网格对复杂起伏地形有更好的灵活性。在本实施例中，震源采用主频为20Hz的雷克子波。
数值模拟过程中均采用完全匹配层技术。图5是利用本发明方法进行数值模拟的波场结果，(a)-(d)分别表示0.3秒，0.4秒，0.5秒和0.6秒时刻的波场快照结果。图中的黑色实线表示起伏地表边界，白色实线表示区域分解线。图6是在0.7秒时刻不同方法的波场快照对比图；其中，图6(a)是使用本发明方法的波场快照结果，图6(b)是使用交错网格方法的波场快照结果，图6(c)是使用间断有限元方法的波场快照结果。
由于交错网格方法自身的局限性，对地形的适应性差。交错网格方法利用阶梯状网格去近似复杂地形，在速度差异大时，在网格点附近会出现明显的绕射波，该绕射波在实际的复杂地层介质中是不存在的，是由于交错网格方法对复杂构造适应性差产生的。如图6(b)箭头所表示的就是绕射波，而图6(a)是用本发明方法进行数值模拟的结果，和使用间断有限元方法的图6(c)进行比较可知，两种方法的波场精度是一致的，箭头所指的位置并没有绕射波。而使用间断有限元方法的效率几乎是本发明方法的两倍，在一定程度上验证了本发明具有较高的效率。结果对比证明了本发明的正确性。
为了进一步验证本发明在实际地质模型当中的有效性，将本发明应用到某油田实际的地质模型当中。图7是本发明在实施例中用到的某油田地区的二维速度模型构造图，其中黑色虚线框所表示的地区用间断有限元方法，其它地区用交错网格方法。对该地质模型进行全波场模拟，震源设置在地表的中心位置，同时在地表进行接收。
对该地质模型应用本发明方法进行数值模拟，图8是模拟得到的全波场地震记，从合成地震记录中可以看出，地震记录上的同相轴连续，同相轴的能量均衡，没有异常地震波出现。结果验证了本发明利用拉格朗日高阶插值方法可以避免不同方法之间由于算法的差异导致波场传递过程中在分解边界出现反射效应，使得两种方法在时间迭代过程中波场平稳地传播，也证明了本发明采用完全匹配层技术可以将两种数值方法很好地耦合在一起。同时，为了进一步验证本发明方法具有高精度性。在图8地震记录上的1和2虚线框区域内进行局部放大。局部放大图如图9所示，同相轴之间连续，能量均衡。反射波、直达波等没有出现频散现象。由于速度差异导致在不同检波器上到时不同，同相轴之间交错也不相互影响。
此外，为了更进一步验证本发明对地质模型高精度的数值模拟，将利用本发明模拟的全波场结果与实际的地震数据进行了对比。图10是本发明的模拟结果与实际地震记录的对比结果图。图10(a)是实际地震数据局部放大图，图10(b)是模拟的地震记录局部放大图，图中的黑色箭头表示的是对应地层记录。数值模拟结果与实际地震记录结果一致，从而验证本发明对实际地质模型模拟的正确性。本发明充分利用不同方法的优势，将不同的方法结合起来，取其精华去其糟粕。不仅可以实现对复杂构造的高精度数值模拟，克服交错网格自身的缺陷，还可以兼顾效率。
需要注意的是，公布实施例的目的在于帮助进一步理解本发明，但是本领域的技术人员可以理解：在不脱离本发明及所附权利要求的范围内，各种替换和修改都是可能的。因此，本发明不应局限于实施例所公开的内容，本发明要求保护的范围以权利要求书界定的范围为准。
Claims (11)
1.一种基于间断有限元和交错网格的高精度地震波场模拟方法，将地质模型分为构造简单区域模型和构造复杂区域模型；在构造简单区域采用交错网格方法进行地震波场模拟，在构造复杂区域采用间断有限元方法进行地震波场模拟；模型输入数据包括交错网格地质模型数据、间断有限元地质模型数据和模型区域分解参数；在模型计算过程中，通过交错网格计算、分解边界波场交换、间断有限元计算、数值模拟阶数设置、数值流设置和人工边界反射压制技术，实现在地质模型复杂构造中地震波场的高精度传播模拟，并获得高精度的勘探地震波传播数据；包括如下步骤：
A.构建地震波场模拟所需的地质模型，通过区域分解方法，将地质模型分为构造简单区域和构造复杂区域的不同情形；获取交错网格地质模型数据、间断有限元地质模型数据和模型区域分解参数；获取地震子波数据；接收地震波场数据记录；
B.进行数据前处理和混合算法稳定性分析；
C.通过进行人工边界反射压制，利用完全匹配层技术使得间断有限元和交错网格方法的吸收边界兼容，在不同的计算区域应用完全匹配层边界条件；
D.进行波场模拟演化计算，包括：
D1.进行初始化，将所有的变量参数初始化为0；变量参数包括波场和地震记录；
D2.在相应的地质模型中输入人工震源，激发地震波场，对模型边界产生的人工反射波进行吸收，在表面对地震波进行接收，获取表面地震观测数据；
D3.计算从n时刻到n+1时刻波场演化：
在n时刻，交错网格区域即构造简单区域内的变量值为(P)n+1/2、(Vx)n、(Vz)n；采用高阶交错网格方法计算n+1时刻的波场变量值(P)n+3/2、(Vx)n+1、(Vz)n+1；其中，P表示正应力，Vx和Vz分别是在x方向和z方向上的质点速度；计算波场演化包括：
首先，利用时间转换公式和高精度插值方法；得到n+1时刻间断有限元区域即构造复杂区域边界处的变量值(P)n+1、(Vx)n+1、(Vz)n+1；
然后，根据得到的间断有限元区域n+1时刻的边界信息，计算构造复杂区域的变量值，得到n+1时刻的(P)n+1、(Vx)n+1、(Vz)n+1；
最后，根据得到的间断有限元区域n+1时刻的信息，再次通过时间转换公式和高精度插值公式反向得到交错网格区域边界处n+1时刻的变量值(P)n+3/2、(Vx)n+1、(Vz)n+1，完成时间迭代更新；
此时区域所有的变量值即从n时刻更新到n+1时刻的波场；
D4.完成波场更新，得到交错网格波场和间断有限元波场系数；将间断有限元波场系数转化为相应区域对应的波场，获得地质模型某一时刻的波场快照；
E.生成地震波场：
生成某一时刻的地震波场切片，当达到最大模拟时间时输出地震记录结果，即实现基于间断有限元法和交错网格法的高精度地震波场的模拟。
2.如权利要求1所述基于间断有限元和交错网格的高精度地震波场模拟方法，其特征是，步骤A具体包括：
A1.初始地质模型包括速度模型和密度模型；根据初始地质模型，获取交错网格地质模型数据、间断有限元地质模型数据和模型区域分解参数；
A2.根据时间采样率、最大模拟时间、震源子波的类型及震源位置，获取地震子波数据；
A3.观测和接收地震波场，建立地面接收地震数据记录，确定炮点和检波点的位置关系。
3.如权利要求1所述基于间断有限元和交错网格的高精度地震波场模拟方法，其特征是，步骤B具体包括：
B1.根据模型区域分解参数，对分解边界进行处理；包括：间断有限元方法和交错网格方法边界确定、间断有限元方法和交错网格方法之间参数变量的转化关系以及时间插值；
B2.对间断有限元地质模型数据进行预处理，包括计算质量矩阵、寻找单元的相邻单元、以及求取每个单元边的法向量；
B3.选择数值流设置和时间演化格式；包括：选择在间断有限元区域使用数值通量的种类，和时间离散格式；
B4.混合算法精度和稳定性分析；
对间断有限元区域和交错网格区域分别进行精度和稳定性分析，确定两种方法对应计算区域的精度，选择各自对应的阶数，使得计算子域的精度保持一致性；确定采样时间大小，保持数值模拟的稳定性。
4.如权利要求1所述基于间断有限元和交错网格的高精度地震波场模拟方法，其特征是，步骤D中，利用间断有限元方法推导模拟计算方程，应用梯度算子标量应力-速度方程表示为如下向量形式：
其中，U＝(P，Vx，Vz)T表示在声波方程中的未知数向量，H(AU，BU)是数值流向量，A和B表示雅克比矩阵；
式(2)对应的积分方程表示为：
其中，n表示单元边界的外法向量，其值为(nx，nz)；/>表示Ωm的边界；Uint表示单元Ωm内部趋近于边界/>的值，Uext表示单元Ωm外部趋近于边界/>的值；H(Uint，Uext，n)是数值流函数。
5.如权利要求4所述基于间断有限元和交错网格的高精度地震波场模拟方法，其特征是，雅克比矩阵表示为：
其中，κ＝ρc2，表示体积模量；c表示波速；ρ表示介质密度。
6.如权利要求4所述基于间断有限元和交错网格的高精度地震波场模拟方法，其特征是，数值流函数H(Uint，Uext，n)表示为如下形式：
对于中心数值流：
对于局部Lax-Friedrichs数值流：
对于Godunov数值流：
其中，是常数，|An|是由对雅克比矩阵An的特征值矩阵取绝对值后获得。
7.如权利要求6所述基于间断有限元和交错网格的高精度地震波场模拟方法，其特征是，步骤B中，通量形式的间断有限元前处理表示为：
在同一时间，对间断有限元模型数据中的每个单元执行式(12)时，只有含变量Uext的边界积分项需要外部波场完成与相邻元素的信息交换；进一步在间断有限元和交错网格的接口中采用三点高斯积分方法计算式(12)中含有Uext的边界积分项。
8.如权利要求6所述基于间断有限元和交错网格的高精度地震波场模拟方法，其特征是，插值具体是通过拉格朗日插值法将交错网格区域规则网格点上的波场值插值到间断有限元边界上；插值方法表示为：
其中，P为在构造简单区域网格点上的正应力波场；Pj-1、Pj、Pj+1分别为在网格点j-1，j，j+1上的波场值。
9.如权利要求6所述基于间断有限元和交错网格的高精度地震波场模拟方法，其特征是，具体使用3阶龙格库塔时间格式，表示为：
式中，C表示波场向量，包括标量应力-速度方程的应力和速度分量；L表示离散空间计算算子，n表示时刻。
10.如权利要求9所述基于间断有限元和交错网格的高精度地震波场模拟方法，其特征是，采用人工边界反射压制技术，在混合计算区域外面设置一个人工区域，当地震波传播到吸收区域时对波进行吸收衰减；通过边界波场衰减达到吸收边界波场，减少人工反射的干扰；生成边界条件包括如下步骤：
首先设置衰减系数得到衰减因子，表示如下：
其中，L为匹配层宽度，R为理想边界层反射系数，取0～10-12，vmax为地质模型中最大的纵波速度，x表示相应位置到吸收边界的距离，d(x)和d(z)分别表示横向和纵向上的衰减因子；
将式(14)中的衰减因子加载在交错网格与间断有限元方法的边界进行边界发射压制；对交错网格方法与间断有限元方法需分别加载衰减因子；
间断有限元的吸收方程变形为：
交错网格方法的吸收方程为：
稳定性分析包括：
对于交错网格方法：
其中，ΔtSSM为交错网格方法时间步长；h为交错网格的网格尺度；vmax为最大速度；
对于间断有限元方法：
式中，ΔtDGM为间断有限元对应的时间步长；αSSM和αDGM分别对应交错网格方法和间断有限元方法的最大库朗数；
建立混合方法，对应的时间步长Δt必须满足：
Δt≤min{ΔtSSM，ΔtDGM} 式(19)
其中，Δt为所述基于间断有限元和交错网格的高精度地震波场模拟方法采用的时间步长。
11.实现权利要求1～10任一项所述基于间断有限元和交错网格的高精度地震波场模拟方法的系统，其特征是，包括：震源子波模块、观测系统模块、数据前处理模块、混合算法稳定性分析模块、人工边界反射压制模块、波场演化计算、地震波场结果生成与输出模块；
其中：
震源子波模块用于根据时间采样率、最大模拟时间、震源子波的类型及震源位置，获取地震子波数据；
观测系统模块，用于接收地震波场，建立地面接收地震数据记录，确定炮点和检波点的位置关系；
数据前处理模块，用于确定间断有限元方法和交错网格方法边界，间断有限元方法和交错网格方法之间参数变量的转化关系以及时间插值、B2.对间断有限元地质模型数据进行预处理、选择数值流设置和时间演化格式；
混合算法稳定性分析模块，用于对间断有限元区域和交错网格区域分别进行精度和稳定性分析，确定两种方法对应计算区域的精度，选择各自对应的阶数，使得数值模拟保持稳定性；
人工边界反射压制模块，用于通过完全匹配层技术使得间断有限元和交错网格方法的吸收边界兼容，在不同的计算区域应用完全匹配层边界条件；
波场演化计算模块，用于计算和更新波场，得到交错网格波场和间断有限元波场系数，并获得地质模型的波场快照；
地震波场结果生成与输出模块，用于生成某一时刻的地震波场切片，输出波场时间切片结果和地面地震记录结果。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202110831354.0A CN113536638B (zh) | 2021-07-22 | 2021-07-22 | 基于间断有限元和交错网格的高精度地震波场模拟方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202110831354.0A CN113536638B (zh) | 2021-07-22 | 2021-07-22 | 基于间断有限元和交错网格的高精度地震波场模拟方法 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN113536638A CN113536638A (zh) | 2021-10-22 |
| CN113536638B true CN113536638B (zh) | 2023-09-22 |
Family
ID=78088615
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202110831354.0A Active CN113536638B (zh) | 2021-07-22 | 2021-07-22 | 基于间断有限元和交错网格的高精度地震波场模拟方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN113536638B (zh) |
Families Citing this family (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN114545505A (zh) * | 2022-01-12 | 2022-05-27 | 吉林大学 | 一种基于质点网格思想的动态海洋波场信息采集方法 |
| CN116559942B (zh) * | 2023-04-07 | 2024-09-17 | 中国科学院地质与地球物理研究所 | 一种流-固耦合的波场模拟方法、装置、设备及存储介质 |
Citations (6)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US5136550A (en) * | 1992-02-25 | 1992-08-04 | Western Atlas International, Inc. | Method for estimating the residual source of receiver coordinates from CMP gathers |
| CN104459773A (zh) * | 2014-08-08 | 2015-03-25 | 中国石油天然气集团公司 | 基于交错网格Lowrank分解的无条件稳定地震波场延拓方法 |
| CN106842320A (zh) * | 2017-01-19 | 2017-06-13 | 北京大学 | Gpu并行三维地震波场生成方法和系统 |
| CN107561585A (zh) * | 2017-09-19 | 2018-01-09 | 北京大学 | 一种多核多节点并行三维地震波场生成方法和系统 |
| CN109725346A (zh) * | 2018-12-17 | 2019-05-07 | 中国石油天然气集团有限公司 | 一种时间-空间高阶vti矩形网格有限差分方法及装置 |
| CN110988993A (zh) * | 2019-11-27 | 2020-04-10 | 清华大学 | 一种偏移成像方法、装置及电子设备 |
-
2021
- 2021-07-22 CN CN202110831354.0A patent/CN113536638B/zh active Active
Patent Citations (6)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US5136550A (en) * | 1992-02-25 | 1992-08-04 | Western Atlas International, Inc. | Method for estimating the residual source of receiver coordinates from CMP gathers |
| CN104459773A (zh) * | 2014-08-08 | 2015-03-25 | 中国石油天然气集团公司 | 基于交错网格Lowrank分解的无条件稳定地震波场延拓方法 |
| CN106842320A (zh) * | 2017-01-19 | 2017-06-13 | 北京大学 | Gpu并行三维地震波场生成方法和系统 |
| CN107561585A (zh) * | 2017-09-19 | 2018-01-09 | 北京大学 | 一种多核多节点并行三维地震波场生成方法和系统 |
| CN109725346A (zh) * | 2018-12-17 | 2019-05-07 | 中国石油天然气集团有限公司 | 一种时间-空间高阶vti矩形网格有限差分方法及装置 |
| CN110988993A (zh) * | 2019-11-27 | 2020-04-10 | 清华大学 | 一种偏移成像方法、装置及电子设备 |
Non-Patent Citations (3)
| Title |
|---|
| 复杂介质中地震波模拟的无单元法;贾晓峰;胡天跃;王润秋;;石油地球物理勘探(第01期);全文 * |
| 弹性波变阶数旋转交错网格数值模拟(英文);王为中;胡天跃;吕雪梅;秦臻;李艳东;张研;应用地球物理(英文版)(第003期);全文 * |
| 王为中 ; 胡天跃.基于GPU的三维地震波变阶数交错网格数值模拟.中国石油学会2017年物探技术研讨会.2017,全文. * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN113536638A (zh) | 2021-10-22 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| KR102020759B1 (ko) | Q-보상된 전 파동장 반전 | |
| US6125330A (en) | Method of determining the response caused by model alterations in seismic simulations | |
| CN113536638B (zh) | 基于间断有限元和交错网格的高精度地震波场模拟方法 | |
| US11493658B2 (en) | Computer-implemented method and system employing nonlinear direct prestack seismic inversion for poisson impedance | |
| US10345466B2 (en) | Memory efficient Q-RTM computer method and apparatus for imaging seismic data | |
| CN111025387B (zh) | 一种页岩储层的叠前地震多参数反演方法 | |
| CN108845351A (zh) | 一种vsp地震资料转换波全波形反演方法 | |
| CN111766628A (zh) | 一种预条件的时间域弹性介质多参数全波形反演方法 | |
| Qu et al. | Fluid-solid coupled full-waveform inversion in the curvilinear coordinates for ocean-bottom cable data | |
| MX2011003850A (es) | Estimado de señal de dominio de imagen a interferencia. | |
| CN101369024A (zh) | 一种地震波动方程生成方法及系统 | |
| Golubev et al. | Simulation of seismic wave propagation in a multicomponent oil deposit model | |
| KR20190075176A (ko) | 다중반사파 없는 데이터 세트를 생성하는 다단식 전 파동장 역산 프로세스 | |
| CN113031068A (zh) | 一种基于反射系数精确式的基追踪叠前地震反演方法 | |
| Nunes et al. | Fast geostatistical seismic inversion coupling machine learning and Fourier decomposition | |
| CN114895351A (zh) | 模拟地震波在任意不连续界面处传播的介质模型化的方法及装置 | |
| NO20190489A1 (en) | Seismic modeling | |
| US9928315B2 (en) | Re-ordered interpolation and convolution for faster staggered-grid processing | |
| CN109239776A (zh) | 一种地震波传播正演模拟方法和装置 | |
| CN105182414A (zh) | 一种基于波动方程正演去除直达波的方法 | |
| CN115236730B (zh) | 一种层间多次波傅里叶有限差分的地震波场偏移成像方法 | |
| US8437219B2 (en) | Correcting an acoustic simulation for elastic effects | |
| CN117388944A (zh) | 地质模型约束的多物性参数反演方法 | |
| CN115373022A (zh) | 一种基于振幅相位校正的弹性波场Helmholtz分解方法 | |
| Roberts et al. | Investigation into the use of 2D elastic waveform inversion from look‐ahead walk‐away VSP surveys |
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 |