Iterative inversion of seismic, ultrasonic, and other wave data by local gradient-based optimization of mean-square data prediction error (Full Waveform Inversion or FWI) can fail to converge to useful model estimates if started from an initial model predicting wave arrival times in error by more than half a wavelength (a phenomenon known as cycle skipping). Matched Source Waveform Inversion (MSWI) extends the wave propagation model by a filter that shifts predicted waves to fit observed data. The MSWI objective adds a penalty for deviation of this filter from the identity to the mean-square data misfit . The extension allows the inversion to make large model adjustments while maintaining data fit and so reduces the chances of local optimization iterates stagnating at non-informative model estimates. Theory suggests that MSWI applied to acoustic transmission data with single-arrival wavefronts may produce an estimate of refractive index similar to the result of travel time inversion, but without requiring explicit identification of travel times. Numerical experiments conform to this expectation, in that MSWI applied to single arrival transmission data gives reasonable model estimates in cases where FWI fails. This MSWI model can then be used to jumpstart FWI for further refinement of the model. The addition of moderate amounts of noise (30%) does not negatively impact MSWI's ability to converge. However, MSWI applied to data with multiple arrivals is no longer theoretically equivalent to travel-time tomography and exhibits the same tendency to cycle-skip as does FWI.

Professor Pierre Sabatier contributed much to the study of inverse problems in theory and practice. Two of these contributions were a focus on theory that actually supports practice, and the identification of well-posed aspects of inverse problems that may be quite ill-posed. This paper illustrates these two themes in the context of electrical impedance tomography, which is both very ill-posed and very practical. We show that for a highly constrained version of this inverse problem, in which a small elliptical inclusion in a homogeneous background is to be identified, applying optimal experiment design principles to choose electrode locations vastly improves the stability of the solution.

Adaptive Waveform Inversion applied to transient transmitted wave data yields estimates of index of refraction (or wave velocity) similar to those obtained by travel time inversion, provided that the data contain a single smooth wavefront.

Corrigendum: Solution of an acoustic transmission inverse problem by extended inversion (2022 Inverse Problems 38 115003) William W Symes, Huiyi Chen2,∗ and Susan E Minkoff 1 Department of Computational and Applied Mathematics, Rice University, Houston, TX 77005, United States of America 2 Department of Mathematical Sciences, FO 35, University of Texas at Dallas, Richardson, TX 75080, United States of America

Source extension is a reformulation of inverse problems in wave propagation, that at least in some cases leads to computationally tractable iterative solution methods. The core subproblem in all source extension methods is the solution of a linear inverse problem for a source (right hand side in a system of wave equations) through minimization of data error in the least squares sense with soft imposition of physical constraints on the source via an additive quadratic penalty. For an acoustic formulation for sources supported on a surface, with a soft contraint enforcing concentration at a point, a variant of the time reversal method from photoacoustic tomography provides an approximate solution. This approximate inverse can be used to precondition Krylov space iteration for rapid convergence to the solution of the core subproblem in this setting. Numerical examples illustrate the effectiveness of this preconditioner.

Study of a simple single-trace transmission example shows how an extended source formulation of full-waveform inversion can produce an optimization problem without spurious local minima (‘cycle skipping’), hence efficiently solvable via Newton-like local optimization methods. The data consist of a single trace extracted from a causal pressure field, propagating in a homogeneous fluid according to linear acoustics, and recorded at a given distance from a transient point energy source. The source intensity (‘wavelet’) is presumed quasi-impulsive, with zero energy for time lags greater than a specified maximum lag. The inverse problem is: from the recorded trace, recover both the sound velocity or slowness and source wavelet with specified support, so that the data is fit with prescribed RMS relative error. The least-squares objective function has multiple large residual minimizers. The extended inverse problem permits source energy to spread in time, and replaces the maximum lag constraint by a weighted quadratic penalty. A companion paper shows that for proper choice of weight operator, any stationary point of the extended objective produces a good approximation of the global minimizer of the least squares objective, with slowness error bounded by a multiple of the maximum lag and the assumed noise level. This paper summarizes the theory developed in the companion paper and presents numerical experiments demonstrating the accuracy of the predictions in concrete instances. We also show how to dynamically adjust the penalty scale during iterative optimization to improve the accuracy of the slowness estimate.

