# An exact inversion formula for determining a planar fault from boundary measurements

Hui Duong Bui, Andrei Constantinescu, H. Maigre

► To cite this version:

Hui Duong Bui, Andrei Constantinescu, H. Maigre. An exact inversion formula for determining a planar fault from boundary measurements. Journal of Inverse and Ill-posed Problems, 2005, 13 (6), pp.553-565. 10.1515/156939405775199514. hal-00167753

HAL Id: hal-00167753

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

Submitted on 23 Nov 2020

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 CC BY 4.0 - Attribution - International License
# An exact inversion formula from determining a planar fault from boundary measurements

H. D. BUI*, A. CONSTANTINESCU*, and H. MAIGRE*

**Abstract** — This paper considers the inverse problem of determining a time-dependent slip fault which releases shear stress elastically. The input data are the accelerations measured on the free external surface. A new formulae for determining explicitly the geometry of the planar fault is proposed. Potential applications to real earthquakes are discussed.

## 1. INTRODUCTION

Practical applications of acoustics and elastodynamics in geophysics are classically based on several techniques: i) mechanical set-ups such as the transducers and the accelerometers which collect acceleration histories data on the ground; ii) numerical methods such as Finite Elements Method, Boundary Integral Equations, Optimal Control theory, which search the best fitting of models predictions to experimental data.

Exact solutions of inverse scattering problems exist in acoustics under some assumptions: incident plane waves, Born or Kirchhoff approximations, Achenbach et al (1987), Bui (1994). Analytical solutions of the inverse acoustic scattering of a rigid inclusion has been given by Bojarski (1981). The classical “Elastic Wave Tomography” using the Kirchhoff approximation, the far-field analysis and the assumption of smallness of the defects, leads to the so-called POFFIS method (Physical Optics Far Field Inverse Scattering).

Recently, a different approach has been proposed in a series of papers giving the exact solution of the inverse problems for determining a planar crack from boundary measurements. The first paper is devoted to planar crack for

<footnote>*Laboratoire de Mécanique des Solides (CNRS UMR 7649), Ecole Polytechnique, 91128 Palaiseau cedex, France. E-mails: bui@lms.polytechnique.fr, constant@lms.polytechnique.fr, maigre@lms.polytechnique.fr</footnote>
the harmonic equation, Andrieux and Ben Abda (1992), Bannour et al (1997). Later, results are presented for quasi-static elasticity, Andrieux et al (1999), heat diffusion equation, Ben Abda and Bui (2003) and transient scalar acoustics, Bui et al (1999). A review of these solutions can be found in the paper by Bui et al (2003). These works are based on the use of the reciprocity gap functional (RGF). Mathematically, the RGF provides information on the crack geometry in the same way as an unknown function is determined by its components with respect to a complete basis of dual or adjoint functions. The key point is that a suitable choice of the adjoint basis of functions, allows these components to be determined by experiments since they depend explicitly on boundary data. The RGF provides therefore a simple method for deriving exact solutions. Similar techniques can be found in dynamic fracture mechanics, for determining the stress-intensity factors directly from boundary data, without solving the forward initial boundary value problem, Bui and Maigre (1988).

Mechanically, the reciprocity gap functional is simply the principle of virtual work known as the reciprocal relationship in linear materials. Recently a new class of adjoint functions conducting to an instantaneous RGF, is introduced by Bui et al (2004) for determining the “arrival time” of plane waves touching the planar crack. The “instantaneous” reciprocity gap functional RGF determines the crack plane and also the crack geometry in the case of a stationary convex shaped crack. The “arrival time” method determines the crack plane and geometry as the smallest possible convex hull containing the crack.

In this paper, we assume that the crack plane has already been determined, for example by the instantaneous RGF, as proposed in Bui et al (2004) and we will show how the time dependent geometry of the crack (or fault) $\Sigma(t)$ can be recovered in a closed form from boundary measurements. Explicit dependence on time of the geometry $\Sigma(t)$ is particularly important for potential applications in seismology since it is precisely the case for real earthquake due to release of shear stress in the propagation of rupture. In classical approaches by optimization, it is difficult to deal with defect propagation. Regularization techniques require many constraints to be imposed on the numerical solution, considered as penalties to be added to the cost functional. These constraints are as follows: positive slip rate, constant slip rate, bounded velocity of the rupture front, etc. Otherwise, there is no stable numerical solution, Das and Suhadolc (1996), Vallée (2003). One important conclusion in the first reference is that, even if the fitting of data seems to be quite good “the faulting process is poorly reproduced, so that in the real case, it would be difficult to know when one has obtained the correct solution”. The reasons for that may be twofold: flat minimum of the functional to be minimized, existence of local minima. Therefore, it is of great interest to derive a closed form formulae, which determines clearly the history of the rupture from boundary measurements. Of course, the quality of the solution using the proposed formulae depends always on the quality and richness of data. The aim of this paper is also to point out the necessary and sufficient conditions for obtaining the exact solution of the seismic inverse problem for exact data.
The paper is organized as follows. In Section 2.1 we set the equations for the inverse problem in consideration. For self-consistency of the paper, we recall in Section 2.2 the method of instantaneous reciprocity gap functional introduced in Bui et al (2004). Section 3 presents the more general RGF for solving a class of inverse problems. Section 4 presents two methods of solutions, the first one based on the P-wave adjoint fields, the second one based on the S-wave adjoint fields. We do not consider numerical analysis of the proposed method which are beyond the scope of this paper.

