Viscoacoustic Wave Equation Research Articles

Viscoacoustic reverse time migration with stable and effective two-way attenuation compensation

Intrinsic attenuation in seismic wave propagation leads to amplitude dissipation and phase dispersion in seismic data. The attenuation (∼1/ Q) effect is a common cause for inaccurate image locations and dimmed image energy in reverse time migration (RTM). Conventional viscoacoustic or viscoelastic RTM methods implement attenuation compensation by reversing the sign of the amplitude loss term and keeping the dispersion term unchanged for forward and backward wavefield extrapolation, which requires the decoupled viscoacoustic wave equation and gives rise to the numerical instability issue. To address these problems, we develop a stable and effective two-way attenuation-compensated viscoacoustic RTM method. Our method uses a new two-way normalized crosscorrelation imaging condition to compensate the attenuation effect. In the new imaging condition, only the attenuated forward wavefield of sources and the viscoacoustic Green’s functions of receivers are required, which eliminates the instability source arising from the operation of reversing the sign of the dissipation term in the viscoacoustic wave equation. Because no modifications are made to the viscoacoustic wave equation, the new imaging condition can be calculated without any additional computational cost and can be applied to decoupled and coupled viscoacoustic wave equation-based RTMs, which is more flexible than other attenuation compensation strategies. Two synthetic tests and one field data example are presented to demonstrate the stability and effectiveness of our method.

Discontinuous Galerkin method for solving viscoacoustic wave equations with amplitude dissipation and phase dispersion separation in isotropic and anisotropic media

SUMMARY The standard-linear-solid (SLS) theory works well for viscoelasticity. However, the coupling of amplitude dissipation and phase dispersion makes it impossible to investigate their effects separately by the discontinuous Galerkin method (DGM). In this paper, we have derived newly viscoacoustic wave equations with amplitude dissipation and phase dispersion separation in isotropic and anisotropic media, based on a Fourier method, which is suitable for using a time–space-domain DGM on unstructured meshes. The basic framework of DGM is constructed and the amplitude-dissipation effect and the phase-dispersion effect in viscoacoustic wave equations are investigated. The original equation is first transformed into the frequency–wavenumber domain, where the amplitude dissipation and phase dispersion are separated effectively, and then the decoupled formulation is converted back to the time–space domain. The new equivalent and approximate equations can be obtained. Compared with the original equation, the newly approximated equation enables us to separate the amplitude-loss and phase-delay terms, respectively, and experiences four kinds of effects, namely acoustic effect, only amplitude-dissipation effect, only phase-dispersion effect and both amplitude-dissipation and phase-dispersion effect. Moreover, the stability condition and numerical dispersion for using DGM to solve the new and old equations are presented. Several numerical examples are used to verify the correctness and effectiveness of the modified approximate equations in viscoacoustic isotropic and anisotropic media. The numerical results in a cave and SEG/EAGE salt models demonstrate that the new equations combined with DGM have performances on viscoacoustic media with complex geological structures.

Attenuation-compensated reverse time migration based on the stereo-modeling operator

The viscosity in the subsurface is ubiquitous. Prestack reverse time migration (RTM) based on the viscoacoustic wave equation is an accurate imaging method to compensate attenuation. The resolution of an imaging result highly depends on the accuracy of the wave equation solution. We extend the nearly analytic center difference (NACD) method to a viscoacoustic wave equation to compensate seismic attenuation. The NACD method can achieve fourth-order accuracy in the time and space domain. We use the pressure and its gradients to approximate the spatial derivatives based on stereo modeling. Numerical simulations and dispersion analysis demonstrate that our NACD has less numerical dispersion and is more accurate than the Lax-Wendroff correction (LWC) method. Thus, the NACD is less expensive in terms of computational time and storage space. Numerical simulations indicate that the forward modeling of the viscoacoustic wavefield via the NACD is able to obtain more accurate wavefields compared with the conventional LWC method and can accurately simulate the effect of the viscosity on waves in attenuative media. When applied to RTM, the viscous effects of the subsurface can be compensated. Synthetic examples demonstrate that the attenuation-compensated RTM using NACD can obtain an image with higher resolution in attenuative media compared with that using the acoustic RTM. Furthermore, the image determined by the NACD is clearer than that of the LWC method in RTM.

