Maya De Buhan
email: debuhan@ann.jussieu.fr

Axel Osses
email: axosses@dim.uchile.cl

Logarithmic stability in determination of a 3D viscoelastic coefficient and a numerical example

Keywords: numbers: 00.00, 20.00, 42.10

We prove a Carleman estimate and a logarithmic stability estimate for an inverse problem in three dimensional viscoelasticity. More precisely, we obtain logarithmic stability for the inverse problem of recovering the spatial part of a viscoelastic coefficient of the form p(x)h(t) from a unique measurement on an arbitrary part of the boundary. The main assumptions are: h (0) = 0, h(0) = 0, p is known in a neighborhood of the boundary and regularity and sensitivity of the reference trajectory. We propose a method to solve the problem numerically and illustrate the theoretical result by a numerical example.

Introduction and main results

This article is concerned with the inverse problem of determining an unknown coefficient in the 3D viscoelasticity system. We set Ω an open bounded domain of R 3 with a sufficiently smooth boundary ∂Ω. The viscoelasticity system, endowed with initial and boundary conditions, is then the following:

         Pu(x, t) = f (x, t), ∀(x, t) ∈ Ω × (0, +∞), u(x, 0) = ū0 (x), ∀x ∈ Ω, ∂ t u(x, 0) = ū1 (x), ∀x ∈ Ω, u(x, t) = 0, ∀(x, t) ∈ ∂Ω × (0, +∞), ( 1
)
where P is the integro-differential hyperbolic operator defined as Pu(x, t) = ∂ 2 t u(x, t) -∇ • µ(x)(∇u(x, t) + ∇u(x, t) T ) + λ(x)(∇ • u)(x, t)I + t 0 ∇ • μ(x, s)(∇u(x, t -s) + ∇u(x, t -s) T ) + λ(x, s)(∇ • u)(x, t -s)I ds.

(2)

This system models the dynamics of a 3D viscoelastic material subjected to a load f , u being the displacement vector, ū0 and ū1 the initial displacement and velocity, (λ, µ) the Lamé coefficients and ( λ, μ) the viscosity coefficients, respectively. The well-posedness nature of the direct problem (1)-( 2) is guaranteed by the following proposition:

Proposition 1 (Theorem 4.2 in [START_REF] Bellout | Existence of global weak solutions for a class of quasilinear hyperbolic integro-differential equations describing visco-elastic materials[END_REF]) Under regularity assumptions on the coefficients (λ, µ, λ, μ), if ū0 ∈ H 1 0 (Ω) 3 , ū1 ∈ L 2 (Ω) 3 and f ∈ L ∞ (0, +∞; L 2 (Ω) 3 ) ∩ L 2 (0, +∞; L 2 (Ω) 3 ), then there exists a unique weak solution u ∈ L ∞ (0, +∞; H 1 0 (Ω) 3 )∩ W 1,∞ (0, +∞; L 2 (Ω) 3 ) ∩ W 2,∞ (0, +∞; H -1 (Ω) 3 ) to the problem (1)- [START_REF] Grasselli | Determining the relaxation tensor in linear viscoelasticity of integral type[END_REF].

If we assume now that the coefficient μ can be decomposed as follows: μ(x, t) = p(x)h(t), ∀(x, t) ∈ Ω × (0, +∞).

(3)

Then, the inverse problem we are interested in is formulated accordingly:

Definition 1 (Inverse problem) Given (λ, µ, λ, h, u 0 , u 1 , f ), recover p(x), for all x ∈ Ω, from measurements of u(x, t), ∀(x, t) ∈ Γ × (0, T ),

where Γ is a part of ∂Ω and T > 0.

Several authors have dealt with the problem of recovering the coefficients of the viscoelasticity system (1)- [START_REF] Grasselli | Determining the relaxation tensor in linear viscoelasticity of integral type[END_REF]. Some of them, e.g. Grasselli [START_REF] Grasselli | Determining the relaxation tensor in linear viscoelasticity of integral type[END_REF], Janno et al. [START_REF] Janno | Inverse problems for identification of memory kernels in viscoelasticity[END_REF], von Wolfersdorf [START_REF] Wolfersdorf | On identification of memory kernels in linear viscoelasticity[END_REF], Cavaterra et al. [START_REF] Cavaterra | On an inverse problem for a model of linear viscoelastic Kirchhoff plate[END_REF], recovered the time dependance h of the coefficient μ by reducing the problem to a non-linear Volterra integral equation by using Fourier's method to solve the direct problem and by applying the contraction principle. Others, e.g. Lorenzi [START_REF] Lorenzi | A multidimensional identification problem related to a hyperbolic integrodifferential equation[END_REF], Lorenzi and Romanov [START_REF] Lorenzi | Stability estimates for an inverse problem related to viscoelastic media[END_REF], recovered the space dependance p of the coefficient μ by using the method of Bukhgeim and Klibanov [START_REF] Bukhgeim | Global uniqueness of class of multidimensional inverse problems[END_REF] based on Carleman estimates [START_REF] Carleman | Sur un problème d'unicité pour les systèmes d'équations aux dérivées partielles à deux variables indépendantes[END_REF]. Thus, in 2007, Lorenzi et al. [START_REF] Lorenzi | Recovering a Lamé kernel in a viscoelastic system[END_REF] recovered the coefficient p from three measurements, relying on the assumption that all the coefficients in the operator [START_REF] Grasselli | Determining the relaxation tensor in linear viscoelasticity of integral type[END_REF] are independent of the third space variable. Here, unlike the latter, we recover p from a unique measurement, assuming that it is known in a neighborhood of the boundary. Furthermore, Theorem 1 establishes the logarithmic stability with respect to a unique measurement on an arbitrary part of the boundary for the inverse problem of recovering the coefficient p. Its proof is given in Section 2. It relies on a Carleman estimate (Theorem 2) and a sharp unique continuation result (Theorem 3) for the operator [START_REF] Grasselli | Determining the relaxation tensor in linear viscoelasticity of integral type[END_REF]. We prove theses results in Section 3 and 4 respectively. In Section 5, we propose an adaptive spectral method to approximate the solution of the inverse problem numerically.

Stability estimate

Inverse problems are ill-posed in the classical sense [START_REF] Lavrent'ev | Some Ill-Posed Problems of Mathematics Physics[END_REF]. Stability estimates play thus a special role in the theory. Bukhgeim and Klibanov [START_REF] Bukhgeim | Global uniqueness of class of multidimensional inverse problems[END_REF] developed a remarkable method based on Carleman estimates [START_REF] Carleman | Sur un problème d'unicité pour les systèmes d'équations aux dérivées partielles à deux variables indépendantes[END_REF] to prove the uniqueness and stability for inverse problems associated to partial differential equations. In [START_REF] De Buhan | Un résultat de stabilité pour la récupération d'un paramètre du système de la viscoélasticité 3[END_REF], using this method, we proved a Hölder stability result with a unique internal measurement for the recovery of p in the system (1)- [START_REF] Janno | Inverse problems for identification of memory kernels in viscoelasticity[END_REF]. In this article, we extend this result using a unique continuation estimate and obtain a logarithmic stability result. Although this estimate is weaker, it is nonetheless related to a measurement on an arbitrary part of the boundary. A similar result was proved in Bellassoued [START_REF] Bellassoued | Global logarithmic stability in inverse hyperbolic problem by arbitrary boundary observation[END_REF] for the hyperbolic scalar equation and in Bellassoued et al. [START_REF] Bellassoued | Inverse problem of determining the density and the Lamé coefficients by boundary data[END_REF] for the Lamé system. We follow their method but adapt it to take into account the additive integral term of operator [START_REF] Grasselli | Determining the relaxation tensor in linear viscoelasticity of integral type[END_REF]. Thus, our result can be seen as an extension of the paper [START_REF] Bellassoued | Inverse problem of determining the density and the Lamé coefficients by boundary data[END_REF]. Nevertheless, the presence of the integral term modifies the structure of the equation. Firstly, the latter is not reversible anymore because of the energy dissipation associated with the viscous phenomenon.

In practice, this prevents the extension of the solution to [-T, 0] as usually performed. This explains why we need using a point-wise Carleman estimate (see Section 1.2).

Secondly, the coefficients of the equation now depend on time and not only on space as for the Lamé system. This implies that the FBI transform classically used to prove the unique continuation result (see Section 1.3) is useless. Hence, we have to introduce a new transformation in order to deal with convolution terms. We need now to define the following condition on a scalar function q:

Condition 1 The scalar function q is said to satisfy Condition 1 if (i) there exists K > 0 such that ∀x ∈ Ω, q(x) ≥ K, (ii) there exists x 0 ∈ R 3 \ Ω such that ∀x ∈ Ω, 1 2 q(x) -∇q(x) • (x -x 0 ) ≥ 0.

We can now state the main result:

Theorem 1 (Logarithmic stability) Let u (respectively ū) be the solution of the system (1)-( 3), associated to the coefficient p (respectively p). We assume that (H1) (λ, µ) ∈ C 2 (Ω) 2 and ( λ, μ) ∈ C 2 (Ω × (0, +∞)) 2 are such that the solutions u and ū ∈ W 8,∞ (Ω × (0, +∞)) 3 , (H2) µ and λ + 2µ satisfy Condition 1 with a same x 0 , (H3) p = p is known in a neighborhood ω of ∂Ω, (H4) h(0) = 0, h (0) = 0, (H5) there exists M > 0 such that, ∀x ∈ Ω \ ω, |(∇ū 0 + ∇ū T 0 )(x) • (x -x 0 )| ≥ M . Then, for Γ ⊂ ∂Ω arbitrarily small, there exist κ ∈ (0, 1) and T 0 > 0 such that, for all T > T 0 , the following estimate holds:

p -p H 2 (Ω) ≤ C      log      2 + C 1≤|α|≤2 ∂ α x (u -ū) 2 L 2 (Γ×(0,T ))           -κ
where C > 0 depends on the C 2 (Ω)-norm of p and p and on the W 8,∞ (Ω × (0, T )) 3 norm of u and ū.

An example of an initial datum ū0 verifying (H5) corresponds to the choice ū0 ∈

H 1 0 (Ω) 3 such that ū0 (x) = x, ∀x ∈ Ω \ ω. Thus, |(∇ū 0 + ∇ū T 0 )(x) • (x -x 0 )| = |x -x 0 | ≥ d 0 > 0, ∀x ∈ Ω \ ω.
From Theorem 1, we immediately deduce the Corollary 1.

Corollary 1 (Uniqueness) Under all the hypothesis of Theorem 1, we have then

∂ α x u(x, t) = ∂ α x ū(x, t), ∀(x, t) ∈ Γ × (0, T ), |α| = 1, 2 =⇒ p(x) = p(x), ∀x ∈ Ω.

A Carleman estimate

In a pioneering work [START_REF] Carleman | Sur un problème d'unicité pour les systèmes d'équations aux dérivées partielles à deux variables indépendantes[END_REF], Carleman introduced what is now commonly known as a Carleman estimate, in the context of proving the uniqueness in the Cauchy elliptic problem. Since then, the theory of Carleman estimates has been extensively studied.

As for a general treatment of Carleman estimates, see Hörmander [START_REF] Hörmander | Linear Partial Differential Operators[END_REF], Isakov [START_REF] Isakov | Inverse Problems for Partial Differential Equations[END_REF], Tataru [START_REF] Tataru | Carleman estimates and unique continuation for solutions to boundary value problems[END_REF], Imanuvilov and Yamamoto [START_REF] Imanuvilov | Global uniqueness and stability in determining coefficients of wave equations[END_REF]. In 2006, Cavaterra et al. [START_REF] Cavaterra | A stability result via Carleman estimates for an inverse source problem related to a hyperbolic integro-differential equation[END_REF] modified the point-wise Carleman inequality of Klibanov and Timonov [START_REF] Klibanov | Carleman Estimates for Coefficient Inverse Problems and Numerical Applications[END_REF] for a hyperbolic scalar equation and integrated, thanks to a change of variable, the integral term of (2). In 2008, Imanuvilov and Yamamoto [START_REF] Imanuvilov | Carleman estimates for the three-dimensional nonstationary Lamé system and applications to an inverse problem[END_REF] proved a Carleman estimate for the 3D Lamé system -operator (2) without the integral term. Their approach consisted in decoupling the system by writing the equations satisfied by u, ∇ ∧ u and ∇ • u in view of applying well known results for hyperbolic scalar equations [START_REF] Imanuvilov | On Carleman estimate for hyperbolic equations[END_REF]. In this article, we prove a Carleman estimate for the complete integro-differential operator (2) by combination of these two techniques. We start by introducing the following Carleman weight function:

Definition 2 (a Carleman weight function) For x 0 ∈ R 3 \ Ω and β > 0, we introduce a function ϕ in the following way: Let us introduce now some notations used in the sequel:

ϕ(x, t) = |x -x 0 | 2 -βt 2 , ∀(x, t) ∈ Ω × (0, +∞). (4)
• Distances (see Figure 2)

For x 0 ∈ R 3 \ Ω, we note

d 0 = inf x∈Ω |x -x 0 | and d = sup x∈Ω |x -x 0 | 2 .
• Domains (see Figure 3) For ε > 0 and δ > 0, we note

Ω(ε) = {x ∈ Ω, dist(x, ∂Ω) > ε}, Ω(0) = Ω, Q(ε, δ) = {(x, t) ∈ Ω(ε) × (0, +∞), ϕ(x, t) > δ} , Q(0, 0) = Q. We notice that Q ⊂ Ω × (0, T ) if T > d √ β . Ω d d 0 x 0 Figure 2.
Ω Ω(ε) ε 0 x x 0 √ δ d Ω(ε) ε t ϕ(x, t) = 0 d √ β Q Q(ε, δ) d 0 Figure 3.
Definitions of the domains Ω(ε) and Q(ε, δ).

• Norms For Q ⊂ R 3 × R, σ > 0 and k ∈ N, we introduce the following norms:

2 H 1,σ (Q) = ∇ x,t 2
L 2 (Q) + σ 2 2 L 2 (Q) , 2 H k,σ x (Q) = |α|≤k σ 2(k-|α|) ∂ α x 2 L 2 (Q) .
• Constants C and C j , 1 ≤ j ≤ 6 are generic positive constants.

Equipped with these definitions, we can now state the following result:

Theorem 2 (Carleman estimate) Let P be the operator defined by [START_REF] Bellout | Existence of global weak solutions for a class of quasilinear hyperbolic integro-differential equations describing visco-elastic materials[END_REF]. We assume that

(H1) (λ, µ) ∈ C 2 (Ω) 2 and ( λ, μ) ∈ C 2 (Ω × (0, +∞)) 2 ,
(H2) µ and λ + 2µ satisfy Condition 1 with a same x 0 .

Then, there exist β 0 > 0 and σ 0 > 0 such that, for all β ≤ β 0 , for all σ ≥ σ 0 , for all ε > 0, for all 0 < δ < d 2 0 and for all u ∈ H 2 (Q) 3 satisfying u(x, 0) = 0 or ∂ t u(x, 0) = 0, ∀x ∈ Ω, the following estimate holds:

C σ ue σϕ 2 H 2,σ x (Q(ε,δ)) ≤ (Pu)e σϕ 2 L 2 (Q) + ∇(Pu)e σϕ 2 L 2 (Q) +e Cσ u 2 H 2 (Q(ε,δ)\Q(2ε,δ)) + σ 3 e 2σδ u 2 H 2 (Q)
, where ϕ is defined by (4) and C > 0 depends on the C 2 (Q)-norm of the coefficient μ but is independent of σ.

