Section: New Results

New schemes for time-domain simulations

Solving the Isotropic Linear Elastodynamics Equations Using Potentials

Participant : Patrick Joly.

This work is done in collaboration with Sébastien Impériale (EPI M3DISIM) and Jorge Albella and Jeronimo Rodríguez from the University of Santiago de Compostela.

We pursue our research on the numerical solution of 2D elastodynamic equations in piecewise homogeneous media using the decomposition of the displacement fields into the sum of the gradient and the rotational (respectively) of two scalar potentials potentials. This allows us to obtain an automatic decomposition of the wave field into the sum of pressure and shear waves (respectively). The approach is expected to be efficient when the velocity of shear waves is much smaller than the velocity of pressure waves, since one can adapt the discretization to each type of waves. This appears as a challenge for finite element methods , the most delicate issue being the treatment of boundary and transmission conditions, where the two potentials are coupled..

A stable (mixed) variational formulation of the evolution problem based on a clever choice of Lagrange multipliers has been proposed as well as various finite element approximations which have been successfully implemented. The analysis of the continuous problem has been published in a long paper in the journal of Scientific computing.The numerical analysis of the discretized problem is in progress.

Time domain Half-Space Matching method

Participants : Sonia Fliss, Hajer Methenni.

This work is done in the framework of the PhD of Hajer Methenni (funded by CEA-LIST) and in collaboration with Sebastien Imperiale (EPI M3DISIM) and Alexandre Imperiale (CEA-LIST).

The objective of this work is to propose a numerical method to solve the elastodynamics equations in a locally perturbed unbounded anisotropic media. Let us mention that all the classical methods to restrict the computation around the perturbations are unstable in anisotropic elastic media (PMLs for instance) or really costly (Integral equations). The idea is to extend the method already developed for the corresponding time harmonic problem, called the Halfspace Matching Method. We have considered, for now, the 2D scalar wave equation but the method is constructed in order to be applied to the elastodynamic problem. The method consists in coupling several representations of the solution in half-planes surrounding the defect with a FE representation in a bounded domain including the defect. In order to ensure the stability of the method, we first semi-discretize in time the equations and apply the method to the semi-discrete problem. Thus, for each time step, by ensuring that all the representations of the solution match, in particular in the intersection of the half-planes, we end up, at each time step, with a system of equations which couples, via integral operators, the solution at this time step in the bounded domain and its traces on the edge of the half-planes, the right hand side being a convolution operator involving the solution at the previous time steps. The method has been implemented and validated with Xlife++.

We are now looking to make the method more efficient by implementing methods of acceleration. Finally, we will also seek to develop another version of the method based on the Convolution quadrature.

Time domain modelling for wave propagation in fractal trees

Participants : Patrick Joly, Maryna Kachanovska.

In order to simulate wave propagation in fractal trees (see section 7.4.3), which have infinite structure, it is necessary to be able to truncate the computations to a finite subtree. This was done using Dirichlet-to-Neumann (DtN) operators in our previous work in collaboration with A. Semin (TU Darmstadt). In this case a DtN operator is a convolution operator, whose kernel is not known in a closed form. Based on the results of this previous work, in 2017 we had proposed two methods for approximating these convolution operators:

  • constructing an exact DtN operator for a semi-discretized system (in the spirit of convolution quadrature methods).

  • truncating meromorphic expansion for the symbol (Fourier transform of the convolution kernel) of the DtN operator, which allows to approximate the DtN operator by local operators.

This year we have performed a complete convergence and stability analysis of these methods, based on the energy techniques.

In particular, for the convolution quadrature methods, we were able to obtain all the estimates using time-domain analysis, by avoiding passage to the Laplace domain.

As for the method based on the meromorphic expansion of the symbol of the DtN operator, we have shown that the error induced by truncating the expansion to L terms can be controlled by a remainder of a series, which, in particular, depends on the eigenvalues of the weighted Laplacian on the fractal trees. To obtain an explicit dependence of the error on L, we have computed Weyl bounds for the eigenvalues, based on a refinement of the ideas of [Kigami, Lapidus, Comm. Math. Phys. 158 (1993)].

Additionally, we have addressed some computational aspects of the two methods, in particular, efficient evaluation of the symbol of the DtN operator (we have an algorithm that allows to evaluate it at the frequency ω in O(logk|ω|) time), as well as a method for efficient computation of the poles of the symbol (based on Möbius transform and polynomial interpolation).