Modeled and recorded seismic data is contaminated by noise due to a variety of factors including geophones recording ambient noise unrelated to seismic exploration, equipment malfunction and limitations, and restrictive modeling assumptions. The level of noise in seismic data is not known a priori, hence being able to estimate the noise level in the data is a valuable tool for assessing the quality of mechanical Earth parameter estimates, especially resulting from inversion. Full waveform inversion (FWI) is a promising tool for estimating subsurface parameters such as wave velocity, but users of FWI encounter their own challenges. Specifically, commonly-used gradient based local optimization techniques are well known to stall in geologically uninformative models if the starting guess for the optimization is not close enough to the desired global optimum. Extension-based methods relax physical constraints on model parameters to enlarge the search space of acceptable solutions, helping to reduce the impact of a poor initial guess and potentially convexifying the objective function. These extended inversion methods involve a penalty term that is added to the FWI least squares misfit. Then the challenge of performing the optimization shifts to adjusting the penalty weight to balance the reduction of data misfit with driving the penalty term towards physically-meaningful solutions. The source extended objective function minimized using the discrepancy algorithm requires an estimate of the noise level in the data to proceed. In this work we illustrate an automated algorithm for simultaneously updating this noise estimate so that we bypass the cycle-skipping problem experienced by FWI and converge towards geologically meaningful parameter estimates.

A single-trace transmission inverse problem for the wave equation seeks to determine both the wave velocity in a homogenous acoustic medium and the transient waveform of an isotropic point source. The duration (support) of the source waveform and the source-to-receiver distance are assumed known. A least squares formulation of this problem exhibits the ``cycle-skipping'' behaviour observed in field scale problems of this type, with many local minima differing greatly from the global minimizer. This behaviour is eliminated by dropping the hard support constraint on the source waveform, replacing it by a soft penalty implemented as a weighted mean-square of the source waveform. For properly chosen weight function, penalizing nonzero values away from $t=0$, the velocity component of any stationary point differs from the global minimizer of the constrained least-squares formulation by a linear combination of the source waveform support radius and data noise-to-signal ratio. Given an estimate of data noise, the penalty weight can be dynamically adjusted during iterative optimization to maximize predicted data accuracy and closely approximate a support-constrained source waveform.

Accurate representation and estimation of seismic sources is crucial to the joint medium-source full waveform inversion problem. We focus on the source estimation subproblem and its difficulties, where seismic sources are modeled by truncated series of multipoles. The source full waveform inversion formulation results in a highly ill-conditioned (potentially ill-posed) linear least squares problem which we attempt to solve iteratively via conjugate gradient. Our main contribution lies in developing a preconditioner to accelerate the performance of conjugate gradient on multipole source inversion. The proposed preconditioner consists of (fractional) time derivative/integral operators based on analytical solutions to the wave equation with multipole sources in an unbounded, homogeneous medium. Numerical results in 2D demonstrate that the conjugate gradient iterations are accelerated when incorporating the proposed preconditioning scheme.

An extremely simple single-trace transmission example shows how an extended source formulation of full waveform inversion can produce an optimization problem without spurious local minima ("cycle skipping"). The data consist of a single trace recorded at a given distance from a point source. The velocity or slowness is presumed homogeneous, and the target source wavelet is presumed quasi-impulsive or focused at zero time lag. The source is extended by permitting energy to spread in time, and the spread is controlled by adding a weighted mean square of the extended source wavelet to the data misfit, to produce the extended inversion objective. The objective function and its gradient can be computed explicitly, and it is easily seen that all local minimizers must be within a wavelength of the correct slowness. The derivation shows several important features of all similar extended source algorithms. For example, nested optimization, with the source estimation in the inner optimization (variable projection method), is essential. The choice of the weight operator, controlling the extended source degrees of freedom, is critical: the choice presented here is a differential operator, and that property is crucial for production of an objective immune from cycle-skipping.

Nonlinear least squares data-fitting driven by physical process simulation is a classic and widely successful technique for the solution of inverse problems in science and engineering. Known as ‘full waveform inversion (FWI)’ in application to seismology, it can extract detailed maps of earth structure from near-surface seismic observations, but also suffers from a defect not always encountered in other applications: the least squares error function at the heart of this method tends to develop a high degree of nonconvexity, so that local optimization methods (the only numerical methods feasible for field-scale problems) may fail to produce geophysically useful final estimates of earth structure, unless provided with initial estimates of a quality not always available. A number of alternative optimization principles have been advanced that promise some degree of release from the multimodality of FWI, amongst them wavefield reconstruction inversion (WRI), the focus of this paper. Applied to a simple 1D acoustic transmission problem, both full waveform and WRI methods reduce to minimization of explicitly computable functions, in an asymptotic sense. The analysis presented here shows explicitly how multiple local minima arise in FWI, and that WRI can be vulnerable to the same ‘cycle-skipping’ failure mode.

