One-way Wave Equation Research Articles

A new variable-order fractional momentum operator for wave absorption when solving Schrödinger equations

Absorbing boundary conditions are often needed for solving Schrödinger equations when the analytical solution is unknown. This article presents a new wave-absorption method by unifying the Schrödinger equations with associated one-way wave equations via a fractional-order momentum operator pˆα+1 where 0≤α≤1. By gradually varying the fractional order α, the computational boundary is made to be transparent. This method allows for incident waves to propagate freely and forces reflected waves to decay exponentially with time. The unified equation is then solved using the Finite-Difference Time-Domain (FDTD) method so that the scheme is iteratively explicit and can be parallelized easily. Finally, the obtained fractional FDTD scheme is tested to simulate soliton and particle propagations and is compared with the existing Perfectly Matched Layer (PML) method. Results show that the new method is promising.

Truncated pseudo-differential operator √−▽2 and its applications in viscoacoustic reverse-time migration

SUMMARY The pseudo-differential operator with symbol |k|α has been widely used in seismic modelling and imaging when involving attenuation, anisotropy and one-way wave equation, which is usually calculated using the pseudo-spectral method. For large-scale problems, applying high-dimensional Fourier transforms to solve the wave equation that includes pseudo-differential operators is much more expensive than finite-difference approaches, and it is not suitable for parallel computing with domain decomposition. To mitigate this difficulty, we present a truncated space-domain convolution method to efficiently compute the pseudo-differential operator $\sqrt{-\nabla ^2}$, and then apply it to viscoacoustic reverse-time migration. Although $\sqrt{-\nabla ^2}$ is theoretically non-local in the space domain, we take the limited frequency band of seismic data into account, and constrain the approximated convolution stencil to a finite length. The convolution coefficients are computed by solving a least-squares inverse problem in the wavenumber domain. In addition, we exploit the symmetry of the resulting convolution stencil and develop a fast spatial convolution algorithm. The applications of the proposed method in Q-compensated reverse-time migration demonstrate that it is a good alternative to the pseudo-spectral method for computing the pseudo-differential operator $\sqrt{-\nabla ^2}$, with almost the same accuracy but much higher efficiency.

Seismic imaging of medical ultrasound data: Towards in vivo applications

This study presents the application to medical ultrasonography of a high value-added imaging method developed in the context of geophysical exploration and based on the undulatory description of the physical process. Currently, the most serious limitation of conventional medical ultrasound systems is the implicit assumption of constant-velocity media stemming from the widespread use of geometrical ray-based imaging algorithms and leading to phase aberration phenomena in the case typical of the human body of two or more tissues with different velocities. As a result, imaging of targets underlying tissue layers can be severely degraded. To address this problem, the proposed algorithm is designed expressly for variable-velocity media, providing high-resolution images by implementing the solution of the one-way wave equation. From a macro-velocity model, the illuminated structures are imaged by correlating, in the frequency domain and for each source element, the simulated emitted wavefield with the back-propagated echoes sensed by the ultrasound transducer. The pointwise sum of the correlations reveals, constructively, all diffraction points and reflectors in the image space. The experiments conducted on in vitro laboratory data reveal remarkable spatial resolution and very accurate target positioning, even in the presence of an aberrant layer. A preliminary in vivo test also demonstrates the method's effectiveness in imaging echoes recorded on human tissues, giving a glimpse of its potentiality for clinical purposes.

Improved least-squares migration through double-sweeping solver

Based on a recently developed approximate wave-equation solver, we have developed a methodology to reduce the computational cost of seismic migration in the frequency domain. This approach divides the domain of interest into smaller subdomains, and the wavefield is computed using a sequential process to determine the downward- and upward-propagating wavefields — hence called a double-sweeping solver. A sequential process becomes possible using a special approximation of the interface conditions between subdomains. This method is incorporated into the least-squares migration framework as an approximate solver. The associated computational effort is comparable to one-way wave-equation approaches, yet, as illustrated by the numerical examples, the accuracy and convergence behavior are comparable to that of the full-wave equation.

On the accuracy of one-way approximate models for nonlinear waves in soft solids.

