Elie Bretin
email: bretin@cmap.polytechnique.fr

Abdul Wahab
email: wahab@cmap.polytechnique.fr

Some Anisotropic Viscoelastic Green Functions *

Keywords: Green Function, Viscoelastic Media, Anisotropic Media MSC 2000: Primary 35A08, 74D99; Secondary 92C55, 74L05

] to express unknown Green function in terms of three scalar functions φi, by using the spectral decomposition of the Christoffel tensor associated with the medium. The problem of computing Green function is, thus reduced to the resolution of three scalar wave equations satisfied by φi, and subsequent equations with φi as source terms. To describe viscosity effects, we choose an empirical power law model which becomes well known Voigt model for quadratic frequency losses.

Introduction

Numerous applications in biomedical imaging [START_REF] Ammari | An introduction to Mathematics of Emerging Biomedical Imaging[END_REF][START_REF] Ammari | Reconstruction of Small Inhomogeneities from Boundary Measurements[END_REF], seismology [START_REF] Aki | Quantitative Seismology[END_REF][START_REF] Carcione | Wave Field in the Real Media[END_REF], exploration geophysics [START_REF] Helbig | Foundations of Anisotropy for Exploration Seismics[END_REF][START_REF] Helbig | 75-plus years of anisotropy in exploration and reservoir seismics: A historical review of concepts and methods[END_REF], material sciences [START_REF] Allaire | Shape Optimization by the Homogenization Method[END_REF][START_REF] Ammari | Polarization and Moment Tensors[END_REF] and engineering sciences [START_REF] Achenbach | Wave Propagation in Elastic Solids[END_REF][START_REF] Ben-Menahem | Seismic waves and sources[END_REF][START_REF] Lekhnitskii | Theory of Elasticity of an Anisotropic Body[END_REF] have fueled research and development in theory of elasticity. Elastic properties and attributes have gained interest in the recent decades as a diagnostic tool for non-invasive imaging [START_REF] Greenleaf | Selected methods for imaging elastic properties of biological tissues[END_REF][START_REF] Sarvazyan | Shear wave elasticity imaging: a new ultrasonic technology of medical diagnostics[END_REF]. Their high correlation with the pathology and the underlying structure of soft tissues has inspired many investigations in biomedical imaging and led to many interesting mathematical problems [START_REF] Ammari | Direct elastic imaging of a small inclusion[END_REF][START_REF] Ammari | A method of biological tissues elasticity reconstruction using magnetic resonance elastography measurements[END_REF][START_REF] Ammari | A new optimal control approach for the reconstruction of extended inclusions[END_REF][START_REF] Ammari | Effective viscosity properties of dilute suspensions of arbitrarily shaped particles[END_REF][START_REF] Ammari | Separation of scales in elasticity imaging: A numerical study[END_REF][START_REF] Bercoff | The role of viscosity in the impulse diffraction field of elastic waves induced by the acoustic radiation force[END_REF][START_REF] Sinkus | High resolution tensor MR Elastography for breast tumor detection[END_REF][START_REF] Sinkus | Imaging anisotropic and viscous properties of breast tissue by magnetic resonance-elastography[END_REF].

Biological materials are often assumed to be isotropic and inviscid with respect to elastic deformation. However, several recent studies indicate that many soft tissues exhibit anisotropic and viscoelastic behavior [START_REF] Gennisson | Transient elastography in anisotropic medium: Application to the measurement of slow and fast shear wave speeds in muscles[END_REF][START_REF] Oida | The measurement of anisotropic elasticity in skeletal muscle using MR Elastography[END_REF][START_REF] Sinkus | High resolution tensor MR Elastography for breast tumor detection[END_REF][START_REF] Sinkus | Imaging anisotropic and viscous properties of breast tissue by magnetic resonance-elastography[END_REF][START_REF] Namani | Shear wave propagation in anisotropic soft tissues and gels[END_REF][START_REF] Yoon | Ultrasonic wave propagation in human cortical bone-I. Theoretical considerations for hexagonal symmetry[END_REF]. Sinkus et al. have inferred in [START_REF] Sinkus | High resolution tensor MR Elastography for breast tumor detection[END_REF] that breast tumor tends to be anisotropic, while Weaver et al. [START_REF] Weaver | Evidence of the anisotropic nature of the mechanical properties of breast tissue[END_REF] have provided an evidence that even non cancerous breast tissue is anisotropic. White matter in brain [START_REF] Namani | Shear wave propagation in anisotropic soft tissues and gels[END_REF] and cortical bones [START_REF] Yoon | Ultrasonic wave propagation in human cortical bone-I. Theoretical considerations for hexagonal symmetry[END_REF] 1 also exhibit similar behavior. Moreover, it has been observed that the shear velocities parallel and orthogonal to the fiber direction in forearm [START_REF] Oida | The measurement of anisotropic elasticity in skeletal muscle using MR Elastography[END_REF] and biceps [START_REF] Gennisson | Transient elastography in anisotropic medium: Application to the measurement of slow and fast shear wave speeds in muscles[END_REF] are different. This indicates that the skeletal muscles with directional structure are actually anisotropic. Thus, an assumption of isotropy can lead to erroneous forward-modelled wave synthetics, while an estimation of viscosity effects can be very useful in characterization and identification of anomaly [START_REF] Bercoff | The role of viscosity in the impulse diffraction field of elastic waves induced by the acoustic radiation force[END_REF].

A possible approach to handle viscosity effects on image reconstruction has been proposed in [START_REF] Bretin | On the Green function in visco-elastic media obeying a frequency power-law[END_REF] using stationary phase theorem. It is shown that the ideal Green function (in an inviscid regime) can be approximated from the viscous one by solving an ordinary differential equation. Once the ideal Green function is known one can identify a possible anomaly using imaging algorithms such as time reversal, back-propagation, Kirchhoff migration or MUSIC [START_REF] Ammari | Direct elastic imaging of a small inclusion[END_REF][START_REF] Ammari | Transient elasticity imaging and time reversal[END_REF][START_REF] Ammari | Reconstruction of Small Inhomogeneities from Boundary Measurements[END_REF][START_REF] Ammari | An introduction to Mathematics of Emerging Biomedical Imaging[END_REF]. One can also find the elastic moduli of the anomaly using the asymptotic formalism and reconstructing a certain polarization tensor in the far field [START_REF] Ammari | A method of biological tissues elasticity reconstruction using magnetic resonance elastography measurements[END_REF][START_REF] Ammari | Transient elasticity imaging and time reversal[END_REF][START_REF] Ammari | Polarization and Moment Tensors[END_REF][START_REF] Ammari | Expansion methods, Handbook of Mathematical Methods in Imaging[END_REF].