Extended modeling is one of a number of modifications suggested to enhance the reliability of iterative FullWaveform Inversion (FWI), by making it less prone to stagnation away from useful model estimates ("cycle skipping"). All extended modeling methods add parameters to wave modeling beyond those suggested by basic wave physics, in some cases entailing substantial computational expense. Source extension methods add parameters to the description of the physical energy source; many of these cost little more computationally than standard FWI. The first purpose of this paper is to present a (necessarily incomplete) taxonomy of proposed source extensions and related inversion methods. The second purpose is to introduce a particular such method, the Surface Source Extension, and give a couple of examples of its use. We will also explain the theoretical justification for some of these methods: for pure transmission (including diving wave) data, a close link to travel time tomography explains the absence of cycle skipping observed in numerical experiments. No such theoretical link is known to exist in the case of pure reflection data, but some positive numerical results suggest that some source extension methods may ameliorate cycle-skipping in that case as well. Presentation Date: Monday, September 16, 2019 Session Start Time: 1:50 PM Presentation Start Time: 4:20 PM Location: 301B Presentation Type: Oral

The subsurface offset is a possible link connecting wave-equation migration methods, such as reverse-time migration, with the angle-domain. However, it describes an artificial reflection geometry that splits the reflection point between incident and scattered waves. The split configuration contradicts basic principles of continuum mechanics, as it represents action-at-a-distance, and has kinematics distinct from those of ordinary physical reflection. The angle-domain image, computed as a Radon transform of the subsurface offset image, implicitly inherits the split configuration, which distinguishes it from conventional angle-domain decomposition. As a consequence, it shows a dissimilar moveout response to erroneous migration velocity, but still carries valuable information about the migration velocity error. Conventional linearized traveltime inversion techniques, applied for velocity optimization, are most likely to fail while inverting angle-domain images produced from subsurface offset split reflections. These images should be labelled differently to highlight the split mechanism, and input to traveltime inversion algorithms only after a proper generalization. We present a modified linearized traveltime inversion formulation, appropriate for split reflections in the image space. The key ingredient in the algorithm is a subsurface offset depended extension of the traveltime change along migrated split reflection rays. A generalized reflection tomography scheme is proposed accordingly to tie imaging errors, in either the subsurface offset domain or its equivalent angle-domain map, to the migration velocity error.

Seismic sources are commonly idealized as point-sources due to their small spatial extent relative to seismic wavelengths. The acoustic isotropic point-radiator is inadequate as a model of seismic wave generation for seismic sources that are known to exhibit directivity. Therefore, accurate modeling of seismic wavefields must include source representations generating anisotropic radiation patterns. Such seismic sources can be modeled as multipoles, that is, a time-dependent linear combination of spatial derivatives of the spatial delta function. Since the solutions of linear hyperbolic systems with point-source right hand sides are necessarily singular, standard results on convergence of grid-based numerical methods (finite difference or finite element) do not imply convergence of numerical solutions. We present a method for discretizing multipole sources in a finite difference setting, an extension of the moment matching conditions developed for the Dirac delta function in other applications, along with numerical evidence demonstrating the accuracy of these approximations. Using this analysis, we develop a weak convergence theory for the discretization of a family of symmetric hyperbolic systems of first-order partial differential equations, with singular source terms, solved via staggered-grid finite difference methods: we show that grid-independent space-time averages of the numerical solutions converge to the same averages of the continuum solution, and provide an estimate for the error in terms of moment matching and truncation error conditions. Numerical experiments confirm this result, but also suggest a stronger one: optimal convergence rates appear to be achieved point-wise in space away from the source.

Seismic migration in the angle-domain generates multiple images of the earth's interior in which reflection takes place at different scattering-angles. Mechanically, the angle-dependent reflection is restricted to happen instantaneously and at a fixed point in space: Incident wave hits a discontinuity in the subsurface media and instantly generates a scattered wave at the same common point of interaction. Alternatively, the angle-domain image may be associated with space-shift (regarded as subsurface offset) extended migration that artificially splits the reflection geometry. Meaning that, incident and scattered waves interact at some offset distance. The geometric differences between the two approaches amount to a contradictory angle-domain behaviour, and unlike kinematic description. We present a phase space depiction of migration methods extended by the peculiar subsurface offset split and stress its profound dissimilarity. In spite of being in radical contradiction with the general physics, the subsurface offset reveals a link to some valuable angle-domain quantities, via post-migration transformations. The angle quantities are indicated by the direction normal to the subsurface offset extended image. They specifically define the local dip and scattering angles if the velocity at the split reflection coordinates is the same for incident and scattered wave pairs. Otherwise, the reflector normal is not a bisector of the opening angle, but of the corresponding slowness vectors. This evidence, together with the distinguished geometry configuration, fundamentally differentiates the angle-domain decomposition based on the subsurface offset split from the conventional decomposition at a common reflection point. An asymptotic simulation of angle-domain moveout curves in layered media exposes the notion of split versus common reflection point geometry. Traveltime inversion methods that involve the subsurface offset extended migration must accommodate the split geometry in the inversion scheme for a robust and successful convergence at the optimal velocity model.

