# High-Resolution Seismic Imaging for Dam-Rock Interface using Full-Waveform Inversion

Mohamed Aziz Boukraa, Lorenzo Audibert, Marcella Bonazzoli, Houssem Haddar, Denis Vautrin

► To cite this version:

Mohamed Aziz Boukraa, Lorenzo Audibert, Marcella Bonazzoli, Houssem Haddar, Denis Vautrin. High-Resolution Seismic Imaging for Dam-Rock Interface using Full-Waveform Inversion. WAVES 2024 - 16th International Conference on Mathematical and Numerical Aspects of Wave Propagation, Jun 2024, Berlin, Germany. 10.17617/3.MBE4AA. hal-04661907

HAL Id: hal-04661907

https://hal.science/hal-04661907v1

Submitted on 25 Jul 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# High-Resolution Seismic Imaging for Dam-Rock Interface using Full-Waveform Inversion

Mohamed Aziz Boukraa¹*, Lorenzo Audibert², Marcella Bonazzoli¹, Houssem Haddar¹, Denis Vautrin²

¹UMA, INRIA, ENSTA Paris, Institut Polytechnique de Paris, 91120 Palaiseau, France

²EDF R&D, département PRISME, bâtiment U, 6 quai Watier, BP 49, 78401 Chatou Cedex, France

*Email: mohamed.boukraa@inria.fr

## Abstract

We are interested in imaging the interface between the concrete structure of the hydroelectric dam and the underlying rock using non-destructive seismic waves. Our method, combining "Full Waveform Inversion" with shape optimization, produces high-resolution images. Numerical results from realistic synthetic data demonstrate the method ability to accurately recover the interface with a limited number of measurements.

**Keywords:** Geophysics, subsurface imaging, seismic waves, Full Waveform Inversion.

## 1 Introduction

Seismic imaging consists in creating detailed images of subsurface structures using seismic waves. In our case, these waves are generated within the dam structure, either on its surface or through an immersed source for non-empty reservoirs, and then propagate through the ground. Receivers, including geophones for emerged locations or hydrophones submerged underwater, capture the reflections and refractions of these waves as they encounter various materials.

Processing these measurements to reconstruct the dam-rock interface is the objective of this study. Previous works employing geotechnical techniques, such as "Seismic Reflection" and "Vertical Seismic Profile," have not yielded satisfactory results. In this study, we propose leveraging the geophysical properties of these measurements by simulating the propagation of seismic waves to achieve a high-resolution image.

This study concerns gravity dams known for having a nearly invariant shape along their longitudinal direction. Hence, our study is limited to two dimensions. The computational domain consists of a vertical section of the dam and the surrounding media (refer to Fig. 1). Both dam and rock materials are treated as isotropic and homogeneous, with dissipation.

<figure><img src="image_2.png" /><figcaption>Figure 1: 2D schematic diagram for the considered excitations of the dam structure, the water and the rock foundation.</figcaption></figure>

## 2 Problem setting

The responses to excitation are described as solutions to frequency-domain wave propagation problems, involving the linear elasticity problem in the solid medium and the acoustics problem in the water medium, with appropriate transmission boundary conditions at their interface.

