Section: New Results

Seismic Imaging and Inverse Problems

Shape-reconstruction and parameter identification of an elastic object immersed in a fluid

Participants : Izar Azpiroz Iragorri, Hélène Barucq, Julien Diaz.

We have developed a procedure to reconstruct the shape and material parameters of an elastic obstacle immersed in a fluid medium from some external measurements given by the so called far-field pattern. It is a nonlinear and ill-posed problem which is solved by applying a Newton-like iterative method involving the Fréchet derivatives of the scattered field. These derivatives express the sensitivity of the scattered field with respect to the parameters of interest. They are defined as the solution of boundary value problems which differ from the direct one only at the right-hand sides level. We have been able to establish the well-posedness of each problem in the case of a regular obstacle and it would be interesting in the near future to extend those results to the case of scatterers with polygonal boundaries. It requires to work with less regular Sobolev spaces for which the definition of traces is not obvious. We have also provided an analytical representation of the Fréchet derivatives in the case of a circle.

Next, we have introduced a series of numerical experiments that have been performed by applying two algorithms which propose two strategies of full reconstruction regarding the material parameters are retrieved simultaneously with the shape or not. It turns out that both work similarly delivering the same level of accuracy but the simultaneous reconstruction requires less iterations. We have thus opted for retrieving all the parameters simultaneously. Since realistic configurations include noisy data, we have performed some simulations for the reconstruction of the shape along with the Lamé coefficients for different noise levels. Other interesting experiments have been carried out using a multistage procedure where the parameters of interest are the density of the solid interior, the shape of the obstacle and its position. We have considered the case of Limited Aperture Data in back-scattering configurations, using multiple incident plane waves, mimicing a physical disposal of non-destructive testing.

We extended the solution methodology to the case of anisotropic media. Since the impact of some of the anisotropic parameters on the FFP is even weaker than the Lamé coefficients, the reconstruction of these parameters together with the shape parameters requires several frequencies and carefully adapted regularization parameters. It is in particular difficult to retrieve the Thomsen parameters $ϵ$ and $\delta$ because their reconstruction requires to have an accurate adjustment on the rest of material and shape parameters. The recovery process is thus computationally intensive and some efforts should be done in the near future to decrease the computational costs. We were able to recover all the anisotropic parameters when the shape were assumed to be known. However, when trying to recover both shape and material parameters, we could only recover the shape and some of the physical parameters (namely the three most important ones : the density and the two velocities $V_p$ and $V_s$ ).

These results have been obtained in collaboration with Rabia Djellouli (California State University at Northridge, USA) and are presented in Izar Azpiroz Ph.D thesis [1]

Time-harmonic seismic inverse problem with dual-sensors data

Participants : Hélène Barucq, Florian Faucher.

We study the inverse problem for the time-harmonic acoustic wave equation. The seismic context implies restrictive set of measurements: it consists of reflection data (resulting from an artificial source) acquired from the near surface area only. The inverse problem aims at recovering the subsurface medium parameters and we use the Full Waveform Inversion (FWI) method, which defines an iterative minimization algorithm of the difference between the measurement and simulation.

We investigate the use of new devices that have been introduced in the acoustic setting. They are able to capture both the pressure field and the vertical velocity of the waves and are called dual-sensors. For solving the inverse problem of interest, we define a new cost function, adapted to these two-components data. We first note that the stability of the problem can be shown to be Lipschitz, assuming the parameters to be piecewise linear.

The usefulness of the cost function is to allow a separation between the observational and numerical sources. Therefore, the numerical sources do not have to coincide with the observational ones, offering new possibilities to create adapted computational acquisitions, and possibilities to reduce the numerical burden. We illustrate our approach with three-dimensional medium reconstructions, where we start with minimal information on the target models.

This work is a collaboration with Giovanni Alessandrini (Università di Trieste), Maarten V. de Hoop (Rice University), Romina Gaburro (University of Limerick) and Eva Sincich (Università di Trieste). It has been presented in the GDR-Meca Wave conference [31].

Stability and convergence analysis for seismic depth imaging using Full Waveform Inversion

Participants : Hélène Barucq, Florian Faucher.

We study the convergence of the inverse problem associated with the time-harmonic wave equations. In the context of seismic, the inverse problem uses reflection data which can only be obtained from the near surface area. We consider the propagation of waves in a domain $\Omega$ and the forward problem is defined from the displacement vector field $u$, solution to

where $g$ stands for the source, $\rho$ is the density and $\sigma$ the stress tensor. The inverse problem aims at the recovery of the medium parameters (contained in the stress tensor) and can be solved with an iterative minimization algorithm: this is the Full Waveform Inversion (FWI) method. We study the convergence of the minimization by introducing the framework of Finite Curvature/Limited Deflection (FC/LD) problems. The idea is to obtain the FC/LD properties by restricting the model space to guarantee strictly quasiconvex attainable set. It allows us to numerically estimate the size of the basin of attraction depending on characteristics of the inverse problem such as the frequency or the geometry of the target. In particular, it allows a quantitative comprehension of frequency progression during the iterative scheme, which is an aspect that appeared mostly intuitive. It also allows a comparison of methods from a convergence point of view. This analysis is to relate with stability estimates in order to provide a consistent scheme where frequency progression is justified from the quantitative estimates. Eventually, we illustrate our approach with elastic medium reconstructions, starting from minimal information on the initial models; this also serves to illustrate the numerical requirement of the large scale optimization seismic experiments.

