Reducing computation cost by Lax-Wendroff methods with fourth-order temporal accuracy

Geophysics ◽  
2019 ◽  
Vol 84 (3) ◽  
pp. T109-T119 ◽  
Author(s):  
Hongwei Liu ◽  
Houzhu Zhang

Explicit time-marching finite-difference stencils have been extensively used for simulating seismic wave propagation, and they are the most computationally intensive part of seismic forward modeling, reverse time migration, and full-waveform inversion. The time-marching step, determined by both the stability condition and numerical dispersion, is a key factor in the computational cost. In contrast with the widely used second-order temporal stencil, the Lax-Wendroff stencil is more cost effective because the time-marching step can be much larger. It can be proved, using theory and numerical tests, that the Lax-Wendroff stencil does enable larger time steps. In terms of numerical dispersion, the time steps for second-order and Lax-Wendroff stencils are functions of the number of shortest wavelengths away from the source location. These functions are derived by evaluating the relative L2-norm of the differences between the analytical and numerical solutions. The method for determining the time-marching step is model adaptive and easy to implement. We use the pseudo spectral method for the computation of spatial derivatives, and the wave equations that we solved are for isotropic media only, but the described principles can be easily implemented for more complicated types of media.

Geophysics ◽  
2011 ◽  
Vol 76 (1) ◽  
pp. A1-A6 ◽  
Author(s):  
Robin P. Fletcher ◽  
Johan O. Robertsson

We propose two new boundary conditions to regulate coherent reflections from the model boundaries in numerical solutions of wave equations. Both boundary conditions have the common feature that the boundary condition is varied with respect to time. The first boundary condition expands or contracts the computational model during a modeling simulation. The effect is to cause a Doppler shift in the reflected wavefield that can be used to shift energy outside a frequency band of interest. In addition, when the computational domain is expanding, the range of possible incidence angles on the boundary is restricted. This can be used to increase the effectiveness of many existing absorbing boundary conditions that are more effective for incidence angles close to normal. The second boundary condition is an extension of random boundaries. By carefully changing the realization of a random boundary over time, a more diffusive wavefield can be simulated. We show results with 2D numerical simulations of the scalar wave equation for both these boundary conditions. The first boundary condition has application to modeling, but both these boundary conditions have potential application within algorithms that rely upon modeling kernels, such as reverse-time migration and full-waveform inversion.


Geophysics ◽  
2020 ◽  
Vol 85 (6) ◽  
pp. T359-T373
Author(s):  
Jeffrey Shragge ◽  
Tugrul Konuk

Numerical solutions of 3D isotropic elastodynamics form the key computational kernel for many isotropic elastic reverse time migration and full-waveform inversion applications. However, real-life scenarios often require computing solutions for computational domains characterized by non-Cartesian geometry (e.g., free-surface topography). One solution strategy is to compute the elastodynamic response on vertically deformed meshes designed to incorporate irregular topology. Using a tensorial formulation, we have developed and validated a novel system of semianalytic equations governing 3D elastodynamics in a stress-velocity formulation for a family of vertically deformed meshes defined by Bézier interpolation functions between two (or more) nonintersecting surfaces. The analytic coordinate definition also leads to a corresponding analytic free-surface boundary condition (FSBC) as well as expressions for wavefield injection and extraction. Theoretical examples illustrate the utility of the tensorial approach in generating analytic equations of 3D elastodynamics and the corresponding FSBCs for scenarios involving free-surface topography. Numerical examples developed using a fully staggered grid with a mimetic finite-difference formulation demonstrate the ability to model the expected full-wavefield behavior, including complex free-surface interactions.


Geophysics ◽  
2018 ◽  
Vol 83 (5) ◽  
pp. S399-S408 ◽  
Author(s):  
Yunyue Elita Li ◽  
Yue Du ◽  
Jizhong Yang ◽  
Arthur Cheng ◽  
Xinding Fang