$$ 
\left\{
\begin{array}{ll}
\nabla \cdot \underline{\underline{\sigma}}(\vec{u}) + \rho\omega^2\vec{u} = 0 & \text{in } \Omega_s = \Omega_{s+} \cup \Omega_{s-} \\
\Delta p + k^2 p = g_f & \text{in } \Omega_f \\
\rho_f\omega^2\vec{u} \cdot \vec{n} = \nabla p \cdot \vec{n} & \text{on } \Gamma_{sf} \\
\underline{\underline{\sigma}}(\vec{u})\vec{n} = -p\vec{n} & \text{on } \Gamma_{sf} \\
\underline{\underline{\sigma}}(\vec{u})\vec{n} = \vec{g}_s & \text{on } \Gamma_s \\
\underline{\underline{\sigma}}(\vec{u})\vec{n} = 0 & \text{on } \Gamma_{s0} \\
\nabla p \cdot \vec{n} = 0 & \text{on } \Gamma_{f0} \\
\vec{u}|_{\Omega_{s-}} & \text{is outgoing} \\
p|_{\Omega_f} & \text{is outgoing}
\end{array}
\right.
\tag{1}
$$

The term $\vec{g}_s$ represents an emerged source on the wall $\Gamma_s$ and $g_f$ an immersed source in the water $\Omega_f$. On-site, only one source is active at a time. The outgoing condition is modeled using *Perfectly Matched Layers* technique.
### 3 The inverse problem

The proposed imaging strategy is based on Full Waveform Inversion (FWI), a quantitative reconstruction method that involves iteratively adjusting an initial model by aligning corresponding synthetic measurements obtained using a model similar to (1) with the seismic data.

Given $M$ source terms $g_1, \dots, g_M$, we can formulate the inverse problem of determining the interface $\Sigma_{ex}$ from (potentially noisy) measurements:

$$ \vec{u}_{\text{mes}}^{1} \approx \vec{u}_{1}(\Sigma_{\text{ex}})|_{\Gamma_{\text{ms}}}, \dots, \vec{u}_{\text{mes}}^{M} \approx \vec{u}_{M}(\Sigma_{\text{ex}})|_{\Gamma_{\text{ms}}}, $$

$$ p_{\text{mes}}^{1} \approx p_{1}(\Sigma_{\text{ex}})|_{\Gamma_{\text{mf}}}, \dots, p_{\text{mes}}^{M} \approx p_{M}(\Sigma_{\text{ex}})|_{\Gamma_{\text{mf}}}, $$

where $\Gamma_{ms} \subset \Gamma_s$ and $\Gamma_{mf} \subset \Gamma_{sf}$ are the measurements locations and $\mathcal{U}_{admissible}$ is a set of admissible shapes, i.e.,

$$ \left\{ \begin{array}{l} \text{Find } \Sigma_{\text{opt}} \in \mathcal{U}_{\text{admissible}} \text{ such that} \\ \Sigma_{\text{opt}} = \underset{\Sigma \in \mathcal{U}_{\text{admissible}}}{\operatorname{argmin}} J(\Sigma) \end{array} \right. $$

where the cost functional incorporates misfit terms denoted by $J_d$ in a least squares sense, along with a perimeter regularization term weighted by a regularization parameter $\eta \ge 0$:

$$ J(\Sigma) = \frac{1}{2} \sum_{i=0}^{M} \| \vec{u}^{i}(\Sigma) |_{\Gamma_{ms}} - u_{\text{mes}}^{i} \|_{L^2(\Gamma_{ms})}^2 \\ + \frac{1}{2} \sum_{i=0}^{M} \| p^{i}(\Sigma) |_{\Gamma_{mf}} - p_{\text{mes}}^{i} \|_{L^2(\Gamma_{mf})}^2 + \eta |\Sigma|. $$

The minimization process is conducted using the gradient descent algorithm. To compute the gradient, we employ a shape derivative approach, enabling us to express the derivative of the misfit terms $J_d$ with respect to the interface $\Sigma$ by solving both direct and adjoint problems. This derivative can be expressed by the following integral:

$$ J'_{d}(\Sigma)(\vec{h}) = \int_{\Sigma} \mathrm{Re}(j(\Sigma))\vec{h} \cdot \vec{n} \, ds, $$

for a given direction $\vec{h}$ and $j(\Sigma)$ is a function of the direct and adjoint states and their derivatives.

### 4 Inversion scheme

In our inversion scheme, we assume that the interface can be characterized by its height, i.e.,

the coordinates of a two-dimensional curve. This assumption is supported by engineers and specialists due to its realistic nature. The set of admissible shapes can be expressed as following:

$$ \mathcal{U}_{\text{admissible}} = \left\{ \Sigma = \{(x, f(x)), \ x \in I, \ f = \sum_{i=1}^{n_{\Sigma}} \xi_i \hat{\phi}_i \}, \right\} $$

for some interpolation functions $\hat{\phi} = \{\hat{\phi}_i\}_{1 \le i \le n_{\Sigma}}$.

Hence, $J(\Sigma) \equiv J(\Xi)$ with $\Xi = (\xi_1, \dots, \xi_{n_{\Sigma}})$. For $\Xi^0 \in \mathbb{R}^{n_{\Sigma}}$, the gradient descent algorithm can be written as:

$$ \Xi^{k+1} = \Xi^k - \tau \left( \nabla_{\Xi} J(\Xi^k) \right), $$

where $\tau > 0$ is the descent step. The components of the gradient of the misfit terms are $\frac{\partial}{\partial \xi_i} J_d(\Xi) = \int_{\Sigma} \mathrm{Re}(j(\Sigma)) \hat{\phi}_i (\vec{n} \cdot \vec{e}_y) ds$, where $\vec{n}$ is the normal vector to the interface and $\vec{e}_y$ represents the vertical component of $\mathbb{R}^2$.

### 5 Numerical results

We used synthetic measurements generated by Code\_Aster, a reliable engineering software developed by EDF, to validate our imaging method. We show the achievement of high-resolution imaging for two sets of data with different interfaces.

<figure><img src="image_3.png" /><figcaption>Figure 2: Illustration of sources and receivers configuration, the initialization and the considered interface geometries.</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>Figure 3: Reconstruction results for the two interface configurations, using frequencies: 100 Hz - 1000 Hz.</figcaption></figure>