The last term of this estimate is rather uncommon since it is global but it will disappear thanks to a clever choice of σ.

A unique continuation result

In [START_REF] Bellassoued | Unicité et contrôle pour le système de Lamé[END_REF], Bellassoued proved a sharp unique continuation result for the Lamé system. He applied a method initially developed by Robbiano [START_REF] Robbiano | Fonction de coût et contrôle des solutions des équations hyperboliques[END_REF] which uses the Fourier-Bros-Iagolnitzer (FBI) transform [START_REF] Bros | Tuboïdes et structure analytique des distributions. ii. Support essentiel et structure analytique des distributions[END_REF] to change the problem near the boundary into a problem for which elliptic estimates can be applied. In this article, we show a similar result for the system (1)-( 2) by adapting these techniques. In particular, we propose a new transformation inspired from the FBI transform but which is able to deal with the additive convolution term of the operator [START_REF] Grasselli | Determining the relaxation tensor in linear viscoelasticity of integral type[END_REF]. The result we prove is the following: Theorem 3 (Unique continuation) Let u be the vector solution of

         Pu(x, t) = R(x, t), ∀(x, t) ∈ Ω × (0, +∞), u(x, 0) = 0, ∀x ∈ Ω, ∂ t u(x, 0) = 0, ∀x ∈ Ω, u(x, t) = 0, ∀(x, t) ∈ ∂Ω × (0, +∞).
We assume that (H1) (λ, µ) ∈ C 2 (Ω) 2 and ( λ, μ) ∈ C 2 (Ω × (0, +∞)) 2 are such that the solution u ∈ W 4,∞ (Ω × (0, +∞)) 3 , (H6) R(x, t) = 0, ∀(x, t) ∈ ω × (0, +∞) where ω is a neighborhood of ∂Ω.

Then, for Γ ⊂ ∂Ω arbitrarily small, there exists T 0 > 0 such that, for all T ≥ T 0 and for all ρ > 0, the following estimate holds:

u 2 H 2 (ω×(0, T 2 -ρ)) ≤ C      log      2 + C 1≤|α|≤2 ∂ α x u 2 L 2 (Γ×(0,3T ))           -1
where C > 0 depends on the C 2 (Ω × (0, T ))-norm of the coefficient μ and on the norm W 4,∞ (Ω × (0, T )) 3 of u.

Numerical results

Section 5 presents a numerical example to illustrate one of our theoretical result [START_REF] De Buhan | Un résultat de stabilité pour la récupération d'un paramètre du système de la viscoélasticité 3[END_REF].

Here, we propose a non-quadratic functional to solve the inverse problem. Moreover, as regularization method [START_REF] Kirsch | An introduction to the mathematical theory of inverse problems[END_REF], we use a spectral basis, adapted in space and in frequency to the solution. This idea, combined with mesh adaptation, allows to improve the accuracy of the method by minimizing the numerical error. The unknown parameter is successfully recovered at each vertex of the discretized domain.

Proof of Theorem 1

This section presents the proof of the stability result (Theorem 1). The idea of the proof is the following. Firstly, we bring the unknown parameter p to the source by writing the equation satisfied by v = ∂ t (u -ū). Then, we use the method of Bukhgeim and Klibanov, i.e. we differentiate the previous equation to bring the parameter in the initial condition and we apply the Carleman estimate of Theorem 2 to the new variable ∂ t v in order to bound the initial energy. Thanks to a Carleman estimate for a first order operator (Lemma 2), we come back to the coefficient in the estimate. Thus, we obtain a first Hölder stability result with an internal measurement on a neighborhood of the boundary. Finally, we combine this result with the unique continuation relation of Theorem 3 to conclude. Let u (respectively ū) be the solution of the system (1)-( 3), associated with the coefficient p (respectively p). We suppose that hypothesis (H1)-(H5) hold and we fix Γ ⊂ ∂Ω arbitrarily small. By linearity and without loss of generality, we can suppose that ū1 = 0.

Bring the unknown parameter p to the source

We suppose (p, ū) known and we introduce p = p -p and û = u -ū which satisfies the following equation:

P û(x, t) = ∂ 2 t û(x, t) -∇ • µ(x)(∇û(x, t) + ∇û(x, t) T ) + λ(x)(∇ • û)(x, t)I + t 0 ∇ • h(s)p(x)(∇û(x, t -s) + ∇û(x, t -s) T ) + λ(x, s)(∇ • û)(x, t -s)I ds = - t 0 h(s)∇ • p(x)(∇ū(x, t -s) + ∇ū(x, t -s) T ) ds, ∀(x, t) ∈ Ω × (0, +∞),
with null initial and boundary conditions. We notice that the coefficient p we wish to recover appears in the source term. The problem is that this source term vanishes at the initial time t = 0, so it does not satisfy the hypothesis of the method of Bukhgeim and Klibanov [START_REF] Bukhgeim | Global uniqueness of class of multidimensional inverse problems[END_REF]. To overcome this problem, we derive the equation in time and set v(x, t) = ∂ t û(x, t), ∀(x, t) ∈ Ω × (0, +∞).

Then, v satisfies the equation

P v(x, t) = - t 0 h(s)∇ • p(x)(∇v(x, t -s) + ∇v(x, t -s) T ) ds -h(t)∇ • p(x)(∇ū 0 (x) + ∇ū 0 (x) T ) , ∀(x, t) ∈ Ω × (0, +∞), (5)
with null initial and boundary conditions. As a consequence, the unknown coefficient is still in the source term but it is not vanishing at the initial time anymore.

Use the method of Bukhgeim and Klibanov

We derive (5) and set

ŵ(x, t) = ∂ t v(x, t) = ∂ 2 t û(x, t), ∀(x, t) ∈ Ω × (0, +∞),
which satisfies then

P ŵ(x, t) = - t 0 h(s)∇ • p(x)(∇ w(x, t -s) + ∇ w(x, t -s) T ) ds -h (t)∇ • p(x)(∇ū 0 (x) + ∇ū 0 (x) T ) , ∀(x, t) ∈ Ω × (0, +∞), (6)
with the following initial and boundary conditions:

     ŵ(x, 0) = 0, ∀x ∈ Ω, ∂ t ŵ(x, 0) = -h(0)∇ • p(x)(∇ū 0 (x) + ∇ū 0 (x) T ) , ŵ(x, t) = 0, ∀(x, t) ∈ ∂Ω × (0, +∞). (7)
Therefore, the coefficient p that we want to recover appears now in the initial conditions. We are going to write an inequality by bounding the energy of the system (7) by the source and the observations. To do that, we apply the Carleman estimate of Theorem 2 to ŵ and to its time derivatives.

Apply the Carleman estimate

Theorem 2 states that for ϕ defined by equation ( 4) with β sufficiently small, for σ sufficiently large, for ε > 0 and 0 < δ < d 2 0 , the following estimate holds:

1 σ ŵe σϕ 2 H 2,σ x (Q(ε,δ)) ≤ C e Cσ ŵ 2 H 2 (Q(ε,δ)\Q(2ε,δ)) + σ 3 e 2σδ ŵ 2 H 2 (Q) + Q t 0 h(s)∇ • p(x)(∇ w(x, t -s) + ∇ w(x, t -s) T ) ds 2 e 2σϕ(x,t) dx dt + Q t 0 h(s)∇ ∇ • p(x)(∇ w(x, t -s) + ∇ w(x, t -s) T ) ds 2 e 2σϕ(x,t) dx dt + h ∇ • p(∇ū 0 + ∇ū T 0 ) e σϕ 2 L 2 (Q) + h ∇ ∇ • p(∇ū 0 + ∇ū T 0 ) e σϕ 2 L 2 (Q) . (8)
We now need a result to bound the integral terms. This is given by the following lemma: Lemma 1 (Lemma 3.1.1 in [START_REF] Klibanov | Carleman Estimates for Coefficient Inverse Problems and Numerical Applications[END_REF]) Let ϕ be defined by [START_REF] Wolfersdorf | On identification of memory kernels in linear viscoelasticity[END_REF]. There exists C > 0 such that, for all σ > 0 and u ∈ L 2 (Q),

Q t 0 |u(x, s)|ds 2 e 2σϕ(x,t) dx dt ≤ C σ Q |u(x, t)| 2 e 2σϕ(x,t) dx dt.
We use the result of Lemma 1 and the fact that (p, h) is bounded in Ω × (0, +∞) (H1) to write:

C σ ŵe σϕ 2 H 2,σ x (Q(ε,δ)) ≤ e Cσ ŵ 2 H 2 (Q(ε,δ)\Q(2ε,δ)) + σ 3 e 2σδ ŵ 2 H 2 (Q) + ∇ • p(∇ū 0 + ∇ū T 0 ) e σϕ 2 L 2 (Q) + ∇ ∇ • p(∇ū 0 + ∇ū T 0 ) e σϕ 2 L 2 (Q) .
Thus, taking into account the fact that ū0 ∈ W 8,∞ (Ω × (0, +∞)) (H1), we obtain

C σ ŵe σϕ 2 H 2,σ x (Q(ε,δ)) ≤ e Cσ ∂ 2 t û 2 H 2 (Q(ε,δ)\Q(2ε,δ)) + σ 3 e 2σδ ∂ 2 t û 2 H 2 (Q) + |α|≤2 (∂ α x p)e σϕ 2 L 2 (Q) .
Similarly, we obtain,

for i = 1, 2, C σ (∂ i t ŵ)e σϕ 2 H 2,σ x (Q(ε,δ)) ≤ e Cσ ∂ 2+i t û 2 H 2 (Q(ε,δ)\Q(2ε,δ)) + σ 3 e 2σδ ∂ 2+i t û 2 H 2 (Q) + |α|≤2 (∂ α x p)e σϕ 2 L 2 (Q) . (9)
In order to apply the Carleman estimate, we verified that the functions ∂ t ŵ and ∂ 2 t ŵ have one of their initial conditions null. This requirement is true since we supposed h (0) = 0 in (H4). We will now bound the initial energy of the system (7).

Bound the initial energy

We define the cut-off function

χ 1 ∈ C ∞ (R 3 × R) by 0 ≤ χ 1 ≤ 1 and χ 1 (x, t) = 1, if (x, t) ∈ Q(0, 2δ), 0, if (x, t) ∈ Q \ Q(0, δ),
and we set

w * (x, t) = χ 1 (x, t) ŵ(x, t), ∀(x, t) ∈ Ω × (0, +∞).
This function is illustrated in Figure 4. Since the time derivative of w * is null outside Q(0, δ), i.e. for

00000000
0 x t ϕ(x, t) = 0 Q √ δ √ 2δ χ 1 = 1 Q(0, 2δ) χ 1 = 0
t ≥ t δ x = |x -x 0 | 2 -δ β
which is well defined for all x ∈ Ω because δ < d 2 0 , we can write:

Ω(ε) |∂ t w * (x, 0)| 2 e 2σϕ(x,0) dx = - Ω(ε) t δ x 0 ∂ t |∂ t w * (x, t)| 2 e 2σϕ(x,t) dx dt = Q(ε,δ) 4βtσ|∂ t w * (x, t)| 2 e 2σϕ(x,t) -2∂ t w * (x, t) • ∂ 2 t w * (x, t)e 2σϕ(x,t) dx dt.
In addition, we have

∂ t w * (x, t) = χ 1 (x, t) ∂ t ŵ(x, t) + ∂ t χ 1 (x, t) ŵ(x, t), ∀(x, t) ∈ Ω × (0, +∞), ∂ 2 t w * (x, t) = χ 1 (x, t) ∂ 2 t ŵ(x, t) + 2∂ t χ 1 (x, t) ∂ t ŵ(x, t) + ∂ 2 t χ 1 (x, t) ŵ(x, t). Therefore, we obtain Ω(ε) |∂ t w * (x, 0)| 2 e 2σϕ(x,0) dx ≤ C Q(ε,δ) σ| ŵ(x, t)| 2 + σ|∂ t ŵ(x, t)| 2 + |∂ 2 t ŵ(x, t)| 2 e 2σϕ(x,t) dx dt .
We do the same work for the derivatives in x for |α| ≤ 2 and thus we obtain:

|α|≤2 σ 2(1-|α|) Ω(ε) |∂ α x (∂ t ŵ)(x, 0)| 2 e 2σϕ(x,0) dx ≤ C σ ŵe σϕ 2 H 2,σ x (Q(ε,δ)) + (∂ t ŵ)e σϕ 2 H 2,σ x (Q(ε,δ)) + (∂ 2 t ŵ)e σϕ 2 H 2,σ x (Q(ε,δ)) ≤ C   |α|≤2 (∂ α x p)e σϕ 2 L 2 (Q) + e Cσ û 2 H 6 (Q(ε,2δ)\Q(2ε,δ)) + σ 3 e 2σδ û 2 H 6 (Q)   .
where we used [START_REF] Carleman | Sur un problème d'unicité pour les systèmes d'équations aux dérivées partielles à deux variables indépendantes[END_REF] in the last inequality. We have to go back to the coefficient in the left hand side.

Return to the coefficient in the estimate

To this end, we need a Carleman estimate for a first order operator. It is given by the following lemma:

Lemma 2 (Lemma 3.2 in [START_REF] Imanuvilov | Global uniqueness and stability in determining coefficients of wave equations[END_REF]) We consider the following first order partial differential operator:

R(x) = a(x) • ∇ +a 0 (x) , ∀x ∈ Ω, with a 0 ∈ C 0 (Ω), a 0 C 2 (Ω) ≤ M 1 and a ∈ C 1 (Ω) n , a C 2 (Ω) ≤ M 2 , ∃x 0 , |a(x) • (x -x 0 )| ≥ M 3 > 0, ϕ 0 (x) = |x -x 0 | 2 , ∀x ∈ Ω.
Then, there exist σ 0 > 0 and C > 0 such that, for all σ ≥ σ 0 and for all q ∈ C 2 0 (Ω),

σ 2 |α|≤2 Ω |∂ α x q(x)| 2 e 2σϕ0(x) dx ≤ C |α|≤2 Ω |∂ α x (R(x)q(x))| 2 e 2σϕ0(x) dx.
Since the coefficient p satisfies the following first order system:

h(0)∇ • p(x)(∇ū 0 (x) + ∇ū 0 (x) T ) = -∂ t w * (x, 0),
we can apply Lemma 2 to each equation of the system with

q(x) = p(x), R(x)q(x) = -∂ t w * (x, 0), ϕ 0 (x) = ϕ(x, 0), a 0 (x) = h(0)∇ • (∇ū 0 (x) + ∇ū 0 (x) T ), a(x) = h(0)(∇ū 0 (x) + ∇ū 0 (x) T ),
which satisfy (H5) and we obtain then

|α|≤2 Ω(ε) |∂ α x p(x)| 2 e 2σϕ(x,0) dx ≤ C |α|≤2 σ 2(1-|α|) Ω(ε) |∂ α x (∂ t w * )(x, 0)| 2 e 2σϕ(x,0) dx.
We supposed in (H3) that we know p in a neighborhood ω of the boundary ∂Ω, i.e. p = 0 in ω and so we can assume ε to be sufficiently small such that Ω \ Ω(ε) ⊂ ω. This allows to integrate on Ω in the left hand side. Then, we have

