Lowrank seismic wave extrapolation on a staggered grid

Error analysis using the Method of Manufactured Solutions

Next we use the method of manufactured solutions (MMS) (Salari and Knupp, 2000) to analyze numerical error of the proposed methods. MMS provides an approach to designing exact reference solutions for wave equations in heterogeneous media. We couple the use of manufactured solutions with mesh refinement to plot the error curve for each method. To begin our investigation, we consider the system of 1D first order wave equations,

$$\begin{array}{l} \displaystyle \frac{\partial u(x,t)}{\partia... ...)v^2(x)\frac{\partial u(x,t)}{\partial x}+s_p(x,t), \end{array}$$ (21)

which are subjected to appropriate boundary and initial conditions, where $s_u$ and $s_p$ represent inject source terms for the particle velocity $u$ and pressure $p$ , respectively. For this study, we are interested in the governing equation. We choose a Gaussian pulse as the solution of the pressure and particle velocity in variable velocity media (Ober et al., 2009),

$$\displaystyle u(x,t) = p(x,t) = e^{-\lambda^2(x-x_0-v(x)t)^2},$$ (22)

where $\lambda$ is the wavelength of the Gaussian pulse. $x_0$ is the source location and $v(x)$ is the variable velocity. It is easy to derive the corresponding source terms from equations 21 and 22,

$$\begin{array}{l} \displaystyle s_u(x,t) = 2\lambda^2[x-x_0-v(... ...2[x-x_0-v(x)t][v(x)+\rho(x)v^2(x)(v'(x)t-1)]p(x,t). \end{array}$$ (23)

For numerical solution, accuracy is often affected by various factors. In this study, we examine the numerical error of the proposed methods when using different temporal and spatial discretization and wavelets of different dominant frequencies. The velocity we used here increases with $x$ , defined as $v(x)=2.1+0.1x^2(km/s)$ . Its gradient is $v'(x)=0.2x$ correspondingly. We use constant density in this experiment. Figure 6a shows the RMS errors of the wavefields for different time interval. For all these calculations, we keep space interval $\Delta x = 25m$ as constant. The blue dashed line and green dotted line indicate the errors of SGL method and SGLFD method respectively, while red solid line plots the errors of the traditional SGFD method. Compared with the traditional SGFD method, the proposed lowrank methods correct the distortion caused by increasing $\Delta t$ . They exhibit a high accuracy in time. Figure 6b shows the recorded trace at the $x$ of $4km$ for $\Delta t=4 ms$ .The black solid line corresponds to the exact solution generated by MMS. The colors of remaining lines has the same meaning as in Figure 6a. Numerical dispersion is more visible when $\Delta t$ is increased; this effect is much less significant in our methods.

et,trec19
Figure 6. (a) RMS errors as function of time interval. (b) Recorded data for $\Delta t=4 ms$ . Blue dashed line: SGL method. Green dotted line: SGLFD method. Red solid line: SGFD method. Black solid line: MMS.

Figure 7a shows the RMS errors of the wavefields for different space intervals. We use a small time interval $\Delta t = 0.4ms$ to keep time error small. Figure 7b shows the recorded data for $\Delta x=40m$ . The first-order $kx$ -space operators (equation 9) have a high accuracy in space which make the error of SGL method increase slowly when $\Delta x$ increase. The coefficients of the SGLFD are obtained by applying least-squares fitting to the SGL method. The error from the least-squares fitting make SGLFD method less accurate than SGL method. On the other hand, the coefficients of lowrank FD are optimized and auto-adapted to variations in velocity, which makes it significantly more accurate than SGFD method.

em,mrec7
Figure 7. (a) RMS errors as function of space interval. (b) Recorded data for $\Delta x=40m$ . Blue dashed line: SGL method. Green dotted line: SGLFD method. Red solid line: SGFD method. Black solid line: MMS.

In applications of seismic wave extrapolation, temporal and spatial intervals are often increased at a certain ratio for saving computational cost. Figure 8 show the error curves for a temporal and spatial refinement study where both $\Delta x$ and $\Delta t$ are increased simultaneously. We define the Courant-Friedrichs-Lewy (CFL) number as $\alpha=v_{max}\Delta t/\Delta x$ to specify the stability condition, where $v_{max}$ is the maximum velocity. Figure 8a shows the RMS errors for $\alpha =0.2$ . Figure 8b shows the RMS error of $\alpha =0.4$ , while figure 8c for $\alpha =0.8$ . The color of the lines has the same meaning as indicated before. We can see that the errors increase with the CFL number $\alpha$ . This increase is especially significant for SGFD method. Both SGL and SGLFD methods keep high accuracy for a larger scale of temporal and spatial intervals than SGFD.

ea2,ea4,ea8
Figure 8. RMS errors for CFL number (a) $\alpha =0.2$ , (b) $\alpha =0.4$ , (c) $\alpha =0.8$ . Blue dashed line: SGL method. Green dotted line: SGLFD method. Red solid line: SGFD method.

Seismic exploration techniques, such as prestack depth migration or FWI, need the modeling engine which remains highly accurate for a wide band seismic wavelet. The numerical dispersion can be serious for a high frequency source. Figure 9a compares the accuracy of different methods for different source frequencies. The meanings of different line colors are the same as before. Figure 9b shows the recorded data with source frequency of $40Hz$ . It can be seen that the SGFD method has a visible dispersion, while SGL and SGLFD methods remain almost dispersion free.

fre,frec15
Figure 9. (a) RMS errors as function of dominant frequency of source wavelet. (b) Recorded data for frequency of $40Hz$ . Blue dashed line: SGL method. Green dotted line: SGLFD method. Red solid line: SGFD method. Black solid line: MMS.

Lowrank seismic wave extrapolation on a staggered grid