Stable attenuation-compensated reverse time migration and its application to land seismic data

Intrinsic attenuation of the earth causes energy loss and phase distortion in seismic wave propagation. To obtain high-resolution imaging results, these negative effects must be considered during reverse time migration (RTM). We can easily implement attenuation-compensated RTM using the constant Q viscoacoustic wave equation with decoupled amplitude attenuation and phase dispersion terms. However, the nonphysical amplitude-compensation process will inevitably amplify the high-frequency noise in the wavefield in an exponential form, causing the numerical simulation to become unstable. This is due to the fact that the amplitude of the compensation grows exponentially with frequency. In order to achieve stable attenuation-compensated RTM, we modify the analytic expression of the attenuation compensation extrapolation operator and make it only compensate for amplitude loss within the effective frequency band. Based on this modified analytic formula, we then derive an explicit time-space domain attenuation compensation extrapolation operator. Finally, the implementation procedure of stable attenuation-compensated RTM is presented. In addition to being simple to implement, the newly proposed attenuation-compensated extrapolation operator is superior to the conventional low-pass filter in suppressing random noise, which will further improve the imaging resolution. We use two synthetic and one land seismic datasets to verify the stability and effectiveness of the proposed attenuation-compensated RTM in improving imaging resolution in viscous media.

GPU performance analysis for viscoacoustic wave equations using fast stencil computation from the symbolic specification

Seismic forward modeling is a computationally and data-intensive stage in the seismic processing workflow. By profiling the kernels of seismic forward modeling algorithms, it was observed that they need to access a wide variety of memory locations, in addition to the computational cost of performing floating-point operations for the numerical solution of wave equations. In this context, the Roofline model was used to analyze six representative computing kernels in seismic modeling on GPU environment to indicate bottlenecks in the performance and suggest improvements of these wave equation propagators. Based on this, six viscoacoustic equations were implemented using the Devito tool. Experimental data have shown that optimizations in increasing data reuse and decreasing off-chip memory traffic can significantly improve performance.

An explicit stable Q-compensated reverse time migration scheme for complex heterogeneous attenuation media

Prestack reverse-time migration (RTM) is a popular imaging technique for complex geological conditions, since the amplitude attenuation and velocity dispersion are common in seismic recordings. To image attenuated seismic recordings accurately, a robust migration algorithm with a stable attenuation compensation approach should be considered. In the context of the Q-compensated RTM approach based on the decoupled fractional Laplacians (DFLs) viscoacoustic wave equation, amplitude compensation can be implemented by flipping the sign of the dissipation term. However, the non-physical magnification of image amplitude could lead to a well-known numerical instability problem. The explicit stabilization operator can rectify the amplitude attenuation and suppress the numerical instability. However, limited by the inconvenient mixed-domain operator, the average Q value rather than the real Q value is often used in the compensation operator, lowering the compensated accuracy of the migration image. To overcome this problem, we propose a novel explicit Q-compensation scheme. The main advantage of the proposed compensation operator is that its order is space-invariant, making it more suitable for handling complex heterogeneous attenuation media. Several two-dimensional (2D) and three-dimensional (3D) synthetic models are used to verify the superiority of the proposed approach in terms of amplitude fidelity and image resolution. Field data further demonstrates that our approach has potential applications and can greatly enhance the resolution of seismic images.

Modeling seismic wave propagation in the Loess Plateau using a viscoacoustic wave equation with explicitly expressed quality factor

The thick Quaternary loess on the Loess Plateau of China produces strong seismic attenuation, resulting in weak reflections from subsurface exploration targets. Accurately simulating seismic wavefield in the Loess Plateau is important for guiding subsequent data processing and interpretation. We present a 2D/3D wavefield simulation method for the Loess Plateau using a viscoacoustic wave equation with explicitly expressed quality factor. To take into account the effect of irregular surface, we utilize a vertically deformed grid to represent the topography, and solve the viscoacoustic wave equation in a regular computational domain that conforms to topographic surface. Grid deformation introduces the partial derivatives such as ∂vx/∂z and ∂vy/∂z in the wave equation, which is difficult to be accurately computed using traditional staggered-grid finite-difference method. To mitigate this issue, a finite-difference scheme based on a fully staggered-grid is adopted to solve the viscoacoustic wave equation. Numerical experiments for a simple layer model and 2D/3D realistic Loess Plateau models demonstrate the feasibility and adaptability of the proposed method. The 3D modeling results show comparable amplitude and waveform characteristics to the field data acquired from the Chinese Loess Plateau, suggesting a good performance of the proposed modeling method.

