Christian Boehm

A Semismooth Newton-CG Method for Constrained Parameter Identification in Seismic Tomography

Seismic tomography is a technique to determine the material properties of the Earths subsurface based on the observation of seismograms. This can be stated as a PDE-constrained optimization problem governed by the elastic wave equation. We present a semismooth Newton-PCG method with a trust-region globalization for full-waveform seismic inversion that uses a Moreau--Yosida regularization to handle additional constraints on the material parameters. We establish results on the differentiability of the parameter-to-state operator and analyze the proposed optimization method in a function space setting. The elastic wave equation is discretized by a high-order continuous Galerkin method in space and an explicit Newmark time-stepping scheme. The matrix-free implementation relies on the adjoint-based computation of the gradient and Hessian-vector products and on an MPI-based parallelization. Numerical results are shown for an application in geophysical exploration at reservoir scale.

A Newton-CG Method for Full-Waveform Inversion in a Coupled Solid-Fluid System

We present a Newton-CG method for full-waveform seismic inversion. Our method comprises the adjoint-based computation of the gradient and Hessian-vector products of the reduced problem and a preconditioned conjugate gradient method to solve the Newton system in matrix-free fashion. A trust-region globalization strategy and a multi-frequency inversion approach are applied. The governing equations are given by a coupled system of the acoustic and the elastic wave equation for the numerical simulation of wave propagation in solid and fluid media. We show numerical results for the application of our method to marine geophysical exploration.

Ambient seismic source inversion in a heterogeneous Earth - Theory and application to the Earth's hum†

The sources of ambient seismic noise are extensively studied both to better understand their influence on ambient noise tomography and related techniques, and to infer constraints on their excitation mechanisms. Here we develop a gradient-based inversion method to infer the space-dependent and time-varying source power spectral density of the Earths hum from cross correlations of continuous seismic data. The precomputation of wavefields using spectral elements allows us to account for both finite-frequency sensitivity and for three-dimensional Earth structure. Although similar methods have been proposed previously, they have not yet been applied to data to the best of our knowledge. We apply this method to image the seasonally varying sources of Earths hum during North and South Hemisphere winter. The resulting models suggest that hum sources are localized, persistent features that occur at Pacific coasts or shelves and in the North Atlantic during North Hemisphere winter, as well as South Pacific coasts and several distinct locations in the Southern Ocean in South Hemisphere winter. The contribution of pelagic sources from the central North Pacific cannot be constrained. Besides improving the accuracy of noise source locations through the incorporation of finite-frequency effects and 3-D Earth structure, this method may be used in future cross-correlation waveform inversion studies to provide initial source models and source model updates.

Optimal experimental design to position transducers in ultrasound breast imaging

We present methods to optimize the setup of a 3D ultrasound tomography scanner for breast cancer detection. This approach provides a systematic and quantitative tool to evaluate different designs and to optimize the con- figuration with respect to predefined design parameters. We consider both, time-of-flight inversion using straight rays and time-domain waveform inversion governed by the acoustic wave equation for imaging the sound speed. In order to compare different designs, we measure their quality by extracting properties from the Hessian operator of the time-of-flight or waveform differences defined in the inverse problem, i.e., the second derivatives with respect to the sound speed. Spatial uncertainties and resolution can be related to the eigenvalues of the Hessian, which provide a good indication of the information contained in the data that is acquired with a given design. However, the complete spectrum is often prohibitively expensive to compute, thus suitable approximations have to be developed and analyzed. We use the trace of the Hessian operator as design criterion, which is equivalent to the sum of all eigenvalues and requires less computational effort. In addition, we suggest to take advantage of the spatial symmetry to extrapolate the 3D experimental design from a set of 2D configurations. In order to maximize the quality criterion, we use a genetic algorithm to explore the space of possible design configurations. Numerical results show that the proposed strategies are capable of improving an initial configuration with uniformly distributed transducers, clustering them around regions with poor illumination and improving the ray coverage of the domain of interest.

Optimized transducer configuration for ultrasound waveform tomography in breast cancer detection

Waveform inversion is a promising method for ultrasound computed tomography able to produce high-resolution images of human breast tissue. However, the computational complexity of waveform inversion remains a considerable challenge, and the costs per iteration are proportional to the number of emitting transducers. We propose a twofold strategy to accelerate the time-to-solution by identifying the optimal number and location of emitters using sequential optimal experimental design (SOED). SOED is a powerful tool to iteratively add the most informative transducer or remove redundant measurements, respectively. This approach simultaneously provides optimized transducer configurations and a cost-benefit curve that quantifies the information gain versus the computational cost. First, we propose a method to identify the emitters that provide reconstructions with minimal expected uncertainties. Using a Bayesian approach, model uncertainties and resolution can be quantified with the trace of the posterior covariance. By linearizing the wave equation, we can compute the posterior covariance using the inverse of the Gauss-Newton approximation of the Hessian. Furthermore, this posterior is independent of the breast model and the experimental data, thus enabling pre-acquisition experimental optimization. Then, for the post-acquisition inversion, we present an approach to select a subsample of sources that accurately approximates the full gradient direction in each iteration. We control the convergence of the angular differences between consecutive gradient directions by randomly adding new emitters into the subsample. We present synthetic studies in 2D and 3D that consider a ring-shaped and a semi-ellipsoidal scanning device, respectively. Numerical results suggest that the provided methods have the potential to identify redundancies from the corresponding cost-benefit curves. Furthermore, the gradient direction rapidly converges to the direction of the full gradient, which appears to be independent of the model and the emitter locations.