## 2. THE FAULT IDENTIFICATION PROBLEM: THE MODEL PROBLEM

Let $\Omega$ be a bounded solid subjected to an internal stress distribution $\sigma^0(x)$ in self-equilibrium in the absence external forces (see Figure 1). At $t \ge 0$, a planar fault $\Sigma(t)$ appears suddenly and starts propagating.

Let us denote by $u(x, t)$ and $\sigma(x, t)$ the displacement and respectively the stress variation between the initial time $t = 0$ and the actual instant $t$.

The appearance of the planar fault $\Sigma(t)$ for $t > 0$, releases an amount of stress characterized by the traction vector $\sigma(x, t) \cdot n$ on the faces of the fault. Where $n = n^\pm$ is the unit normal to the faces of the fault.

If we suppose that the fault does not open, then the movement at the faces of the fault corresponds to an amount of tangential slip $[u(x, t)] \cdot \tau$, where $\tau$ is an unit vector parallel to the fault plane $\Pi$ and $[u(x, t)]$ denotes the displacement jump at the crack. We recall that the displacement $u(x, t)$ field is actually a variation between the initial and actual state and this variation is purely elastic.

In the inverse problem defined next, we shall assume that the vector field of accelerations $\ddot{u}(x, t)$ is measured on the external boundary of the body $S = \partial\Omega(t) \setminus \Sigma(t)$, for any time $t > 0$.

The inverse problem consists in finding the amount of slip, the geometry of the propagating fault $\Sigma(t)$, the stress release etc.

The displacement field $u(x, t)$ satisfies the system of elastodynamic equations:

$$ \mathrm{div}(\sigma[u]) - \rho\ddot{u} = 0 \quad \text{in } (\Omega \setminus \Sigma) \times [0, \infty), \quad (1) $$

$$ \sigma[u] = A : \nabla u \quad \text{in } (\Omega \setminus \Sigma) \times [0, \infty), \quad (2) $$

$$ \sigma[u] \cdot n = 0 \quad \text{on } S \times [0, \infty), \quad (3) $$

$$ [u(x, t)] \cdot n = 0 \quad \text{on } S \times [0, \infty), \quad (4) $$

where $A$ is the positive definite symmetric tensor of the elastic moduli.

It is important to notice that the main unknown of the problem is the tangential slip $[u(x, t)] \cdot \tau$. Once the slip has been determined, one can remark from that the system of equations (1), (2) together with the boundary conditions (3), (4) and the complete set of initial conditions:

$$ u(x, 0) = 0, \quad \dot{u}(x, 0) = 0 \quad \text{in } \Omega \quad (5) $$
<figure><img src="image_2.png" /><figcaption>Figure 1. The planar slip fault $\Sigma(t)$ in a bounded solid $\Omega$. The fault plane $\Pi$ is assumed to be known. Measured acceleration data on the free surface $S$</figcaption></figure>

that the problem transforms in a well-posed problem. The solution provides
then information on the stress release on the fault.

Let us recall the main input data of the inverse problem discussed next:
i) $\ddot{u}(\boldsymbol{x}, t)$ the acceleration field on the free surface (3) and ii) the a priori knowledge of a zero normal jump (4) on the fault.

Moreover we have to assume that the pair of boundary values on $S$, the null
stress (3) and the acceleration field $\ddot{u}(\boldsymbol{x}, t)$ are a compatible data pair, in the
sense that they correspond to the same slip history on the fault.

This is an inverse problem, similar to the ones formulated in geophysics
for the determination of the slip on at the epicenter of the earthquakes. The
problem is generally ill-posed and therefore it is of interest to derive an closed
form solution which links the time-dependent fault geometry $\Sigma(t)$ to the given
data.

3. THE ARRIVAL TIME METHOD
FOR THE PLANE DETERMINATION

