Hélène Barucq

Henri Calandra
email: henri.calandra@total.com

Julien Diaz

Elvira Shishenina
email: elvira.shishenina@inria.fr

Trefftz -Discontinuous Galerkin Approach for Solving Elastodynamic Problem

Methods based on Discontinuous Finite Element approximation (DG FEM) are basically well-adapted to specifics of wave propagation problems in complex media, due to their numerical accuracy and flexibility. However, they still lack of computational efficiency, by reason of the high number of degrees of freedom required for simulations.

The Trefftz-DG solution methodology investigated in this work is based on a formulation which is set only at the boundaries of the mesh. It is a consequence of the choice of test functions that are local solutions of the problem. It owns the important feature of involving a space-time approximation which requires using elements defined in the space-time domain.

Herein, we address the Trefftz-DG solution of the Elastodynamic System. We establish its well-posedness which is based on mesh-dependent norms. It is worth noting that we employ basis functions which are space-time polynomial. Some numerical experiments illustrate the proper functioning of the method.

Introduction

Among the different possible approaches to solve partial differential equations there exists a distinct family of methods based on the use of trial functions in the form of exact solution of the governing equations (but not the boundary conditions). The idea was first proposed by Trefftz in 1926 [START_REF] Trefftz | Ein Gegenstuck zum Ritzschen Verfahren[END_REF], and since then it has been largely developed and generalized. The main step for its implementation as an efficient computational tool was achieved in 1978 when Jirousek and his collaborators proposed the Hybrid-Trefftz (HT) finite element model [START_REF] Jirousek | Basis for development of large finite elements locally satisfying all field equations[END_REF]. The results of their work allowed solving different boundary value problems, thus giving roots to multiple applications in different fields such as potential problems, plane elasticity, plate bending (thin, thick, post-buckling), heat conduction as well as advective-diffusive transport (see [START_REF] Zienkiewic | Trefftz type approximation and the generalized finite element method -history and development[END_REF] and the references therein).