Elastic wave imaging has been a significant challenge in the exploration industry due to the complexities in wave physics and numerical implementation. We have separated the governing equations for P- and S-wave propagation without the assumptions of homogeneous Lamé parameters to capture the mode conversion between the two body waves in an isotropic, constant-density medium. The resulting set of two coupled second-order equations for P- and S-potentials clearly demonstrates that mode conversion only occurs at the discontinuities of the shear modulus. Applying the Born approximation to the new equations, we derive the PP, PS, SP, and SS imaging conditions from the first gradients of waveform matching objective functions. The resulting images are consistent with the physical perturbations of the elastic parameters, and, hence, they are automatically free of the polarity reversal artifacts in the converted images. When implementing elastic reverse time migration (RTM), we find that scalar wave equations can be used to back propagate the recorded P-potential, as well as individual components in the vector field of the S-potential. Compared with conventional elastic RTM, the proposed elastic RTM implementation using acoustic propagators not only simplifies the imaging condition, it but also reduces the computational cost and the artifacts in the images. We have determined the accuracy of our method using 2D and 3D numerical examples.


Geophysics ◽  
2019 ◽  
Vol 84 (6) ◽  
pp. R827-R844 ◽  
Author(s):  
Zongcai Feng ◽  
Gerard Schuster

We present a quasi-elastic wave equation as a function of the pressure variable, which can accurately model PP reflections with elastic amplitude variation with offset effects under the first-order Born approximation. The kinematic part of the quasi-elastic wave equation accurately models the propagation of P waves, whereas the virtual-source part, which models the amplitudes of reflections, is a function of the perturbations of density and Lamé parameters [Formula: see text] and [Formula: see text]. The quasi-elastic wave equation generates a scattering radiation pattern that is exactly the same as that for the elastic wave equation, and only requires the solution of two acoustic wave equations for each shot gather. This means that the quasi-elastic wave equation can be used for true-amplitude linearized waveform inversion (also known as least-squares reverse time migration) of elastic PP reflections, in which the corresponding misfit gradients are with respect to the perturbations of density and the P- and S-wave impedances. The perturbations of elastic parameters are iteratively updated by minimizing the [Formula: see text]-norm of the difference between the recorded PP reflections and the predicted pressure data modeled from the quasi-elastic wave equation. Numerical tests on synthetic and field data indicate that true-amplitude linearized waveform inversion using the quasi-elastic wave equation can account for the elastic PP amplitudes and provide a robust estimate of the perturbations of P- and S-wave impedances and, in some cases, the density. In addition, true-amplitude linearized waveform inversion provides images with a wider bandwidth and fewer artifacts because the PP amplitudes are accurately explained. We also determine the 2D scalar quasi-elastic wave equation for P-SV reflections and the 3D vector equation for PS reflections.


Geophysics ◽  
2009 ◽  
Vol 74 (6) ◽  
pp. WCA153-WCA158 ◽  
Author(s):  
Faqi Liu ◽  
Guanquan Zhang ◽  
Scott A. Morton ◽  
Jacques P. Leveille

The acoustic wave equation has been widely used for the modeling and reverse time migration of seismic data. Numerical implementation of this equation via finite-difference techniques has established itself as a valuable approach and has long been a favored choice in the industry. To ensure quality results, accurate approximations are required for spatial and time derivatives. Traditionally, they are achieved numerically by using either relatively very fine computation grids or very long finite-difference operators. Otherwise, the numerical error, known as numerical dispersion, is present in the data and contaminates the signals. However, either approach will result in a considerable increase in the computational cost. A simple and computationally low-cost modification to the standard acoustic wave equation is presented to suppress numerical dispersion. This dispersion attenuator is one analogy of the antialiasing operator widely applied in Kirchhoff migration. When the new wave equation is solved numerically using finite-difference schemes, numerical dispersion in the original wave equation is attenuated significantly, leading to a much more accurate finite-difference scheme with little additional computational cost. Numerical tests on both synthetic and field data sets in both two and three dimensions demonstrate that the optimized wave equation dramatically improves the image quality by successfully attenuating dispersive noise. The adaptive application of this new wave equation only increases the computational cost slightly.


Geophysics ◽  
2017 ◽  
Vol 82 (5) ◽  
pp. T225-T235 ◽  
Author(s):  
Leandro Di Bartolo ◽  
Leandro Lopes ◽  
Luis Juracy Rangel Lemos

Pseudoacoustic algorithms are very fast in comparison with full elastic ones for vertical transversely isotropic (VTI) modeling, so they are suitable for many applications, especially reverse time migration. Finite differences using simple grids are commonly used to solve pseudoacoustic equations. We have developed and implemented general high-order 3D pseudoacoustic transversely isotropic formulations. The focus is the development of staggered-grid finite-difference algorithms, known for their superior numerical properties. The staggered-grid schemes based on first-order velocity-stress wave equations are developed in detail as well as schemes based on direct application to second-order stress equations. This last case uses the recently presented equivalent staggered-grid theory, resulting in a staggered-grid scheme that overcomes the problem of large memory requirement. Two examples are presented: a 3D simulation and a prestack reverse time migration application, and we perform a numerical analysis regarding computational cost and precision. The errors of the new schemes are smaller than the existing nonstaggered-grid schemes. In comparison with existing staggered-grid schemes, they require 25% less memory and only have slightly greater computational cost.