Before explaining the identification $\Pi$ of the plane of the fault by the arrival time
method, let us indicate briefly how this technique is used in the identification
of defects in classical acoustic scatterings.

Let us consider a convex defect $D$ situate inside the solid. At a point on the
solid boundary, a device called “emitter-receiver” sends a traveling acoustic
signal inside $\Omega$. If we denote the arrival time at a point on the defect by $T$
then the reflected signal is recorded at the receiver at time $2T$. If the shortest
delay time between emission and recording is $2T$, then $T$ corresponds to the
travel time from the point on the boundary to the point on the defect, the path
being an orthogonal projection on the surface of the defect. The complete defect
can be reconstructed from a series of measurements obtained by changing the
position of the “emitter-receiver” over the boundary.
Let us now consider, that the defect is capable of emitting a wave, for example the seismic wave produced by the slip on the fault at the epicenter. Using a similar method, one can reconstruct the position of the epicenter using different receivers on S determine a convex hull containing the defect.

Next we will use a similar technique, but we shall enhance the efficiency of the method by employing the instantaneous reciprocity gap functional (IRGF) introduced in Bui et al (2004). The IRGF is defined as:

$$ \mathcal{R}[u, v](t) = \int_{\Sigma(t)} n \cdot \sigma[v] \cdot [u] ds \quad (6) $$

$$ = \int_S \{n \cdot \sigma[v] \cdot u - n \cdot \sigma[u] \cdot v\} ds \quad (7) $$

where $u(x, t)$ is the current elastodynamic field and $v$ is the adjoint wave in the uncracked body. In (6), the second integral over S denoted by $\mathcal{R}[u, v]$ is a known linear functional of $v$.

As adjoined waves, one can choose plane waves of the form:

$$ v(x, t) = kY(t - \frac{1}{c_2} x \cdot p - \tau) \quad (8) $$

which propagates in the $p$ direction, with the shear wave velocity $c_2$ and $\tau$ a parameter chosen in such a way that waves are outside the solid $\Omega$ at $t = 0$. Plane waves can exist in a bounded solid if the corresponding traction field $\sigma[v] \cdot n$ is applied to its boundary, i. e. if the boundary tractions are equal to the restriction of the planar 3D wave to boundary. In (8), the Heaviside function is denoted by $Y(\cdot)$.

For the $p$ wave, the corresponding stress field $\sigma[v]$ becomes a travelling Dirac delta impulse which does not produce any virtual work with the crack displacement $[u]$ until the adjoint field $v(x, t)$ interacts with the crack $\Sigma(t)$.

The IRGF has two expressions, one depending on the unknown field $[u]$ and the chosen adjoint field $v$, defined on the fault surface $\Sigma$ and another one depending on the measured boundary data and $v$ on the external surface S.

As a function of $v$, the RGF depends implicitly of the parameter $p$: $R(t) = R(t; p)$. It vanishes identically for $t < T(p)$ and begins to take a non zero value after the arrival time $T(p)$. The arrival times $T(p)$ and $T(-p)$, for all directions $p$, determine the crack if the apparatus is sensitive enough to detect the true starting time $T$ such that $R(T^-) = 0$, $|\partial_t R(T^+)| > 0$ or $|R(T^+)| > 0$.

A numerical computation of $R(t; p)$ for a numerical simulation of 2D crack which suddenly releases shear stress in shown in Figure 2, where the arrival times are $T(p)$ and $T(-p)$ for waves of opposite directions. In general, in real experiments or numerical simulations, a small value of $R(t)$ is observed for the arrival time, $|R(t)| = \varepsilon$. For small $\varepsilon$, the plot of different plane waves at their arrival times shows a flat convex hull containing the crack (see Figure 3). The mathematical instantaneous RGF acts as the physical emitter-receiver device for measuring time arrival and therefore the instantaneous RGF method does not work for determining defects with a concave shape.
<figure><img src="image_3.png" /><figcaption>Figure 2. The instantaneous RGF for the identification of a 2D crack. The time arrivals of waves $T(p)$ and $T(-p)$</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>Figure 3. The convex hull of a 2D crack determined by the instantaneous RGF. The bold straight line is the real crack</figcaption></figure>

The determination of the fault plane in 3D is very simple. Assume that all adjoint waves at $t = 0$ are at the same distance $\Delta$ from the origin.

The function

$$F(p) = 2\Delta - c_2(T(p) + T(-p)) \quad (9)$$

measures the distance between the locations of the fronts of the adjoint waves given by $p$ and $-p$ at their respective arrival times: $T(p)$ and $T(-p)$. The solution $p_0$ of $F(p) = 0$, determines the normal to the fault plane $n = \pm p_0$ and consequently the location of the plane fault itself.