Target-oriented Q-compensated reverse-time migration by using optimized pure-mode wave equation in anisotropic media

Research on seismic anisotropy and attenuation plays a significant role in exploration geophysics. To enhance the imaging quality for complicated structures, we develop several effective improvements for anisotropic attenuation effects in reverse-time migration (Q-RTM) on surface and vertical seismic profiling (VSP) acquisition geometries. First, to suppress pseudo-shear wave artifact and numerical instability of the commonly used anisotropic pseudo-acoustic wave equations, an optimized pure P-wave dispersion relation is derived and the corresponding pure-mode wave equation is solved by combining the finite-difference and Possion methods. Second, a simplified anisotropic pure-mode visco-acoustic wave equation (PVAWE) based on standard linear solid model is established. Third, a time-dispersion correlation strategy is applied to improve the modeling accuracy. Fourth, we extend a target-oriented scheme to anisotropic attenuated modeling and imaging. Instead of the conventional wavefield modeling and RTM, the proposed approach can extract available wavefield information near the target regions and produce high imaging resolution for target structures. Last, both anisotropic surface and VSP Q-RTMs are executed by combining optimized PVAWE, time-dispersion correlation and target-oriented algorithm. Modeling examples demonstrate the advantages of our schemes. Moreover, our modified Q-compensated imaging workflow can be regarded as a supplement to the classical anisotropic RTM.

Frequency-Domain Reverse-Time Migration With Attenuation Compensation

Seismic wave suffers from amplitude attenuation and phase distortion when propagating in the attenuating media, thus reverse time migration (RTM) for viscous media should take the attenuation effects into consideration. Compensating for the attenuation effects in RTM may occur numerical instability because of the exponential amplification of the extrapolated wavefields. To obtain stable imaging results, we have developed a stabilized Q-compensated reverse time migration in the frequency domain. This algorithm is implemented by the following steps: first, we use the Kolsky-Futterman model to derive a frequency-domain viscoacoustic wave equation, which can simulate the amplitude loss and phase dispersion effects separately. Then, we calculate the source wavefields in the viscoacoustic media. Next, treating the recorded (viscoacoustic) data as the receiver sources, we can obtain the phase-dispersion-only and viscoacoustic receiver wavefields, which can be used to construct the stabilized Q-compensated receiver wavefields. Finally, we apply the deconvolution imaging condition for obtaining a Q-compensated image. A simple anticline model and gas chimney model are used to verify the effectiveness of the proposed approach. The Q-compensated images for the noise-free data indicate the algorithmic stability and compensation accuracy of the proposed scheme. The noisy data tests for the gas chimney model demonstrate the good anti-noise property of our method. The field data applications further prove its feasibility and practicability.

Approximating Constant-Q Seismic Wave Propagations in Acoustic and Elastic Media Using a Cole–Cole Model

ABSTRACT Seismic wave propagation in the Earth’s interior inevitably encounters attenuation and dispersion effects, which usually can be represented by a constant-Q model. However, solving the constant-Q wave equations formulated by fractional Laplacians is computationally intensive. Alternatively, the Cole–Cole model provides an optimal description of seismic attenuation. Because of the fractional time derivatives of both stress and strain in the expression, this method exhibits good adaptability and flexibility. In this article, we investigate the performance of the Cole–Cole model to approximate constant-Q behaviors with different fractional orders in acoustic and elastic media. The phase velocity and quality factor are compared to determine an optimal fractional order. After that the Cole–Cole model can be easily represented by the other three parameters (i.e., relaxed modulus, minimum angular frequency, and minimum quality factor), which are similar to the constant-Q theory (reference modulus, reference angular frequency, and frequency-independent quality factor). The first-order viscoacoustic and viscoelastic wave equations are derived to implement seismic wavefield simulations by combining conservation equation and stress–displacement relation. Furthermore, a time-domain algorithm is developed to solve the wave equations based on the Grunwald–Letnikov approximation and finite-difference scheme. Numerical modeling results calculated by our proposed method have good consistencies with the reference solutions from the constant-Q theory, suggesting that a small fractional order can well approximate the constant-Q outputs in a broad frequency band.