Time-domain spectral-element ultrasound waveform tomography using a stochastic quasi-Newton method

Waveform inversion for ultrasound computed tomography (USCT) is a promising imaging technique for breast cancer screening. However, the improved spatial resolution and the ability to constrain multiple parameters simultaneously demand substantial computational resources for the recurring simulations of the wave equation. Hence, it is crucial to use fast and accurate methods for numerical wave propagation, on the one hand, and to keep the number of required simulations as small as possible, on the other hand. We present an efficient strategy for acoustic waveform inversion that combines (i) a spectral-element continuous Galerkin method for solving the wave equation, (ii) conforming hexahedral mesh generation to discretize the scanning device, (iii) a randomized descent method based on mini-batches to reduce the computational cost for misfit and gradient computations, and (iv) a trust-region method using a quasi-Newton approximation of the Hessian to iteratively solve the inverse problem. This approach combines ideas and state-of-the-art methods from global-scale seismology, large-scale nonlinear optimization, and machine learning. Numerical examples for a synthetic phantom demonstrate the efficiency of the discretization, the effectiveness of the mini-batch approximation and the robustness of the trust-region method to reconstruct the acoustic properties of breast tissue with partial information.

Automated Large‐Scale Full Seismic Waveform Inversion for North America and the North Atlantic

We present a new anisotropic seismic tomography model based on a multiscale full seismic waveform inversion for crustal and upper-mantle structure from the western edge of North America across the North Atlantic and into Europe. The gradient-based inversion strategy utilizes the adjoint state method coupled with an L-BFGS quasi-Newton optimization scheme. To improve the handling of large data quantities in the context of full seismic waveform inversions, we developed a workflow framework automating the procedure across all stages, enabling us to confidently invert for waveforms from 72 events recorded at 7,737 unique stations, resulting in a total of 144,693 raypaths, most of them with three-component recordings. The final model after 20 iterations is able to explain complete waveforms including body as well as surface waves of earthquakes that were not used in the inversion down to periods of around 30s. The model is complemented by a detailed resolution analysis in the form of 3-D distributions of direction-dependent resolution lengths.

Towards Full Waveform Ambient Noise Inversion

We develop a method for the joint inversion of noise correlation functions for the distribution of noise sources and for Earth structure. The forward problem is free of assumptions required to equate noise correlations with Green functions and allows us to compute inter-station correlations for arbitrary distributions of noise sources in space and time. Using adjoint techniques, we design an iterative inversion scheme for noise sources and Earth structure based on waveform and energy differences as misfit functional. Starting from an initial model from a wave equation traveltime inversion, we recover the target velocity model with high accuracy. A key prerequisite is a good inference of the noise source distribution.

Automatic Global Multiscale Seismic Inversion: Insights into Model, Data, and Workflow Management

Modern global seismic waveform tomography is formulated as a PDE-constrained nonlinear optimization problem, where the optimization variables are Earths visco-elastic parameters. This particular problem has several defining characteristics. First, the solution to the forward problem, which involves the numerical solution of the elastic wave equation over continental to global scales, is computationally expensive. Second, the determinedness of the inverse problem varies dramatically as a function of data coverage. This is chiefly due to the uneven distribution of earthquake sources and seismometers, which in turn results in an uneven sampling of the parameter space. Third, the seismic wavefield depends nonlinearly on the Earths structure. Sections of a seismogram which are close in time may be sensitive to structure greatly separated in space. In addition to these theoretical difficulties, the seismic imaging community faces additional issues which are common across HPC applications. These include the storage of massive checkpoint files, the recovery from generic system failures, and the management of complex workflows, among others. While the community has access to solvers which can harness modern heterogeneous computing architectures, the computational bottleneck has fallen to these memory- and manpower-bounded issues. We present a two-tiered solution to the above problems. To deal with the problems relating to computational expense, data coverage, and the increasing nonlinearity of waveform tomography with scale, we present the Collaborative Seismic Earth Model (CSEM). This model, and its associated framework, takes an open-source approach to global-scale seismic inversion. Instead of attempting to monolithically invert all available seismic data, the CSEM approach focuses on the inversion of specific geographic subregions, and then consistently integrates these subregions via a common computational framework. To deal with the workflow and storage issues, we present a suite of workflow management software, along with a custom designed optimization and data compression library. It is the goal of this paper to synthesize these above concepts, originally developed in isolation, into components of an automatic global-scale seismic inversion.