### 3.1. The elastodynamic RGF and the observation equation

Let us assume that the current elastodynamic field $u(x, t)$ satisfies conditions (1)–(5) defined in Section 1.2. The adjoint field $v(x, t)$ is assumed to satisfy the same wave equation in infinite medium containing no fault inside:

$$\mathrm{div}(\sigma[v]) - \rho\ddot{v} = 0 \quad \text{in } \mathbb{R}^3 \times [0, \infty) \quad (10)$$

The boundary values of $\sigma[v] \cdot n$ on $S$ and $\Sigma$ are simply computed from the field $v(x, t)$.
The reciprocal relationship takes the form

$$ \int_0^\infty \int_\Sigma \llbracket \mathbf{u} \rrbracket \cdot \boldsymbol{\sigma}[\mathbf{v}] \cdot \mathbf{n} \, ds \, dt = \int_0^\infty \int_S \{ \mathbf{n} \cdot \boldsymbol{\sigma}[\mathbf{v}] \cdot \mathbf{u} - \mathbf{n} \cdot \boldsymbol{\sigma}[\mathbf{u}] \cdot \mathbf{v} \} \, ds \, dt \\ + \int_\Omega \{ \mathbf{v} \cdot \ddot{\mathbf{u}} - \mathbf{u} \cdot \ddot{\mathbf{v}} \}^\infty_0 \, dv. \quad (11) $$

The last term vanishes if we choose the adjoint field such that: $v(t = \infty) = 0$, $v(t = 0) = 0$:

$$ \int_{\Omega} [\mathbf{v} \cdot \ddot{\mathbf{u}} - \mathbf{u} \cdot \ddot{\mathbf{v}}]^{\infty}_{0} \, dv. $$

Therefore, the elastodynamic RGF defined by the right hand side of (6) is equal to the mechanical work done by the adjoint stress field $\boldsymbol{\sigma}[\mathbf{v}]$ on the jump of the crack opening displacement $\llbracket \mathbf{u} \rrbracket$:

$$ \int_0^\infty \int_\Sigma \llbracket \mathbf{u} \rrbracket \cdot \boldsymbol{\sigma}[\mathbf{v}] \cdot \mathbf{n} \, ds \, dt = \int_0^\infty \int_S \{ \mathbf{n} \cdot \boldsymbol{\sigma}[\mathbf{v}] \cdot \mathbf{u} - \mathbf{n} \cdot \boldsymbol{\sigma}[\mathbf{u}] \cdot \mathbf{v} \} \, ds \, dt. \quad (12) $$

This equation can be interpreted as the “observation” equation for the unknown $\llbracket \mathbf{u} \rrbracket$.

Let us remark that the last equation is of the general form $a(\mathbf{u}, \mathbf{v}) = b(\mathbf{v})$ for a linear system, with a bilinear form $a(\mathbf{u}, \mathbf{v})$ and with $\mathbf{u}$ satisfying conditions of Section 1.2, for any admissible fields $\mathbf{v}$ as stated above.

## 4. EXACT FORMULAS FOR THE FAULT POSITION

The key step of the solution method is to choose suitable adjoint elastodynamic fields which reveal the geometry of the fault in the reciprocity gap functional. In the following we present two possible choices based on the representation of adjoint fields either by $\mathbf{v} = \nabla \mathbf{u}(\mathbf{x}, t)$ or by $\mathbf{v} = \text{rot} \mathbf{u}(\mathbf{x}, t)$. In order to simplify the notation, we shall suppose that the coordinates have been chosen such that the equation of plane $\Pi$ is $x_3 = 0$.

For the sake of simplicity, we restrain the presentation to the case of material isotropy, i. e. the general form of the tensor of elastic moduli $\mathbf{A}$ is:

$$ A_{ijkl} = \lambda \delta_{ij} \delta_{kl} + \mu (\delta_{ik} \delta_{jl} + \delta_{il} \delta_{jk}) \quad i, j, k, l = 1, 3. $$

### 4.1. Adjoint P-wave field method

To identify a vector function $\llbracket \mathbf{u} \rrbracket(\hat{\mathbf{x}}, t)$ of three variables $\hat{\mathbf{x}} = (x_1, x_2, 0)$ and time $t$, we need to introduce an adjoint vector field $\mathbf{v}$, which depends on three parameters, denoted hereafter by $\hat{\mathbf{s}} = (s_1, s_2, 0)$ and $q$. The $\mathbf{p}$-wave adjoint field satisfying (10) is defined as follows:
$$v^{(\hat{s},q)}(\boldsymbol{x}, t) = \nabla_x \Phi^{(\hat{s},q)}(\boldsymbol{x}, t) \qquad (13)$$

