Seismic imaging and inversion based on spectral-element and adjoint methods

Abstract

One of the most important topics in seismology is to construct detailed tomographic images beneath the surface, which can be interpreted geologically and geochemically to understand geodynamic processes happening in the interior of the Earth. Classically, these images are usually produced based upon linearized traveltime anomalies involving several particular seismic phases, whereas nonlinear inversion fitting synthetic seismograms and recorded signals based upon the adjoint method becomes more and more favorable. The adjoint tomography, also referred to as waveform inversion, is advantageous over classical techniques in several aspects, such as better resolution, while it also has several drawbacks, e.g., slow convergence and lack of quantitative resolution analysis.

In this dissertation, we focus on solving these remaining issues in adjoint tomography, from a theoretical perspective and based upon synthetic examples. To make the thesis complete by itself and easy to follow, we start from development of the spectral-element method, a wave equation solver that enables access to accurate synthetic seismograms for an arbitrary Earth model, and the adjoint method, which provides Frechet derivatives, also named as sensitivity kernels, of a given misfit function. Then, the sensitivity kernels for waveform misfit functions are illustrated, using examples from exploration seismology, in other words, for migration purposes. Next, we show step by step how these gradient derivatives may be utilized in minimizing the misfit function, which leads to iterative refinements on the Earth model. Strategies needed to speed up the inversion, ensure convergence and improve resolution, e.g., preconditioning, quasi-Newton methods, multi-scale measurements and combination of traveltime and waveform misfit functions, are discussed. Through comparisons between the adjoint tomography and classical tomography, we address the resolution issue by calculating the point-spread function, the action of the Hessian on an arbitrarily-chosen model perturbation, and the resolution function, the action of the resolution matrix on the arbitrarily-chosen model perturbation. Inner products between the two functions and the chosen model perturbation (properly normalized) are two scalars -- the point-spread parameter and the resolution parameter. The two functions serve as trade-off maps between the chosen model perturbation and all other model parameters, whereas the two parameters indicate whether the chosen model perturbation is well resolved in the inversion. While the point-spread function and the point-spread parameter work in relative sense, the resolution function and the resolution parameter are absolute quantities, regardless of the misfit function used in the inversion. Besides the optimization point of view, we also treat inverse problems from Tarantola's perspective -- the Bayesian inference, where each Earth model is associated with certain probability, preferably obeying multivariate normal distribution by choosing Cartesian model parameters, such as the logarithm of wavespeed. With a new limit-memory square root variable metric algorithm, we may sample the a posteriori distribution of model parameters, which allows statistical analysis on the inversion, e.g., addressing uncertainty and non-uniqueness of the inversion. Although, due to limit of time, seismic examples are to be added, analytical examples involving 20,000 model parameters validate our theory and algorithm, and it is promising that they can be easily adapted to real seismic applications. After solving both resolution and non-uniqueness issues, we finally extend capability of seismic inversions to consider noise simulations, i.e., by cross correlating noisy seismograms between pairs of seismic stations, without help of natural earthquakes and man-made explosions. At the end, we talk about implications of our studies on the model parameterization, in terms of both types of model parameters, partially mentioned throughout all chapters, and (spatial) basis functions for each type of model parameters, where wavelet/curvelet bases or kernel-driven bases might be used.