|α|≤2 Ω |∂ α x p(x)| 2 e 2σϕ(x,0) dx ≤ C   |α|≤2 (∂ α x p)e σϕ 2 L 2 (Q) + e Cσ û 2 H 6 (Q(ε,δ)\Q(2ε,δ)) + σ 3 e 2σδ û 2 H 6 (Q)   . (10)
We can absorb the first term of the right hand side of [START_REF] Lorenzi | Recovering a Lamé kernel in a viscoelastic system[END_REF] in the left hand side, thanks to the Carleman weights. Indeed,

|α|≤2 (∂ α x p)e σϕ 2 L 2 (Q) = |α|≤2 Ω |∂ α x p(x)| 2 e 2σϕ(x,0) t 0
x 0 e 2σ(ϕ(x,t)-ϕ(x,0)) dt dx, with, ∀x ∈ Ω:

t 0 x 0 e 2σ(ϕ(x,t)-ϕ(x,0)) dt ≤ +∞ 0 e -2σβt 2 dt = 1 √ 2σβ +∞ 0 e -z 2 dz = C √ σ .
Then,

|α|≤2 (∂ α x p)e σϕ 2 L 2 (Q) ≤ C √ σ |α|≤2 Ω |∂ α x p(x)| 2 e 2σϕ(x,0) dx.
Finally, we obtain

e 2σd 2 0 p 2 H 2 (Ω) ≤ |α|≤2 Ω |∂ α x p(x)| 2 e 2σϕ(x,0) dx ≤ C e Cσ û 2 H 6 (Q(ε,δ)\Q(2ε,δ)) + σ 3 e 2σδ û 2 H 6 (Q)
, and thus, we conclude that

p 2 H 2 (Ω) ≤ C e Cσ û 2 H 6 (Q(ε,δ)\Q(2ε,δ)) + e -σ(d 2 0 -δ) , with l = d 2 0 -δ > 0, since sup σ (σ 3 e -σl
) < +∞ and since we assumed in (H1) that û is in W 8,∞ (Q). In order to have

e Cσ û 2 H 6 (Q(ε,δ)\Q(2ε,δ)) = e -σl , we propose to set σ = - 1 l + C log û 2 H 6 (Q(ε,δ)\Q(2ε,δ)) . (11)
We have σ > 0 if we suppose that û 2 H 6 (Q(ε,δ)\Q(2ε,δ)) < 1. Thus, the right hand side of the equation can be written as follows:

e Cσ û 2 H 6 (Q(ε,δ)\Q(2ε,δ)) + e -σl = 2 û 2 H 6 (Q(ε,δ)\Q(2ε,δ)) l l+C .
We use now the following interpolation result:

Lemma 3 (Proposition 4 in [START_REF] Dautray | Analyse mathématique et calcul numérique pour les sciences et les techniques[END_REF]) For all m ∈ N * , there exists a constant C > 0 such that, for all r satisfying 0 ≤ r ≤ m and for all u ∈ H m (Q), we have

u H r (Q) ≤ C u 1-r m L 2 (Q) u r m H m (Q)
We apply Lemma 3 to the second derivative of û ∈ H 8 (Q) with m = 6 and r = 2 in order to write:

û H 6 (Q(ε,δ)\Q(2ε,δ)) ≤ C û 1/3 H 2 (Q(ε,δ)\Q(2ε,δ)) û 2/3 H 8 (Q(ε,δ)\Q(2ε,δ)) . Hence, p H 2 (Ω) ≤ C û 1/3 H 2 (Q(ε,δ)\Q(2ε,δ)) l l+C . ( 12
)
So we already have a stability result with observation in

Q(ε, δ) \ Q(2ε, δ) which is included in (Ω(ε) \ Ω(2ε)) × (0, T ) if T > d √ β .

Use the unique continuation result

We can now apply Theorem 3 to û. We verify that hypothesis (H6) holds, i.e. that

R(x, t) = - t 0 h(s)∇ • p(x)(∇ū(x, t -s) + ∇ū(x, t -s) T ) ds
vanishes in Ω \ Ω(2ε) ⊂ ω and choose ρ < δ. We deduce that, if T is large enough, we have

û 2 H 2 (Q(ε,δ)\Q(2ε,δ)) ≤ C log 2 + C 2 |α|=1 ∂ α x û 2 L 2 (Γ×(0,6T )) -1 . Therefore, p H 2 (Ω) ≤ C log 2 + C 2 |α|=1 ∂ α x û 2 L 2 (Γ×(0,6T )) - l 6(l+C)
.

We set κ = l 6(l + C) ∈ (0, 1) and change 6T by T . This achieves to prove Theorem 1.

Proof of Theorem 2

This section is devoted to the proof of the Carleman estimate (Theorem 2). The proof consists in decoupling the system of equations as in [START_REF] Imanuvilov | Carleman estimates for the three-dimensional nonstationary Lamé system and applications to an inverse problem[END_REF] by writing the equations satisfied by u, ∇ ∧ u and ∇ • u. Then, the idea is to introduce a change of variable as in [START_REF] Cavaterra | A stability result via Carleman estimates for an inverse source problem related to a hyperbolic integro-differential equation[END_REF] to reduce the problem to a scalar hyperbolic equation for which the Carleman estimate is well known (Corollary 2). Finally, we come back to the initial variable by a series of inequalities assuming the regularity of the coefficients and using the fundamental Lemma 1.

Let u be the solution of the system

     Pu(x, t) = F (x, t), ∀(x, t) ∈ Ω × (0, +∞), u(x, 0) = 0 or ∂ t u(x, 0) = 0, ∀x ∈ Ω, u(x, t) = 0, ∀(x, t) ∈ ∂Ω × (0, +∞), (13)
with

Pu(x, t) = ∂ 2 t u(x, t) -∇ • µ(x)(∇u(x, t) + ∇u(x, t) T ) + λ(x)(∇ • u)(x, t)I + t 0 ∇ • μ(x, s)(∇u(x, t -s) + ∇u(x, t -s) T ) + λ(x, s)(∇ • u)(x, t -s)I ds.
Let assume that hypothesis (H1) and (H2) are satisfied. We have first to decouple the equations.

Decouple the system of equations

To this end, we develop the system (13) as follows:

∂ 2 t u(x, t) -µ(x)∆u(x, t) -(µ + λ)(x)∇(∇ • u)(x, t) -(∇ • u)(x, t)∇λ(x) -∇u(x, t) + ∇u(x, t) T • ∇µ(x) + t 0 μ(x, t -s)∆u(x, s)ds + t 0 (∇ • u)(x, s)∇ λ(x, t -s)ds + t 0 (μ + λ)(x, t -s)∇(∇ • u)(x, s)ds + t 0 ∇u(x, s) + ∇u(x, s) T • ∇μ(x, t -s)ds = F (x, t) ∀(x, t) ∈ Ω × (0, +∞).
Notice the change of variable in the convolution terms. Then, we take the curl and the divergence of this system. We introduce the vectors u 1 = u, u 3 = ∇ ∧ u and the scalar u 2 = ∇ • u which satisfy then the following system of seven equations.

∂ 2 t u i (x, t) -q i (x)∆u i (x, t) + t 0 qi (x, t -s)∆u i (x, s)ds = F i (x, t) + A i (u 1 , u 2 , u 3 )(x, t), ∀i ∈ {1, 2, 3}, (14)
where we set

q 1 = q 2 = µ, q 3 = λ + 2µ, q1 = q2 = μ, q3 = λ + 2μ, F 1 = F, F 2 = ∇ ∧ F, and
F 3 = ∇ • F.
Here, the coupling terms A i are first order integro-differential operators in x and t, with coefficients bounded in Ω × (0, +∞), according to hypothesis (H1). These equations are coupled only at order 1, therefore we can apply the results known for the scalar equations. However, the problem is now that the trace of the functions ∇ • u and ∇ ∧ u on ∂Ω are not defined anymore. We need now to change the integrodifferential hyperbolic equation ( 14) into a hyperbolic equation in order to apply the classical results.

Use a change of variable

We introduce the following change of variable:

ũi (x, t) = q i (x)u i (x, t) - t 0 qi (x, t -s)u i (x, s)ds, ∀i ∈ {1, 2, 3}. (15)
Then, we have, for all i ∈ {1, 2, 3} and for all (x, t) ∈ Ω × (0, +∞),

∂ 2 t ũi (x, t) = q i (x)∂ 2 t u i (x, t) - t 0 ∂ 2 t qi (x, t -s)u i (x, s)ds +∂ t qi (x, 0)u i (x, t) + qi (x, 0)∂ t u i (x, t),
and

∆ũ i (x, t) = q i (x)∆u i (x, t) - t 0 qi (x, t -s)∆u i (x, s)ds + 2∇q i (x) • ∇u i (x, t) +∆q i (x)u i (x, t) -2 t 0 ∇q i (x, t -s) • ∇u i (x, s)ds - t 0 ∆q i (x, t -s)u i (x, s)ds.
Therefore, for all i ∈ {1, 2, 3} and for all (x, t) ∈ Ω × (0, +∞), the ũi satisfy an hyperbolic system of the type:

∂ 2 t ũi (x, t) -q i (x)∆ũ i (x, t) = q i (x) (F i (x, t) + A i (u 1 , u 2 , u 3 )(x, t)) + L i (u i )(x, t),
where the L i are first order integro-differential operators, the coefficients of which are bounded in Ω × (0, +∞). We now need a Carleman estimate for a scalar hyperbolic equation.

Use a Carleman estimate for an hyperbolic scalar equation

The equations we consider are valid in (0, +∞) and because of the presence of the integral term, we can not extend the solution to (-∞, 0). That is why we can not use a classical global Carleman estimate found in the literature. Therefore, we will start with the following point-wise Carleman estimate stated in Q(0, δ):

Lemma 4 (Theorem 2.2.4 in [START_REF] Klibanov | Carleman Estimates for Coefficient Inverse Problems and Numerical Applications[END_REF]) Let q ∈ C 2 (Ω) satisfy Condition 1 and ϕ be defined by (4). Then, there exist β 0 > 0 and σ 0 > 0 such that, for all β ≤ β 0 , for all σ ≥ σ 0 , for all δ > 0 and for all u ∈ H 2 (Q(0, δ)), we have

(σ|∇ x,t u(x, t)| 2 + σ 3 |u(x, t)| 2 )e 2σϕ(x,t) + ∇ • U (x, t) + ∂ t V (x, t) ≤ C|∂ 2 t u(x, t) -q(x)∆u(x, t)| 2 e 2σϕ(x,t) , ∀(x, t) ∈ Q(0, δ). (16)
Here, (U, V ) is a vector-valued function and satisfies

|U (x, t)| + |V (x, t)| ≤ C(σ|∇ x,t u(x, t)| 2 + σ 3 |u(x, t)| 2 )e 2σϕ(x,t) , ∀(x, t) ∈ Q(0, δ). Moreover, V (x, 0) = 0, ∀x ∈ Ω if u(x, 0) = 0 or ∂ t u(x, 0) = 0, ∀x ∈ Ω.
From this Lemma, we deduce the following Carleman estimate:

Corollary 2 Let q ∈ C 2 (Ω)
satisfy Condition 1 and ϕ be defined by ( 4), with β > 0 sufficiently small. Then, there exist β 0 > 0 and σ 0 > 0 such that, for all β ≤ β 0 , for all σ ≥ σ 0 , for all δ > 0 and for all u

∈ H 2 (Q(0, δ)) satisfying u(x, 0) = 0 or ∂ t u(x, 0) = 0, ∀x ∈ Ω, we have σ ue σϕ 2 H 1,σ (Q(0,δ)) ≤ C (∂ 2 t u -q∆u)e σϕ 2 L 2 (Q(0,δ)) + σ ue σϕ 2 H 1,σ (∂Q(0,δ)\(Ω×{0})) .
Proof We integrate ( 16) over Q(0, δ):

Q(0,δ) (σ|∇ x,t u(x, t)| 2 + σ 3 |u(x, t)| 2 )e 2σϕ(x,t) dx dt ≤ - Q(0,δ) ∇ • U (x, t) -∂ t V (x, t) + C|∂ 2 t u(x, t) -µ(x)δu(x, t)| 2 e 2σϕ(x,t) dx dt.
We define now (see Figure 5)

∀t ≥ 0, Ω t (ε, δ) = Q(ε, δ) ∩ (Ω × {t}) and ∀x ∈ Ω, t δ x =      |x -x 0 | 2 -δ β if δ < |x -x 0 | 2 , 0 else.
We notice that

Q(0,δ) ∇ • U (x, t)dx dt = ∞ 0 Ωt(0,δ) ∇ • U (x, t) dx dt = ∞ 0 ∂Ωt(0,δ) U (x, t) • n dx dt = ∂Q(0,δ)\(Ω×{0}) U (x, t) • n dx dt
and that Thus,