$$\Phi^{(\hat{s},q)}(\boldsymbol{x}, t) = \exp(iq t - \varepsilon t) \exp(i\hat{s} \cdot \boldsymbol{x}) \exp\left(x_3\left(\|\hat{s}\|^2 + \frac{(iq - \varepsilon)^2}{c_1^2}\right)^{1/2}\right),$$

where $c_1 = ((\lambda + 2\mu)/\rho)^{1/2}$ is the velocity of the pressure wave (P-wave) and $\varepsilon$ is a small positive number tending to zero, i.e $\varepsilon = 0+$.

The proposed adjoint field $v^{(\hat{s},q)}$ together with the solution $u$ satisfies (10), because at $t = 0$:

$$u = 0, \quad \ddot{u} = 0$$

and at $t = \infty$:

$$v = 0, \quad \ddot{v} = 0.$$

In order to simplify notation we shall denote the instantaneous reciprocity gap functional given in (6) as:

$$\mathcal{R}[d; \hat{s}, q] = \mathcal{R}[u, v^{(\hat{s},q)}],$$

where $d$ denotes the measured data on the exterior boundary $S$, i. e. the displacement $u|_S$ and traction $\sigma \cdot n|_S$ fields.

Taking into account that the normal of the fault is $n = e_3$ and that the slip $[u]$ has components only in the $\Pi$ plane along $e_1$ and $e_2$ and computing the stress components $\sigma_{13}[v]$ and $\sigma_{13}[v]$ we obtain:

$$\mathcal{R}[d; \hat{s}, q] = \int_{0}^{\infty} \int_{\Sigma(t)} \boldsymbol{n} \cdot \boldsymbol{\sigma}[v] \cdot [\![u]] \, ds \, dt = 2\mu \left[ \|\hat{s}\|^2 + \frac{(iq - \varepsilon)^2}{c_1^2} \right]^{1/2} \\ \times \int_{0}^{\infty} \int_{\Sigma(t)} (i\hat{s} \cdot [\![u]]) \exp(i\hat{s} \cdot \boldsymbol{x}) \exp(iqt - \varepsilon t) \, dx_1 dx_2 dt. \quad (14)$$

Let us denote $q' = q + i0^+$ and by $G(\boldsymbol{x}; q') = \mathcal{F}_t(g(\boldsymbol{x}, t))$ the integral Fourier transform of function $g(\boldsymbol{x}, t)$ with respect to time. The Fourier transform of a function $g_+$ vanishing for $t < 0$:

$$g_+(\boldsymbol{x}, t) = 0 \quad \text{if} \quad t < 0 \quad \text{and} \quad g_+(\boldsymbol{x}, t) = g(\boldsymbol{x}, t) \quad \text{if} \quad t > 0$$

has the following expression:

$$G(\boldsymbol{x}; q') = \mathcal{F}_t(g(\boldsymbol{x}, t)) = \int_{-\infty}^{\infty} \exp(iq't)g_+(\boldsymbol{x}, t) \, dt \\ = \int_{0}^{\infty} \exp(iq't)g(\boldsymbol{x}, t) \, dt. \quad (15)$$

Let us now introduce $\mathcal{F}_{\hat{x}}$ the Fourier transform with respect to spatial variable $\hat{x} = (x_1, x_2, 0)$ of a function $h(\boldsymbol{x}; q')$:

$$H(\hat{s}; q') = \mathcal{F}_{\hat{x}}(h(\boldsymbol{x}; q')) = \int_{-\infty}^{\infty} \int_{-\infty}^{\infty} \exp(i\hat{s} \cdot \boldsymbol{x})h(\boldsymbol{x}; q') \, dx_1 dx_2. \quad (16)$$
In order to express the reciprocity gap functional as a spatial and temporal Fourier transform, we extend the function $[\mathbf{u}]$ on the whole plane $\Pi$ by $[\mathbf{u}](\hat{\mathbf{x}}, t) = 0$ outside the fault, i. e. $\hat{\mathbf{x}} \notin \Sigma(t)$.

Equation (14) becomes:

$$2\mu \left( \|\hat{\mathbf{s}}\|^2 + \frac{(iq')^2}{c_1^2} \right)^{1/2} \mathcal{F}_{\hat{x}}(i\hat{\mathbf{s}} \cdot \mathcal{F}_t(\llbracket \mathbf{u} \rrbracket)). \quad (17)$$

Remarking that

$$\mathcal{F}_{\hat{x}}(i\hat{\mathbf{s}} \cdot \mathcal{F}_t(\llbracket \mathbf{u} \rrbracket)) = \mathcal{F}_{\hat{x}}(\mathcal{F}_t(\mathrm{div}[\mathbf{u}]) .$$

It follows from (17) that:

$$\mathcal{F}_{\hat{x}}(\mathcal{F}_t(\mathrm{div}[\mathbf{u}])) = \frac{1}{2\mu} \mathcal{R}^P[\mathbf{d}; \hat{s}, q] \left( \|\hat{\mathbf{s}}\|^2 + \frac{(iq')^2}{c_1^2} \right)^{-1/2}.$$

As a consequence, the solution for $\mathrm{div}[\mathbf{u}]$ can be expressed by means of inverse space and time Fourier transforms:

$$\mathrm{div}[\mathbf{u}] = \mathcal{F}_t^{-1} \mathcal{F}_{\hat{x}}^{-1} \left( \frac{1}{2\mu} \mathcal{R}^P[\mathbf{d}; \hat{s}, q] \left( \|\hat{\mathbf{s}}\|^2 + \frac{(iq')^2}{c_1^2} \right)^{-1/2} \right). \quad (18)$$

This is the desired formula we are expecting for, the function $\mathrm{div}[\mathbf{u}]$ has the same spatial support as the slip itself $[\mathbf{u}]$. In order to complete the proof we shall later show that the function of $\hat{\mathbf{x}}$ in the right hand side of the last expression is a function with a compact support.

## 4.2. Adjoint S-wave field method

The second method is based on the adjoint function defined as follows:

$$\mathbf{v}^{(\hat{s}, q)}(\mathbf{x}, t) = \mathrm{rot} (\Phi^{(\hat{s}, q)}(\mathbf{x}, t) e_3), \quad (19)$$

$$\Phi^{(\hat{s}, q)}(\mathbf{x}, t) = \exp(iqt - \varepsilon t) \exp(i\hat{\mathbf{s}} \cdot \mathbf{x}) \exp\left(x_3\left(\|\hat{\mathbf{s}}\|^2 + \frac{(iq - \varepsilon)^2}{c_2^2}\right)^{1/2}\right),$$

where $e_3$ is the normal to the fault plane $\Pi$, $c_2 = (\mu/\rho)^{1/2}$ is the velocity of the shear wave (S-wave). The method of derivation is the same as is given for the P-wave in the previous section.

In order to simplify the writing, we define a vector $[\mathbf{u}^{\perp}]$ orthogonal to $[\mathbf{u}]$ by setting:

$$[\mathbf{u}^{\perp}] = (-[u_2], [u_1], 0). \quad (20)$$

Equation (9) becomes:

$$2\mu \left( \|\hat{\mathbf{s}}\|^2 + \frac{(iq - \varepsilon)^2}{c_1^2} \right)^{1/2} \mathcal{F}_{\hat{x}}(\mathrm{div}\mathcal{F}_t(\llbracket \mathbf{u}^{\perp} \rrbracket)) = \mathcal{R}^S[\mathbf{d}; \hat{s}, q] \quad (21)$$

or by interchanging the operators $\mathcal{F}_t$ and $\mathrm{div}$:

$$2\mu \left( \|\hat{\mathbf{s}}\|^2 + \frac{(iq - \varepsilon)^2}{c_1^2} \right)^{1/2} \mathcal{F}_{\hat{x}}(\mathcal{F}_t(\mathrm{div}[\mathbf{u}^{\perp}])) = \mathcal{R}^S[\mathbf{d}; \hat{s}, q]. \quad (22)$$
As a consequence we obtain:

$$ \mathrm{d}iv[\mathbf{u}^{\perp}] = \mathcal{F}_t^{-1} \mathcal{F}_{\hat{x}}^{-1} \left( \frac{1}{2\mu} \mathcal{R}^S[\mathbf{d}; \hat{s}, q] \left( \|\hat{s}\|^2 + \frac{(iq')^2}{c_2^2} \right)^{-1/2} \right) \quad (23) $$

with $\mathcal{R}^S[\mathbf{u}; \hat{s}, q]$ and $c_2$ instead of $\mathcal{R}^P[\mathbf{u}; \hat{s}, q]$ and $c_1$. The support of the function $\mathrm{d}iv[\mathbf{u}^{\perp}]$ is the time-dependent fault $\Sigma(t)$ itself.

### 4.3. Compactness of the support functions

Let us give some mathematical arguments for proving that the support functions right hand side of (18), considered as a function of $\hat{x}$, is a compact support function.

Consider equation (17):

$$ 2\mu \left( \|\hat{s}\|^2 + \frac{(iq')^2}{c_1^2} \right)^{1/2} \mathcal{F}_{\hat{x}}(i\hat{s} \cdot \mathcal{F}_t(\llbracket \mathbf{u} \rrbracket)) = \mathcal{R}^P[\mathbf{u}; \hat{s}, q]. \quad (24) $$

Let

$$ \mathcal{G}(\mathbf{x}, q') = \mathcal{F}_t(\llbracket \mathbf{u} \rrbracket). $$

If the support $\Sigma$ of $\llbracket \mathbf{u} \rrbracket$ is time independent, the support of $\mathcal{G}(\mathbf{x}, q') = \mathcal{F}_t(\llbracket \mathbf{u} \rrbracket)$ is the same as for $\Sigma$ and we have to check that $\mathrm{supp} \mathcal{F}_t(\llbracket \mathbf{u} \rrbracket)$ is compact.

For the time dependent case, $\mathrm{supp}\llbracket \mathbf{u} \rrbracket$ is necessarily bounded by the solid $\Omega$. The compactness of $\mathrm{supp} \mathcal{F}_t(\llbracket \mathbf{u} \rrbracket)$ implies the compactness of $\mathrm{supp}\llbracket \mathbf{u} \rrbracket = \mathrm{supp}(\mathrm{div}\llbracket \mathbf{u} \rrbracket)$.

Let us check that $\mathrm{supp} \mathcal{F}_t(\llbracket \mathbf{u} \rrbracket)$ is compact. Starting from equation (17) we have:

$$ \int_{\Pi} \exp(i\hat{s} \cdot \mathbf{x}) (\mathrm{div} \mathcal{G}(\mathbf{x}, q')) dx_1 dx_2 = \mathcal{R}^P[\mathbf{u}, \hat{s}, q] \frac{1}{2\mu} \left( \|\hat{s}\|^2 + \frac{(iq')^2}{c_1^2} \right)^{-1/2} . \quad (25) $$

Owing to the term $q' = q + i0^+$, the right hand side of the last formula is infinitely differentiable in $\hat{s}$. Consider now the analytical extension of functions in the complex plane:

$$ \begin{gathered} f(\hat{s}) \to f(z), \quad z \in \mathbb{C}^2, \quad z = \hat{s} + i\mathbf{w}, \\ (z_1, z_2) = (s_1 + i\mathbf{w}_1, s_2 + i\mathbf{w}_2), \quad \mathbf{w} = (\mathbf{w}_1, \mathbf{w}_2) \in \mathbb{R}^2. \end{gathered} $$

Extended functions in (25) satisfy:

$$ \int_{\Pi} \exp(i\mathbf{z} \cdot \mathbf{x}) (\mathrm{div} \mathcal{G}(\mathbf{x}, q')) dx_1 dx_2 = \mathcal{R}^P[\mathbf{u}, \hat{z}, q] \frac{1}{2\mu} \left( \|\mathbf{z}\|^2 + \frac{(iq')^2}{c_1^2} \right)^{-1/2} . \quad (26) $$

If the solid $\Omega$ is bounded by the cube $[-a, a]^3$, then the extension of the reciprocity gap functional: $\mathcal{R}^P[\mathbf{u}, \hat{z}, q]$ satisfies the following inequalities at infinity, i.e. for large $\|\mathbf{z}\|$:

$$ \mathcal{R}^P[\mathbf{u}, \hat{z}, q] \leq C \exp(a(\|z_1\| + \|z_2\|)) \exp(a(\|z_1\|^2 + \|z_2\|^2)) \quad (27) $$

$$ \leq C \exp(2a(\|z_1\| + \|z_2\|)). \quad (28) $$
Therefore the right hand side of (26) satisfies a similar inequality with another positive constant $C'$

$$Z(z, q) \leq \mathcal{R}^P[\mathbf{u}, \hat{\mathbf{z}}, q] \frac{1}{2\mu} \left( \|z\|^2 + \frac{(iq')^2}{c_1^2} \right)^{-1/2} \leq C' \exp(2a(\|z_1\| + \|z_2\|)). \quad (29)$$

This last inequality implies that the extension $Z(z, q)$ is of the “exponential type” with coefficient $2a$, and thus, according to the Paley—Wiener theorem (Schwartz (1978), p. 271), its original $Z(\hat{\mathbf{z}}, q)$ is a tempered distribution, belonging to the space of Fourier transform of compact support functions of $\hat{\mathbf{x}} = (x_1, x_2, 0)$ on the plane $\Pi$.

The function of $\hat{\mathbf{x}}$ with compact support on the plane $\Pi$, considered here is exactly the function:

$$\mathrm{div} \, \mathbf{G}(\mathbf{x}, q') = \mathcal{F}_t(\mathrm{div} \, [\mathbf{u}]).$$

This proves the compactness of $\mathrm{supp}(\mathrm{div} \, [\mathbf{u}])$. The proof for the compactness of $\mathrm{supp}(\mathrm{div} \, [\mathbf{u}^\perp])$ defined in (23) follows from the same arguments.

Finally, we remark that the proof does is independent of the given boundary data and would be the same when the boundary data involves a non zero stress vector.

## 5. CONCLUDING REMARKS

In this paper, we have proposed an explicit and exact formula for solving inverse seismic problem. In contrast with different solutions previously known in the literature for unbounded media and for the frequency domain, (Born’s approximation of far-fields, Kirchhoff’s approximation of high-frequencies), we have considered the time domain inversion for a bounded solid which releases elastically its internal stress. The solution involves time and space Fourier transforms of quantities related to boundary measurements. It makes use of suitable adjoint waves fields, considered as dual bases of functions in the Reciprocity Gap Functional. For a stress free solid which releases its internal or residual stresses during the slip rupture process, the solution for the time-dependent fault history is exact. The method applies to a partially stress free solid, such as the ground surface, provided that the boundary conditions on remaining part of the external boundary, the underground part, can be estimated approximately.

## REFERENCES

1. J. D. Achenbach, D. A. Sotitopoulos, and H. Zhu, Characterisation of cracks from ultrasonic scattering data. *J. Appl. Mechanics* (1987) **54**, p. 754.

2. K. Aki and P. G. Richards, *Quantitative Seismology: Theory and Methods*. W.H. Freeman, New York.
3. S. Andrieux and A. Ben Abda, Identification de fissures planes par une donnée de bord unique : un procédé direct de localisation et d'identification. *Comptes Rendus Académie des Sciences, I* (1992) **315**, 1323-1328.

4. S. Andrieux and A. Ben Abda, Identification of planar cracks by complete overdetermined data inversion formulae. *Inverse Problems* (1996) **12**, 553-563.

5. S. Andrieux, H. D. Bui, and A. Ben Abda, Reciprocity and crack identification. *Inverse Problems* (1999) **15**, 59-65.

6. C. J. S. Alves and T. Ha-Duong, On inverse scattering by screen. *Inverse Problems* (1999) **15**, 91-97.

7. T. Bannour, A. Ben Abda, and M. Jaoua, A semi-explicit algorithm for the reconstruction of 3D planar cracks. *Inverse Problems* (1997) **13**, 899-917.

8. A. Ben Abda and H. D. Bui, Planar cracks identification for the transient heat equation. *J. Inv. Ill-Posed Problems* (2003) **11**, No. 1, 67-86.

9. N. N. Bojarski, Exact inverse scattering theory. *Radio Sciences* (1981) **16**, p. 1025.

10. H. D. Bui, *Inverse Problems in the Mechanics of Materials: An Introduction*. CRC Press, Boca Raton, 1994.

11. H. D. Bui, A. Constantinescu, and H. Maigre, Diffraction acoustique inverse de fissure plane: solution explicite pour un solide borné. *Comptes Rendus Académie des Sciences, II* (1999) **327**, 971-976.

12. H. D. Bui, A. Constantinescu, and H. Maigre, *The Reciprocity Gap Functional for Identifying Defects and Cracks*. CISM Course and Lectures. Udine, Springer, Wien, New York (in press).

13. H. D. Bui, A. Constantinescu, and H. Maigre, Numerical identification of planar cracks in elastodynamics using the instantaneous reciprocity gap. *Inverse Problems* (2004) **20**, 993-1001.

14. H. D. Bui and H. Maigre, Facteur d'intensité des contraintes tire des grandeurs mécaniques globales. *Comptes Rendus Acad. Sci. Paris, II* (1988) 306, p. 1213.

15. S. Das and P. Suhadolc, On the inverse problem for earthquake rupture: The Haskell-type source model. *J. Geophysical Research* (1996) **101** (B3), 5725-5738.

16. S. Das and B. V. Kostrov, Inversion for seismic slip rate and distribution with stabilizing constraints: Application to the 1986 Andeanof Islands earthquake. *J. Geophysical Research* (1990) **95**, 6899-6913.
17. P. E. Gill and W. Murray, *Numerical Methods for Constraints Optimization*. Academic Press, New York, 1974.

18. L. M. Guelfand and G. E. Chilov, *Les Distributions*. Dunod, Paris, 1962.

19. L. Schwartz, *Théorie des Distributions*. Hermann, Paris, 1978.

20. A. Tikhonov and V. Arsenine, *Méthode de Résolution de Problèmes mal Posés*. Mir, Moscou, 1976.

21. M. Vallée, *Etude Cinématique de la Rupture Sismique en Champ Lointain: Méthodes et Résolution*. Thèse Doctorat, Université de Grenoble, 2003.