Simple strain-rate viscoelasticity models of isotropic soft solid are introduced. The constitutive equations account for finite strain, incompressibility, material frame-indifference, nonlinear elasticity, and viscous dissipation. A nonlinear viscous wave equation for the shear strain is obtained exactly and corresponding one-way Burgers-type equations are derived by making standard approximations. Analysis of the travelling wave solutions shows that these partial differential equations produce distinct solutions, and deviations are exacerbated when wave amplitudes are not arbitrarily small. In the elastic limit, the one-way approximate wave equation can be linked to simple wave theory and shock wave theory, thus, allowing direct error measurements.

One-Way Wave Operator

The second-order partial differential wave Equation (Cauchy’s first equation of motion), derived from Newton’s force equilibrium, describes a standing wave field consisting of two waves propagating in opposite directions, and is, therefore, a “two-way wave equation”. Due to the second order differentials analytical solutions only exist in a few cases. The “binomial factorization” of the linear second-order two-way wave operator into two first-order one-way wave operators has been known for decades and used in geophysics. When the binomial factorization approach is applied to the spatial second-order wave operator, this results in complex mathematical terms containing the so-called “Dirac operator” for which only particular solutions exist. In 2014, a hypothetical “impulse flow equilibrium” led to a spatial first-order “one-way wave equation” which, due to its first order differentials, can be more easily solved than the spatial two-way wave equation. To date the conversion of the spatial two-way wave operator into spatial one-way wave operators is unsolved. By considering the one-way wave operator containing a vector wave velocity, a “synthesis” approach leads to a “general vector two-way wave operator” and the “general one-way/two-way equivalence”. For a constant vector wave velocity the equivalence with the d’Alembert operator can be achieved. The findings are transferred to commonly used mechanical and electromagnetic wave types. The one-way wave theory and the spatial one-way wave operators offer new opportunities in science and engineering for advanced wave and wave field calculations.

One-Way Vibration Absorber

A vibration absorber consisting of a one-dimensional waveguide with a reflectionless termination extracts vibrational energy from a structure that is to be damped. An optimum energy dissipation occurs for the so-called power adjustment, i.e, the same level of resistance and the opposite reactance of structure and absorber. The dimensioning of these impedance parameters on the base of the classic second order “two-way” wave equation provides analytical solutions for a few simple waveguide shapes; solutions for all other waveguides are only accessible via numerical finite-element computation. However, the competing first order “one-way” wave equation allows for an analytical conception of both the known broadband vibration absorber and the “Acoustic Black Hole” absorber. For example, for an exponential waveguide, the two-way calculation shows no resistance (and hence no real wave propagation) below a cut-off frequency, while the one-way wave equation predicts absorption in the whole frequency range.

Finite Difference Pre-Stack Depth Migration of One-Way Wave Equation in Isotropic Visco-Acoustic Media of a Ray-Centered Coordinate System

In seismic exploration, obtaining accurate bandwidth of the reflected signal is essential for determining seismic resolution. A portion of the stored energy is attenuated during signal propagation, narrowing the received seismic signal bandwidth. Therefore, compensating the energy attenuation is important for improving the seismic resolution. The current method of compensating absorption and attenuation based on a single channel can only compensate the post-stack data (self-exciting and self-receiving), whereas in practice, seismic waves do not propagate along the self-exciting and self-receiving seismic wave path; the propagation path is complex. The absorption and attenuation depend on the propagation path. The primary methods used for Q-compensation along the propagation paths are one-way wave extrapolation in the Cartesian coordinate system and Gaussian beam Q-compensation migration in the ray-centered coordinate system. However, the large angle limits the one-way wave method, and the Gaussian beam method refers to the high-frequency approximation solution of the two-way wave equation. Therefore, a 15-degree equation in the ray-centered coordinate system is proposed. Seismic waves extrapolate along the ray, which compensates the absorption and attenuation along the real propagation path. The 15-degree equation in the ray-centered coordinate system does not perform high-frequency approximation in the ray beam and has no large angle limit, facilitating the accurate description of local wavefields in the ray beam.

Theoretical and quantitative evaluation of hybrid PML-ABCs for seismic wave simulation