This work is a collaboration with Guy Chavent (Inria Rocquencourt) and Henri Calandra (TOTAL). The results have been presented in the conference “Reconstruction Methods for Inverse Problems” [15].

Quantitative localization of small obstacles with single-layer potential fast solvers

Participants : Hélène Barucq, Florian Faucher, Ha Howard Faucher.

In this work, we numerically study the inverse problem of locating small circular obstacles in a homogeneous medium using noisy backscattered data collected at several frequencies. The main novelty of our work is the implementation of a single-layer potential based fast solver (called FSSL) in a Full-Waveform inversion procedure, to give high quality reconstruction with low-time cost. The efficiency of FSSL was studied in our previous works. We show reconstruction results with up to 12 obstacles in structured or random configurations with several initial guesses, all allowed to be far and different in nature from the target. This last assumption is not expected in results using nonlinear optimization schemes in general. For results with 6 obstacles, we also investigate several optimization methods, comparing between nonlinear gradient descent and quasi-Newton, as well as their convergence with different line search algorithms.

The work is published in Journal of Computational Physics [9]. This work has been presented at GDR-Meca Wave conference in Fréjus cf. [21].

Time Domain Full Waveform Inversion Adjoint Studies

Participants : Hélène Barucq, Julien Diaz, Pierre Jacquet.

Full Waveform Inversion (FWI) allows retrieving the physical parameters (e.g. the velocity, the density) from an iterative procedure underlying a global optimization technique. The recovering of the medium corresponds to the minimum of a cost function quantifying the difference between experimental and numerical data. In this study we have considered the adjoint state method to compute the gradient of this cost function. At each iteration the parameters are updated with the solution of an adjoint equation which can be defined either as the adjoint of the continuous equation or the discrete problem. Some studies have addressed the question of establishing what the best strategy is. The answer is still unclear and turns out to be strongly dependent on the problem under study.

The purpose of this study was to investigate several computations of the adjoint state as a preamble of a FWI method applied to the time-dependent acoustic wave approximated in a Discontinuous Galerkin framework involving Bernstein elements. We have considered different time schemes to feature the inherited properties of the computed adjoint state. By comparing the different discrete adjoint operators both from a mathematical and numerical point of view, we aim at defining the best option for computing the adjoint state with accuracy at least cost.

This work is a collaboration with Henri Calandra (TOTAL). It was presented at Total MATHIAS conference in Paris [28].

Seismic imaging of remote targets buried in an unknown media

Participant : Yder Masson.

Box Tomography:first application to the imaging of upper-mantle shear velocity and radial anisotropy structure beneath the North American continent: The EarthScope Transpotable Array (TA) deployment provides dense array coverage through- out the continental United States and with it, the opportunity for high-resolution 3-D seismic velocity imaging of the stable part of the North American (NA) upper mantle. Building upon our previous long-period waveform tomographic modeling, we present a higher resolution 3-D isotropic and radially anisotropic shear wave velocity model of the NA lithosphere and asthenosphere. The model is constructed using a combination of teleseismic and regional waveforms down to 40 s period and wavefield computations are performed using the spectral element method both for regional and teleseismic data. Our study is the first tomographic ap- plication of ‘Box Tomography’, which allows us to include teleseismic events in our inversion, while computing the teleseismic wavefield only once, thus significantly reducing the numerical computational cost of several iterations of the regional inversion. We confirm the presence of high-velocity roots beneath the Archean part of the continent, reaching 200–250 km in some areas, however the thickness of these roots is not everywhere correlated to the crustal age of the corresponding cratonic province. In particular, the lithosphere is thick ( 250 km) in the western part of the Superior craton, while it is much thinner ( 150 km) in its eastern part. This may be related to a thermomechanical erosion of the cratonic root due to the passage of the NA plate over the Great Meteor hotspot during the opening of the Atlantic ocean 200–110 Ma. Below the lithosphere, an upper-mantle low-velocity zone (LVZ) is present everywhere under the NA continent, even under the thickest parts of the craton, although it is less developed there. The depth of the minimum in shear velocity has strong lateral variations, whereas the bottom of the LVZ is everywhere relatively flat around 270–300 km depth, with minor undulations of maximum 30 km that show upwarping under the thickest lithosphere and downwarping under tectonic regions, likely reflecting residual temperature anomalies. The radial anisotropy structure is less well resolved, but shows distinct signatures in highly deformed regions of the lithosphere.

This is the first application to a real case study of a novel imaging method called "Box Tomography". These results were obtained through collaborations with Barbara Romanowicz (Berkeley Seimological Laboratory, UC Berkeley; Collège de France) and Pierre Clouzet (Institut de Physique du Globe de Paris). The results have been published in the Geophysical Journal International [11].

Additional developments are conducted in collaboration with Sevan Adourian and Barbara Romanowicz at the Berkeley Seismological Laboratory, UC Berkeley, in particular, to efficiently account for receivers located outside the imaged region. These new results have been presented in different international conferences [23], [19].

To strengthen existing collaborations, a proposal has been submitted to the France-Berkeley Fund (13000$ for travelling and living expenses). We propose a joint effort to further develop and apply a novel seismic tomographic approach, Box-Tomography, to image and characterize small scale structures of interest in the deep Earth, such as the roots of mantle plumes, ultra-low velocity zones, or the edges of large low shear velocity provinces. Our objective is to forge a long-term collaboration between applied mathematicians at Magique3D developing wave propagation modeling methods and the seismologists at the Berkeley Seismological Laboratory (UC Berkeley) using these methods to investigate the Earth's internal structure.