Q(0,δ) ∂ t V (x, t)dx dt = Ω t δ x 0 ∂ t V (x, t)dt dx = Ω V (x, t δ x ) -V (x, 0) dx = ∂Q(0,δ)\(Ω×{0}) V (x,
Ω t (0, 0) Ω t (ε, δ) t ϕ(x, •) ϕ(x, 0) 0 d √ β δ t δ x
- Q(0,δ) ∇ • U (x, t)dx dt - Q(0,δ) ∂ t V (x, t)dx dt ≤ ∂Q(0,δ)\(Ω×{0}) (|U (x, t)| + |V (x, t)|)dx dt ≤ C ∂Q(0,δ)\(Ω×{0}) (σ|∇ x,t u(x, t)| 2 + σ 3 |u(x, t)| 2 )e 2σϕ(x,t) dx dt.
This achieves to prove the expected result.

For ε > 0 and δ > 0, we introduce then the cut-off function χ 2 ∈ C ∞ (R 3 × R) which satisfies 0 ≤ χ 2 ≤ 1 and such that

χ 2 (x, t) = 1, if (x, t) ∈ Q(2ε, 2δ), 0, if (x, t) ∈ Q \ Q(ε, δ). (17
0 x t ϕ(x, t) = 0 Q √ δ √ 2δ Q(2ε, 2δ) χ 2 = 1 χ 2 = 0 ε ε
u * i (x, t) = χ 2 (x, t)ũ i (x, t), ∀i ∈ {1, 2, 3}, ∀(x, t) ∈ Ω × (0, +∞). Thus, u * i satisfies the equation ∂ 2 t u * i (x, t) -q i (x)∆u * i (x, t) = χ(x, t) ∂ 2 t ũi (x, t) -µ(x)∆ũ i (x, t) + 2∂ t χ(x, t)∂ t ũi (x, t) +ũ i (x, t) ∂ 2 t χ(x, t) -µ(x)∆χ(x, t) -2µ(x)∇χ(x, t) • ∇ũ i (x, t) = χ(x, t) (q i (x) (F i (x, t) + A i (u 1 , u 2 , u 3 )(x, t)) + L i (u i )(x, t)) + Li (ũ i )(x, t),
where Li is a first order integro-differential operator with coefficients bounded in Ω × (0, +∞). Since we supposed u(x, 0) = 0 or ∂ t u(x, 0) = 0, ∀x ∈ Ω, and that, according to (H2), the coefficient q i satisfies Condition 1, we can apply the Corollary 2 to u * i in Q(0, δ/2), i.e. there exists β > 0 sufficiently small such that, for σ > 0 sufficiently large, we have • all the coefficients λ, µ, λ, μ and χ are bounded in Ω×(0, +∞) according to (H1), • all the operators A i , L i and Li are of order inferior or equal to 1,

σ u * i e σϕ 2 H 1,σ (Q(ε,δ)) ≤ C (∂ 2 t u * i -q i ∆u * i )e σϕ 2 L 2 (Q(ε,δ)) , without boundary term since u * i is null in Q \ Q(ε, δ).
• Lemma 1 holds, to calculate:

(∂ 2 t u * i -q i ∆u * i )e σϕ 2 L 2 (Q(ε,δ)) = χ (µ (F i + A i (u 1 , u 2 , u 3 )) + L i (u i )) + Li (ũ i ) e σϕ 2 L 2 (Q(ε,δ)) . Then, σ u * i e σϕ 2 H 1,σ (Q(ε,δ)) ≤ C F i e σϕ 2 L 2 (Q) + C 1≤j≤3 u j e σϕ 2 H 1,σ (Q(ε,δ)) . (18)
3.4.2. Lower bound for the left hand side We use again the change of variable [START_REF] Hörmander | Linear Partial Differential Operators[END_REF], for all (x, t) ∈ Ω × (0, +∞):

u i (x, t) = 1 q i (x) ũi (x, t) + t 0 qi (x, t -s) q i (x) u i (x, s)ds, ∀i ∈ {1, 2, 3},
Then, we can write, taking into account (H2) and thanks to Lemma 1, ∀i ∈ {1, 2, 3}:

Q(ε,δ) |u i (x, t)| 2 e 2σϕ(x,t) dx dt ≤ C Q(ε,δ) |ũ i (x, t)| 2 e 2σϕ(x,t) dx dt + C σ Q(ε,δ) |u i (x, t)| 2 e 2σϕ(x,t) dx dt.
For σ sufficiently large, we obtain

Q(ε,δ) |u i (x, t)| 2 e 2σϕ(x,t) dx dt ≤ C Q(ε,δ) |ũ i (x, t)| 2 e 2σϕ(x,t) dx dt.
In the same way,

Q(ε,δ) |∇ x,t u i (x, t)| 2 e 2σϕ(x,t) dx dt ≤ C Q(ε,δ) (|ũ i (x, t)| 2 + |∇ x,t ũi (x, t)|
2 )e 2σϕ(x,t) dx dt which allows to write, ∀i ∈ {1, 2, 3}:

u i e σϕ 2 H 1,σ (Q(ε,δ)) ≤ C ũi e σϕ 2 H 1,σ (Q(ε,δ)) . (19)
In addition,

χ 2 = 1 sur Q(2ε, 2δ), so we have ũi = u * i . Hence ∀i ∈ {1, 2, 3}, ũi e σϕ 2 H 1,σ (Q(ε,δ)) = u * i e σϕ 2 H 1,σ (Q(2ε,2δ)) + ũi e σϕ 2 H 1,σ (Q(ε,δ)\Q(2ε,2δ)) . (20)
Finally, using ( 18), ( 19) and ( 20) and again Lemma 1, we obtain

1≤i≤3 σ u i e σϕ 2 H 1,σ (Q(ε,δ)) ≤ C 1≤i≤3 F i e σϕ 2 L 2 (Q) + u i .e. σϕ 2 H 1,σ (Q(ε,δ)) + σ u i e σϕ 2 H 1,σ (Q(ε,δ)\Q(2ε,2δ)) ,
and thanks to the Carleman weights, we absorb the second term of the right hand side in the left hand side. Finally, we need the following lemma:

Lemma 5 Let ϕ defined by (4). There exist C > 0 and σ 0 > 0 such that, for all σ ≥ σ 0 and for all u ∈ H 2 (Q), we have

1 σ ue σϕ 2 H 2,σ x (Q(2ε,2δ) ≤ C σ ue σϕ 2 H 1,σ x (Q(ε,δ)) + ∇(∇ • u)e σϕ 2 L 2 (Q(ε,δ)) + ∇(∇ ∧ u)e σϕ 2 L 2 (Q(ε,δ)) .
Proof With χ 2 defined in [START_REF] Tataru | Carleman estimates and unique continuation for solutions to boundary value problems[END_REF], we compute, for all (x, t) ∈ Q, ∆ χ 2 (x, t)u(x, t)e σϕ(x,t) = 2(σχ 2 (x, t)∇ϕ(x, t) + ∇χ 2 (x, t)) • ∇u(x, t)e σϕ(x,t)

+ σ 2 χ 2 (x,
= 0 in O \ Ω t (ε, δ), we obtain χ 2 ue σϕ 2 H 2,σ x (Q(ε,δ)) ≤ ue σϕ 2 L 2 (0,+∞;H 2 (Ωt(ε,δ))) + σ 2 ∇(ue σϕ ) 2 L 2 (Q(ε,δ)) + σ 4 ue σϕ 2 L 2 (Q(ε,δ)) ≤ C e σϕ ∆u 2 L 2 (Q(ε,δ)) + σ 2 e σϕ ∇u 2 L 2 (Q(ε,δ)) + σ 4 ue σϕ 2 L 2 (Q(ε,δ)) . We use the classical formula ∆u = ∇ (∇ • u) -∇ ∧ (∇ ∧ u) , and the fact that χ 2 = 1 in Q(2ε, 2δ) to write ue σϕ 2 H 2,σ x (Q(2ε,2δ)) ≤ χ 2 ue σϕ 2 H 2,σ x (Q(ε,δ)) ≤ C e σϕ ∇(∇ • u) 2 L 2 (Q(ε,δ)) + e σϕ ∇ (∇ ∧ u) 2 L 2 (Q(ε,δ)) + σ 2 ue σϕ 2 H 1,σ x (Q(ε,δ))
. This is the expected result.

Using this lemma, we can finally write:

1 σ ue σϕ 2 H 2,σ x (Q(2ε,2δ)) ≤ C σ u 1 e σϕ 2 H 1,σ x (Q(ε,δ)) + ∇u 2 e σϕ 2 L 2 (Q(ε,δ)) + ∇u 3 e σϕ 2 L 2 (Q(ε,δ)) ≤ C 1≤i≤3 σ u i e σϕ 2 H 1,σ (Q(ε,δ)) ≤ C 1≤i≤3 F i e σϕ 2 L 2 (Q) + σ u i e σϕ 2 H 1,σ (Q(ε,δ)\Q(2ε,2δ)) , then 1 σ ue σϕ 2 H 2,σ x (Q(ε,δ)) = 1 σ ue σϕ 2 H 2,σ x (Q(2ε,2δ)) + 1 σ ue σϕ 2 H 2,σ x (Q(ε,δ)\Q(2ε,2δ)) ≤ C 1≤i≤3 F i e σϕ 2 L 2 (Q) + σ u i e σϕ 2 H 1,σ (Q(ε,δ)\Q(2ε,2δ)) .
Moreover, we have

1≤i≤3 σ u i e σϕ 2 H 1,σ (Q(ε,δ)\Q(2ε,2δ)) = 1≤i≤3 σ u i e σϕ 2 H 1,σ (Q(ε,δ)\Q(2ε,δ) + σ u i e σϕ 2 H 1,σ (Q(ε,δ)\Q(ε,2δ)) ≤ e Cσ u 2 H 2 (Q(ε,δ)\Q(2ε,δ)) + σ 3 e 2σδ u 2 H 2 (Q) , since ϕ(x, t) ≤ δ, ∀(x, t) ∈ Q(ε, δ) \ Q(ε, 2δ). Thus 1 σ ue σϕ 2 H 2,σ x (Q(ε,δ)) ≤ C F e σϕ 2 L 2 (Q) + (∇F )e σϕ 2 L 2 (Q) + e Cσ u 2 H 2 (Q(ε,δ)\Q(2ε,δ)) + σ 3 e 2σδ u 2 H 2 (Q)
, and this concludes the proof of Theorem 2.

Proof of Theorem 3

In this section, we prove the unique continuation result of Theorem 3. It consists in transforming the integro-differential hyperbolic system (1)-( 2) into an elliptic one thanks to a FBI type transform. Then, we show a Carleman estimate (Theorem 4) for the resulting integro-differential elliptic operator, using the same techniques as in the proof of Theorem 2. And we use this Carleman estimate to obtain interpolation inequalities which link the value of the solution in an interior domain near the boundary to the one on the boundary. Finally, we come back to the solution of the initial problem by a series of inequalities. Let u be the solution of

Pu(x, t) = R(x, t), ∀(x, t) ∈ Ω × (0, +∞), (22)
with null initial and boundary conditions and let assume that (H1) and (H6) are satisfied.

Transform the hyperbolic system into an elliptic one

Let us choose ε > 0 such that Ω \ Ω(3ε) ⊂ ω.

000000000000 000000000000 000000000000 000000000000 000000000000 000000000000 000000000000 000000000000 000000000000 000000000000 We introduce the cut-off function χ 3 ∈ C ∞ 0 (R 3 ) (see Figure 8) which satisfies 0 ≤ χ 3 ≤ 1 and is such that

χ 3 = 1 χ 3 = 0 Ω Ω \ Ω(3ε) ε 3ε
χ 3 (x) = 0, if x ∈ Ω(4ε), 1, if x ∈ Ω \ Ω(3ε),
and we set

u * (x, t) = χ 3 (x) u(x, t), ∀(x, t) ∈ Ω × (0, +∞).
Thus, the new variable u * satisfies the following equation:

Pu * (x, t) = χ 3 (x) Pu(x, t) + [P, χ 3 ]u(x, t) =⇒ ∂ 2 t u * (x, t) -L(x)u * (x, t) + t 0 L(x, s)u * (x, t -s)ds = [P(x, t), χ 3 ]u(x, t), (23)
since, according to (H6), R(x, t) is equal to zero in ω. We have introduced the following operators:

L(x) = ∇ • µ(x)(∇ +∇ T ) + λ(x)(∇ • )I , L(x, t) = ∇ • p(x)h(t)(∇ +∇ T ) + λ(x, t)(∇ • )I .
We define a special transformation which is inspired from the classical Fourier-Bros-Iagolnitzer transform [START_REF] Bros | Tuboïdes et structure analytique des distributions. ii. Support essentiel et structure analytique des distributions[END_REF]:

Definition 3 Let T > 1 and 0 < η < 1 6 . We introduce Q = Ω×(-T 2 , 5T
2 )×(-3η, 3η). Then, ∀(x, t, r) ∈ Q, we define

(F γ u * )(x, t, r) = γ 2π +∞ 0 e -γ 2 (t+ir-y) 2 θ(t -y)u * (x, y)dy, where θ ∈ C ∞ (R), 0 ≤ θ ≤ 1, θ(z) =        1, if z > - T 2 + 3η, 0, if z < - T 2 + 3η 2 .
This function is plotted in Figure 9. Thanks to this transform, we will convert locally the hyperbolic system (23) into an elliptic system. Firstly, we notice that, ∀(x, t, r

0 z - T 2 + 3η θ 1 - T 2 - T 2 + 3 η 2
) ∈ Q, ∂ r (F γ u * )(x, t, r) = γ 2π +∞ 0 -iγ(t + ir -y)e -γ 2 (t+ir-y) 2 θ(t -y)u * (x, y)dy = γ 2π +∞ 0 -i∂ y (e -γ 2 (t+ir-y) 2 )θ(t -y)u * (x, y)dy = γ 2π +∞ 0 ie -γ 2 (t+ir-y) 2 ∂ y (θ(t -y)u * (x, y))dy.
The boundary terms of the integration by parts vanish because u * (x, 0) = 0 and e -γ 2 (t+ir-y)

2 -→ 0 if y -→ ∞. Moreover θ(t -y) = 0 if y ≥ T . Then, ∂ 2 r (F γ u * )(x, t, r) = - γ 2π +∞ 0 e -γ 2 (t+ir-y) 2 ∂ 2 y (θ(t -y)u * (x, y)) dy = -(F γ (∂ 2 y u * ))(x, t, r) -F γ (x, t, r), ∀(x, t, r) ∈ Q,
where

F γ (x, t, r) = γ 2π +∞ 0 e -γ 2 (t+ir-y) 2
(θ (t -y)u * (x, y) -2θ (t -y)∂ y u * (x, y)) dy.

Similarly, we can write:

∂ 2 t (F γ u * )(x, t, r) = (F γ (∂ 2 y u * ))(x, t, r), ∀(x, t, r) ∈ Q,
and it is easy to see that

L(x)(F γ u * )(x, t, r) = (F γ (L(x)u * ))(x, t, r), ∀(x, t, r) ∈ Q.
We have now to look at the integral term

F γ y 0 L(x, s)u * (x, y -s)ds (x, t, r) = γ 2π +∞ 0 e -γ 2 (t+ir-y) 2 θ(t -y) y 0 L(x, s)u * (x, y -s)ds dy = γ 2π +∞ 0 L(x, s) +∞ s e -γ 2 (t+ir-y) 2 θ(t -y)u * (x, y -s)dy ds = +∞ 0 L(x, s) γ 2π +∞ 0 e -γ 2 (t-s+ir-z) 2 θ(t -s -z)u * (x, z)dz ds = +∞ 0 L(x, s)(F γ u * )(x, t -s, r)ds = t -T /2 L(x, t -s)(F γ u * )(x, s, r)ds.
The main advantage of this transformation with respect to the FBI transform is that it transforms the convolution product of two functions into the convolution of the first function by the transform of the second function. The function u γ = F γ u * , function of three variables (x, t, r), satisfies then a system of elliptic integro-differential equations, ∀(x, t, r) ∈ Q:

Qu γ (x, t, r) = -2∂ 2 r u γ (x, t, r) -∂ 2 t u γ (x, t, r) -L(x)u γ (x, t, r) + t -T /2 L(x, t -s)u γ (x, s, r)ds = 2F γ (x, t, s) + G γ (x, t, s), (24)
where

G γ (x, t, r) = γ 2π +∞ 0 e -γ 2 (t+ir-y) 2 θ(t -y)[P(x, t), χ 3 ]û(x, y)dy,
with the boundary conditions:

     u γ (x, t, r) = 0, ∀(x, t, r) ∈ ∂Ω × (- T 2 , 5T 2
) × (-3η, 3η),

u γ (x, t, r) = 0, ∀(x, t, r) ∈ Ω × {- T 2 } × (-3η, 3η).
We notice then that

G γ (x, t, r) = 0, ∀x ∈ Ω \ Ω(3ε), (25)
because [P(x, t), χ 3 ] only involves the derivatives of χ 3 and its support is in Ω(3ε) \ Ω(4ε). We also have

F γ H 1 ( e Q) ≤ Ce -mγT u * H 2 (Ω×(0,3T )) (26)
as well as

u γ H 2 ( e Q) ≤ Ce Cγ u * H 2 (Ω×(0,3T )) , (27)
where C depends of η, T, Ω but not of γ and where m is independent of T . Indeed,

F γ 2 L 2 ( e Q) = γ 2π e Q +∞ 0 e -γ 2 (t+ir-y) 2 (θ (t -y)u * (x, y) -2θ (t -y)∂ y u * (x, y)) dy 2 dr dt dx. with +∞ 0 e -γ 2 (t+ir-y) 2 (θ (t -y)u * (x, y) -2θ (t -y)∂ y u * (x, y)) dy 2 ≤ y≥0, -T /2+3η/2≤t-y≤-T /2+3η |e -γ 2 (t+ir-y) 2 | 2 dy × t+T /2 0 |θ (t -y)u * (x, y) -2θ (t -y)∂ y u * (x, y)| 2 dy ≤ y≥0, -T /2+3η/2≤t-y≤-T /2+3η e -γ((t-y) 2 -r 2 ) dy C u * (x, •) 2 H 1 (0,3T ) ≤ Ce -mγT u * (x, •) 2 H 1 (0,3T ) , ∀(x, t, r) ∈ Q. Indeed, (t -y) 2 ≥ -T 2 + 3η 2 and r 2 ≤ 9η 2 therefore (t -y) 2 -r 2 ≥ T 4 (T -6η
). Thus, if we suppose T > 1, since η < 1/6, then there exists m > 0 independent of T such that (t -y) 2 -r 2 ≥ mT . We can do the same work for the derivatives of F γ . Likewise, we have

u γ 2 L 2 ( e Q) = Ω 5T /2 -T /2 3η -3η γ 2π +∞ 0 e -γ 2 (t+ir-y) 2 θ(t -y)u * (x, y)dy 2 dr dt dx with, ∀(x, t, r) ∈ Q, +∞ 0 e -γ 2 (t+ir-y) 2 θ(t -y)u * (x, y)dy 2 ≤ t+T /2 0 |e -γ 2 (t+ir-y) 2 | 2 dy t+T /2 0 |θ(t -y)u * (x, y)| 2 dy ≤ 3T 0 e -γ((t-y) 2 -r 2 ) dy C u * (x, •) 2 L 2 (0,3T ) ≤ Ce Cγ u * (x, •) 2 L 2 (0,3T ) .
And we use the same procedure for the derivatives of u γ to obtain the desired result.

Prove a Carleman estimate

We write a Carleman estimate for the system of elliptic integro-differential equations [START_REF] Robbiano | Fonction de coût et contrôle des solutions des équations hyperboliques[END_REF], applying the same methods we used in the hyperbolic case. That is, we combine the decoupling of the equations proposed by [START_REF] Bellassoued | Inverse problem of determining the density and the Lamé coefficients by boundary data[END_REF] with the change of variables of [START_REF] Lorenzi | Recovering a Lamé kernel in a viscoelastic system[END_REF], and we use a classical Carleman estimate for a scalar elliptic equation (Lemma 2). The difference is that here we have 3 variables (x, t, r), so the weight function has to be modified. Let us introduce first some notations:

• Domains For T > 0 and η > 0, we note

Q = Ω × (- T 2 , 5T 2 ) × (-3η, 3η), Σ = ∂Ω × (- T 2 , 5T 2
) × (-3η, 3η).

• Norms

2 H k,σ ( e Q) = |α|≤k σ 2(k-|α|) ∂ α 2 L 2 ( e Q) , 2 H k,σ ( e Σ) = |α|≤k σ 2(k-|α|) ∂ α 2 L 2 ( e Σ) .
Definition 4 (a Carleman weight function) Let x 0 ∈ R 3 \ Ω and ξ > 0. We introduce a function ψ in the following way:

ψ(x, t, r) = |x -x 0 | 2 + (t + T 2 ) 2 + r 2 , ∀(x, t, r) ∈ Q,
and we set ϕ(x, t, r) = e -ξψ(x,t,r) , ∀(x, t, r) ∈ Q.

Because we introduced this new weight function, we must now check that we still have a result, similar to to one of Lemma 1, to bound the integral terms in the proof of the Carleman estimate. It is given by the following lemma: Lemma 6 Let ϕ be defined by ( 28), σ > 0 and u ∈ L 2 ( Q). Then, there exists a constant C > 0, independent of σ, such that

e Q t -T /2 |u(x, s, r)|ds 2 e 2σϕ(x,t,r) dx dt dr ≤ C σ e Q |u(x, t, r)| 2 e 2σϕ(x,t) dx dt dr.
Proof Using the Cauchy-Schwarz inequality, we can write:

e Q t -T /2
|u(x, s, r)|ds 2 e 2σϕ(x,t,r) dx dt dr

≤ e Q t -T /2 |u(x, s, r)| 2 ds (t + T
2 )e 2σϕ(x,t,r) dx dt dr.

We notice then that, ∀(x, t, r) ∈ Q,

(t + T 2 )e 2σϕ(x,t,r) ≤ - e ξψ(x,t,r) 4ξσ ∂ t e 2σϕ(x,t,r) ≤ - C σ ∂ t e 2σϕ(x,t,r) .
Therefore,

e Q t -T /2
|u(x, s, r)|ds 2 e 2σϕ(x,t,r) dx dt dr

≤ Ω 3η -3η 5T /2 -T /2 - C σ t -T /2
|u(x, s, r)| 2 ds ∂ t (e 2σϕ(x,t,r) )dt dx dr

≤ C σ Ω 3η -3η 5T /2 -T /2
|u(x, t, r)| 2 e 2σϕ(x,t,r) dt dx dr -e 2σϕ(x,5T /2,r)

5T /2 -T /2
|u(x, s, r)| 2 ds dx dr .

Thus, we deduce the result.

We are now ready to state the Carleman estimate for the operator Q:

Theorem 4 (Carleman estimate) Let Q be the operator defined by [START_REF] Robbiano | Fonction de coût et contrôle des solutions des équations hyperboliques[END_REF]. Let K be a compact set in Ω × (-T 2 , 5T 2 ) × (-3η, 3η). We assume that (H1) (λ, µ) ∈ C 2 (Ω) 2 and ( λ, μ) ∈ C 2 (Ω × (0, +∞)) 2 .

Then, there exists ξ 0 > 0 such that, for all ξ ≥ ξ 0 , there exists σ 0 > 0 such that, for all σ ≥ σ 0 and for all u ∈ C ∞ 0 (K) 3 , we have the following estimate:

1 σ ue σϕ 2 H 2,σ ( e Q) ≤ C (Qu)e σϕ 2 H 1 ( e Q) + ue σϕ 2 H 2,σ ( e Σ) ,
where ϕ is defined by [START_REF] Imanuvilov | Global Carleman estimates for weak solutions of elliptic nonhomogeneous Dirichlet problems[END_REF], and C > 0 depends on the C 2 (Q)-norm of the coefficient μ but is independent of σ.

Proof Let K be a compact set in Ω × (-T 2 , 5T 2 ) × (-3η, 3η), let u ∈ C ∞ 0 (K) 3 be the solution of: Qu(x, t, r) = S(x, t, r), ∀(x, t, r) ∈ Q, (29)
and let assume (H1). We must first decouple the equations. To this end, we take the divergence of ( 29) and we introduce the scalar v = ∇ • u and the vector w = ∇ ∧ u. We obtain the following system of seven scalar equations only coupled at the first order:

-2∂ 2 r u(x, t, r) -∂ 2 t u(x, t, r) -µ(x)∆u(x, t, r) + t -T /2 μ(x, t -s)∆u(x, s, r)ds = S(x, t, r) + A 1 (u, v)(x, t, r), ∀(x, t, r) ∈ Q, -2∂ 2 r v(x, t, r) -∂ 2 t v(x, t, r) -(λ + 2µ)(x)∆v(x, t, r) + t -T /2
( λ + 2μ)(x, t -s)∆v(x, s, r)ds

= (∇ • S)(x, t, r) + A 2 (u, v, w)(x, t, r), ∀(x, t, r) ∈ Q, -2∂ 2 r w(x, t, r) -∂ 2 t w(x, t, r) -µ(x)∆w(x, t, r) + t -T /2 μ(x, t -s)∆w(x, s, r)ds = (∇ ∧ S)(x, t, r) + A 3 (u, v, w)(x, t, r), ∀(x, t, r) ∈ Q, (30)
where A 1 , A 2 and A 3 are first order integrodifferential operators. We only treat the first equation of ( 30), the others equations can be treated accordingly. We introduce the following change of variable:

ũ(x, t, r) = µ(x)u(x, t, r) - t -T /2 μ(x, t -s, r)u(x, s, r)ds, ∀(x, t, r) ∈ Q.
Thus, ũ satisfies the classical elliptic equation

-2∂ 2 r ũ(x, t, r) -∂ 2 t ũ(x, t, r) -µ(x)∆ũ(x, t, r) = µ(x) (S(x, t, r) + A 1 (u, v)(x, t, r)) + L 1 (u)(x, t, r),
where L 1 is a first order integro-differential operator with bounded coefficients. We now need a Carleman estimate for a scalar elliptic equation. It is given by the following lemma:

Lemma 7 ( [START_REF] Imanuvilov | Global Carleman estimates for weak solutions of elliptic nonhomogeneous Dirichlet problems[END_REF] and [START_REF] Lebeau | Contrôle exact de l'équation de la chaleur[END_REF])

Let Q be an open domain in R n , K be a compact set in Q and ψ be a C ∞ (Q) function satisfying ∇ψ(x) = 0, ∀x ∈ K. Let ϕ(x) = e -ξψ(x) , ∀x ∈ Q,
where ξ > 0 is sufficiently large. We consider the scalar second-order elliptic operator

R(x) = a(x) : ∇ 2 +b(x) • ∇ +c , ∀x ∈ Q,
where all the coefficients are C 2 (Q). Then, there exist σ 0 > 0 and C > 0 such that, for all σ ≥ σ 0 and for all u ∈ C ∞ 0 (K), the following Carleman estimates hold true:

C σ ue σϕ 2 H 2,σ (Q) ≤ (Ru)e σϕ 2 L 2 (Q) + ue σϕ 2 H 2,σ (∂Q) , Cσ ue σϕ 2 H 1,σ (Q) ≤ (Ru)e σϕ 2 L 2 (Q) + σ ue σϕ 2 H 1,σ (∂Q)
. As the weight function ϕ defined in [START_REF] Imanuvilov | Global Carleman estimates for weak solutions of elliptic nonhomogeneous Dirichlet problems[END_REF] satisfies the hypothesis of 7, we can apply the first inequality to ũ:

1 σ ũe σϕ 2 H 2,σ ( e Q) ≤ C (-2∂ 2 r ũ -∂ 2 t ũ -µ∆ũ)e σϕ 2 L 2 ( e Q) + ũe σϕ 2 H 2,σ ( e Σ) .
We take into account the fact that the coefficients λ, λ, µ, μ are bounded in Q, that the Lemma 6 holds and that all the terms of the operators A 1 and L 1 are of order less or equal than 1. Then, we obtain

1 σ ũe σϕ 2 H 2,σ ( e Q) ≤ C Se σϕ 2 L 2 ( e Q) + ue σϕ 2 H 1,σ ( e Q) + ve σϕ 2 H 1,σ ( e Q) + ue σϕ 2 H 2,σ ( e Σ)
. We use again the change of variable

u(x, t, r) = 1 µ(x) ũ(x, t, r) + t -T /2 μ(x, t -s) µ(x) u(x, s, r)ds, ∀(x, t, r) ∈ Q.
Hence,

e Q |u(x, t, r)| 2 e 2σϕ(x,t,r) dx dt dr ≤ C e Q |ũ(x, t, r)| 2 e 2σϕ(x,t,r) dx dt dr + C σ e Q
|u(x, t, r)| 2 e 2σϕ(x,t,r) dx dt dr.

For σ sufficiently large, the second term is absorbed. We do the same for the derivatives of u, which leads to

1 σ ue σϕ 2 H 2,σ ( e Q) ≤ C 1 σ ũe σϕ 2 H 2,σ ( e Q) . Finally, 1 σ ue σϕ 2 H 2,σ ( e Q) ≤ C Se σϕ 2 L 2 ( e Q) + ue σϕ 2 H 1,σ ( e Q) + ve σϕ 2 H 1,σ ( e Q) + ue σϕ 2 H 2,σ ( e Σ) .
Thanks to the Carleman weights, we absorb the second term of the right hand side in the left hand side. We make the same work for v and w but using the second Carleman inequality of the Lemma 7 in order to obtain

σ ve σϕ 2 H 1,σ ( e Q) ≤ C (∇ • S)e σϕ 2 L 2 ( e Q) + ue σϕ 2 H 1,σ ( e Q) + we σϕ 2 H 1,σ ( e Q) + σ ve σϕ 2 H 1,σ ( e Σ)
σ we σϕ 2

H 1,σ ( e Q) ≤ C (∇ ∧ S)e σϕ 2 L 2 ( e Q) + ue σϕ 2 H 1,σ ( e Q) + ve σϕ 2 H 1,σ ( e Q) + σ we σϕ 2 H 1,σ ( e Σ) . Hence, σ ve σϕ 2 H 1,σ ( e Q) ≤ C (∇S)e σϕ 2 L 2 ( e Q) + ue σϕ 2 H 1,σ ( e Q) + we σϕ 2 H 1,σ ( e Σ) + σ ve σϕ 2 H 1,σ ( e Σ) . Finally, 1 σ ue σϕ 2 H 2,σ ( e Q) ≤ C Que σϕ 2 L 2 ( e Q) + 1 σ (∇Qu)e σϕ 2 L 2 ( e Q) + 1 σ ue σϕ 2 H 1,σ ( e Q) + ue σϕ 2 H 2,σ ( e Σ) .
And thanks to the Carleman weights, we absorb the third term of the right hand side in the left hand side to obtain the result:

1 σ ue σϕ 2 H 2,σ ( e Q) ≤ C Que σϕ 2 H 1 ( e Q) + ue σϕ 2 H 2,σ ( e Σ)
and this concludes the proof of Theorem 4.

Obtain local estimations

We are now going to apply the Carleman estimate we proved in Theorem 4 to the function u γ . To ensure that its support is in a compact set of Ω×(-T 2 , 5T 2 )×(-3η, 3η), we multiply it by a cut-off function. A good choice of this cut-off function leads to two local estimations. Summing these estimations, we obtain the desired interpolation result. This work is inspired from [START_REF] Bellassoued | Unicité et contrôle pour le système de Lamé[END_REF] which uses some results of [START_REF] Robbiano | Fonction de coût et contrôle des solutions des équations hyperboliques[END_REF] but noticing that we have here an additional variable.

First estimation

Let Γ ⊂ ∂Ω be arbitrarily small. We begin by estimating u γ in a ball, which is close to Γ. Let us choose 0 < η < ε and

x (0) ∈ R 3 \ Ω such that B(x (0) , η) ∩ Ω = ∅, B(x (0) , 2η) ∩ Ω = ∅, B(x (0) , 4η) ∩ ∂Ω ⊂ Γ.
Figure 10 illustrates this construction. We define then

ψ (0) (x, t, r) = |x -x (0) | 2 + η 2 T 2 (t + T 2 ) 2 + r 2 , ∀(x, t, r) ∈ Q,
and we set

ϕ (0) (x, t, r) = e -ξ η 2 ψ (0) (x,t,r) , ∀(x, t, r) ∈ Q. Ω 4η x (0) 2η Γ Figure 10
. Definition of x (0) and η.

As in Figure 11, we introduce the cut-off function

χ 4 ∈ C ∞ 0 (R) such that χ 4 (z) =      0, if z < 1 2 , z > 8, 1, if 3 4 < z < 7,
and we set Here,

u * γ (x, t, r) = χ 4 ψ (0) η 2 u γ (x, t, r), ∀(x, t, r) ∈ Q.
u * γ has a compact support in Ω × (-T 2 , 5T 2 ) × (-3η, 3η) since ψ (0) (x, t, ±3η) > 10η 2 =⇒ u * γ (x, t, r) = 0 if |r| ≥ 3η, ψ (0) (x, 5T 2 , r) > 10η 2 =⇒ u * γ (x, t, r) = 0 if t ≥ 5T 2 ,
and

u γ (x, t, r) = 0 if t < - T 2 + 3η 2 .
Therefore, we can apply the Carleman estimate we showed in Theorem 4 to this function:

1 σ u * γ e σϕ (0) 2 H 2,σ ( e Q) ≤ C (Qu * γ )e σϕ (0) 2 H 1 ( e Q) + u * γ e σϕ (0) 2 H 2,σ ( e Σ) .
The last norm is actually in Γ = Γ × (-T 2 , 5T 2 ) × (-3η, 3η) ⊂ Σ. Indeed, we supposed B(x (0) , 4η) ∩ ∂Ω ⊂ Γ, therefore, if x ∈ ∂Ω \ Γ, then u * γ (x, t, r) = 0. In addition, ∀(x, t, r) ∈ Q, Qu * γ (x, t, r)

= χ 4 ψ (0) η 2 Qu γ (x, t, r) + Q, χ 4 ψ (0) η 2 u γ (x, t, r) = χ 4 ψ (0) η 2 (F γ (x, t, r) + G γ (x, t, r)) + Q, χ 4 ψ (0) η 2 u γ (x, t, r) = χ 4 ψ (0) η 2 F γ (x, t, r) + Q, χ 4 ψ (0) η 2 u γ (x, t, r),
since we saw in [START_REF] Bros | Tuboïdes et structure analytique des distributions. ii. Support essentiel et structure analytique des distributions[END_REF] that G γ (x, t, r) = 0, ∀x ∈ Ω \ Ω(3ε). And here, χ 4

ψ (0) η 2 is different from zero if |x -x (0) | 2 ≤ 8η 2 , thus, in particular, if |x -x (0) | ≤ 3η < 3ε.
Then, we can write:

1 σ u * γ e σϕ (0) 2 H 2,σ ( e Q) ≤ C (Qu * γ )e σϕ (0) 2 H 1 ( e Q) + u * γ e σϕ (0) 2 H 2,σ ( e Γ) =⇒ 1 σ e 2σe -6ξ u γ 2 H 2,σ ({η 2 ≤ψ (0) ≤6η 2 }∩ e Q) ≤ C σ 2 e 2σe -7ξ u γ 2 H 2 ( e Q) + σ 2 e 2σe -ξ 2 F γ 2 H 1 ( e Q) + σ 4 e 2σe -ξ 2 u γ 2 H 2 ( e Γ) ,
We choose then ρ and x (1) such that dist(x (1) , ∂Ω) ≥ 4ρ, B 1 = B(x (1) , ρ) × (-

T 2 , T 2 ) × (-ρ, ρ) ⊂ {η 2 ≤ ψ (0) ≤ 6η 2 }. 2η x (0)
Γ Ω B(x (1) , ρ)

Figure 12. Definition of x (1) and ρ.

This choice is valid since, in B 1 , we have

ψ (0) (x, t, r) ≤ |x -x (1) | 2 + |x (1) -x (0) | 2 + η 2 T 2 (t + T 2 ) 2 + r 2 ≤ ρ 2 + |x (1) -x (0) | 2 + η 2 + ρ 2
and since B(x (0) , 2η) ∩ Ω = ∅, we can choose x (1) such that |x (1) -x (0) | 2 < 4η 2 and 4ρ ≤ η so that ψ (0) (x, t, r) ≤ 6η 2 . Then, if σ is sufficiently large, we have

u γ 2 H 2 ( e B1) ≤ Ce -C1σ u γ 2 H 2 ( e Q) + Ce C2σ F γ 2 H 1 ( e Q) + u γ 2 H 2 ( e Γ) .
We minimize with respect to σ (cf. [START_REF] Robbiano | Fonction de coût et contrôle des solutions des équations hyperboliques[END_REF]) to obtain, with

ν 0 = C 2 C 1 + C 2 : u γ H 2 ( e B1) ≤ C u γ H 2 ( e Q) 1-ν0 F γ H 1 ( e Q) + u γ H 2 ( e Γ) ν0 . (31)

Second estimation

Now, we extend the estimation in B(x (1) , ρ) into Ω(ε) \ Ω(2ε). Let B(x (j) , ρ), 2 ≤ j ≤ N be a cover of Ω(ε)\Ω(2ε). We suppose that x (j) is such that for 2 ≤ j ≤ N , dist(x (j) , ∂Ω) ≥ 4ρ, and B(x (j+1) , ρ) ⊂ B(x (j) , 2ρ).

00000000 00000000 00000000 00000000 00000000 00000000 00000000 00000000 11111111 11111111 11111111 11111111 11111111 We introduce

11111111 11111111 11111111 ε B j x (0) 2η B 1 Ω(ε) \ Ω(2ε) B(x (j) , 2ρ) B(x (j+1) , ρ) B(x (j) , ρ)
B j = B(x (j) , ρ) × (- T 2 , T 2 ) × (-ρ, ρ).
This construction is illustrated in Figure 13. We define

ψ (j) (x, t, r) = |x -x (j) | 2 + ρ 2 T 2 (t + T 2 ) 2 + r 2 , ∀(x, t, r) ∈ Q, and
ϕ (j) (x, t, r) = e -ξ ρ 2 ψ (j) (x,t,r) , ∀(x, t, r) ∈ Q.
We set

u * γ (x, t, r) = χ 4 ψ (j) ρ 2 u γ (x, t, r), ∀(x, t, r) ∈ Q,
to which we apply Theorem 4, without the boundary term, since the support of u * γ is an interior domain:

1 σ u * γ e σϕ (j) 2 H 2,σ ( e Q) ≤ C (Qu * γ )e σϕ (j) 2 H 1 ( e Q) .
In addition, ∀(x, t, r) ∈ Q,

Qu * γ (x, t, r) = χ 4 ψ (j) ρ 2 F γ (x, t, r) + Q, χ 4 ψ (j) ρ 2 u γ (x, t, r).
Thus, we have

C σ e 2σe -6ξ u γ 2 H 2,σ ({ρ 2 ≤ψ (j) ≤6ρ 2 }∩ e Q) ≤ σ 2 e 2σe -ξ 2 u γ 2 H 2 (ψ (j) ≤ρ 2 ) + σ 2 e 2σe -7ξ u γ 2 H 2 (ψ (j) ≤8ρ 2 ) + σ 2 e 2σe -ξ 2 F γ 2 H 1 ( e Q)
, and, if we choose σ large, we can write:

e 2σe -6ξ u γ 2 H 2,σ ({ψ (j) ≤6ρ 2 }∩ e Q) ≤ Ce 2σe -ξ 3 u γ 2 H 2 (ψ (j) ≤ρ 2 ) + Ce 2σe -13ξ 2 u γ 2 H 2 (ψ (j) ≤8ρ 2 ) + Ce 2σe -ξ 3 F γ 2 H 1 ( e Q)
. This allows to conclude that

u γ 2 H 2 ( e Bj+1) ≤ Ce C3σ u γ 2 H 2 ( e Bj ) + F γ 2 H 1 ( e Q) + Ce -C4σ u γ 2 H 2 ( e Q) ,
since B j+1 ⊂ {ψ (j) ≤ 6ρ 2 } and {ψ (j) ≤ ρ 2 } ⊂ B j . We minimize with respect to σ to obtain, with

ν 1 = C 4 C 3 + C 4 : u γ H 2 ( e Bj+1) ≤ C u γ H 2 ( e Q) 1-ν1 F γ H 1 ( e Q) + u γ H 2 ( e Bj ) ν1
We use then the recurrence result of the following lemma:

Lemma 8 (Lemma 4 in [START_REF] Lebeau | Contrôle exact de l'équation de la chaleur[END_REF]) Let α j > 0 satisfying, for all j ≥ 0,

α j ≤ B 1-ν (α j-1 + A) ν and α j ≤ B
where A > 0, B > 0 and ν ∈]0, 1[. Then, for all µ ∈]0, ν N [, we have:

α N ≤ 2 1/(1-ν) B 1-µ (α 0 + A) µ .
Here,

α j = u γ H 2 ( e Bj+1) , A = F γ H 1 ( e Q) , B = u γ 2 H 2 ( e Q)
. Therefore, we obtain

u γ H 2 ( e Bn) ≤ C u γ H 2 ( e Q) 1-ν F γ H 1 ( e Q) + u γ H 2 ( e B1) ν
.

We apply the Young inequality and we obtain

u γ H 2 ( e Bn) ≤ q u γ H 2 ( e Q) + -q F γ H 1 ( e Q) + u γ H 2 ( e B1) ,
where q = 1 1 -ν and q = 1 ν . We use the estimations ( 26) and ( 27) on F γ and u γ to write:

u γ H 2 ( e Bn)
≤ q e M γ u * H 2 (Ω×(0,3T )) + -q e -mγT u * H 2 (Ω×(0,3T )) + u γ H 2 ( e B1) .

We choose then = e -2M γ/q , so that

u γ H 2 ( e Bn) ≤ e -M γ u * H 2 (Ω×(0,3T ))
+e -(mT -2M q q )γ u * H 2 (Ω×(0,3T )) + e 2M γ q q u γ H 2 ( e B1) . If we fix T > T n with mT n -2M q q = M, we obtain, with κ = 2M γ q q , u γ H 2 ( e Bn) ≤ e -M γ u * H 2 (Ω×(0,3T )) + e κγ u γ H 2 ( e B1) .

We can do the same work on the first estimation [START_REF] Modecai | Nonlinear Programming: Analysis and Methods[END_REF], the one in B 1 . We apply the Young inequality and the estimations ( 26) and ( 27) to obtain

u γ H 2 ( e B1)
≤ q0 e M γ u * H 2 (Ω×(0,3T )) + -q 0 e -mγT u * H 2 (Ω×(0,3T )) + u γ H 2 ( e Γ) .

where q 0 = 1 1 -ν 0 and q 0 = 1 ν 0 . We choose then = e -(2M +κ)γ/q0 , so that

u γ H 2 ( e B1) ≤ e -(M +κ)γ u * H 2 (Ω×(0,3T )) +e -(mT -(2M +κ) q 0 q 0 )γ u * H 2 (Ω×(0,3T )) + e (2M +κ)γ q 0 q 0 u γ H 2 ( e Γ) . If we fix T > T 0 with mT 0 -(2M + κ) q 0 q 0 = M + κ, we obtain u γ H 2 ( e B1) ≤ e -(M +κ)γ u * H 2 (Ω×(0,3T )) + e Cγ u γ H 2 ( e Γ) .

And we regroup the two estimations in

u γ H 2 ( e Bn) ≤ e -M γ u * H 2 (Ω×(0,3T )) + e Cγ u * 2 H 2 (Γ×(0,3T ))
. Finally, we take T > max n T n and we sum the inequalities to obtain

u γ 2 H 2 (Ω(ε)\Ω(2ε)×(-T 2 , T
2 )×(-ρ,ρ)) ≤ e -Cγ u * 2 H 2 (Ω×(0,3T )) + e Cγ u * 2 H 2 (Γ×(0,3T )) . We can now go back to the variable u in the estimations.

Return to the variable u

We introduce the new variable

w γ (x, t) = u γ (x, t, r = 0) = γ 2π ∞ 0 e -γ 2 (t-y) 2 θ(t -y)u * (x, y)dy = (K γ * u * )(x, t), ∀(x, t) ∈ Ω × R, with K γ (t) = γ 2π e -γ 2 t 2 θ(t), ∀t ∈ R.
We see that w γ converges to u * when γ tends to infinity. We have

w γ 2 H 2 (Ω(ε)\Ω(2ε)×(0, T 2 -ρ)) ≤ C u γ 2 H 2 (Ω(ε)\Ω(2ε)×(-T 2 , T 2 )×(-ρ,ρ)) ≤ e -Cγ u * 2 H 2 (Ω×(0,3T )) + e Cγ u * 2 H 2 (Γ×(0,3T ))
. Indeed, the Cauchy formula says that, for 0 < d < ρ, we have

w γ (x, a) = 1 2iπ |w-a|=d w γ (x, w) w -a dw,
which, using polar coordinates, implies that

|w γ (x, a)| 2 ≤ C 2π 0 |w γ (x, a + de iθ )| 2 dθ.
We integrate then, for 0 < d < ρ, into

|w γ (x, a)| 2 ≤ C ρ ρ 0 2π 0 |w γ (x, a + de iθ )| 2 dθ dρ,
that we can write:

|w γ (x, a)| 2 ≤ C ρ -ρ |t-a|≤ρ |w γ (x, t + ir)| 2 dt dr = C ρ -ρ |t-a|≤ρ |u γ (x, t, r)| 2 dt dr,
according to the definition of w γ . Then, we integrate for x ∈ Ω(ε) \ Ω(2ε):

w γ (., a) 2 L 2 (Ω(ε)\Ω(2ε)) ≤ C ρ -ρ |t-a|≤ρ u γ (., t, r) 2 L 2 (Ω(ε,2ε)) dt dr ≤ C u γ 2 L 2 (Ω(ε,2ε)×(-T 2 , T 2 )×(-ρ,ρ)) . We integrate for a ∈ (0, T 2 -ρ): w γ 2 L 2 (Ω(ε)\Ω(2ε)×(0, T 2 -ρ)) ≤ C u γ 2 L 2 (Ω(ε,2ε)×(-T 2 , T 2 )×(-ρ,ρ))
and we do the same for the derivatives of w γ . We can now come back to u * (and u) using the classical Fourier transform (noted with a big hat symbol) because we notice that

w γ (x, t) = (K γ * u * )(x, t) =⇒ w γ (x, τ ) = K γ (τ ) u * (x, τ ),
where

K γ (τ ) = γ 2π +∞ -∞ e -iτ t e -γ 2 t 2 θ(t)dt = γ 2π +∞ -∞ e -iτ t e -γ 2 t 2 dt + γ 2π +∞ -∞ e -iτ t e -γ 2 t 2 (θ(t) -1) dt. We have γ 2π +∞ -∞ e -iτ t e -γ 2 t 2 dt = γ 2π +∞ -∞ e - " √ γ 2 t+ q 1 2γ iτ " 2 e -τ 2 2γ dt = 1 √ π e -τ 2 2γ +∞ -∞ e -z 2 dz = e -τ 2 2γ .
Hence,

K γ (τ ) = e -τ 2 2γ + γ 2π -T 2 +3η -∞ e -iτ t e -γ 2 t 2 (θ(t) -1) dt.
Therefore, we have

|(1 -K γ )(τ )| ≤ 1 -e -τ 2 2γ + γ 2π -T 2 +3η -∞ e -iτ t e -γ 2 t 2 |θ(t) -1| dt ≤ τ 2 2γ + γ 2π -T 2 +3η -∞ e -γ 2 t 2 dt,
where we used for the first term the fact that the function 1 -e -t 2 t 2 is bounded in R.

We win one order of convergence since we have 1 γ and no 1 √ γ as in [START_REF] Bellassoued | Inverse problem of determining the density and the Lamé coefficients by boundary data[END_REF]. Then, we can introduce in the second term the change of variable z = e -γ 2 t 2 , that is t = -2 γ ln(z), and dz = -γ t e -γ 2 t 2 dt, and so we have:

γ 2π -T 2 +3η -∞ e -γ 2 t 2 dt ≤ 1 2 √ π e - γ 2 ( -T 2 +3η ) 2 0 dz -ln(z) ≤ 1 2 √ π e -γ 2 (-T 2 +3η) 2 γ 2 -T 2 + 3η 2 ≤ C γ .
Then,

u * -w γ 2 L 2 (Ω(ε)\Ω(2ε)×(0, T 2 -ρ)) ≤ u * -w γ 2 L 2 (Ω×R) = u * -w γ 2 L 2 (Ω×R) ≤ Ω +∞ -∞ τ 2 + C γ 2 | u * (x, τ )| 2 dτ ≤ C γ 2 Ω +∞ -∞ (| u * (x, τ )| 2 + | ∂ 2 t u * (x, τ )| 2 )dτ ≤ C γ 2 u * 2 H 2 (Ω×R) .
Finally, if we suppose that u * is prolonged by zero outside (0, 3T ), we have

u * -w γ L 2 (Ω(ε)\Ω(2ε)×(0, T 2 -ρ)) ≤ C γ u * H 2 (Ω×(0,3T )) . Thus, u * H 2 (Ω(ε)\Ω(2ε)×(0, T 2 -ρ)) ≤ C u * -w γ H 2 (Ω(ε)\Ω(2ε)×(0, T 2 -ρ)) + w γ H 2 (Ω(ε)\Ω(2ε)×(0, T 2 -ρ)) ≤ C γ u * H 4 (Ω×(0,3T )) + C w γ H 2 (Ω(ε)\Ω(2ε)×(0, T 2 -ρ))
. Coming back to u and taking into account that u * (x, t) = χ 3 (x)u(x, t), so that u = u * in Ω(ε) \ Ω(2ε), we have

u H 2 (Ω(ε)\Ω(2ε)×(0, T 2 -ρ)) ≤ C γ u H 4 (Ω×(0,3T )) + e Cγ u H 2 (Γ×(0,3T )) .
We use the fact that u H 4 (Ω×(0,3T )) is bounded and that

u 2 H 2 (Γ×(0,3T )) = 2 |α|=1 ∂ α x u 2 L 2 (Γ×(0,3T )) ,
since û and its time derivatives vanish on Γ, to write:

u 2 H 2 (Ω(ε)\Ω(2ε)×(0, T 2 -ρ)) ≤ C 5 γ 2 + e C6γ 2 |α|=1 ∂ α x u 2 L 2 (Γ×(0,3T )) .
Then, we choose γ such that the first term dominate the second. We can for example choose

γ = 1 2C 6 log 2 + C 5 2 |α|=1 ∂ α x u 2 L 2 (Γ×(0,3T ))
.

We achieve the proof of Theorem 3.

Numerical results

To illustrate the theoretical result [START_REF] De Buhan | Un résultat de stabilité pour la récupération d'un paramètre du système de la viscoélasticité 3[END_REF], we present here one numerical example in 2D, even if the result is also valid in 3D. The numerical resolution that we propose leads to recover the values of the unknown coefficient p at each vertex of the discretized domain.

Since the linear viscoelastic system ( 1) is a simplified model for the mechanical behavior of brain structures [START_REF] Sarron | Identification of linear viscoelastic constitutive models[END_REF], we could use the method to retrieve the localization of a brain tumor. However, this example can be considered as an academic test case, sufficient to emphasize the motivation of the study, to illustrate the stability result and to propose an efficient numerical method.

Direct problem

We consider now the system (1)-( 3), in two dimensions, with the following coefficients:

• µ(x) = λ(x) = 1200, h(t) = e -t/τ with τ = 1, λ(x, t) = 400 h(t),

• p(x) = 400, in the healthy tissue, > 400, in the tumor (cf. Figure 14),

• f (x, t) = 0, ū1 (x) = 0 and ū0 is the solution of the stationary problem associated to (1):

-∇ • µ(x)(∇ū 0 (x) + ∇ū 0 (x) T ) + λ(x)(∇ • ū0 )(x)I = 1, ∀x ∈ Ω, u 0 (x) = 0, ∀x ∈ ∂Ω.

When dealing with a real experiment, it is difficult to ensure that the initial data ū0 satisfies hypothesis (H5). In practice, the recuperation is carry out with ū0 = 0. Moreover, we assume that p is known in the boundary of Ω. This hypothesis, less restrictive than hypothesis (H3), is more likely to be acknowledged because experimental measurements of the coefficient on the boundary may be possible. We can solve numerically the direct problem (1)-( 3) by discretizing the equations • in space using P 1 Lagrange Finite Elements in the mesh shown on Figure 1, • in time using a θ-scheme with θ = 0.5 (implicit centered scheme) and δt = 1,

• by using the trapezium formula for the integral term.

The time of observation is taken equal to T = 50 whereas the observation zone ω is the one shown on Figure 1. A uniform relative error of 2% corresponding to the experimental error is added to the solution i.e. the observed displacement measured experimentally u obs is obtained by the formula: u obs (x, t) = (1 + r 0.02)ū(x, t), ∀(x, t) ∈ ω × (0, T ), where r is a constant in [-1, 1], randomly chosen (according to a uniform distribution law) at each discretization point in the domain ω × (0; T ).

Inverse problem

We are looking for the minimizer of the non quadratic functional

J(p) = 1 2 T 0 ω |u -u obs | 2 + |∇(u -u obs )| 2 dx dt,
with u = M(p), M being the nonlinear operator from P to U associated to system (1)-( 3). In general, u obs ∈ M(P ). Notice that adding a Tikhonov regularizing term to J does not numerically help in this case since we will filter high frequencies when choosing the parameter space for p (see section 5.4). We solve the minimization problem by a Broyden-Fletcher-Goldfarb-Shanno (BFGS) algorithm [START_REF] Modecai | Nonlinear Programming: Analysis and Methods[END_REF]. Thus, we calculate ∇J(p, δp) = lim (

) 32
We introduce the adjoint operator P * of (2)

P * v(x, t) = ∂ 2 t v(x,

The results

We solve the inverse problem by the method introduced in 5.2-5.4. On Figure 17, we plot the relative error between p * and p in L 2 -norm with respect to the number K of eigenfunctions in the basis and for different steps of the iterative process. On Figure 18, we show the exact coefficient and the numerical result we obtain at different steps.

Perspectives

Notice that we only try to recover the coefficient p, assuming that all the other coefficients are known. Indeed, the aim of this example is to illustrate the main result of the paper, which is the theoretical stability for recovering the spatial part of only one viscoelastic coefficient from only one measurement. The other objective of this example is to show the interesting choice we can do for the spectral basis in an adaptive way. To recover simultaneously all the coefficients, there are various methods (see [START_REF] Pagnacco | Inverse strategies for the identification of elastic and viscoelastic material parameters using full-field measurements[END_REF]). One method consists in adding to the minimal criteria more parameters, but this could be unstable if only a single observation is available, so we would need more independent observations. Another method consists in using an iterative procedure fixing all coefficients except one at each step. In both cases, we should start in a neighborhood of the real coefficients and play with noise levels. All this deep analysis of the numerical simulations is of course very relevant but it is out of the scope of the present publication. We are preparing a forthcoming publication that will include these numerical issues and more details about the choice of the spectral adaptive basis.

Figure 1 .

 1 Figure 1. The Carleman weight function ϕ(x, t) (left) and its isovalues ϕ(x, t) = δ (right) for x 0 = 0 and β = 1.

Figure 4 .

 4 Figure 4. Definition of the cut-off function χ 1 .

Figure 5 .

 5 Figure 5. Definition of the domain Ωt(ε, δ) and the time t δ x .

Figure 6 .

 6 Figure 6. Definition of the cut-off function χ 2 .

Figure 8

 8 Figure 8 illustrates this cut-off function. Then, we set

3. 4 .

 4 Return to the initial variable u 3.4.1. Upper bound for the right hand side We take into account the fact that

Figure 8 .

 8 Figure 8. Definition of the neighborhood ω of ∂Ω and the cut-off function χ 3 .

Figure 9 .

 9 Figure 9. Definition of the cut-off function θ.

4 Figure 11 .

 411 Figure 11. Definition of the cut-off function χ 4 .

Figure 13 .

 13 Figure 13. Definition of the balls B(x (j) , ρ).

Figure 14 .

 14 Figure 14. Unknown coefficient p (on the left), computational mesh and initial data ū0 (in the center), observation zone ω (in red, on the right).

M 0 ∇

 0 p + εδp) -J(p)) p δp • ((u -u obs ) -∆(u -u obs )) dx dt,with M p the linearized operator of M around p, i.e. M p δp = δu satisfies • δp(x)h(t -s)(∇u(x, s) + ∇u(x, s) T ) ds, δu(x, 0) = 0, ∀x ∈ Ω, ∂ t δu(x, 0) = 0, ∀x ∈ Ω, δu(x, t) = 0, ∀(x, t) ∈ ∂Ω × (0, T ).

P 0 Ω

 0 t) -∇ • µ(x)(∇v(x, t) + ∇v(x, t) T ) + λ(x)(∇ • v(x, t))I + T t ∇ • p(x)h(s -t)(∇v(x, s) + ∇v(x, s) T ) + λ(x, s -t)(∇ • v)(x, s)I ds ,and we define δu * , the solution of * δu * (x, t) = (u -u obs )(x, t) -∆(u -u obs )(x, t), ∀x ∈ ω, 0, ∀x ∈ (Ω \ ω), δu * (x, T ) = 0, ∀x ∈ Ω, ∂ t δu * (x, T ) = 0, ∀x ∈ Ω, δu * (x, t) = 0, ∀(x, t) ∈ ∂Ω × (0, T ).Therefore we can write:∇J(p, δp) = T 0 Ω δu(x,t) • P * δu * (x, t)dxdt = T Pδu(x, t) • δu * (x, t)dxdt -s)(∇u(x, s) + ∇u(x, s) T ) : ∇δu * (x, t)dsdt dx.

.

  The initial and adapted meshes and basis functions are shown on Figure15and Figure16

Figure 16 .

 16 Figure 16. Adapted mesh (left), adapted basis functions n o 1, n o 2 and n o 5 (right).

Figure 17 .

 17 Figure 17. Relative error for the coefficient p in L 2 -norm with respect to K.

Figure 18 .

 18 Figure 18. Exact coefficient p (on the left), recovered coefficient p * at step 0 with K 0 = 50 and K 0 = 100 (in the center), recovered coefficient p * at step 2 with K 0 = K 1 = K 2 = 50 (on the right).

  Definitions of the distances d 0 and d.

	0000000000 0000000000 0000000000 0000000000 0000000000 0000000000 0000000000 0000000000 0000000000 0000000000 0000000000 1111111111 1111111111 1111111111 1111111111 1111111111 1111111111 1111111111 1111111111 1111111111 1111111111 1111111111	0000000 0000000 0000000 0000000 0000000 0000000 0000000 0000000 0000000 0000000 1111111 1111111 1111111 1111111 1111111 1111111 1111111 1111111 1111111 1111111

  t) dx dt.

	x 0	0000000 0000000 0000000 0000000 0000000 0000000 0000000 0000000 0000000 0000000 0000000 1111111 1111111 1111111 1111111 1111111 1111111 1111111 ε 1111111 1111111 1111111 1111111

  For all t ≥ 0, we can find an open O with smooth boundary (at least C 1,1 ) such that Ω t (ε, δ) ⊂ O ⊂ Ω t (0, 0). Since, for t fixed, the function χ 2 (•, t)u(•, t)e σϕ(•,t) is in

	H 2 (O) ∩ H 1 0 (O), we can write:
	χ 2 (•, t)u(•, t)e σϕ(•,t)	H 2 (O) ≤ C ∆(χ 2 (•, t)u(•, t)e σϕ(•,t) ) L 2 (O) . (21)
	If we integrate inequality (21) for t ∈ (0, +∞) and notice that χ 1
	x 0	Ω t (0, 0)	0000000 0000000 0000000 0000000 0000000 0000000 0000000 0000000 0000000 0000000 0000000 1111111 1111111 1111111 1111111 1111111 1111111 1111111 ε 1111111 1111111 1111111 1111111 Ω t (ε, δ)

t)|∇ϕ(x, t)| 2 + σχ 2 (x, t)∆ϕ(x, t) + ∆χ 2 (x, t) u(x, t)e σϕ(x,t) + (σ∇ϕ(x, t) • ∇χ 2 (x, t)) u(x, t)e σϕ(x,t) + χ 2 (x, t)∆u(x, t)e σϕ(x,t) . O Figure 7. Definition of the open O.

Acknowledgments

This work has been partially supported by the MathAmSud CIP-EDP project, by doctoral fellowships from CNRS and CONICYT for the first author and by the FONDECYT-CONICYT 1061263 grand for the second author.

Regularization method

We denote by p * the numerical parameter and by p the exact parameter, and we suppose p * close to p. We can write:

If M p, linear from P (of infinite dimension) to U , is compact for . P and . U , then its inverse M -1 p is not bounded. So, even if we find p * such that M pp * is close to ū, we cannot guarantee that p * will be close to p. Thus, as regularization method [START_REF] Kirsch | An introduction to the mathematical theory of inverse problems[END_REF], we introduce the operator

p , where Π K is the projection of P to a finite dimensional space P K . This operator R K is bounded in U and verifies

and we now look for p * ∈ P K .

Mesh and basis adaptation

We choose to look for the unknown coefficient p * in the space P K of the K first eigenfunctions of the mesh, that is

where p| ∂Ω is a raising of the trace of the exact value of p (which is known as we assumed p known in ω) and We propose an adaptive method to solve the problem accurately. After computing a first solution p * 0 on the initial mesh, we use it to refine the mesh and to adapt the spectral basis. Indeed, we consider next as basis functions the solutions of the following problem:

-∇ • (a(x)∇ϕ i (x)) = σ i ϕ i (x), ∀x ∈ Ω, ϕ i (x) = 0, ∀x ∈ ∂Ω,


Zhuang Zhao

Oussama Mhibik

Malik Nafa

Sébastien Chénais

Sébastien Forget
email: sebastien.forget@univ-paris13.fr.

High brightness diode-pumped organic solid-state laser

High-power, diffraction-limited organic solid-state laser operation has been achieved in a vertical external cavity surface-emitting organic laser (VECSOL), pumped by a low-cost compact blue laser diode. The diode-pumped VECSOLs were demonstrated with various dyes in a polymer matrix, leading to laser emissions from 540 nm to 660 nm. Optimization of both the pump pulse duration and output coupling leads to a pump slope efficiency of 11% for a DCM based VECSOLs. We report output pulse energy up to 280 nJ with 100 ns long pump pulses, leading to a peak power of 3.5 W in a circularly symmetric, diffraction-limited beam.

Thin film organic solid-state lasers (OSSLs) have emerged as an important category of coherent sources, offering wavelength tunability over the whole visible spectrum. 1,2 One of the main advantages of organic materials as gain media for lasers is their ability to be easily processed and their intrinsic low cost, making them interesting light sources for many applications such as spectroscopy, 3 bio/chemo-sensing, 4,5 and polymer fiber telecommunications. 6 However, as electrical pumping is still an open challenge, OSSLs are generally optically pumped by frequency-doubled or tripled solid-state lasers, which are bulky and much more expensive than the organic laser structure itself, jeopardizing the promises of a truly low-cost system. Within the last decade, the spectacular progress in the field of high power inorganic blue laser diodes (LD) and light-emitting diodes (LED) enabled several teams to report on diode-pumped and LED-pumped OSSLs, thus demonstrating true low-cost systems. [7][8][9][10][11][12] Nevertheless, as the peak power of laser diodes (and all the more LEDs) is weak compared to the values achieved with pulsed solid-state lasers, very low threshold laser resonators are needed. Diode pumping has been achieved up to now through the use of distributed feedback (DFB) [7][8][9]12 or organic VCSEL (1D microcavity) 10,11 resonators.

As low threshold operation involves a low output coupling, minimizing threshold and maximizing output power are generally incompatible in a single laser device. [START_REF] Siegman | Lasers[END_REF] As a consequence, the reported output powers of such diode-pumped lasers are very low (either given in arbitrary units 8,9 or when measured, in the pJ 10 or nJ 9 range) and the only reported optical-to-optical efficiencies are around 1%. 9,10 Furthermore, DFB low loss compact resonators do not generally provide diffraction-limited beams, 12 which combined with low output powers leads to a weak brightness, which is a strong limitation especially for applications requiring a free-propagating beam, such as e.g. Raman spectroscopy 3 . Brightness is related to spatial coherence and is a figure that is typically orders-of-magnitude higher in lasers than in any incoherent source. The power (or energy, respectively) brightness or radiance is the peak power (energy, resp.) over the optical étendue, and is an indication of how much a given optical power (or energy) can be channeled into a given mode or focused to a tight spot.

In this letter, we report on the first diode-pumped vertical external cavity surface-emitting organic laser (VECSOL). This kind of structure [START_REF] Rabbani-Haghighi | [END_REF] (see figure 1 ) has been demonstrated previously under solid-state laser pumping, and enabled obtaining a diffraction-limited output, a high optical-to-optical efficiency (around 50 % [START_REF] Chenais | Organic Solid State Lasers[END_REF]16 ), and output energies up to several µJ [START_REF] Rabbani-Haghighi | [END_REF] , that is several orders of magnitude higher than what is usually obtained in OSSLs with DFB or microcavity resonators. This is because the output coupling is set to optimize the efficiency (and not minimize the threshold) and because mode area can be made large while keeping the beam quality to the diffraction limit. Reaching threshold in VECSOLs is thus more challenging than in low-loss resonators, but it can now be overcome thanks to the recent availability of high power blue laser diodes (emitting several watts in CW regime) at very reasonable costs. The VECSOL open cavity enables pump-to-cavity mode matching which is the key for high conversion efficiency, but also offers flexibility for output coupling control and addition of intracavity elements (e.g. for wavelength tuning). We studied the influence of the pump pulse duration and output coupler transmission on the laser performance and we achieved a best performance of 3.5 W output peak power. This corresponds to a peak power brightness of 8.9 TW/m 2 /sr (at a wavelength around 620 nm).

The VECSOL design is shown in figure 1. The gain chip consists of an 18-µm-thick film of dye-doped Polymethylmethacrylate (PMMA, Mw=9.5×10 5 g•mol -1 ), spin coated directly onto a highly-reflective plane dielectric mirror (R> 99.5% in the range of 500-670 nm), that is also designed to be transparent (T>94%) at the pump wavelength of 450 nm.

The remote output coupler has a radius-of-curvature of 200 mm and is set 1 mm away from the active mirror.

The pump laser source is a low-cost high power blue InGaN laser diode (maximum output power of 1.7 W, central emission wavelength 450 nm) driven by a pulsed power supply (Picolas GmbH), with controllable pulsewidth. The pump beam is focused to a 90 µm*100 µm -in-diameter spot size onto the dye-doped PMMA layer through a high-numerical aperture longworking-distance microscope objective, reasonably matched to the cavity mode (130 µm in diameter).

The time evolution of organic laser emission is recorded with a reverse biased Si pin photodiode attached to a 1 GHz oscilloscope (Agilent) and is synchronized to the pump pulses. All the fabrication process and measurements are made in ambient air with no specific encapsulation of the samples. To illustrate the versatility of organic laser, we used in this work different lasing dyes to cover a large part of the visible spectrum: Coumarin 510 (C510), Coumarin 540 (C540), Pyrromethene 567 (Py567), Pyrromethene 597 (Py597), 4-Dicyanomethylene-2-methyl-6-p dimethylaminostyryl-4H-pyran (DCM) and Rhodamine 640 (Rh640). All active materials were dispersed into PMMA either alone or in the form of donor/acceptor mixtures (host-guest blends) according to their absorption properties at the pump wavelength. For C540 and DCM, the weight ratio in PMMA is set to 1% which enables the absorption at 450 nm to be 80% in a single pass, while circumventing concentrating quenching. For Py567, Py597 and Rh640, direct pumping at 450 nm is not efficient because of a too low absorption cross-section at this wavelength. We used C510/Py567, C510/Py597 and C540/Rh640 as excitation donor/acceptor pairs, in which the donor dye absorbed the pump light and transferred its excitation to the acceptor dye by a Forster resonant energy transfer (FRET) mechanism. 17 In those cases we blended the "donor" and the "acceptor" with the same weight ratio of 1% in PMMA. Figure 2 shows the lasing spectra recorded using an optical spectrometer (USB 2000+, Ocean optics) with a resolution of ~ 0.3 nm. The evenly spaced peaks in the spectrum correspond to the free spectral range of the Fabry-Perot etalon formed by the 18µm-thick active layer itself. Each of these peaks contains several thin peaks (linewidth below 0.1 nm, limited by our spectrometer resolution) corresponding to the modes of the cavity (see inset in figure 2). Inserting an additional Fabry-Perot filter allowed single-peak operation 18 whereas continuous tunability over 10 nanometers could be achieved upon taking advantage of the thickness variation at the edge of the sample. [START_REF] Rabbani-Haghighi | [END_REF] In the following, the laser characterization and optimization are detailed with the DCM-based VECSOL.

Beside their interest as low cost pumping sources, diode lasers offer the advantage that pulse duration is totally usercontrollable from a few ns to the CW regime. This is of utmost interest for studying organic lasers, where the laser dynamics is complex and strongly linked to the long-lived triplet state population, causing triplet-triplet absorption and triplet-singlet quenching that affect the maximum pulsewidth obtainable from an OSSL. 19 Figure 3 shows measurements of the VECSOL temporal dynamics when operating with pump pulse widths of 50 ns, 100 ns, 200 ns and 500 ns. The clearly visible delay between the onsets of pump and laser pulses is closely related to the oscillation build up time and therefore depends on the cavity length, the transmission of the output coupler and the pump power density.

The clearly visible delay between the onsets of pump and laser pulses is closely related to the oscillation build up time and therefore depends on the cavity length, the transmission of the output coupler and the pump power density. This large delay can be well accounted by solving Statz-DeMars coupled population equations, as detailed in ref 16 in a slightly different case under laser pumping, and is a salient feature of VECSOLs, due to the fact that the active medium thickness represents only a tiny fraction of the cavity length. In figure 3 (a-d), pump and laser pulses are displayed for various pump durations. A typical relaxation oscillation is observed at the beginning of the laser pulse. For long pump pulses (fig. 3 (b-d)), the laser intensity vanishes during the pump pulse, as a consequence of the slow increase of the triplet state population. Those results show that the maximum pulsewidth of a DCM-based VECSOLs is around 100 ns (FWHM) and is limited by triplet population buildup, as also observed in ref. 19. The variation of the laser slope efficiency (slope of the laser pulse energy versus the incident pump pulse energy curve) with pump pulse duration is plotted in figure 4. When a 20-ns pump pulse is used, the efficiency is only 2.5%. as the laser pulse starts only at the very end of the pump pulse. The laser slope efficiency increases with the pump pulse width and reaches a maximum value of 11% (corresponding to a quantum slope efficiency of ~ 16%) with a pump pulse width of 100 ns, and then decays due to the triplet population buildup.

Fig. 1 :

 1 Fig.1: Schematic diagram of the diode-pumped plano-concave VECSOL resonator (MO : microscope objective, LD : Laser Diode). Inset: Long-exposure time photo of the working diode-pumped VECSOL.

Fig. 2 :

 2 Fig.2: Laser spectra measured from LD pumped VECSOLs with five different dyes, from left to right: C540 (green dash), Py567/C510 (yellow dash-dot), Py597/C510 (light orange dash-double-dot), DCM(orange solid line) and Rh640/C540 (red dot). Inset : high-resolution spectrum (0.1 nm) showing the structure of each peak (here for a 450 µm long cavity).

Fig. 3

 3 Fig. 3 Temporal profile of 1%-DCM VECSOL (orange line) pumped with varied pump pulse durations (blue line): (a) 50 ns; (b) 100 ns; (c) 200 ns; (d) 500ns. An output coupler with a transmission of 3% and a peak pump power density of 132 kW/cm 2 were used in the experiment.

Fig. 4

 4 Fig. 4 pump slope efficiency of 1%-DCM VECSOLs at varied pump pulse duration from 20 ns to 500 ns. The maximum pump slope efficiency was obtained with a pump pulse duration of ~ 100 ns.

Fig. 1 :

 1 Fig.1: Schematic diagram of the diode-pumped plano-concave VECSOL resonator (MO : microscope objective, LD : Laser Diode). Inset: Long-exposure time photo of the working diode-pumped VECSOL.

Fig. 2 :

 2 Fig. 2 : Laser spectra measured from LD pumped VECSOLs with five different dyes, from left to right: C540 (green dash), Py567/C510 (yellow dash-dot), Py597/C510 (light orange dash-double-dot), DCM(orange solid line) and Rh640/C540 (red dot). Inset : high-resolution spectrum (0.1 nm) showing the structure of each peak (here for a 450 µm long cavity).

Fig. 3

 3 Fig. 3 Temporal profile of 1%-DCM VECSOL (orange line) pumped with varied pump pulse durations (blue line): (a) 50 ns; (b) 100 ns; (c) 200 ns; (d) 500ns. An output coupler with a transmission of 3% and a peak pump power density of 132 kW/cm 2 were used in the experiment.

Fig. 4

 4 Fig. 4 pump slope efficiency of 1%-DCM VECSOLs at varied pump pulse duration from 20 ns to 500 ns. The maximum pump slope efficiency was obtained with a pump pulse duration of ~ 100 ns.

Fig. 5

 5 Fig.5Laser performance including pump slope efficiency and pump density threshold, were investigated by changing the transmission of output couplers: 0.5%, 2%, 3% and 6%.

Fig. 6

 6 Fig.6Input-output characteristic of 1%-DCM VECSOL, with 1 LD (below 1.6 µJ) and 2 polarization-coupled LDs (above 1.6 µJ.) as the pump. Inset: spatial profile of the output beam and its evolution along the emission direction, yielding M 2 ~1 in both transverse directions (diffraction-limited output.)

Acknowledgements:

The authors wish to acknowledge the Agence Nationale de la Recherche (ANR-12-EMMA-0040 "Vecspresso" project) for funding this work.

Another factor playing a key role for power extraction in lasers is the choice of the transmission of the output coupler. [START_REF] Siegman | Lasers[END_REF] Whereas every laser cavity has an optimal output coupler transmission (resulting from the trade-off between cavity photon lifetime and energy extraction from the resonator), this feature is not easily controlled in typical monolithic 2 nd order DFB structures where the output laser beam comes from the diffraction of the laser wave on the grating itself. In a VECSOL structure, changing the amount of output coupling is straightforward and can be done without affecting the gain properties, providing an additional degree of freedom to optimize the resonator. To this end, the input-output characteristic curves of a 1mm long VECSOL (with optimal pump-pulse duration of 100 ns, corresponding to a laser pulsewidth of 80 ns) are compared (figure 5) for four different output couplers having transmission of 0.5%, 2%, 3% and 6% respectively. As expected, the minimum pump threshold (32 kW/cm 2 ) is obtained using the output coupler with the lowest transmission (0.5%). The optimum coupler transmission in order to achieve the maximum pump slope efficiency of 11% is 3%; in this case, the pump density threshold is 69 kW/cm 2 .

Fig. 5 Laser performance including pump slope efficiency and pump density threshold, were investigated by changing the transmission of output couplers: 0.5%, 2%, 3% and 6%.

In figure 6, the emission pulse energy is given as a function of incident pump pulse energy for two output couplers. With a single LD pumping system, the maximum pump pulse energy deposited on the VECSOL chip is 1.6 µJ in a 100 ns pulse, and a maximum emission pulse energy of 130 nJ is achieved with a 2% transmission coupler. Using the 3% transmission coupler, the maximum emission pulse energy is slight lower (115 nJ) while it shows higher pump slope efficiency. To push the output power to its maximum, we drove two LDs with a single LD driver. The two diodes were polarization-coupled through a polarizing beam splitter to provide an incident pump pulse energy on the VECSOL chip of 3.2 µJ in a 100 ns pulse. An energy of 280 nJ is obtained with the 3% transmission coupler, and the peak power (estimated from the pulse energy over FWHM of laser pulse duration) is around 3.5 W. Beam quality was quantified by measuring the M² factor, shown in the inset of figure 6. The LD pumped VECSOL emits a diffraction limited beam with M 2 ~ 1.02. A typical image of the TEM 00 beam profile is also shown. The measured power brightness, defined as P/(M²) 2 •λ², where P is the peak power (in W) and λ the wavelength, is then 9.1 TW/m²/sr at 620 nm. This figure is not easily comparable with other OSSLs, as the M² factor is rarely measured. However, for DFB lasers, its value is always greater that one. 20,21 M² values range from a few units (an astigmatic annular mode with a M² value as low as 2.2 has been reported in a two dimensional DFB resonator pumped by a solid-state microchip laser, 21 ) to several tens (a high divergence of 5° is described for example in ref 22). As a consequence, the highest power brightness (resp. energy brightness) that can be deduced (at 620 nm, with the optimistic hypothesis of a state-of-the art value of 2.2 for M²) from the literature for LD-pumped OSSL is around 54 GW/m²/sr (resp. 1.1 kJ/m²/sr). 9 The herein reported power brightness of 9.1 TW/m²/sr thus corresponds to a two-order of magnitude increase. For the energy brightness, our value of 700 kJ/m²/sr is more than 600 times higher than the previously reported results.

Fig. 6 Input-output characteristic of 1%-DCM VECSOL, with 1 LD (below 1.6 µJ) and 2 polarization-coupled LDs (above 1.6 µJ.) as the pump. Inset: spatial profile of the output beam and its evolution along the emission direction, yielding M 2 ~1 in both transverse directions (diffraction-limited output.)

Finally, photostability was investigated at fixed pump density above threshold (80 kW/cm 2 ). Under 100-ns long pulse pumping, the output emission decreased to half its initial value after ~ 1000 pulses. As VECSOLs are made of planar, uniform organic layers, it is straightforward to scan the sample (perpendicularly to the cavity axis) to enhance the lifetime. As an example, scanning a 1-inch-in-diameter sample with a 100-µm-in-diameter pump spot, allows more than 8000 fresh spots to be used (with no realignment issue), 18 leading to a lifetime of 8×10 6 pulses (equivalent to 13 min operation at 10 kHz).

Furthermore, encapsulation can be used since it has shown a significant improvement in the lifetime of the organic device. 1,2 In summary, we have demonstrated diode-pumped VECSOLs that cover the visible range from 450 nm to 660 nm, by choosing proper dyes or dye mixtures. The lasing performance of diode-pumped VECSOLs was optimized by varying the pump pulse duration and the output coupler transmission. An optimum pump slope efficiency is demonstrated for 100 ns pump pulses, corresponding to a duration long enough to allow for efficient build-up of laser oscillation and short enough to avoid losses and quenching resulting from the buildup of a triplet state population. With DCM-based VECSOLs, a pump slope efficiency of 11% has been obtained with a diode pump pulse width of 100 ns. A maximum energy level of 280 nJ (3.5 W peak power), limited by the available pump pulse peak power (32 W), has been obtained with two polarizationcoupled laser diodes, with a nearly diffraction limited beam. Those results pave the way toward the development of truly low-cost, compact, high-brightness tunable organic thin-film solid-state lasers at the µJ level.