A good artificial boundary treatment in a seismic wave grid-based numerical simulation can reduce the size of the computational region and increase the computational efficiency, which is becoming increasingly important for seismic migration and waveform inversion tasks requiring hundreds or thousands of simulations. Two artificial boundary techniques are commonly used: perfectly matched layers (PMLs), which exhibit the excellent absorption performance but impose a greater computational burden by using finite layers to gradually reduce wave amplitudes; and absorbing boundary conditions (ABCs), which have the high computational efficiency but are less effective in absorption because they employ the one-way wave equation at the exterior boundary. Naturally, PMLs have been combined with ABCs to reduce the number of PMLs, thus improving the computational efficiency; many studies have proposed such hybrid PMLs. Depending on the equations from which the ABCs are derived, there are two hybrid PML variants: the PML+unstretched ABC (UABC), in which the ABC is derived from a physical equation; or the PML+stretched ABC (SABC), in which the ABC is derived from the PML equation. Even though all the previous studies concluded that hybrid PMLs can improve the absorption performance, none of them quantified how many PMLs can be removed by combining the PML with the ABC compared with the pure PML. In this paper, we systematically study the absorption performance of the two hybrid PML variants. We develop a method to distinguish the artificial reflections from the PML-interior interface and those caused by the PML exterior boundary to accurately approximate the additional absorption achieved by using the UABC and the SABC. The reflection coefficients based on a theoretical derivation and numerical tests both show that the UABC amplifies most reflections and is not recommended in any situation; conversely, the SABC can always diminish reflections, but the additional absorption achieved by the SABC is relatively poor and cannot effectively reduce the number of PMLs. In contrast, we find that simply increasing the damping parameter improves absorption better than the PML+SABC. Our results show that the improvement in absorption achieved by combining the PML with either the SABC or the UABC is not better than that obtained by simply adjusting the damping profile of the PML; thus, combining the PML with the ABC is not recommended in practice.

Direct expansion of Fourier extrapolator for one-way wave equation using Chebyshev polynomials of the second kind

The Fourier method for one-way wave propagation is efficient but potentially inaccurate in complex media. The implicit finite-difference (FD) method can handle arbitrarily complex media, but it can be inefficient in 3D and has limited dip bandwidth. We have adopted a new Fourier method based on the Chebyshev expansion of the second kind. Theoretical analyses and numerical experiments indicate that our method is comprehensively superior to a similar method based on the Chebyshev expansion of the first kind in terms of balanced amplitude and error tolerance. Within the dip bandwidth from 0° to 65°, the fourth-order form of our method has an error tolerance of 2%, which is approximately one-third that of the Chebyshev expansion of the first kind. Our method also is superior to the implicit FD method in several important aspects: effective bandwidth, computational efficiency, numerical dispersion, and two-way splitting error. It can be easily extended from 2D to 3D compared with the FD method and from low orders to high orders compared with the optimized Chebyshev-Fourier method. Our method finds better imaging results of the SEG/EAGE model by providing a well-focused salt dome, flank, and bottom as well as the detailed structures beneath the salt body, compared with the implicit FD method and the Chebyshev expansion of the first kind; meanwhile, our method has fewer imaging artifacts because it can better position the reflectors.

A three-dimensional Laguerre one-way wave equation solver

A finite difference algorithm based on the integral Laguerre transform in time for solving a three-dimensional one-way wave equation is proposed. This allows achieving high accuracy of calculation results. In contrast to the Fourier method, the approach does not need to solve systems of linear algebraic equations with indefinite matrices. To filter the unstable components of a wave field, Richardson extrapolation or spline approximation can be used. However, these methods impose additional limitations on the integration step in depth. This problem can be solved if the filtering is performed not in the direction of extrapolation of the wave field, but in a horizontal plane. This approach called for fast methods of converting the Laguerre series coefficients into the Fourier series coefficients and vice versa. The high stability of the new algorithm allows calculations with a large depth step without loss of accuracy and, in combination with Marchuk-Strang splitting, this can significantly reduce the calculation time. Computational experiments are performed. The results have shown that this algorithm is highly accurate and efficient in solving the problems of seismic migration.

Factorized One-Way Wave Equations

The method used to factorize the longitudinal wave equation has been known for many decades. Using this knowledge, the classical 2nd-order partial differential Equation (PDE) established by Cauchy has been split into two 1st-order PDEs, in alignment with D’Alemberts’s theory, to create forward- and backward-traveling wave results. Therefore, the Cauchy equation has to be regarded as a two-way wave equation, whose inherent directional ambiguity leads to irregular phantom effects in the numerical finite element (FE) and finite difference (FD) calculations. For seismic applications, a huge number of methods have been developed to reduce these disturbances, but none of these attempts have prevailed to date. However, a priori factorization of the longitudinal wave equation for inhomogeneous media eliminates the above-mentioned ambiguity, and the resulting one-way equations provide the definition of the wave propagation direction by the geometric position of the transmitter and receiver.