Frequency-Domain Q-Compensated Reverse Time Migration Using a Stabilization Scheme

Seismic attenuation occurs during seismic wave propagation in a viscous medium, which will result in a poor image of subsurface structures. The attenuation compensation by directly amplifying the extrapolated wavefields may suffer from numerical instability because of the exponential compensation for seismic wavefields. To alleviate this issue, we have developed a stabilized frequency-domain Q-compensated reverse time migration (FQ-RTM). In the algorithm, we use a stabilized attenuation compensation operator, which includes both the stabilized amplitude compensation operator and the dispersion correction operator, for the seismic wavefield extrapolation. The dispersion correction operator is calculated based on the frequency-domain dispersion-only viscoacoustic wave equation, while the amplitude compensation operator is derived via a stabilized division of two propagation wavefields (the dispersion-only wavefield and the viscoacoustic wavefield). Benefiting from the stabilization scheme in the amplitude compensation, the amplification of the seismic noises is suppressed. The frequency-domain cross-correlation imaging condition is exploited to obtain the compensated image. The layered model experiments demonstrate the effectiveness and great compensation performance of our method. The BP gas model examples further verify its feasibility and stability. The field data applications indicate the practicability of the proposed method. The comparison between the acoustic and compensated results confirms that the proposed method is able to compensate for the seismic attenuation while suppressing the amplification of the high-frequency seismic noise.

A stable Q reverse time migration method with regularization

Q reverse time migration (Q-RTM) is an effective method to compensate the deep amplitude energy of seismic records and improve the resolution. However, an important problem of the Q-RTM method is that the wavefield extrapolation process will generate high-frequency noise, resulting in unstable imaging results. On the basis of the decoupled viscousacoustic wave equation based on standard linear solid model, a viscoacoustic wave equation with regularization term is derived, and a stable Q-RTM based on this equation is developed. And the detailed derivation process of the viscousacoustic wave equation containing the regularization term is given. Compared with the commonly used low-pass filtering method, the method based on the regularization term avoids the design of the filter and the selection of the cutoff frequency. We test the Q-RTM method on the BP gas layer model, the Marmousi model and real seismic data. The numerical example results show that the Q-RTM method proposed in this paper can effectively solve the instability problem, compensate the amplitude energy, and obtain a higher resolution seismic profile.

Q-compensated least-squares reverse time migration with velocity-anisotropy correction based on the first-order velocity-pressure equations

Anisotropy and Q attenuation bring great challenges to seismic wave migration. On migrated images, anisotropy creates structural and positioning errors, and Q attenuation leads to weak amplitudes and misplacement of reflectors. A 2D Q-compensated least-squares reverse time migration with velocity-anisotropy correction ( QLSRTM-VA) is proposed through the construction of velocity-anisotropic Q-compensated forward modeling, Q-compensated adjoint, and Q-attenuated demigration operators to simultaneously correct velocity-anisotropy and Q-attenuation in the migration process. The preceding operators are derived using first-order velocity-anisotropic viscoacoustic quasi-differential wave equations with variable densities, which are stable, capable of conveniently dealing with variable density media and are easy to transform between velocity-anisotropic Q-compensation and Q-attenuation versions. As exemplified by two synthetic and field data sets, our QLSRTM-VA method increases the imaging resolution, signal-to-noise ratio, and amplitude preservation in deep regions. Our method is capable of producing better images than viscoacoustic isotropic least-squares reverse time migration (LSRTM) and acoustic anisotropic LSRTM.

Two exact first-order k-space formulations for low-rank viscoacoustic wave propagation on staggered grids

