ADMM-based Wavefield Reconstruction inversion

Meeting the challenge of cycle skipping

Imaging large-contrast media containing salt bodies by FWI is a major challenge since it is difficult to satisfy the cycle skipping criterion (fit the data with a kinematic error smaller than half the period from the first FWI iteration). In order to cope with this challenge, H. Aghamiry with Pr. A. Gholami developed the Wavefield Reconstruction Inversion (WRI) method, originally proposed by van Leeuwen & Herrmann (2013) and van Leeuwen & Herrmann (2016), with edge-preserving regularization in the framework of the alternating-direction method of multiplier ( ADMM). The primary goal of the method is to fit data with inaccurate subsurface model to mitigate the ill-famed cycle skipping pathology. To achieve this goal, we simulate the wavefields with a wave equation relaxation generated by a feedback term to the data. Then, we update the subsurface parameters by minimizing the wave equation violation such that the reconstructed wavefields are pushed back toward the wave equation. Wavefield simulation and parameter estimation are performed in alternating mode in the framework of ADMM. This allows us to break down the waveform inversion as a sequence of two linear problems, capitalizing on the bilinearity of the wave equation. Moreover, recasting the parameter estimation problem as a linear problem simplifies the interfacing of nonsmooth regularization such as total variation regularization with the optimization algorithm using the split Bregman recipe developed by Goldstein & Osher and proximal algorithms. Check the publication list of Aghamiry et al. for more details.

Why should WRI be suitable for long offset data?

In WRI based upon alternating directions, subsurface parameters are updated keeping fixed the simulated wavefields, which are assumed to be representative of the true wavefields (those propagating in the true subsurface medium). These simulated wavefields are computed by jointly solving in a least-squares sense the wave equation (physical constraint) and the observation equation (feedback to the data). When long offsets are available, we record energetic wide-angle data (diving waves, post-critical reflections) at the surface of the earth that have deeply penetrate the subsurface. In this long-offset framework, we expect the data-assimilated wavefields to be more representative of the true wavefields than those which would have been reconstructed from short-offset data because the information collected at the surface has more fully illuminated the subsurface. This improved wavefield reconstruction should translate into a better posed, more convex and faster converging inverse problem. This statement is illustrated in Figure 1 with a two-layer model synthetic experiment. The goal of the exercise is to reconstruct the velocities in the two homogeneous half spaces and the position of the reflector with the 5Hz frequency and starting from a homogeneous velocity model whose velocity differs from those of the two layers.The experimental setup is designed such that classical FWI would remain rapidly stuck into a local minimum due to cycle skipping. The first experiment (Figure 1a-d) is performed with a limited-offset experiment. The inversion remains stuck into a local minimum due to the ill-famed velocity-depth ambiguity (the pathology is not cycle skipping here). In the second experiment, the maximum offset is doubled. The inversion converges to the correct solution because the long-offset acquisition allows for the recording of post-critical reflections and refracted waves from the interface, allowing for a more accurate wavefield reconstruction during WRI. In Figure 1, we also clearly see that WRI proceeds from shallow to depth for surface acquisitions, following some kinds of layer stripping approach. This may result because the reconstructed wavefields are more accurate nearby the receivers. Therefore, WRI starts updating the shallow part of the subsurface where receivers are located and energetic waves propagate. When WRI updates the shallow part during the early iterations, it generates artificial boundaries between the updated and non-updated parts. These artificial boundaries feed the data misfit function with a new sets of residuals, which are used to prolongate in depth the subsurface update.

Figure 1: A two-layer synthetic experiment to illustrate how WRI takes advantage of long offsets to perform more accurate wavefield reconstruction and mitigates the ill-posedness of the waveform inversion (a-d) Maximum offset is 8km. (e-h) Maximim offset is 16km. The top row shows the true velocity model and the true 5Hz monochromatic wavefield. The next rows show from top to bottom the results of ADMM-based WRI at iterations 1, 15, 100. From Aghamiry (2019).

After intensive testing of ADMM-based WRI on synthetic tests (Figures 2-3), a central objective of WIND is to assess more precisely the potential and limits of this approach against real data case studies as those illustrated in (1) and (2).

Figure 2: The target of the BP model reconstructed by ADMM-based WRI with different regularizations. (a) Bound constraints (zero-order Tikhonov). (b) Tikhonov regularization. (c) First-order total variation (TV) regularization. (d) Hybrid Tikhonov + TV regularization implemented by convex combination. (e) Hybrid Tikhonov + TV (TT) regularization implemented by infimal convolution. (f) Total generalized variation (TGV), the infimal convolution of 1st and 2nd order total variation. The TT regularization provides the best reconstruction because the subsurface parameters have been broken down explicitly into smooth (m2) and piecewise constant (m1) components (Figure 16) such that the regularization can be tailored to each component (Tikhonov for the smooth component, TV for the blocky one). These two components are tied together (m=m1+m2) by the wave equation and the observation equation. See for Gholami and Hosseini (2013) and Aghamiry et al. (2019) for more details.

fig_m1_m2.png

Figure 3: The blocky (left) and smooth (right) components reconstructed by TV+Tikhonov regularization (top) and Total Generalized Variation (bottom) regularizations based upon infimal convolution.