Constant-Q wave propagation and compensation by pseudo-spectral time-domain methods

Highlights

• An efficient viscoacoustic pseudo-spectral time-domain (PSTD) modeling method.

• A stable amplitude-boosted PSTD wave propagation engine.

• Application of a time-variant filter in amplitude-boosted PSTD simulation.

• A stable attenuation-compensated reverse-time migration method.

Abstract

Seismic attenuation is usually described by a constant-Q (CQ) model that assumes the seismic quality factor (Q) is independent of the frequency. To simulate the attenuation behaviors of seismic waves, we develop a pseudo-spectral time-domain (PSTD) method to solve a CQ viscoacoustic wave equation. This method is nearly fourth-order accurate in time. Compared to the conventional temporal second-order PSTD method, the new PSTD scheme is verified to be more efficient. In some applications such as Q-compensated reverse-time migration (Q-RTM) and time-reversal imaging, one requires to simulate an anti-attenuation process. To realize this purpose, we switch our viscoacoustic PSTD modeling scheme into an amplitude-compensated PSTD modeling scheme by flipping the signs of the operators that dominate the amplitude loss. To control the numerical instability caused by high-frequency overcompensation in the amplitude-boosted modeling, we integrate a time-variant filter to the PSTD modeling scheme to suppress the high-frequency noise. Wavefield simulation examples in homogeneous media verify the temporal accuracy of our nearly fourth-order PSTD modeling scheme. A Q-RTM test of synthetic data is also presented to demonstrate the robustness of our amplitude-compensated PSTD modeling scheme.

Introduction

In the context of seismic applications, the constant-Q (CQ) model assuming the seismic quality factor (Q) is independent of the frequency is a good approximation of the earth’s viscosity in most cases (McDonal et al., 1958, Kjartansson, 1979). Therefore, the techniques of Q-compensated reverse-time migration (Q-RTM) (Deng and McMechan, 2007, Zhu et al., 2014, Sun et al., 2015, Xie et al., 2018, Zhou et al., 2018) and full-waveform inversion (FWI) (Chen et al., 2016b, Chen et al., 2017, Chen and Zhou, 2017, Fabien-Ouellet et al., 2017, Xue et al., 2017, Fang et al., 2020, Zhang et al., 2020, Xing and Zhu, 2020, Chen et al., 2020b) prefer to use a constant-Q (CQ) wave equation as wave propagation engines. The finite-difference time-domain (FDTD) (Wang et al., 2019, Zou et al., 2020, Di et al., 2021) and pseudo-spectral time-domain (PSTD) (Xie et al., 2018) methods are two popular numerical schemes to solve the seismic wave equations. FDTD is usually used to solve the generalized Maxwell or Zener body (GMB, GZB) models-based wave equations (Carcione et al., 1988, Robertsson et al., 1994). These wave equations are flexible in describing both frequency-dependent and frequency-independent Q behaviors. However, to approximate the constant-Q (CQ) behavior in seismic frequency band, the GMB/GZB model-based wave equations are reported to require more than two relaxation mechanisms (Blanch et al., 1995, Zhu et al., 2013), which increases the computational cost, compared to their counterparts of acoustic or elastic FDTD modeling. Another drawback of GMB/GZB viscoacoustic wave equations appears when they are used as wave propagation engines in Q-RTM. The Q compensation requires boosting the amplitudes while preserving phase the same as that in viscoacoustic forward modeling (Zhu, 2014). However, this cannot be realized directly for GMB/GZB models-based wave equations (Guo et al., 2016), since the operators controlling the amplitude loss and phase distortion are coupled together in these wave equations.

PSTD simulations of the decoupled fractional Laplacians (DFL) viscoacoustic wave equations have become popular for Q-RTM recently (Zhu et al., 2014, Sun et al., 2015, Li et al., 2016, Wang et al., 2018, Zhao et al., 2018), due to the high accuracy of these wave equations in matching the CQ model. Additionally, the operators controlling the amplitude loss are separated explicitly from the operators controlling the phase distortion in the DFL wave equations (Zhu and Harris, 2014, Chen et al., 2016a, Li et al., 2019, Shukla et al., 2019, Xing and Zhu, 2019). In this case, one only needs to flip the signs of the operators controlling the amplitude loss in reverse time propagation of seismic data. By this way, both the amplitude loss and phase distortion can be compensated correctly. However, the current PSTD modeling schemes for the DFL wave equations still have two shortcomings.

First, the conventional PSTD scheme applies low-order finite- difference (FD) operators to approximate the time derivatives and the wavefield temporal extrapolation accuracy is limited to second-order. The low-order time discretization could introduce visible time dispersion errors and is easy to cause numerical instability in the case of a large timestepping size. The low-rank approximation is later used to improve the accuracy and stability (Sun et al., 2015, Chen et al., 2016a, Chen et al., 2019a, Wang et al., 2020). However, this method is more complicated to program than PSTD, even though the involved low-rank decomposition is efficient (Fomel et al., 2013). Second, the conventional PSTD schemes used in Q-RTM adopt a time-invariant low-pass filter to suppress the numerical instability (Zhu et al., 2014, Sun et al., 2015, Li et al., 2016, Li et al., 2019) but the filter often sacrifices the accuracy by keeping the stability in the presence of strong noise. Wang et al. (2018) develops a time-variant filter to balance the accuracy and stability better. Chen et al. (2020a) further integrate this filter into the wavefield temporal extrapolation formula to avoid extra filtering, however the low-rank decomposition is still involved. As a different approach, Sun and Zhu (2018) present a new imaging condition to realize Q compensation without boosting the amplitudes, but introducing extra computational cost.

Based on an overall consideration of programming complexity and simulation accuracy, we present two PSTD modeling schemes in this study to simulate wave propagation with amplitude attenuation and compensation, respectively. Compared to the low-rank approximation-based PSTD methods, our PSTD schemes are developed with simpler Taylor approximation and they are easier to program. We notice that a few generalized FDTD methods including a Hermite distributed approximating functional (HDAF) method (Yao et al., 2017) and a matrix-transform method (Chen et al., 2019b) are developed recently to solve the DFL wave equations. However, we do not discuss these methods due to their higher computational cost.

We organize this work as follows. We first present a viscoacoustic PSTD modeling scheme with nearly fourth-order accuracy in time. Then, we analyze the Courant–Friedrichs–Lewy (CFL) stability condition and demonstrate the simulation accuracy. Next, we discuss an amplitude-compensated PSTD modeling scheme with a stabilization strategy involved. These contents are followed by numerical examples and conclusions.