Wave propagation in the viscoacoustic media is physically dispersive and dissipated. Completely excluding the numerical dispersion error from the physical dispersion in the viscoacoustic wave simulation is indispensable to understanding the intrinsic property of the wave propagation in attenuated media for the petroleum exploration geophysics. In recent years, a viscoacoustic wave equation characterized by fractional Laplacian gains wide attention in geophysical community. However, the first-order form of the viscoacoustic wave equation, often solved by the conventional staggered-grid pseudospectral method, suffers from the numerical dispersion error in time due to the low-order finite-difference approximation. It is challenging to completely eliminate the error because the viscoacoustic wave equation contains two temporal derivatives, which stem from the time stepping and the amplitude attenuation terms, respectively. To tackle the issue, we derive two exact first-order k-space viscoacoustic formulations that can fully exclude the numerical error from the physical dispersion. For the homogeneous case, two formulations agree with the viscoacoustic analytical solution very well and have the same efficiency. For the heterogeneous case, our second k-space formulation is more efficient than the first one because the second formulation significantly reduces the number of the wavenumber-space mixed-domain operators, which are the expensive part of the viscoacoustic k-space simulation. Numerical cases validate that the two first-order k-space formulations are effective and efficient alternatives to the current staggered-grid pseudospectral formulation for the viscoacoustic wave simulation.

Viscoelastic and viscoacoustic modelling using the Lie product formula

ABSTRACTWe present efficient and accurate modelling of seismic wavefields in anelastic media. We use a first‐order viscoelastic wave equation based on the generalized Robertsson's formulation. In our work, viscoacoustic and viscoelastic wave equations use the standard linear solid mechanism. To numerically solve the first‐order wave equation, we employed a scheme derived from the Lie product formula, where the time evolution operator of the analytic solution is written as a product of exponential matrices, and each exponential matrix term is approximated by the Taylor series expansion. The accuracy of the proposed scheme is evaluated by comparison with the analytical solution for a homogeneous medium. We also present simulations of some geological models with different structural complexities, whose results confirm the accuracy of the proposed scheme and illustrate the attenuation effect on the seismic energy during its propagation in the medium. Our results demonstrate that the numerical scheme employed can be used to extrapolate wavefields stably for even larger time steps than those usually used by traditional finite‐difference schemes.

Viscoacoustic wave equation for qP-wave in transversely isotropic media

The wave equation acts a pivotal part in seismic data processing. We derived two time-domain pure qP viscoacoustic wave equations suitable for attenuation anisotropic media, providing a choice for pure qP-wave forward modeling and inversion in attenuation anisotropic media. Our new qP-wave equation in vertical transversely isotropic (VTI) media are based on the dispersion relation for qP-wave in VTI media and obtained by the relaxation functions of the standard linear solid model and the Kjartansson's constant Q model respectively. Anelasticity including amplitude attenuation and phase dispersion can be dealt with by the proposed wave equation separately. With the coordinate transform, the corresponding two pure qP-wave viscoacoustic anisotropic wave equations for titled transversely isotropic media also are derived. Different from the common viscoacoustic anisotropic wave equation (a coupled second-order partial differential equations system), the new wave equations are simplified by decomposing the original pseudo-differential operator into a differential and a scalar opera- tors. For numerical calculation, the new wave equations can be solved more efficiently than the coupled wave equations. Meanwhile the new wave equations could simulate pure qP-wave steadily in attenuation anisotropic media without quasi-SV wave artifacts. Some synthetic examples, including attenuation Hess and Marmousi-2 models, demonstrate the feasibility of the derived pure qP-wave viscoacoustic wave equations in transversely isotropic media.

Q-compensated least-squares reverse time migration in TTI media using the visco-acoustic TTI wave equation based on the SLS model

The anelasticity and anisotropy widely exist in real subsurface media. Strong anelasticity will lead to phase dispersion and amplitude attenuation during seismic wave propagation. Anisotropy cause seismic waves to have obviously different kinematic characteristics from that of isotropy. For seismic migration, ignoring the anelasticity and anisotropy of subsurface media will significantly reduce the quality of migration images, even cause imaging failure. We propose a Q-compensated least-squares reverse time migration (Q-LSRTM) in tilted transversely isotropic (TTI) media to correct these effects. According to the Born approximation of seismic wave equation, a linearized visco-acoustic TTI pure qP-wave modeling operator is derived using a new visco-acoustic TTI wave equation for one standard linear solid (SLS) model, which can deal with the anelasticity and can simulate pure qP-wave steadily in attenuating anisotropic media without qSV-wave artifacts. Then, the corresponding adjoint equation is formulated using the adjoint-state method to calculate the gradient sensitivity kernel for the visco-acoustic TTI media. Because of the least-squares inversion, the Q compensation can be achieved during the iterations, so that the over-amplification of noises can be avoided naturally. In addition, compared with conventional LSRTM, the proposed method compensates for the anelasticity and corrects the anisotropy, so as to produce images with better spatial resolution and amplitude fidelity. Numerical examples demonstrate the feasibility and advantages of the proposed method for the data including strong attenuation effects over conventional LSRTM.