Extending subsurface imaging beyond ocean-bottom node coverage using multiples: A case study in deepwater Nigeria

We have evaluated the results of a receiver decimation study in a deepwater context using separated wavefield imaging (SWIM) algorithms to provide extended illumination for imaging without ocean-bottom node (OBN) positioning constraints. We carried out subsurface imaging using the SWIM imaging technique with a reduced OBN layout, and we compared the results with those from conventional one-way wave-equation migration. We found from the results that the SWIM algorithm makes it possible to reduce the OBN layout while obtaining a similar subsurface image with the same shot geometry, which allows a reduced receiver acquisition effort, offers more geometry flexibility without affecting the image quality, with a potentially significant reduction of acquisition cost and 4D processing turnaround time.

Introduction of the Hessian in joint migration inversion and improved recovery of structural information using image-based regularization

Joint migration inversion (JMI) is a method based on one-way wave equations that aims at fitting seismic reflection data to estimate an image and a background velocity. The depth-migrated image describes the high spatial-frequency content of the subsurface and, in principle, is true amplitude. The background velocity model accounts mainly for the large spatial-scale kinematic effects of the wave propagation. Looking for a deeper understanding of the method, we briefly review the continuous equations that compose the forward-modeling engine of JMI for acoustic media and angle-independent scattering. Then, we use these equations together with the first-order adjoint-state method to arrive at a new formulation of the model gradients. To estimate the image, we combine the second-order adjoint-state method with the truncated-Newton method to obtain the image updates. For the model (velocity) estimation, in comparison to the image update, we reduce the computational cost by adopting a diagonal preconditioner for the corresponding gradient in combination with an image-based regularizing function. Based on this formulation, we build our implementation of the JMI algorithm. Our image-based regularization of the model estimate allows us to carry over structural information from the estimated image to the jointly estimated background model. As demonstrated by our numerical experiments, this procedure can help to improve the resolution of the estimated model and make it more consistent with the image.

Full wavefield least-squares reverse time migration

Waveform inversion based on least-squares reverse time migration (LSRTM) usually involves Born modeling, which models the primary-only data. As a result, the inversion process handles only primaries and corresponding multiple elimination preprocessing of the input data is required prior to imaging and inversion. Otherwise, multiples left in the input data are mapped as false reflectors, also known as crosstalk, in the final image. At the same time, the developed full-wavefield migration (FWM) methodology can handle internal multiples in an inversion-based imaging process. However, because it is based on the framework of the one-way wave equation, it cannot image dips close to and beyond 90°. Therefore, we have aimed at upgrading the LSRTM framework by bringing in the functionality of FWM to handle internal multiples. We inject the secondary source term, as used in the original formulation of FWM, to define a wavefield relationship that allows modeling multiple scattering via reflectivity. The secondary source term is based on the estimated reflectivity and can be injected into the pressure component when simulating the two-way wave equation using finite-difference modeling. We use this modified forward model for estimating the reflectivity model in an FWM-type manner and validate the method on synthetic and field data containing visible internal multiples.

3D adaptive internal multiples modeling based on globally optimised Fourier finite-difference method

Abstract Numerical methods have been widely applied to simulate seismic wave propagation. However, few studies have focused on internal multiples modeling. The formation mechanism and response of internal multiples are still unclear. Therefore, we develop a weighted-optimised-based internal multiples simulation method under 3D conditions. Using a one-way wave equation and full-wavefield method, the different-order internal multiples are computed numerically in a recursive manner. The traditional Fourier finite-difference (FFD) method has low numerical accuracy in a horizontal direction. A globally optimised FFD (OFFD) method is used to improve the lateral propagation accuracy of the seismic waves. Meanwhile, we adopt an adaptive variable-step technique to improve computational efficiency. The 3D internal multiples modeling technique is capable of calculating the different-order multiple reflections in complex structures. We use the present method to simulate internal multiples in several models. Theoretical analyses are consistent with the numerical results. Numerical examples demonstrate that the 3D internal multiples modeling technique has superior performance when adapting to lateral velocity changes and steep dip. This also implies that our method is fit for the simulation of internal multiples propagation in a 3D complex medium and can assist in identifying the internal multiples from full-wavefield data.