Full-waveform inversion (FWI) faces the persistent challenge of cycle skipping, which can result in stagnation of the iterative methods at uninformative models with poor data fit. Extended reformulations of FWI avoid cycle skipping through adding auxiliary parameters to the model so that a good data fit can be maintained throughout the inversion process. The volume-based matched source waveform inversion algorithm introduces source parameters by relaxing the location constraint of source energy: It is permitted to spread in space, while being strictly localized at time [Formula: see text]. The extent of source energy spread is penalized by weighting the source energy with distance from the survey source location. For transmission data geometry (crosswell, diving wave, etc.) and transparent (nonreflecting) acoustic models, this penalty function is stable with respect to the data-frequency content, unlike the standard FWI objective. We conjecture that the penalty function is actually convex over much larger region in model space than is the FWI objective. Several synthetic examples support this conjecture and suggest that the theoretical limitation to pure transmission is not necessary: The inversion method can converge to a solution of the inverse problem in the absence of low-frequency data from an inaccurate initial velocity model even when reflections and refractions are present in the data along with transmitted energy.

Optimization-based migration velocity analysis updates long-wavelength velocity information by minimizing an objective function that measures the violation of a focusing criterion, applied to an image volume. Differential semblance optimization forms a smooth objective function in velocity and data, regardless of the data-frequency content. Depending on how the image volume is formed, however, the objective function may not be minimized at a kinematically correct velocity, a phenomenon characterized in the literature (somewhat inaccurately) as "gradient artifacts." We find that the root of this pathology is imperfect image volume formation resulting from reverse time migration (RTM), and that the use of linearized inversion (least-squares migration) more or less eliminates it. A synthetic Marmousi example and a 2D real data example are used to demonstrate that an approximate inverse operator, a little more expensive than RTM, leads to recovery of a kinematically correct velocity.

Summary Subsurface offset extended imaging provides a possible link connecting migrated reflections in CIGs with the migration velocity error. However, it describes an artificial scattering geometry that splits the reflection point between incident and scattered waves. The split configuration contradicts basic principles of continuum mechanics, as it represents action-at-a-distance, and has kinematics distinct from those of ordinary physical reflection. The angle-domain image, computed as a Radon transform of the subsurface offset image, implicitly inherits the split configuration and distinguished by its response to erroneous migration velocity. Conventional traveltime inversion techniques, applied for velocity optimization, are incompatible with the split configuration. Hence, they are most likely to fail while inverting traveltimes from the subsurface offset image or its angle-domain transformation. Nevertheless, those images that are being identified with the split geometry still carry valuable information about the migration velocity error. They should be input to traveltime inversion algorithms only after a proper generalization. We present a modified traveltime inversion formulation, appropriate for split reflections in the image space. The net result is a generalized reflection tomography scheme where the traveltime estimation accounts for a possible subsurface offset split.

Seismic migration may be altered to provide an asymptotic inversion of the seismic data. Proper weight operators scale the migration to become an approximate inverse of the Born modeling operator. We extend wave-equation based techniques developed for approximate inversion of constant density acoustics to variable density acoustics. Accordingly, seismic reflection data may be matched by a proper choice of extended bulk modulus and density perturbations. Estimation for those acoustic model perturbations is recovered from the angle-dependent response of Born’s approximate inverse application. We describe the essential steps for constructing the variable density acoustic approximate inverse, and demonstrate a multi-parameter inversion for bulk modulus and density perturbations. Presentation Date: Thursday, October 18, 2018 Start Time: 8:30:00 AM Location: 206A (Anaheim Convention Center) Presentation Type: Oral

Source signature estimation from seismic data is a crucial ingredient for successful application of seismic migration and full-waveform inversion (FWI). If the starting velocity deviates from the target velocity, FWI method with on-the-fly source estimation may fail due to the cycle-skipping problem. We have developed a source-based extended waveform inversion method, by introducing additional parameters in the source function, to solve the FWI problem without the source signature as a priori. Specifically, we allow the point source function to be dependent on spatial and time variables. In this way, we can easily construct an extended source function to fit the recorded data by solving a source matching subproblem; hence, it is less prone to cycle skipping. A novel source focusing annihilator, defined as the distance function from the real source position, is used for penalizing the defocused energy in the extended source function. A close data fit avoiding the cycle-skipping problem effectively makes the new method less likely to suffer from local minima, which does not require extreme low-frequency signals in the data. Numerical experiments confirm that our method can mitigate cycle skipping in FWI and is robust against random noise.