Viscoacoustic wave equations for the power-law dependence of Q on frequency

The power-law dependence of the quality factor Q on frequency (i.e. Q = Q 0 | f / f 0 | γ ), which is called the power-law frequency-dependent Q for brevity, is a common phenomenological description of seismic wave attenuation in the interior of the Earth. The wave equation in differential form is essential to forward and inverse modelling of dissipative seismic waveforms in an accurate and efficient manner. However, all existing methods are seemingly unable to explicitly incorporate the exponent parameter γ into the wave equation in differential form. This drawback apparently limits the development of gradient-based inverse methods for spatially varying γ via the wave equation. In this paper, we use a weighting function method to derive the viscoacoustic wave equations for the power-law frequency-dependent Q , which explicitly involve the parameters γ and Q 0 . A critical step in this method is to construct a dissipative model, for which the complex modulus is expressed as an N -order series in terms of the weighting function. Numerical examples are used to illustrate the accuracy of the dissipative model, the effect of γ on waveforms, and the application of the newly proposed wave equations in viscoacoustic wavefield modelling for the power-law frequency-dependent Q .

Decoupled Fréchet kernels based on a fractional viscoacoustic wave equation

We have formulated the Fréchet kernel computation using the adjoint-state method based on a fractional viscoacoustic wave equation. We first numerically prove that the 1/2- and the 3/2-order fractional Laplacian operators are self-adjoint. Using this property, we find that the adjoint wave propagator preserves the dispersion and compensates the amplitude, whereas the time-reversed adjoint wave propagator behaves identically to the forward propagator with the same dispersion and dissipation characters. Without introducing rheological mechanisms, this formulation adopts an explicit [Formula: see text] parameterization, which avoids the implicit [Formula: see text] in the conventional viscoacoustic/viscoelastic full-waveform inversion ([Formula: see text]-FWI). In addition, because of the decoupling of operators in the wave equation, the viscoacoustic Fréchet kernel is separated into three distinct contributions with clear physical meanings: lossless propagation, dispersion, and dissipation. We find that the lossless propagation kernel dominates the velocity kernel, whereas the dissipation kernel dominates the attenuation kernel over the dispersion kernel. After validating the Fréchet kernels using the finite-difference method, we conduct a numerical example to demonstrate the capability of the kernels to characterize the velocity and attenuation anomalies. The kernels of different misfit measurements are presented to investigate their different sensitivities. Our results suggest that, rather than the traveltime, the amplitude and the waveform kernels are more suitable to capture attenuation anomalies. These kernels lay the foundation for the multiparameter inversion with the fractional formulation, and the decoupled nature of them promotes our understanding of the significance of different physical processes in [Formula: see text]-FWI.

Q-compensated full waveform inversion for velocity and density

Full waveform inversion (FWI) makes full use of the seismic waveforms to find high-resolution velocity and density models by minimising residuals between the calculated and recorded data. Q attenuation widely exists in the subsurface media, leading to weak amplitude and misplacement of reflectors. However, the commonly used Q-compensated FWI (QFWI) based on the second-order wave equation has difficulties in simultaneously inverting velocity and density fields. A QFWI method based on new first-order viscoacoustic quasi-differential equations is proposed to simultaneously produce velocity and density fields. Based on the adjoint state inversion theory, Q-compensated forward-propagated operators, adjoint operators, and gradient equations are derived using the newly derived first-order viscoacoustic quasi-differential wave equations. The time-domain multi-scale decomposition method is introduced to update the velocity and density models from a low to a high wavenumber. Numerical examples on an actual work area model and a modified attenuating Marmousi model show that the proposed QFWI method produces higher-accuracy velocity and density models with iterations by correcting the Q attenuation than the conventional acoustic FWI. Even when the Q model is extremely inaccurate, the proposed QFWI obtains acceptable inversion results. Compared to the conventional QFWI, our QFWI better inverts velocity field in the case of an inaccurate density model. Finally, we verify the adaptability of our QFWI to field data.