A multi-level parallel algorithm for seismic imaging based on one-way wave equation migration

This article presents a parallel algorithm for seismic imaging based on depth wavefield extrapolation (the solution of a one-way wave equation [OWE]) employing the pseudospectral method. The algorithm is essentially oriented on common-offset vector (COV) gathers seismic migration based on an OWE, includes several parallelism levels, and utilizes message passing interface (MPI), Nested OpenMP, and CUDA technologies. The uppermost level of parallelism involves data splitting to ensure that each dataset can be processed independently to construct parts of the COV images. The algorithm is embarrassingly parallel at this level. Next, each independent run processes several datasets to compute COV images. Each MPI process computes all COV images for a single dataset, then MPI processes exchange data to build up a single COV image for all datasets at a single node. Computation of the COV image at a node requires OpenMP parallelization so that one thread governs the GPU-based calculations and facilitates the construction of images within a thick slab. All other threads perform image interpolation within the slab. Finally, CUDA technology is used for the most computationally intense part of the algorithm—wavefield extrapolation. Such a complex algorithm structure makes it possible to process all COV images for full-azimuth seismic data employing OWE-based migration.

Amplitude-Preserved Wave Equation: An Example to Image the Gas Hydrate System

Natural gas hydrate is an important energy source. Therefore, it is extremely important to provide a clear imaging profile to determine its distribution for energy exploration. In view of the problems existing in conventional migration methods, e.g., the limited imaging angles, we proposed to utilize an amplitude-preserved one-way wave equation migration based on matrix decomposition to deal with primary and multiple waves. With respect to seismic data gathered at the Chilean continental margin, a conventional processing flow to obtain seismic records with a high signal-to-noise ratio is introduced. Then, the imaging results of the conventional and amplitude-preserved one-way wave equation migration methods based on primary waves are compared, to demonstrate the necessity of implementing amplitude-preserving migration. Moreover, a simple two-layer model is imaged by using primary and multiple waves, which proves the superiority of multiple waves in imaging compared with primary waves and lays the foundation for further application. For the real data, the imaging sections of primary and multiple waves are compared. We found that multiple waves are able to provide a wider imaging illumination while primary waves fail to illuminate, especially for the imaging of bottom simulating reflections (BSRs), because multiple waves have a longer travelling path and carry more information. By imaging the actual seismic data, we can make a conclusion that the imaging result generated by multiple waves can be viewed as a supplementary for the imaging result of primary waves, and it has some guiding values for further hydrate and in general shallow gas exploration.

Spherical One-Way Wave Equation

The coordinate-free one-way wave equation is transferred in spherical coordinates. Therefore it is necessary to achieve consistency between gradient, divergence and Laplace operators and to establish, beside the conventional radial Nabla operator ∂Φ/∂r, a new variant ∂rΦ/r∂r. The two Nabla operator variants differ in the near field term Φ/r whereas in the far field r≫0 there is asymptotic approximation. Surprisingly, the more complicated gradient ∂rΦ/r∂r results in unexpected simplifications for – and only for – spherical waves with the 1/r amplitude decrease. Thus the calculation always remains elementary without the wattless imaginary near fields, and the spherical Bessel functions are obsolete.

Mesh-free radial-basis-function-generated finite differences and their application in reverse time migration

Compared with one-way wave equation migration and ray-based migration, reverse time migration (RTM) using two-way wave propagation information can produce accurate imaging results for complex structures. Its computational accuracy and efficiency are mainly determined by numerical methods for wavefield simulation. When using traditional regular grids for seismic modeling, scattering artifacts may occur due to the stepped approximation of layer interfaces and rugged topography. On the other hand, irregular grids require a complex grid generation algorithm, despite having certain geometric flexibility. Mesh-free RTM can effectively reduce the scattered noise under regular grids and avoid the extra computation in the process of irregular grid generation. For implementation of the mesh-free RTM method, an algorithm with fast generation of node distributions is used to discretize the underground velocity model, and radial-basis-function-generated finite difference is used to realize the numerical simulation of wave propagation; the crosscorrelation imaging condition is adopted for imaging. The mesh-free RTM method, which has the flexibility of the simulation region and an abundance of wavefield information, reduces the storage required for RTM, shows the potential of high-accuracy migration in the case of an undulating surface, and provides more accurate migration imaging results for oil and gas exploration under complex geologic conditions.