Geophysics ◽  
2014 ◽  
Vol 79 (3) ◽  
pp. WA3-WA11 ◽  
Author(s):  
Wei Xie ◽  
Dinghui Yang ◽  
Faqi Liu ◽  
Jingshuang Li

With the capability of handling complicated velocity models, reverse-time migration (RTM) has become a powerful imaging method. Improving imaging accuracy and computational efficiency are two significant but challenging tasks in the applications of RTM. Despite being the most popular numerical technique applied in RTM, finite-difference (FD) methods often suffer from undesirable numerical dispersion, leading to a noticeable loss of imaging resolution. A new and effective FD operator, called the high-order stereo operator, has been developed to approximate the partial differential operators in the wave equation, from which a numerical scheme called the three-step stereo method (TSM) has been developed and has shown effectiveness in suppressing numerical dispersion. Numerical results show that compared with the conventional numerical methods such as the Lax-Wendroff correction (LWC) scheme and the staggered-grid (SG) FD method, this new method significantly reduces numerical dispersion and computational cost. Tests on the impulse response and the 2D prestack Hess acoustic VTI model demonstrated that the TSM achieves higher image quality than the LWC and SG methods do, especially when coarse computation grids were used, which indicated that the new method can be a promising algorithm for large-scale anisotropic RTM.


2018 ◽  
Vol 26 (02) ◽  
pp. 1850005 ◽  
Author(s):  
Jian Wang ◽  
Xiaohong Meng ◽  
Hong Liu ◽  
Wanqiu Zheng ◽  
Zhiwei Liu

Staggered-grid finite-difference forward modeling in the time domain has been widely used in reverse time migration and full waveform inversion because of its low memory cost and ease to implementation on GPU, however, high dominant frequency of wavelet and big grid interval could result in significant numerical dispersion. To suppress numerical dispersion, in this paper, we first derive a new weighted binomial window function (WBWF) for staggered-grid finite-difference, and two new parameters are included in this new window function. Then we analyze different characteristics of the main and side lobes of the amplitude response under different parameters and accuracy of the numerical solution between the WBWF method and some other optimum methods which denotes our new method can drive a better finite difference operator. Finally, we perform elastic wave numerical forward modeling which denotes that our method is more efficient than other optimum methods without extra computing costs.


Geophysics ◽  
2011 ◽  
Vol 76 (2) ◽  
pp. S65-S75 ◽  
Author(s):  
Eric Duveneck ◽  
Peter M. Bakker

We present an approach for P-wave modeling in inhomogeneous transversely isotropic media with tilted symmetry axis (TTI media), suitable for anisotropic reverse-time migration. The proposed approach is based on wave equations derived from first principles — the equations of motion and Hooke’s law — under the acoustic TI approximation. Consequently, no assumptions are made about the spatial variation of medium parameters. A rotation of the stress and strain tensors to a local coordinate system, aligned with the TI-symmetry axis, makes it possible to benefit from the simple and sparse form of the TI-elastic tensor in that system. The resulting wave equations can be formulated either as a set of five first-order or as a set of two second-order partial differential equations. For the constant-density case, the second-order TTI wave equations involve mixed and nonmixed second-order spatial derivatives with respect to global, nonrotated coordinates. We propose a numerical implementation of these equations using high-order centered finite differences. To minimize modeling artifacts related to the use of centered first-derivative operators, we use discrete second-derivative operators for the nonmixed second-order spatial derivatives and repeated discrete first-derivative operators for the mixed derivatives. Such a combination of finite-difference operators leads to a stable wave propagator, provided that the operators are designed properly. In practice, stability is achieved by slightly weighting down terms that contain mixed derivatives. This has a minor, practically negligible, effect on the kinematics of wave propagation. The stability of the presented scheme in inhomogeneous TTI models with rapidly varying anisotropic symmetry axis direction is demonstrated with numerical examples.


Sign in / Sign up

Export Citation Format

Share Document