The importance of Green function stems from its role as a tool for the numerical and asymptotic techniques in biomedical imaging. Many inverse problems involving the estimation and acquisition of elastic parameters become tractable once the associated Green function is computed [START_REF] Alves | Boundary integral formulae for the reconstruction of imperfections of small diameter in an elastic medium[END_REF][START_REF] Ammari | Direct elastic imaging of a small inclusion[END_REF][START_REF] Ammari | Transient elasticity imaging and time reversal[END_REF][START_REF] Bretin | On the Green function in visco-elastic media obeying a frequency power-law[END_REF]. Several attempts have been made to compute Green functions in purely elastic and/or isotropic regime. (See e.g. [START_REF] Bretin | On the Green function in visco-elastic media obeying a frequency power-law[END_REF][START_REF] Bercoff | The role of viscosity in the impulse diffraction field of elastic waves induced by the acoustic radiation force[END_REF][START_REF] Burridge | Fundamental elastodynamic solutions for anisotropic media with ellipsoidal slowness surfaces[END_REF][START_REF] Carcione | Wave Field in the Real Media[END_REF][START_REF] Payton | Elastic Wave Propagation in Transversely Isotropic Media[END_REF][START_REF] Vavrycuk | Asymptotic Green's function in homogeneous anisotropic viscoelastic media[END_REF][START_REF] Vavrycuk | Exact elastodynamic Green functions for simple types of anisotropy derived from higher-order ray theory[END_REF][START_REF] Vavrycuk | Elastodynamic and elastostatic Green tensors for homogeneous weak transversely isotropic media[END_REF] and references therein). However, it is not possible to give a closed form expression for general anisotropic Green functions without certain restrictions on the media. In this work, we provide anisotropic viscoelastic Green function in closed form for three particular anisotropic media.

The elastodynamic Green function in isotropic media is calculated by separating wave modes using Helmholtz decomposition of the elastic wavefield [START_REF] Aki | Quantitative Seismology[END_REF][START_REF] Bretin | On the Green function in visco-elastic media obeying a frequency power-law[END_REF][START_REF] Bercoff | The role of viscosity in the impulse diffraction field of elastic waves induced by the acoustic radiation force[END_REF]. Unfortunately, this simple approach does not work in anisotropic media, where three different waves propagate with different phase velocities and polarization directions [START_REF] Carcione | Wave Field in the Real Media[END_REF][START_REF] Ben-Menahem | Seismic waves and sources[END_REF][START_REF] Červenỳ | Seismic Ray Theory[END_REF]. A polarization direction of quasi-longitudinal wave that differs from that of wave vector, impedes Helmholtz decomposition to completely separate wave modes [START_REF] Dellinger | Anisotropic Seismic Wave Propagation[END_REF].

The phase velocities and polarization vectors are the eigenvalues and eigenvectors of the Christoffel tensor Γ associated with the medium. So, the wavefield can always be decomposed using the spectral basis of Γ. Based on this observation, Burridge et al. [START_REF] Burridge | Fundamental elastodynamic solutions for anisotropic media with ellipsoidal slowness surfaces[END_REF] proposed a new approach to calculate elastodynamic Green functions. Their approach consists of finding the eigenvalues and eigenvectors of Christoffel tensor Γ(∇ x ) using the duality between algebraic and differential objects. Therefore it is possible to express the Green function G in terms of three scalar functions φ i satisfying partial differential equations with constant coefficients. Thus the problem of computing G reduces to the resolution of three differential equations for φ i and of three subsequent equations (which may or may not be differential equations) with φ i as source terms. See [START_REF] Burridge | Fundamental elastodynamic solutions for anisotropic media with ellipsoidal slowness surfaces[END_REF] for more details.

Finding the closed form expressions of the eigenvalues of Christoffel tensor Γ is usually not so trivial because its characteristic equation is a polynomial of degree six in the components of its argument vector. However, with some restrictions on the material, roots of the characteristic equation can be given [START_REF] Payton | Elastic Wave Propagation in Transversely Isotropic Media[END_REF]. In this article, we consider three different media for which not only the explicite expressions of the eigenvalues of Γ are known [START_REF] Burridge | Fundamental elastodynamic solutions for anisotropic media with ellipsoidal slowness surfaces[END_REF][START_REF] Vavrycuk | Exact elastodynamic Green functions for simple types of anisotropy derived from higher-order ray theory[END_REF], but they are also quadratic homogeneous forms, in the components of the argument vector. As a consequence, equations satisfied by φ i become scalar wave equations. Following Burridge et al. [START_REF] Burridge | Fundamental elastodynamic solutions for anisotropic media with ellipsoidal slowness surfaces[END_REF], we find the viscoelastic Green functions for each medium. It is important to note that the elastodynamic Green function in a purely elastic regime, for the media under consideration, are well known [START_REF] Vavrycuk | Exact elastodynamic Green functions for simple types of anisotropy derived from higher-order ray theory[END_REF][START_REF] Burridge | Fundamental elastodynamic solutions for anisotropic media with ellipsoidal slowness surfaces[END_REF]. Also, the expression of the Green function for viscoelastic isotropic medium, which is computed as a special case, matches the one provided in [START_REF] Bretin | On the Green function in visco-elastic media obeying a frequency power-law[END_REF].

It has been shown in [START_REF] Catheline | Measuring of viscoelastic properties of homogeneous soft solid using transient elastography: an inverse problem approach[END_REF] that Voigt model is well adopted to describe the viscosity response of many soft tissues to low frequency excitations. In this work, we consider a more general model proposed by Szabo and Wu in [START_REF] Szabo | A model for longitudinal and shear wave propagation in viscoelastic media[END_REF], which describes an empirical power law behavior of many viscoelastic materials including human myocardium. This model is based on a time-domain statement of causality [START_REF] Szabo | Time domain wave equations for lossy media obeying a frequency power law[END_REF][START_REF] Titchmarsh | Introduction to the Theory of Fourier Integrals[END_REF] and reduces to Voigt model for the specific case of quadratic frequency losses.

We provide some mathematical notions, theme and the outlines of the article in the next section.

Mathematical Context and Paper Outlines 2.1 Viscoelastic Wave Equation

Consider an open subset Ω of R 3 , filled with a homogeneous anisotropic viscoelastic material. Let u(x, t) : Ω × R + → R 3 be the displacement field at time t of the material particle at position x ∈ Ω and ∇ x u(x, t) be its gradient.

Under the assumptions of linearity and small perturbations, we define the order two strain tensor by

ε : (x, t) ∈ Ω × R + -→ 1 2 ∇ x u + ∇ x u T (x, t), (2.1)
where the superscript T indicates a transpose operation. Let C ∈ L 2 s (R 3 ) and V ∈ L 2 s (R 3 ) be the stiffness and viscosity tensors of the material respectively. Here L 2 s (R 3 ) is the space of symmetric tensors of order four. These tensors are assumed to be positive definite, i.e. there exists a constant δ > 0 such that

(C : ξ) : ξ ≥ δ|ξ| 2 and (V : ξ) : ξ ≥ δ|ξ| 2 , ∀ξ ∈ L s (R d ),
where L s (R 3 ) denotes the space of symmetric tensors of order two.

The generalized Hooke's Law [START_REF] Szabo | A model for longitudinal and shear wave propagation in viscoelastic media[END_REF] for power law media states that the stress distribution

σ : Ω × R + → L s (R 3 )
produced by deformation ε, satisfies:

σ = C : ε + V : A[ε] (2.2)
where A is a causal operator defined as

A[ϕ] = -(-1) γ/2 ∂ γ-1 ϕ ∂t γ-1 γ is an even integer, 2 π (γ -1)!(-1) (γ+1)/2 H(t) t γ * t ϕ γ is an odd integer, - 2 π Γ(γ) sin(γπ/2) H(t) |t| γ * t ϕ γ is a non integer. (2.3)
Note that by convention,

A[u] i = A[u i ] and A[ε] ij = A[ε ij ] 1 ≤ i, j ≤ 3.
Here H(t) is the Heaviside function, Γ is the gamma function and * t represents convolution with respect to variable t. See [START_REF] Alekseev | Equations of state for viscoelastic biological media[END_REF][START_REF] Caputo | Linear models of dissipation whose Q is almost frequency independent-II[END_REF][START_REF] Szabo | A model for longitudinal and shear wave propagation in viscoelastic media[END_REF][START_REF] Szabo | Time domain wave equations for lossy media obeying a frequency power law[END_REF][START_REF] Titchmarsh | Introduction to the Theory of Fourier Integrals[END_REF] for comprehensive details and discussion on fractional attenuation models, causality and the loss operator A.

The viscoelastic wave equation satisfied by the displacement field u(x, t)

reads now ρ ∂ 2 u ∂t 2 -F = ∇ x • σ = ∇ x • C : ε + V : A[ε] ,
where F(x, t) is the applied force and ρ is the density (supposed to be constant) of the material.

Remark 2.1 For quadratic frequency losses, i.e, when γ = 2, operator A reduces to a first order time derivative. Therefore, power-law attenuation model turns out to be the Voigt model in this case.

Spectral decomposition by Christoffel tensors

We introduce now the Christoffel tensors Γ c , Γ v : R 3 → L s (R 3 ) associated respectively with C and V defined by:

Γ c ij (n) = 3 k,l=1 C kilj n k n j , Γ v ij (n) = 3 k,l=1 V kilj n k n j , ∀n ∈ R 3 , 1 ≤ i, j ≤ 3.
Remark that the viscoelastic wave equation can be rewritten in terms of Christoffel tensors as :

ρ ∂ 2 u ∂t 2 -F = Γ c [∇ x ]u + Γ v [∇ x ]A[u]. (2.4)
Note that Γ c and Γ v are symmetric and positive definite as C and V are already symmetric positive definite. Let L c i be the eigenvalues and D c i be the associated eigenvectors of Γ c for i = 1, 2, 3. We define the quantities M c i and E c i by

M c i = D c i • D c i , and
E c i = (M c i ) -1 D c i D c i . (2.5)
As Γ c is symmetric, the eigenvectors D c i are orthogonal and the spectral decomposition of the Christoffel tensor Γ c can be given as:

Γ c = 3 i=1 L c i E c i with I = 3 i=1 E c i (2.6)
where I ∈ L s (R 3 ) is the identity tensor.

Similarly, consider Γ v the Christoffel tensor associated with V and define the quantities

L v i , D v i , M v i and E v i such as Γ v = 3 i=1 L v i E v i with I = 3 i=1 E v i . (2.7)
We assume that the tensors Γ c and Γ v have the same structure in the sense that the eigenvectors D c i and D v i are equal. (See Remark 3.3). In the sequel we use D instead of D c or D v and similar for E and M , by abuse of notation.

Paper Outline

The aim of this work is to compute the elastodynamic Green function G associated to viscoelastic wave equation (2.4). More precisely, G is the solution of the equation

(Γ c [∇ x ]G(x, t) + Γ v [∇ x ]A[G](x, t)) -ρ ∂ 2 G(x, t) ∂t 2 = δ(t)δ(x)I, (2.8)
The idea is to use the spectral decomposition of G of the form

G = 3 i=1 E i (∇ x )φ i = 3 i=1 (D i ⊗ D i )M -1 i φ i , (2.9)
where φ i are three scalar functions satisfying

(L c i (∇ x )φ i + L v i (∇ x )A[φ i ]) -ρ ∂ 2 φ i ∂t 2 = δ(t)δ(x) (2.10)
(See Appendix A for more details about this decomposition.) Therefore, to obtain an expression of G, we need to:

1-solve three partial differential equations (2.10) in φ i

2-subsequent equations

ψ i = M -1 i φ i (2.11)
3-and calculate second order derivatives of ψ i to compute

(D i ⊗ D i )ψ i
In the following Section, we give simple examples of anisotropic media which satisfy some restrictive properties and assumptions (see Subsection 3.4) defining the limits of our approach. In Section 4, we derive the solutions φ i of equations (2.10). In Section 5, we give an explicite resolution of

ψ i = M -1 i φ i and (D i ⊗ D i )ψ i .
Finally, in the last section, we compute the Green function for three simple anisotropic media.

Some Simple Anisotropic Viscoelastic Media

In this section, we present three viscoelastic media with simple type of anisotropy. We also describe some important properties of the media and our basic assumptions in this article.

Definition 3.1 We will call a tensor c = (c mn ) ∈ L s (R 6 ) the Voigt represen- tation of an order four tensor C ∈ L 2 s (R 3 ) if c mn = c p(i,j)p(k,l) = C ijkl 1 ≤ i, j, k, l ≤ 3 where p(i, i) = i, p(i, j) = p(j, i), p(2, 3) = 4, p(1, 3) = 5, p(1, 2) = 6.
We will use c and v for the Voigt representations of stiffness tensor C and viscosity tensor V respectively.

We will let tensors c and v to have a same structure. For each media, the expressions for Γ c , L c i (∇ x ), D c i (∇ x ) and M c i (∇ x ) are provided [START_REF] Burridge | Fundamental elastodynamic solutions for anisotropic media with ellipsoidal slowness surfaces[END_REF][START_REF] Vavrycuk | Exact elastodynamic Green functions for simple types of anisotropy derived from higher-order ray theory[END_REF]. Throughout this section, µ pq will assume the value c pq for c and v pq for v where the subscripts p, q ∈ {1, 2, • • • , 6}. Moreover, we assume that the axes of material are identical with the Cartesian coordinate axes e 1 , e 2 and e 3 and

∂ i = ∂ ∂x i .

Medium I

The first medium for which we present a closed form elastodynamic Green function is an orthorhombic medium with the tensors c and v of the form:

        µ 11 -µ 66 -µ 55 0 0 0 -µ 66 µ 22 -µ 44 0 0 0 -µ 55 -µ 44 µ 33 0 0 0 0 0 0 µ 44 0 0 0 0 0 0 µ 55 0 0 0 0 0 0 µ 66        
The Christoffel tensor is given by

Γ c =   c 11 ∂ 2 1 + c 66 ∂ 2 2 + c 55 ∂ 2 3 0 0 0 c 66 ∂ 2 1 + c 22 ∂ 2 2 + c 44 ∂ 2 3 0 0 0 c 55 ∂ 2 1 + c 44 ∂ + c 33 ∂ 2 3   Its eigenvalues L c i (∇ x
) and the associated eigenvectors

D c i (∇ x ) are: L c 1 (∇ x ) = c 11 ∂ 2 1 + c 66 ∂ 2 2 + c 55 ∂ 2 3 L c 2 (∇ x ) = c 66 ∂ 2 1 + c 22 ∂ 2 2 + c 44 ∂ 2 3 L c 3 (∇ x ) = c 55 ∂ 2 1 + c 44 ∂ 2 2 + c 33 ∂ 2 3 D c i = e i with M c i = 1 ∀i = 1, 2, 3

Medium II

The second medium which we consider is a transversely isotropic medium having symmetry axis along e 3 and defined by the stiffness and the viscosity tensors c and v of the form:

        µ 11 µ 12 -µ 44 0 0 0 µ 12 µ 11 -µ 44 0 0 0 -µ 44 -µ 44 µ 33 0 0 0 0 0 0 µ 44 0 0 0 0 0 0 µ 44 0 0 0 0 0 0 µ 66         with µ 66 = (µ 11 -µ 12 )/2.
Here

Γ c =   c 11 ∂ 2 1 + c 66 ∂ 2 2 + c 44 ∂ 2 3 (c 11 -c 66 ) ∂ 1 ∂ 2 0 (c 11 -c 66 ) ∂ 1 ∂ 2 c 66 ∂ 2 1 + c 11 ∂ 2 2 + c 44 ∂ 2 3 0 0 0 c 44 ∂ 2 1 + c 44 ∂ + c 33 ∂ 2 3   The eigenvalues L c i (∇ x ) of Γ c (∇ x ) in this case are L c 1 (∇ x ) = c 44 ∂ 2 1 + c 44 ∂ 2 2 + c 33 ∂ 2 3 L c 2 (∇ x ) = c 11 ∂ 2 1 + c 11 ∂ 2 2 + c 44 ∂ 2 3 L c 3 (∇ x ) = c 66 ∂ 2 1 + c 66 ∂ 2 2 + c 44 ∂ 2
and the associated eigenvectors D c i (∇ x ) are:

D c 1 =   0 0 1   , D c 2 =   ∂ 1 ∂ 2 0   , D c 3 =   ∂ 2 -∂ 1 0   . Thus M c 1 = 1, and M c 2 = M c 3 = ∂ 2 1 + ∂ 2 2

Medium III

Finally, we will present the elastodynamic Green function for another transversely isotropic media with the axis of symmetry along e 3 and having c and v of the form:

        µ 11 µ 11 -2µ 66 µ 11 -2µ 44 0 0 0 µ 11 -2µ 66 µ 11 µ 11 -2µ 44 0 0 0 µ 11 -2µ 44 µ 11 -2µ 44 µ 11 0 0 0 0 0 0 µ 44 0 0 0 0 0 0 µ 44 0 0 0 0 0 0 µ 66        
The Christoffel tensor in this case is

Γ c =   c 11 ∂ 2 1 + c 66 ∂ 2 2 + c 44 ∂ 2 3 (c 11 -c 66 ) ∂ 1 ∂ 2 (c 11 -c 44 ) ∂ 1 ∂ 3 (c 11 -c 66 ) ∂ 1 ∂ 2 c 66 ∂ 2 1 + c 11 ∂ 2 2 + c 44 ∂ 2 3 (c 11 -c 44 ) ∂ 2 ∂ 3 (c 11 -c 44 ) ∂ 1 ∂ 3 (c 11 -c 44 ) ∂ 2 ∂ 3 c 44 ∂ 2 1 + c 44 ∂ 2 2 + c 11 ∂ 2 3   Its eigenvalues L c i (∇ x ) are: L c 1 (∇ x ) = c 11 ∂ 2 1 + c 11 ∂ 2 2 + c 11 ∂ 2 3 = c 11 ∆ x L c 2 (∇ x ) = c 66 ∂ 2 1 + c 66 ∂ 2 2 + c 44 ∂ 2 3 L c 3 (∇ x ) = c 44 ∂ 2 1 + c 44 ∂ 2 2 + c 44 ∂ 2 3 = c 44 ∆ x
and the eigenvectors D c i (∇ x ) are:

D c 1 =   ∂ 1 ∂ 2 ∂ 3   , D c 2 =   ∂ 2 -∂ 1 0   , D c 3 =   -∂ 1 ∂ 3 -∂ 2 ∂ 3 ∂ 2 1 + ∂ 2 2   (3.1)
In this case,

M c 1 = ∆ x M c 2 = ∂ 2 1 + ∂ 2 2 and M c 3 = (∂ 2 1 + ∂ 2 2 )∆ x

Properties of the Media and Main Assumptions

In all anisotropic media discussed above, it holds that

• The Christoffel tensors Γ c and Γ v have the same structure in the sense that

D c i = D v i , ∀i = 1, 2, 3.
• The eigenvalues L c i (∇ x ) are homogeneous quadratic forms in the components of the argument vector ∇ x i.e.

L c i [∇ x ] = 3 j a 2 ij ∂ 2 ∂x 2 j ,
and therefore equations (2.10) are actually scalar wave equations.

• In all the concerning cases, the operator M c i (∇ x ) is either constant or has a homogeneous quadratic form

M c i = 3 j m 2 ij ∂ 2 ∂x 2 j .
In addition, we assume that

• the eigenvalues of Γ c and Γ v satisfy

L v i (∇ x ) = β i L c i (∇ x ).
• and the loss per wave length is small, i.e.

β i << 1. Remark 3.2 The expression M c 3 = (∂ 2 1 + ∂ 2 2
)∆ x will be avoided in the construction of the Green function by using the expression

G = φ 3 I + E 1 (∇ x )(φ 1 -φ 3 ) + E 2 (∇ x )(φ 2 -φ 3 )
for the elastodynamic Green function.

Remark 3.3 In general, D c

i and D v i are dependant on the parameters c pq and v pq . Consequently, Γ c and Γ v can not be diagonalized simultaneously. However, in certain restrictive cases where the polarization directions of different wave modes ( i.e. quasi longitudinal (qP) and quasi shear waves (qSH and qSV)) are independent of the stiffness or viscosity parameters, it is possible to diagonalize both Γ c and Γ v simultaneously. Moreover, the assumption on the eigenvalues L v i and L c i , implies that for a given wave mode, the decay rate of its velocity in different directions is uniform, but for different wave modes (qP, qSH and qSV) these decay rates are different.

Solution of the Model Wave Problem

Let us now study the scalar wave problems (2.10). We consider a model problem and drop the subscript for brevity in this section as well as in the next section. Consider

(L c [∇ x ]φ + L v [∇ x ]A[φ]) -ρ ∂ 2 φ ∂t 2 = δ(t)δ(x). (4.1)
Our assumptions on the media imply that L c and L v have the following form;

L c [∇ x ] = 3 j=1 a 2 j ∂ 2 ∂x 2 j and L v [∇ x ] = βL c [∇ x ] = 3 j=1 βa 2 j ∂ 2 ∂x 2 j
Therefore, the model equation (4.1) can be rewritten as:

3 j=1 a 2 j ∂ 2 φ ∂x 2 j + βa 2 j A ∂ 2 φ ∂x 2 j -ρ ∂ 2 φ ∂t 2 = δ(t)δ(x),
By a change of variables x j = a j √ ρ ξ j , we obtain in function φ(ξ) = φ(x) the following transformed equation :

∆ ξ φ + βA ∆ ξ φ - ∂ 2 φ ∂t 2 = √ ρ a δ(t)δ(ξ). (4.2)
where the constant a = a 1 a 2 a 3 . Now, we apply A on both sides of the equation (4.2), and replace the resulting expression for A ∆ ξ φ back into the equation (4.2). This yields:

∆ ξ φ + βA ∂ 2 φ ∂t 2 -β 2 A 2 ∆ ξ φ - ∂ 2 φ ∂t 2 = √ ρ a δ(ξ) {δ(t) -βA[δ(t)]}
Recall that β << 1 and the term in β 2 is negligible. Therefore, it holds

∆ ξ φ + βA ∂ 2 φ ∂t 2 - ∂ 2 φ ∂t 2 ≃ √ ρ a δ(ξ) {δ(t) -βA[δ(t)]} . (4.3)
Finally, taking temporal Fourier transform on both sides of (4.3), we obtain the corresponding Helmholtz equation:

∆ ξ Φ + ω 2 1 -β A(ω) Φ = 1 -β A(ω) √ ρ a δ(ξ) (4.4)
where Φ(ξ, ω) and A(ω) are the Fourier transforms of φ(ξ, t) and the kernel of the convolution operator A respectively. Let

κ(ω) = ω 2 1 -β A(ω) .
Then the solution of the Helmholtz equation (4.4) (see for instance [START_REF] Courant | Methods of Mathematical Physics[END_REF][START_REF] Nédélec | Acoustic and Electromagnetic Equations[END_REF]) is expressed as

Φ(x, ω) = √ ρ 1 -β A(ω) e √ -1κ(ω)τ (x)
4aπτ (x) .

where

τ (x) = √ ρ x 2 1 a 2 1 + x 2 2 a 2 2 + x 2 3 a 2 3
Using density normalized constants b j = a j √ ρ , we have

Φ(x, ω) = 1 -β A(ω) e √ -1κ(ω)τ (x) 4bρπτ (x) . ( 4
τ (x) = x 2 1 b 2 1 + x 2 2 b 2 2 + x 2 3 b 2 3

Solution of the Model Potential Problem

In this section, we find the solution of equation (2.11). We once again proceed with a model problem. Once the solution is obtained, we will aim to calculate, its second order derivatives for the evaluation of D ⊗ Dψ.

Solution of the Potential Problem

Let ψ(x, t), be the solution of equation (2.11) and Ψ(x, ω) be its Fourier transform with respect to variable t. Then Ψ(x, ω) satisfies,

M Ψ(x, ω) = Φ(x, ω) = 1 -β A(ω) e √ -1κ(ω)τ (x)
4bρπτ (x) .

(5.1)

When M is constant, the solution of this equation is directly calculated. As M = (∂ 2 1 + ∂ 2 2 )∆ x will not be used in the construction of Green function, we are only interested in the case where M is a homogeneous quadratic form in the component of

∇ x i.e. M = 3 j=1 m 2 j ∂ 2 ∂x 2 j .
So, the model equation (5.1) can be rewritten as:

3 j=1 m 2 j ∂ 2 Ψ ∂x 2 j = 1 -β A(ω) e √ -1κ(ω)τ (x) 4bρπτ (x) m j = 0, ∀j (5.2)
By a change of variables x j = m j η j , equation (5.2) becomes Poisson equation in Ψ(η, ω) = Ψ(x, ω) i.e.

∆ η Ψ = 1 -β A(ω) e √ -1κ(ω)τ (η) 4bρπτ (η) = Φ(η, ω) (5.3)
where,

τ (η) = m 2 1 η 2 1 b 2 1 + m 2 2 η 2 2 b 2 2 + m 2 3 η 2 3 b 2 3 = τ (x) and Φ(η, ω) = Φ(x, ω)
Notice that the source Φ(η, ω) is symmetric with respect to ellipsoid τ , i.e.

Φ(η, ω) = Φ(τ , ω).

Therefore, the solution Ψ of the Poisson equation ( 5.3) is the potential field of a uniformly charged ellipsoid due to a charge density Φ(τ , ω). The potential field Ψ can be calculated with a classical approach using ellipsoidal coordinates.

(See for example [START_REF] Chandrasekhar | Ellipsoidal Figures of Equilibrium[END_REF][START_REF] Kellogg | Foundations of Potential Theory[END_REF] for the theory of potential problems in ellipsoidal coordinates.) For the solution of the Poisson equation ( 5.3) we recall following result from [START_REF] Kellogg | Foundations of Potential Theory[END_REF]Ch. 7,Sec.6].

Proposition 5.1 Let f (z) = 3 j=1 ζ 2 j (α j h) 2 + z -1 and g(z) = Π 3 j=1 (α j h) 2 + z
and let Z(h, ζ) be the largest algebraic root of f (z)g(z) = 0. Then the solution of the Poisson equation

∆ 2 Y (ζ) = 4πχ ζ 2 1 α 2 1 + ζ 2 2 α 2 2 + ζ 2 2 α 2 1 ζ ∈ R 3
is given by

Y (ζ) = 2πα 1 α 2 α 3 ∞ 0 χ(h)I(h, ζ)dh.

The integrand I(h, ζ) is defined as

I(h, ζ) = h 2 ∞ Z(h,ζ) 1 g(z) dz Z > 0 h 2 ∞ 0 1 g(z) dz Z < 0
Hence, the solution of (5.3) can be given as

Ψ(η, ω) = 2πb m 1 -β A(ω) 1 4π ∞ 0 e √ -1κ(ω)h 4bρπh I(h, η)dh or equivalently, Ψ(x, ω) = 1 8ρπm 1 -β A(ω) ∞ 0 e √ -1κ(ω)h h I(h, x)dh, m = m 1 m 2 m 3 (5.4)
By a change of variable s = h -2 z, we can write I(h, x) as:

I(h, x) = mh ∞ S(h,x) 1 G(s) ds h < τ mh ∞ 0 1 G(s) ds h > τ (5.5)
with S(h, x) = h -2 Z(h, x) being the largest algebraic root of the equation

F (s)G(s) = 0
where

F (s) = h 2 f (h 2 s) = 3 j=1 {V j (s)} -1 x 2 j -h 2 G(s) = m 2 h 6 g(h 2 s) = Π 3 j=1 {V j (s)} with V j (s) = b 2 j + m 2 j s (5.6) Remark 5.2 Note that, F (s) ≡ 0 corresponds to a set of confocal ellipsoids s -→ h 2 (s) = 3 j=1 {V j (s)} -1 x 2 j (5.7)
such that τ (x) = h(0) i.e. S(τ ) = 0. Moreover, S > 0 if the ellipsoid h lies inside τ and S < 0 if the ellipsoid h lies outside τ .

Derivatives of the Potential field

Now we compute the derivatives of the potential Ψ. We note that I(h, x) is constant with respect to x when h > τ. So,

∂I(h, x) ∂x k = -mh ∂S(h, x) ∂x k 1 G(S(h, x)) h < τ 0 h > τ
for k = 1, 2, 3 and by consequence,

∂Ψ ∂x k = - 1 8ρπm 1 -β A(ω) ∞ 0 e √ -1κ(ω)h h ∂I(h, x) ∂x k dh or ∂Ψ ∂x k = - 1 8ρπ 1 -β A(ω) τ 0 e √ -1κ(ω)h ∂S(h, x) ∂x k 1 G(S(h, x)) dh. (5.8)
Now, we apply ∂ ∂x l for l = 1, 2, 3 on (5.8) to obtain the second order derivatives of Ψ:

-8ρπ ∂ 2 Ψ ∂x k x l = 1 -β A(ω) ∂ ∂x l τ 0 e √ -1κ(ω)h ∂S ∂x k 1 G(S) dh = 1 -β A(ω) ∂τ ∂x l e √ -1κ(ω)τ ∂S(τ ) ∂x k 1 G(S(τ )) + 1 -β A(ω) τ 0 e √ -1κ(ω)h 1 G(S) ∂ 2 S ∂x k ∂x l - 1 2 ∂S ∂x k ∂S ∂x l G ′ (S)

G(S) dh

As F (S)G(S) = 0 and G(s) is normally non-zero on S, therefore by differentiating F (S) = 0, we obtain [19, eq. (5.21)-(5.23)]

∂S ∂x k = -2x k V k (S)F ′ (S)
(5.9)

∂ 2 S ∂x k x l = -4x k x l V k (S)V l (S) [F ′ (S)] 2 F ′′ (S) F ′ (S) + m 2 k V k (S) + m 2 l V l (S) - 2δ kl V k (S)F ′ (S) (5.10) where,
F ′ (s) = 3 j=1 -m 2 j x 2 j V 2 j (s) , F ′′ (s) = 3 j=1 2m 4 j x 2 j V 3 j (s) , G ′ (s) = G(s) 3 j=1 m 2 j V j (s)
(5.11) and prime represents a derivative with respect to variable s.

Substituting the values from (5.9) and (5.10), the second order derivative of Ψ becomes 4ρπ

∂ 2 Ψ ∂x k x l = -x k x l 1 -β A(ω) aa 2 k a 2 l F ′ (0) e √ -1κ(ω)τ τ + 1 -β A(ω) τ 0 e √ -1κ(ω)h 1 F ′ (S) G(S) × 2x k x l V k (S)V l (S)F ′ (S) F ′′ (S) F ′ (S) + m 2 k V k (S) + m 2 l V l (S) + 1 2 G ′ (S) G(S) + δ kl V k (S) dh
(5.12)

Remark 5.3 If for some i ∈ {1, 2, 3}, m i → 0 one semi axis of the ellipsoid τ tends to infinity but no singularity occurs. Therefore the results of this section are still valid in this case.

Elastodynamic Green Function

In this section we present the expressions for the elastodynamic Green functions for the media presented in section 3. Throughout this section c p = c pp ρ with

p ∈ {1, 2, • • • , 6}. We recall that κ i (ω) = ω 2 1 -β i A(ω) .

Medium I

All the eigenvectors of Γ are constants in this case i.e. D i = e i , therefore

M i = 1 and E i = e i ⊗ e i .
If G is the Fourier transform of the viscoelastic Green function G with respect to variable t, then:

G = 3 i=1 Φ i (x, ω)e i ⊗e i = 1 4πρ 3 i=1   c i+3 1 -β i A(ω) c i c 4 c 5 c 6 τ i exp( √ -1κ i (ω)τ i )   e i e i (6.1)
where

τ 1 = x 2 1 c 2 1 + x 2 2 c 2 6 + x 2 3 c 2 5 , τ 2 = x 2 1 c 2 6 + x 2 2 c 2 2 + x 2 3 c 2 4 , τ 3 = x 2 1 c 2 5 + x c + x 2 3 c 2 3 , Medium b 1 b 2 b 3 m 1 m 2 m 3 M i I c 1 c 6 c 5 1 0 0 M 1 c 6 c 2 c 4 0 1 0 M 2 c 5 c 4 c 3 0 0 1 M 3 II c 4 c 4 c 3 0 0 1 M 1 c 1 c 1 c 4 1 1 0 M 2 c 6 c 6 c 4 * * * M 3 III c 1 c 1 c 1 1 1 1 M 1 c 6 c 6 c 4 1 1 0 M 2 c 4 c 4 c 4 * * * M 3
Table 1: Values of b i and m i for different media. Here * represents a value which is not used for reconstructing Green function.

Medium II

According to section 4, the functions Φ i have following expressions:

Φ 1 (x, ω) = 1 -β 1 A(ω) e √ -1κ1(ω)τ1(x) 4c 2 4 c 3 ρπτ 1 (x) Φ 2 (x, ω) = 1 -β 2 A(ω) e √ -1κ2(ω)τ2(x) 4c 2 1 c 4 ρπτ 2 (x) Φ 3 (x, ω) = 1 -β 3 A(ω) e √ -1κ3(ω)τ3(x) 4c 2 6 c 4 ρπτ 3 (x)
where

τ 1 (x) = x 2 1 c 2 4 + x 2 2 c 2 4 + x 2 3 c 2 3 , τ 2 (x) = x 2 1 c 2 1 + x 2 2 c 2 1 + x 2 3 c 2 4 , τ 3 (x) = x 2 1 c 2 6 + x c + x 2 3 c 2 4
To calculate Green function, we use the expression

G = Φ 3 I + D 1 ⊗ D 1 M -1 1 (Φ 1 -Φ 3 ) + D 2 ⊗ D 2 M -1 2 (Φ 2 -Φ 3 ) . D 1 = e 3 and M 1 = 1, yield D 1 ⊗ D 1 M -1 1 (Φ 1 -Φ 3 ) = (Φ 1 -Φ 3 ) e 3 ⊗ e 3 To compute D 2 ⊗ D 2 M -1 2 (Φ 2 -Φ 3 ), suppose Ψ 2 = M -1 2 Φ 2 and Ψ 3 = M -1 2 Φ 3
and notice that m 1 = m 2 = 1 and m 3 = 0. Moreover for Φ 2 and Φ 3 , b 1 = b 2 .

(See Table 1). Thus, we have 4ρπ

1 -β 2 A(ω) ∂ 2 Ψ 2 ∂x k x l = R k R l e √ -1κ2(ω)τ2 c 2 1 c 4 τ 2 - 1 c 4 R 2 (δ kl -2 R k R l ) τ2 0 e √ -1κ2(ω)h dh 4ρπ 1 -β 3 A(ω) ∂ 2 Ψ 3 ∂x k x l = R k R l e √ -1κ3(ω)τ3 c 2 6 c 4 τ 3 - 1 c 4 R 2 (δ kl -2 R k R l ) τ3 0 e √ -1κ3(ω)h dh where R k = x k R for k = 1, 2.
See Appendix C for the derivation of this result. By using the second derivatives of Ψ 2 and Ψ 3 and the expression

D 2 ⊗ D 2 M -1 2 (Φ 2 -Φ 3 ) = 2 k,l=1 ∂ k ∂ l (Ψ 2 -Ψ 3 ) e k ⊗ e l
we finally arrive at

G = 1 -β 3 A(ω) e √ -1κ3(ω)τ3(x) 4c 2 6 c 4 ρπτ 3 (x) J + 1 -β 1 A(ω) e √ -1κ1(ω)τ1(x) 4c 2 4 c 3 ρπτ 1 (x) e 3 ⊗ e 3 + 1 -β 2 A(ω) e √ -1κ2(ω)τ2(x) 4c 2 1 c 4 ρπτ 2 (x) -1 -β 3 A(ω) e √ -1κ3(ω)τ3(x) 4c 2 6 c 4 ρπτ 3 (x) R ⊗ R - 1 4ρπc 4 R 2 (J -2 R ⊗ R)× 1 -β 2 A(ω) τ2 0 e √ -1κ2(ω)h dh -1 -β 3 A(ω) τ3 0 e √ -1κ3(ω)h dh Or equivalently, G = Φ 1 e 3 ⊗ e 3 + Φ 2 R ⊗ R + Φ 3 (J -R ⊗ R) - 1 R 2 c 2 1 τ2 0 hΦ 2 (h, ω)dh -c 2 6 τ3 0 hΦ 3 (h, ω)dh (J -2 R ⊗ R)
Here J = Ie 3 ⊗ e 3 and R = R 1 e 1 + R 2 e 2

Medium III

The solutions of the wave equation Φ i in this case are

Φ 1 (x, ω) = 1 -β 1 A(ω) e √ -1κ1(ω)τ1(x) 4c 3 1 ρπτ 1 (x) Φ 2 (x, ω) = 1 -β 2 A(ω) e √ -1κ2(ω)τ2(x) 4c 2 6 c 4 ρπτ 2 (x) Φ 3 (x, ω) = 1 -β 3 A(ω) e √ -1κ3(ω)τ3(x) 4c 3 4 ρπτ 3 (x)
where

τ 1 (x) = 1 c 1 x 2 1 + x 2 2 + x 2 3 = r c 1 , τ 2 (x) = x 2 1 c 2 6 + x 2 2 c 2 6 + x 2 3 c 2 4 , τ 3 (x) = r c 4
To calculate Green function, we once again use the expression

G = Φ 3 I + D 1 ⊗ D 1 M -1 1 (Φ 1 -Φ 3 ) + D 2 ⊗ D 2 M -1 2 (Φ 2 -Φ 3 ) . Suppose Ψ 1 = M -1 1 Φ 1 and Ψ 3 = M -1 1 Φ 3 . Notice that m 1 = m 2 = m 3 = 1
for M 1 and b 1 = b 2 = b 3 for Φ 1 as well as Φ 3 (see Table 1). Thus, 4ρπ 1β 1 A(ω)

∂ 2 Ψ 1 ∂x k x l = r k r l e √ -1κ1(ω)τ1 c 3 1 τ 1 - 1 r 3 (δ kl -3 r i r j ) τ1 0 he √ -1κ1(ω)h dh 4ρπ 1 -β 3 A(ω) ∂ 2 Ψ 3 ∂x k x l = r k r l e √ -1κ3(ω)τ1 c 3 1 τ 3 - 1 r 3 (δ kl -3 r i r j ) τ3 0 he √ -1κ3(ω)h dh
See Appendix B for the derivation of this result. It yields

D 1 ⊗ D 1 M -1 1 (Φ 1 -Φ 3 ) = 1 4ρπ 1 -β 1 A(ω) e √ -1κ1(ω)τ1(x) c 3 1 τ 1 (x) + 1 -β 3 A(ω) e √ -1κ3(ω)τ3(x) c 3 4 τ 3 (x) r ⊗ r -1 -β 1 A(ω) τ1 0 he √ -1κ1(ω)h dh -1 -β 3 A(ω) τ3 0 he √ -1κ3(ω)h dh × 1 4ρπr 3 (I -3 r ⊗ r) = [Φ 1 (x, ω) -Φ 3 (x, ω)] r ⊗ r - 1 r 3 τ1 0 h 2 Φ 1 (h, ω)dh - τ3 0 h 2 Φ 3 (h, ω)dh (I -3 r ⊗ r)
.

where r = r 1 e 1 + r 2 e 2 + r 3 e 3 with r i = x i r for all i = 1, 2, 3.

To compute,

D 2 ⊗ D 2 M -1 2 (Φ 2 -Φ 3 ), suppose Ψ 2 = M -1 2 Φ 2 and Ψ 4 = M -1
2 Φ 3 . By using formula (C.3) with m 1 = m 2 = 1 and m 3 = 0, we obtain:

4ρπ 1 -β 2 A(ω) ∂ 2 Ψ 2 ∂x k x l = R k R l e √ -1κ2(ω)τ2 c 2 6 c 4 τ 2 - 1 c 4 R 2 (δ kl -2 R k R l ) τ2 0 e √ -1κ2(ω)h dh 4ρπ 1 -β 3 A(ω) ∂ 2 Ψ 4 ∂x k x l = R k R l e √ -1κ3(ω)τ3 c 3 4 τ 3 - 1 c 4 R 2 (δ kl -2 R k R l ) τ3 0 e √ -1κ3(ω)h dh with R k = x k R
and k, l ∈ {1, 2}. This allows us to write

D 2 ⊗ D 2 M -1 2 (Φ 2 -Φ 3 ) = 1 4ρπ 1 -β 2 A(ω) e √ -1κ2(ω)τ2(x) c 3 1 τ 2 (x) + 1 -β 3 A(ω) e √ -1κ3(ω)τ3(x) c 3 4 τ 3 (x) × R 2 2 e 1 ⊗ e 1 -R 1 R 2 [e 1 ⊗ e 2 + e 2 ⊗ e 1 ] + R 2 1 e 2 ⊗ e 2 - 1 4c 4 ρπR 2 1 -β 2 A(ω) τ2 0 e √ -1κ2(ω)h dh -1 -β 3 A(ω) τ3 0 e √ -1κ3(ω)h dh × (1 -2 R 2 2 )e 1 ⊗ e 1 -2 R 1 R 2 [e 1 ⊗ e 2 + e 2 ⊗ e 1 ] + (1 -2 R 2 1 )e 2 ⊗ e 2 = [Φ 2 (x, ω) -Φ 3 (x, ω)] R ⊥ ⊗ R ⊥ - 1 R 2 c 2 6 τ2 0 hΦ 2 (h, ω)dh -c 2 4 τ3 0 hΦ 3 (h, ω)dh J -2 R ⊥ ⊗ R ⊥ .
where R ⊥ = R 2 e 1 -R 1 e 2 and J = Ie 3 ⊗ e 3 . Finally, we arrive at

G = Φ 1 r ⊗ r + Φ 2 R ⊥ ⊗ R ⊥ + Φ 3 (I -r ⊗ r -R ⊥ ⊗ R ⊥ ) - 1 r 3 τ1 0 h 2 Φ 1 (h, ω)dh - τ3 0 h 2 Φ 3 (h, ω)dh (I -3 r ⊗ r) - 1 R 2 c 2 1 τ2 0 hΦ 2 (h, ω)dh -c 2 6 τ3 0 hΦ 3 (h, ω)dh (J -2 R ⊥ ⊗ R ⊥ )

Isotropic Medium

When c 66 = c 44 , medium III becomes isotropic. In this case

Φ 2 (x, ω) = Φ 3 (x, ω), β 1 = β 2 , τ 1 (x) = r c 1 , and τ 2 (x) = r c 4 = τ 3 (x)
Thus, the Green function in an isotropic medium with independent elastic parameters c 11 and c 44 can be given in frequency domain as:

G = Φ 2 I + D 1 ⊗ D 1 M -1 1 (Φ 1 -Φ 2 ) = Φ 1 r ⊗ r + Φ 2 (I -r ⊗ r) - 1 r 3 r c 1 0 h 2 Φ 1 (h, ω)dh - r c 4 0 h 2 Φ 2 (h, ω)dh (I -3 r ⊗ r)
where Φ 1 and Φ 2 are the same as in medium III. This expression of the Green function has already been reported in a previous work [START_REF] Bretin | On the Green function in visco-elastic media obeying a frequency power-law[END_REF].

Acknowledgement

We would like to thank Prof. Habib Ammari ( École Normale Supérieure-Paris) for his continuous support and encouragement. We would also like to thank Meisam Sharify ( École Polytechnique-Paris) for his help and remarks. E. Bretin was supported by the foundation Digitio of France, in terms of a Post-doctoral Fellowship. A. Wahab was supported by Higher Education Commission of Pakistan in terms of a doctoral fellowship.

A Decomposition of the Green function

Consider the elastic equation satisfied by G:

(Γ c (∇ x )G(x, t) + Γ v (∇ x )A[G](x, t)) -ρ ∂ 2 G(x, t) ∂t 2 = δ(t)δ(x)I. (A.1) If G is given in the form G = 3 i=1 E i (∇ x )φ i (A.2)
Then substituting (A.2) in (A.1) yields:

δ(t)δ(x)I = (Γ c (∇ x )G(x, t) + Γ v (∇ x )A[G](x, t)) -ρ ∂ 2 G(x, t) ∂t 2 = 3 i,j=1 L c j (∇ x )φ i + L v j (∇ x )A[φ i ] E j (∇ x )E i (∇ x ) -ρ 3 i=1 E i (∇ x ) ∂ 2 φ i (x, t) ∂t 2
By definition E i (∇ x ) is a projection operator which satisfies

E i (∇ x )E j (∇ x ) = δ ij E j (∇ x )
Consequently, we can have

δ(t)δ(x)I = 3 i,j=1 E j (∇ x )δ ij ρ -1 L c j (∇ x )φ i + L v j (∇ x )A[φ i ] -ρ 3 i=1 E i (∇ x ) ∂ 2 φ i (x, t) ∂t 2 = 3 i=1 E i (∇ x ) (L c i (∇ x )φ i + L v i (∇ x )A[φ i ]) -ρ ∂ 2 φ i (x, t) ∂t 2 . Moreover I = 3 i=1 E i (∇ x ), therefore 3 i=1 E i (∇ x ) (L c i (∇ x )φ i + L v i (∇ x )A[φ i ]) -ρ ∂ 2 φ i (x, t)
V 1 (s) = V 2 (s) = V 3 (s) = b 2 1 + m 2 1 s F (s) = 3 j=1 x 2 j V 1 (s) -h 2 = r 2 V 1 (s) -h 2 F ′ (s) = 3 j=1 -m 2 1 x 2 j V 2 1 (s) = -m 2 1 r 2 V 2 1 (s)
and

F ′ (0) = -m 2 1 r 2 b 4 1 F ′′ (s) = 3 j=1 2m 4 1 x 2 j V 3 1 (s) = 2m 4 1 r 2 V 3 1 (s) G(s) = (V 1 (s)) 3 and G ′ (s) = G(s) 3m 2 1 V 1 (s) (B.1)
with r = x 2 1 + x 2 2 + x 2 3 . When F (S) = 0, we have

V 1 (S) =
F ′ (s) = 2 j=1 -m 2 1 x 2 j V 2 1 (s) = -m 2 1 R 2 V 2 1 (s)
and

F ′ (0) = -m 2 1 R 2 b 4 1 F ′′ (s) = 2 j=1 2m 4 1 x 2 j V 3 1 (s) = 2m 4 1 R 2 V 3 1 (s) G(s) = b 2 3 (V 1 (s))

. 5 )

 5 where constant b = b 1 b 2 b 3 and

∂t 2 -

 2 δ(t)δ(x) = 0 Finally, remark that G we can express in the form (2.8) if the functions φ i satisfy equation (2.10). B Derivative of Potential: Case I If b 1 = b 2 = b 3 and m 1 = m 2 = m 3 , we have

r 2 h 2 , 1 VC

 1 k (S)V l (S)F ′ (.1) and (B.2) in (5.12) we finally arrive at:4ρm 2 1 π 1β A(ω) ∂ 2 Ψ ∂x k x l = r k r l e √ -1κ(ω)τ bτ -1 r 3 (δ kl -3 r i r j ) Derivative of Potential: Case II If b 1 = b 2 , m 1 = m 2 and m 3 = 0, we have V 1 (s) = V 2 (s) = b 2 1 + m 2 1 s andV 3 (s) = b 2 3

2 and 1 V 1 x 2 2 . 1 -

 21121 G ′ (s) = G(s) 2m 2 For all l, k ∈ {1, 2}, we have 1 V k (S)V l (S)F ′ (β A(ω) ∂ 2 Ψ ∂x k x l = R k R l e √ -1κ(ω)τ bτ -1 b 3 R 2 (δ kl -2 R k R l ) ω)h dh (C.3) where R k = x k R for k = 1, 2.