Trefftz type methods have been widely used with time-harmonic formulations by Farhat, Tezaur, Harari, Hetmaniuk (2003[START_REF] Tezaur | Three-dimensional discontinuous Galerkin elements with plane waves and Lagrange multipliers for the solution of midfrequency Helmholtz problems[END_REF] (see [START_REF] Farhat | A discontinuous Galerkin method with Lagrange multipliers for the solution of Helmholtz problems in the midfrequency regime[END_REF][START_REF] Tezaur | Three-dimensional discontinuous Galerkin elements with plane waves and Lagrange multipliers for the solution of midfrequency Helmholtz problems[END_REF]), [START_REF] Gabard | Discontinuous Galerkin methods with plane waves for timeharmonic problems[END_REF] (see [START_REF] Gabard | Discontinuous Galerkin methods with plane waves for timeharmonic problems[END_REF]), [START_REF] Badics | Trefftz-discontinuous Galerkin and finite element multi-solver technique for modeling time-harmonic EM problems with high-conductivity regions[END_REF] (see [START_REF] Badics | Trefftz-discontinuous Galerkin and finite element multi-solver technique for modeling time-harmonic EM problems with high-conductivity regions[END_REF]), Hiptmair, Moiola, Perugia (2011 -2013) (see [START_REF] Hiptmair | Plane wave discontinuous Galerkin methods for the 2D Helmholtz equation: analysis of the p-version[END_REF][START_REF] Moiola | Plane wave approximation of homogeneous Helmholtz solutions[END_REF][START_REF] Hiptmair | Error analysis of Trefftzdiscontinuous Galerkin methods for the time-harmonic Maxwell equations[END_REF]) and others, while studies are still limited for reproducing temporal phenomena. Only few papers are interested in Maxwell equations in time [START_REF] Petersen | A space-time discontinuous Galerkin method for the solution of the wave equation in the time domain[END_REF][START_REF] Kretzschmar | Discontinuous Galerkin methods with Trefftz approximations[END_REF][START_REF] Egger | A Space-Time Discontinuous Galerkin Trefftz Method for Time Dependent Maxwell's Equations[END_REF][START_REF] Kretzschmar | A priori error analysis of space-time Trefftz discontinuous Galerkin methods for wave problems[END_REF]. They are mostly devoted to a theoretical analysis of the method, showing the convergence and stability. To the best of our knowledge, numerical tests involving plane wave approximation are restricted to 1D + time dimensional case. Space-time Trefftz approximation by Lagrange multipliers for the second order formulation of the transient wave equation was explored in [START_REF] Wang | A hybrid discontinuous in space and time Galerkin method for wave propagation problems[END_REF][START_REF] Banjai | A Trefftz polynomial spacetime discontinuous Galerkin method for the second order wave equation[END_REF]. A Trefftz-DG method for the first-order transient acoustic wave equations in arbitrary space dimensions has been introduced in a recent paper of Moiola and Perugia (2017) [START_REF] Moiola | A space-time Trefftz discontinuous Galerkin method for the acoustic wave equation in first-order formulation[END_REF]. It is an extension of the one-dimensional scheme of Kretzschmar et al. [START_REF] Kretzschmar | A priori error analysis of space-time Trefftz discontinuous Galerkin methods for wave problems[END_REF]. The authors provide a complete a priori error analysis in both mesh-dependent and mesh-independent norms.

In this work we develop a theory for Trefftz-DG method applied to the Elastodynamic System (ES) of wave propagation. We confirm well-posedness of the variational problem based on the estimates in mesh-dependent norms. We consider a space-time Trefftz polynomial basis and provide some numerical results for 2D ES. We give a short conclusion and discuss the perspectives in the end of this paper.

Trefftz-DG formulation for elastodynamics

The Elastodynamic System is based on three fundamental laws of continuum mechanics: movement equations, constitutive equations (Hooke's law), and geometric equations (infinitesimal strain tensor definition) [START_REF] Tallec | Modélisation et calcul des milieux continus[END_REF].

We consider a global space-time domain Q = Ω × I, where Ω ⊂ R d is a bounded Lipschitz physical space domain and I = (0, T ) is a time interval.

The Lamé coefficients λ ≡ λ(x), µ ≡ µ(x) and solid density ρ ≡ ρ(x) are the solid parameters, assumed to be piecewise constant and positive.

We consider the first order ES in terms of velocity v ≡ v(x, t) and stress σ ≡ σ(x, t) fields:

       ∂tσ -C ε(v) = 0 in Q, ρ∂tv -divσ = 0 in Q, v(•, 0) = v0, σ(•, 0) = σ0 in Ω, σ = gD in ∂Ω × I.
Here C is the elastic coefficient tensor, ε(v) = (∇v+∇v T )/2 is the infinitesimal strain tensor. The boundary conditions g D ≡ g D (x, t), the velocity v 0 and the stress σ 0 are the initial data. By symmetry and positivity of the tensor C, the application ε -→ C ε is an isomorphism in the symmetrical tensor space [START_REF] Tallec | Modélisation et calcul des milieux continus[END_REF]. Thus, we may consider the corresponding inverse application A, verifying the same properties of symmetry and positivity:

         A∂tσ -ε(v) = 0 in Q, ρ∂tv -divσ = 0 in Q, v(•, 0) = v0, σ(•, 0) = σ0 in Ω, σ = gD in ∂Ω × I.
(1)

Space-time Trefftz-DG formulation

We choose a Lipschitz sub-domain K ⊂ Q such that λ, µ and ρ are constant in K. We define n K ≡ (n x K , n t K ) as the outward pointing unit normal vector on ∂K and (v,

σ) ∈ H 1 (K) d × H 1 (K) d 2
, where d is the dimension of the physical space Ω.

Multiplying both equations of (1) by the test functions (ω, ξ) ∈ H 1 (K) d × H 1 (K) d 2 respectively, and integrating by part in space and time twice, we obtain:

- K σ(A∂tξ -ε(ω)) + v • (ρ∂tω -divξ) dv + ∂K A σ : ξn t K -v • n x K ξ + ρvn t K • ω -σn x K • ω ds = 0.
(

) 2
Without losing generality with respect to the classical space DG methods, we introduce a non-overlapping mesh T h on Q, whose elements are right prisms, with vertical sides parallel to the time axis. All solid parameters inside the elements K are supposed to be constant, so that all discontinuities lie on the inter-element boundaries.

The mesh skeleton F h = ∪ K∈T h ∂K can be decomposed into the families of element faces:

F I h internal I-faces (x -fixed) F Ω h internal Ω-faces (t -fixed) F D h external Dirichlet boundary faces (∂Ω × [0, T ]) F 0 h external initial time faces (Ω × {0}) F T h external final time faces (Ω × {T })
The space-time DG formulation for (1) consists in searching

(v hp , σ hp ) ∈ V(T h ) ⊂ H 1 (T h ) d × H 1 (T h ) d 2 such
that, for all K ∈ T h and for all (ω, ξ) ⊂ V(T h ) the following identity holds true:

- K σ hp : (A∂tξ -ε(ω)) + v hp • (ρ∂tω -divξ) dv + ∂K A σhp : ξn t K -vhp • n x K ξ + ρv hp n t K • ω -σhp n x K • ω ds = 0. (3)
The numerical fluxes vhp , σhp are defined in the standard DG notations [START_REF] Hesthaven | Nodal discontinuous Galerkin methods: algorithms, analysis, and applications[END_REF] on the mesh skeleton F h as follows:

Here δ ∈ L ∞ (F I h ∪ F D h ) and γ ∈ L ∞ (F I h ) are positive penalty parameters. This choice is recommended in order to improve the numerical stability of the scheme.

We define the Trefftz approximation space, such that the chosen test functions (ω, ξ) satisfy the initial elastodynamic system in the homogeneous

F I h vhp σhp ≡ { {v hp } } -δ[[σ hp ]]x { {σ hp } } -γ[[v hp ]]x F Ω h vhp σhp ≡ v hp - σ hp - F D h vhp • n x K σhp n x K ≡ v hp • n x K -δ(σ hp -gD)n x K gDn x K F T h vhp σhp ≡ v hp σ hp F 0 h vhp σhp ≡ v0 σ0
sense (without source term and boundary conditions):

T(T h ) ≡ (ω, ξ) ∈ V(T h ) s. t. ρ∂tω -divξ = A∂tξ -ε(ω) = 0 in all K ∈ T h .
Thanks to this choice of discrete space we remove a volume integration term in [START_REF] Zienkiewic | Trefftz type approximation and the generalized finite element method -history and development[END_REF]. Summing over all elements K ∈ T h , we obtain a space-time Trefftz-DG formulation for (1):

Seek (v hp , σ hp ) ∈ V(T h ) such that, for all (ω, ξ) ∈ T(T h ), it holds true:

-

F I h { {σ hp } }[[ω]]x + { {v hp } }[[ξ]]x -γ[[v hp ]]x[[ω]]x -δ[[σ hp ]]x[[ξ]]x ds + F Ω h A σ hp -: [[ξ]]t + ρv hp -[[ω]]t ds - F D h v hp • n x K ξ -δσ hp : ξ ds + F T h A σ hp : ξ + ρv hp • ω ds - 1 2 F 0 h A σ hp : ξ + ρv hp • ω ds = 1 2 F 0 h A σ hp : ξ + ρv hp • ω ds + F D h gDn x K • ω + δgD : ξ) ds,
or, by introducing bilinear A T DG (• ; •) and linear T DG (•) operators: Seek (v hp , σ hp ) ∈ V(T h ) such that, for all (ω, ξ) ∈ T(T h ), it holds true:

AT DG((vhp, σ hp ); (ω, ξ hp )) = T DG(ω, ξ hp ). (4)
It is worth mentioning that in addition to its setting at the boundaries of the elements only, the Trefftz-DG formulation does not involve differential operators. It is thus straightforward to implement.

Well-posedness of Trefftz-DG formulation

In this section we show the coercivity and continuity estimates proving wellposedness of the obtained Trefftz-DG method for ES in mesh-dependent norms. We refer to the Appendix B in [START_REF] Barucq | Space-Time Trefftz -Discontinuous Galerkin Approximation for Elasto-Acoustics[END_REF] for more details. The analysis is carried out inside the framework developed in [START_REF] Kretzschmar | A priori error analysis of space-time Trefftz discontinuous Galerkin methods for wave problems[END_REF] for the time-dependent Maxwell problem.

We introduce two mesh-dependent norms in T(T h ):

|||(ω, ξ)||| 2 T DG ≡ 1 2 (A) 1 / 2 [[ξ]]t 2 L 2 (F Ω h ) + 1 2 ρ 1 / 2 [[ω]]t 2 L 2 (F Ω h ) + γ 1 / 2 [[ω]]x 2 L 2 (F I h ) + δ 1 / 2 [[ξ]]x 2 L 2 (F I h ) + 1 2 (A) 1 / 2 ξ 2 L 2 (F T h ) + 1 2 ρ 1 / 2 ω 2 L 2 (F T h ) + δ 1 / 2 ξ 2 L 2 (F D h ) , |||(ω, ξ)||| 2 T DG * ≡ |||(ω, ξ)||| 2 T DG + ρ 1 / 2 ω -2 L 2 (F T h )F Ω h + (A) 1 / 2 ξ -2 L 2 (F T h )F Ω h + δ -1 / 2 { {ω} } 2 L 2 (F I h ) + γ -1 / 2 { {ξ} } 2 L 2 (F I h ) + δ -1 / 2 ξ 2 L 2 (F T h )F D h .
Thus, for the bilinear A T DG (•, •) and linear T DG (•) forms we obtain the following coercivity

AT DG((ω, ξ); (ω, ξ)) = |||(ω, ξ)||| 2 T DG , ∀(ω, ξ) ∈ T(T h ),
and continuity properties with respect to the chosen norms.

|AT DG((v, σ); (ω, ξ))| ≤ 2 |||(v, σ)|||T DG * |||(ω, ξ)|||T DG, | T DG(ω, ξ)| ≤ √ 2 ρ 1 / 2 v0 2 L 2 (F 0 h ) + A 1 / 2 σ0 2 L 2 (F 0 h ) 1 / 2 (gD ≡ 0).
The above estimates confirm the well-posedness of the Trefftz-DG problem for ES, moreover:

|||(v -v hp , σ -σ hp )|||T DG ≤ 3 inf (ω,ξ)∈V(T h ) |||(v -ω, σ -ξ)|||T DG * .

Numerical implementation

In this Section we discuss the choice of the discrete approximation space, and provide some numerical tests for 2D + time elastic model.

Polynomial basis

The flexibility in the choice of basis functions is one of the advantages of Trefftz-type methods. The main condition is to satisfy the governing equations inside each element. The natural choice in the case of harmonic problems can be the plane wave trigonometric basis. However, when applied to the spacetime formulations, it demands a high number of trigonometric functions of different frequencies, in order to provide a better approximation order. Thus, it increases the number of degrees of freedom, and as a result -the global numerical cost of the algorithm.

We have computed a space-time polynomial basis, using generating exponential functions -local solutions of the initial systems of equations (see [START_REF] Maciag | Solution of the two-dimensional wave equation by using wave polynomials[END_REF][START_REF] Barucq | Space-Time Trefftz -Discontinuous Galerkin Approximation for Elasto-Acoustics[END_REF] for more details). Numerically speaking, this basis generates a lower computational bound than the standard trigonometric ones. It contains the couples of polynomial functions for velocity -stress, of degrees less or equal to n (n = 0, 1, 2, 3), satisfying the initial ES, to provide an approximation of order n.

Numerical results

In order to explore the method and its algorithm in general, and to perform some basic numerical tests for its validation, we have developed a prototype MATLAB R code, which is, technically speaking, quite limited.

We consider the elastic medium, represented by a unit square with the Dirichlet conditions at the boundaries. The final time of propagation is T = 1. The medium parameters are ρ S = 1, λ = 1, µ = 2. All model parameters are dimensionless quantities. We introduce a source term in the center of the medium. The source signal is represented by the Gaussian function, so that it takes approximately 5 elements per wavelength. Zero initial conditions are imposed for the tests.

Figure 1 shows some results of convergence of the (a) P -velocity and (b) S-velocity. The convergence curves have been computed for different approximation orders (n = 0, 1, 2, 3), and they represent the numerical error as a function of cell size in logarithmic scale. Even though the initial model is limited (large mesh of 30 × 30 × 30 elements, Dirichlet boundaries, which causes many reflections), the numerics reproduce the expected propagation characteristics quite well, and the convergence in both cases is of order higher than the corresponding approximation order.

Conclusion

We have applied the theory for Trefftz-DG method to the Elastodynamic System, and we have studied the well-posedness of the problem. The new polynomial basis has been computed for numerical implementation of the method. The computed numerical solutions have been validated by the analytical ones. The convergence results are of higher order, compared with the classical DG methods.

Even though the obtained results are very promising, from the optimisation point of view, it seems necessary to study the alternative to a global matrix inversion, which brings the main computational cost. We have also in perspective to pass from simple rectangular meshes to more complicate forms (in space domain) -which is one of main advantages of Trefftz-DG methods. It gives the possibility of developing a hybrid method, based on numerical coupling of Trefftz-DG method in Elastics, with less expensive Finite Volume Method (FVM) in Acoustics, in order to create a software able to solve more realistic problems.

Fig. 1 .

 1 Fig. 1. Convergence of velocities vP and vS in function of cell size h = ∆x. (a) velocity vP (b) velocity vS

Acknowledgments

This work is supported by the Inria -Total S.A. strategic action "Depth Imaging Partnership" (http://dip.inria.fr).