Part of my research program has been devoted to developing new, comprehensive quantitative probes of our planet’s interior which include novel convergent nonlinear iterative techniques (surpassing so-called standard ‘adjoint state’ methods). One of our most recent breakthroughs concerns the understanding of the fundamental role of conditional Hölder (or Lipschitz) stability estimates of the nonlinear seismic inverse problem.
Novel iterative methods, conditional Lipschitz stability and convergence
Forward modelling and optimization to match waveforms is an approach that was introduced in seismology more than four decades ago. Even though this approach, commonly using implementations based on the ‘adjoint state’ method, is being rather widely applied today, no results on the recovery of coefficients of models, however, were obtained. The first fundamental question to be addressed is whether the inverse boundary value problem corresponding with idealized seismic data has a unique solution. The second question is whether explicit reconstructions can be obtained, and the third question is what one can say about the case of partial data, also subject to sampling. The fourth question concerns the incorporation of inaccurate data and uncertainty quantification. All these questions are intimately interconnected. The second question typically involves iterative methods (‘full waveform inversion’) and leads to a local result.
We succeeded in obtaining conditional convergence of iterative methods results; we obtained explicit expressions for the radii of convergence . The convergence relies on a Hölder, or Lipschitz, stability estimate for (essentially well-posedness of) the inverse problem. These insights have a significant impact on our understanding and develop- ment of optimization techniques very different from optimization with traditional regularization strategies.
- Lipschitz stability estimates: scalar waves. The time-harmonic or ‘frequency-domain’ formulation lends itself to obtain conditional Lipschitz stability estimates for the nonlinear seismic inverse boundary value problem, even assuming partial boundary data . The data are the Dirichlet-to-Neumann map. Uniqueness of recovery for this inverse problem in the case of scalar waves has been known since the late 1980s under minimal regularity assumptions. We established these estimates in the case of piecewise constant models with underlying domain partitionings containing discontinuities.
The Lipschitz stability constant grows exponentially with the number of subdomains in the domain partitioning. The stability constant appears in the expressions for the radii of convergence, in a way such that the radius tends to zero as this constant becomes large. We obtained a sharp estimate of the stability constant in terms of the number of subdomains in the domain partitioning [preprint 7]. This estimate plays a critical role in the development of a multi-level (multi-scale) scheme.
Moreover, we obtained the Lipschitz stable recovery of polyhedral interfaces and – if the domain partitioning in the above is an (unstructured) tetrahedral mesh – the Lipschitz stable recovery of the number of elements in and the shape of the mesh [preprint 8]. Thus we fundamentally relaxed knowledge about the parametrization.
The Dirichlet-to-Neumann map is related to a single layer potential operator, the kernel of which is directly related to seismic data.
- Time-harmonic elastic waves. We obtained the unique recovery of Lame´ parameters and density that are piece- wise constant on an unstructured tetrahedral mesh from the Neumann-to-Dirichlet map (which corresponds with vibroseis data). We also obtained Lipschitz stability.
- Multi-level, multi-scale schemes. Based on the conditional stability estimates we introduce convex subsets, fol- lowing a refinement of domain partitionings. We designed an iterative method which converges to the point in each subset closest to the unique solution. The distance, or approximation error, depends on frequency, and an estimate can be obtained. Mitigating the growth of the stability constants with the approximation errors controlled by frequency, we obtained a multi-level (nonlinear projected steepest descent iteration) scheme with precise conditions for convergence  defining a nonlinear hierarchical compressed reconstruction. Through the compression and corresponding approximation we incorporate inaccuracies in the data. Other more ad- vanced methods, in particular the Gauss-Newton method, can also be adapted for this purpose.
- Attenuation. The Lipschitz stability estimates hold for complex (frequency dependent) coefficients and, hence, the time-harmonic formulation allows in principle to obtain localized information about the (reciprocal) Q factor, and spatial variability in physical dispersion. This might open the way to the further studying of attenuation in localized regions in the mantle, such as in the asthenosphere.
Hyperbolic inverse boundary value problem
The time-domain formulation of the inverse boundary value problem for scalar and elastic waves involves different techniques and strategies. Naturally, in this formulation, one can introduce time windows identifying structure in and emphasizing parts of the data. We have been studying this problem for coefficients with structure. On the one hand, in the case of piecewise smooth coefficients with interfaces of unknown shapes, we use connections with geometric inverse problems with corresponding conditions; on the other hand, we establish explicit connections with the time- harmonic inverse boundary value problem, a multi-level strategy, and conditional high-frequency stability estimates.