## Section: New Results

### Mathematical modeling of multi-physics involving wave equations

#### Atmospheric Radiation Boundary Conditions for the Helmholtz Equation

Participants : Hélène Barucq, Juliette Chabassier, Marc Duruflé.

An article is to be published in M2AN, see [14]. This work offers some contributions to the numerical study of acoustic waves propagating in the Sun and its atmosphere. The main goal is to provide boundary conditions for outgoing waves in the solar atmosphere where it is assumed that the sound speed is constant and the density decays exponentially with radius. Outgoing waves are governed by a Dirichlet-to-Neumann map which is obtained from the factorization of the Helmholtz equation expressed in spherical coordinates. For the purpose of extending the outgoing wave equation to axisymmetric or 3D cases, different approximations are implemented by using the frequency and/or the angle of incidence as parameters of interest. This results in boundary conditions called Atmospheric Radiation Boundary Conditions (ARBC) which are tested in ideal and realistic configurations. These ARBCs deliver accurate results and reduce the computational burden by a factor of two in helioseismology applications. This work has been done in collaboration with Laurent Gizon and Michael Leguèbe (Max-Planck-Institut für Sonnensystemforschung, Gottingen, Germany).

#### Atmospheric radiation boundary conditions for high frequency waves in time-distance helioseismology

Participants : Hélène Barucq, Juliette Chabassier, Marc Duruflé.

An article has been published in Astronomy and Astrophysics [22]. The temporal covariance between seismic waves measured at two locations on the solar surface is the fundamental observable in time-distance helioseismology. Above the acoustic cutoff frequency ( 5.3 mHz), waves are not trapped in the solar interior and the covariance function can be used to probe the upper atmosphere. We wish to implement appropriate radiative boundary conditions for computing the propagation of high-frequency waves in the solar atmosphere. We consider the radiative boundary conditions recently developed by Barucq et al. (2017) for atmospheres in which sound-speed is constant and density decreases exponentially with radius. We compute the cross-covariance function using a finite element method in spherical geometry and in the frequency domain. The ratio between first-and second-skip amplitudes in the time-distance diagram is used as a diagnostic to compare boundary conditions and to compare with observations. We find that a boundary condition applied 500 km above the photosphere and derived under the approximation of small angles of incidence accurately reproduces the 'infinite atmosphere' solution for high-frequency waves. When the radiative boundary condition is applied 2 Mm above the photosphere, we find that the choice of atmospheric model affects the time-distance diagram. In particular, the time-distance diagram exhibits double-ridge structure when using a VAL atmospheric model. This is a collaboration with Damien Fournier, Laurent Gizon, Chris Hanson and Michael Leguèbe (Max-Planck-Institut für Sonnensystemforschung, Gottingen, Germany).

#### Computational helioseismology in the frequency domain: acoustic waves in axisymmetric solar models with flows

Participants : Hélène Barucq, Juliette Chabassier, Marc Duruflé.

An article has been published in Astronomy and Astrophysics [23]. Context. Local helioseismology has so far relied on semi-analytical methods to compute the spatial sensitivity of wave travel times to perturbations in the solar interior. These methods are cumbersome and lack flexibility. Aims. Here we propose a convenient framework for numerically solving the forward problem of time-distance helioseismology in the frequency domain. The fundamental quantity to be computed is the cross-covariance of the seismic wavefield. Methods. We choose sources of wave excitation that enable us to relate the cross-covariance of the oscillations to the Green's function in a straightforward manner. We illustrate the method by considering the 3D acoustic wave equation in an axisymmetric reference solar model, ignoring the effects of gravity on the waves. The symmetry of the background model around the rotation axis implies that the Green's function can be written as a sum of longitudinal Fourier modes, leading to a set of independent 2D problems. We use a high-order finite-element method to solve the 2D wave equation in frequency space. The computation is embarrassingly parallel, with each frequency and each azimuthal order solved independently on a computer cluster. Results. We compute travel-time sensitivity kernels in spherical geometry for flows, sound speed, and density perturbations under the first Born approximation. Convergence tests show that travel times can be computed with a numerical precision better than one millisecond, as required by the most precise travel-time measurements. Conclusions. The method presented here is computationally efficient and will be used to interpret travel-time measurements in order to infer, e.g., the large-scale meridional flow in the solar convection zone. It allows the implementation of (full-waveform) iterative inversions, whereby the axisymmetric background model is updated at each iteration. This work is a collaboration with Aaron Birch, Damien Fournier, Laurent Gizon, Chris Hanson, Michael Leguèbe and Emanuele Papini(Max-Planck-Institut für Sonnensystemforschung, Gottingen, Germany) and with Thorsten Hohage (Gôttingen University, Germany).

#### The virtual workshop : towards versatile optimal design of musical wind instruments for the makers

Participants : Juliette Chabassier, Robin Tournemenne.

Our project aims at proposing optimization solutions for wind instrument making. Our approach is based on a strong interaction with makers and players, aiming at defining interesting criteria to optimize from their point of view. After having quantified those criteria under the form of a cost function and a design parameters space, we wish to implement state-of-the-art numerical methods (finite elements, full waveform inversion, neuronal networks, diverse optimization techniques...) that are versatile (in terms of models, formulations, couplings...) in order to solve the optimization problem. More precisely, we wish to take advantage of the fact that sound waves in musical instruments satisfy the laws of acoustics in pipes (PDE), which gives us access to the full waveform inversion technique, usable in harmonic or temporal regime. The methods that we want to use are attractive because the weekly depend on the chosen criterion, and they are easily adaptable to various physical situations (multimodal decomposition in the pipe, coupling with the embouchure, ...), which can therefore be modified a posteriori. The goal is to proceed iteratively between instrument making and optimal design (the virtual workshop) in order to get close to tone quality related and playability criteria.

#### Energy based model and simulation in the time domain of linear acoustic waves in a radiating pipe

Participants : Juliette Chabassier, Robin Tournemenne.

We model in the time domain linear acoustic waves in a radiating pipe without daming. The acoustic equations system in formulated in flow and pressure, which leads to a first order space time equations system. The radiation condition is also written as a first order in time equation, and is parametrized by two real coefficients. Moreover, an auxiliary variable is introduced at the radiating boundary. The choice of this variable is adapted to the considered source type in order to ensure the model stability by energy techniques, under some conditions on the radiating condition. We then propose a stable space time explicit discretization, which ensures the dissipation of a discrete energy. The novelty of the discretization lies, on the one hand, in the variational nature of the space approximation ( which leads to artibrary order finite elements with no required matrix inversion), and on the other hand, on the definition of the auxiliary variable for any acoustic source type (which leads to the decay of a well defined energy). Finaly, we quantify the frequential domain of validity of the used radiation condition by comparison with theoretical and experimental models of the litterature. This is a collaboration with Morgane Bergot (Université Claude Bernard, Lyon 1).

#### Computation of the entry impedance of a dissipative radiating pipe

Participants : Juliette Chabassier, Robin Tournemenne.

Modeling the entry impedance of wind instruments pipes is essential for sound synthesis or instrument qualification. We study this modeling with the finite elements method in one dimension (FEM1D) and with the more classically used transfer matrix method (TMM). The TMM gives an analytical formula of the entry impedance depending on the bore (intern geometry of the instrument) defined as a concatenation of simple elements (cylinders, cones, etc). The FEM1D gives the entry impedance for any instrument geometry. The main goals of this work are to assess the viability of the FEM1D and to study the approximations necessary for the TMM in dissipative pipes. First, lossless Weber's equation in one dimension is studied with arbitrary radiation conditions. In this context and for cylinders or cones, the TMM is exact. We verify that the error made with FEM1D for fine enough elements is as small as desired. When we consider viscothermal losses, the TMM does not solve the classical Kirchhoff model because two terms are supposed constant. In order to overcome this model approximation, simple elements, on which are based the TMM, are decomposed into much smaller elements. The FEM1D does not necessitate any model approximation, and it is possible to show that it solves the dissipative equation with any arbitrarily small error. With this in hand, we can quantify the TMM model approximation error.

#### Hybrid discontinuous finite element approximation for the elasto-acoustics.

Participants : Hélène Barucq, Julien Diaz, Elvira Shishenina.

Discontinuous Finite Element Methods (DG FEM) have proven their numerical accuracy and flexibility. However, numerically speaking, the high number of degrees of freedom required for computation makes them more expensive, compared to the standard techniques with continuous approximation. Among the different variational approaches to solve boundary value problems there exists a distinct family of methods, based on the use of trial functions in the form of exact solutions of the governing equations. The idea was first proposed by Trefftz in 1926, and since then it has been largely developed and generalized. By its definition, Trefftz-DG methods reduce numerical cost, since the vari- ational formulation contains the surface integrals only. Thus, it makes possible exploration of the meshes with different geometry, in order to create more realistic application. Trefftz-type approaches have been widely used for time-harmonic problems, while their implementation is still limited in time domain. The particularity of Trefftz-DG methods applied to the time-dependent formulations consists in the use of space-time meshes. Even though it creates another computational difficulty, due to a dense form of the matrix, which represents the global linear system, the inversion of the full ”space-time” matrix can be reduced to the inversion of one block-diagonal matrix, which corresponds to the interactions in time. In the present work, we develop a theory for solving the coupled elasto-acoustic wave propagation system. We study well-posedness of the problem, based on the error estimates in mesh-dependent norms. We consider a space-time polynomial basis for numerical discretization. The obtained numerical results are validated with analytical solutions. Regarding the advantages of the method, following properties have been proven by the numerical tests: high flexibility in the choice of basis functions, better order of convergence, low dispersion. These results have been obtained in collaboration with Henri Calandra (TOTAL) and have been published in a research report [43]. A paper has been submitted and a second one is being prepared.

#### Construction of stabilized high-order hybrid Galerkin schemes.

Participants : Hélène Barucq, Aurélien Citrain, Julien Diaz.

We have compared the perfomance of Discontinuous Galerkin Methods and Spectral Element Methods on academic benchmark and on realistic geophysical model in two dimensions. We have shown that, for a given accuracy, SEm on quadrilateral meshes could be 10 times faster than DGm, which justifies our strategy to consider SEm wherever it is possible to use quadrilateral/hexahedral cells. These first results have been presented in Matthias conference. Then, we have considered the SEM/DG coupling proposed for electromagnetics in [78] and we have implemented it in our acoustics code. We are now analyzing the performance of this strategy and we are extending it to deal with elastodynamic and elasto-acoustic coupling. The following steps will be the extension of the analysis to 3D dimensional problems and the application to realistic test case. The main bottleneck is obviously to the definition of an efficient strategy to couple tetrahedra and hexahedra. Indeed, if in the 2D case, the edges of both triangles and quadrilaterals are all segments the faces of tetrahedra are triangle while the faces of hexahedra are quadrilaterals. Hence, in 2D it sufficed to define integration on segment, while in 3D it will be necessary to consider integration of various polygon resulting of the intersection of triangle and quadrilaterals. Once this strategy is defined and implemented, we expect to be able to reduce the computational of the plateform that we develop jointly with Total by a factor between 5 and 10. These results have been obtained in collaboration with Henri Calandra (TOTAL) and Christian Gout (INSA Rouen).

#### Modeling of dissipative porous media.

Participants : Juliette Chabassier, Julien Diaz, Fatima Jabiri.

In this work we have considered the modeling of 1D acoustic wave propagation coupled with visco-thermical losses that occur in porous media. We have proposed a family of dissipative models from which we have been able to obtain a quasi-constant quality factor (which is an indicator of the dissipation as a function of the frequency). We have derived stability conditions on the parameters of the model thanks to an energy analysis and we have rewritten the problem of designing a quasi-constant quality factor as a contrained least-square optimization problem. The parameters to optimize are the parameters of the family of dissipative models and the constraints are the stability of the final model. We are now considering the extension of the family to more general formulations and to heterogeneous media, before tackling multidimensional problems. These results have been obtained during the Master internship of Fatima Jabiri, in collaboration with Sébastien Imperiale (Inria Project-Team M3DISIM)

#### Asymptotic models for the electric potential across a highly conductive casing. Application to the field of resitivity measurements.

Participants : Hélène Barucq, Aralar Erdozain, Victor Péron.

A configuration that involves a steel-cased borehole is analyzed, where the casing that covers the borehole is considered as a highly conductive thin layer. Asymptotic techniques are presented as the suitable tool for deriving reduced problems capable of dealing with the numerical issues caused by the casing when applying the traditional numerical methods. The derivation of several reduced models is detailed by employing two different approaches, each of them leading to different classes of models. The stability and convergence of these models is studied and uniform estimates are proved. The theoretical orders of convergence are supported by numerical results obtained with the finite element method. We develop an application to the field of resistivity measurements. The second derivative of the potential which solves a reduced model has been employed to recover the resistivity of rock formations. These results are in accordance with an experiment of Kaufmann for the reference solution and have been obtained in collaboration with David Pardo (UPV/EHU).

#### Semi-Analytical Solutions for the Electric Potential across a Highly Conductive Casing.

Participants : Hélène Barucq, Aralar Erdozain, Victor Péron.

A transmission problem for the electric potential is considered, where one part of the domain is a high-conductive casing. Semi-analytical solutions are derived for several asymptotic models. These asymptotic models are designed to replace the casing by appropriate impedance conditions in order to avoid numerical instabilities. A decomposition in Fourier series of the solution to these asymptotic models is characterized. As an application we reproduced successfully the experiment of Kaufmann, using his same parameters, but computing with a fourth order asymptotic model. This experiment allows to recover the resitivity of rock formations employing a second derivative of the potential along the vertical direction. These results have been obtained in collaboration with Ignacio Muga (Pontificia Universidad Catolica of Valparaiso).

#### Asymptotic modeling of the electromagnetic wave scattering problem by a small sphere perfectly conducting.

Participants : Justine Labat, Victor Péron, Sébastien Tordeux.

In the context of non-destructive testing in medical imaging or civil engineering, the detection of small heterogeneities can be a difficult task in three dimensional domains. The complexity for solving numerically the direct problem both in terms of computation time and memory cost is due to the small size of obstacles in comparison with the incident wavelength and the large size of the domain of interest. Then the fine mesh size makes unsuitable or too expensive the use of classical numerical methods type continuous and discontinuous finite element methods or boundary element methods.
The use of reduced models allows to get an approximation of the exact solution at a certain accuracy with a lower cost. We develop a Matched Asymptotic Expansions method to solve a time-harmonic electromagnetic scattering problem by a small sphere. This method allows to replace the scatterer by an equivalent asymptotic point source. In practice, it consists in defining an approximate solution using multi-scale expansions over far and near fields, related in a matching area. When the scatterer is a sphere, we make explicit the asymptotic expansions until the second order of approximation, relatively to the sphere radius. Numerical results make evident the convergence rate with respect to the sphere radius. Reference solutions are analytical solutions computed thanks to Montjoie Code. This work has been presented in the *Caleta Numerica* seminar, Pontificia Universidad Catolica of Valparaiso, Chili [48].

#### Asymptotic Models and Impedance Conditions for Highly Conductive Sheets in the Time-Harmonic Eddy Current Model.

Participant : Victor Péron.

This work is concerned with the time-harmonic eddy current problem for a medium with a highly conductive thin sheet. We present asymptotic models and impedance conditions up to the second order of approximation for the electromagnetic field. The conditions are derived asymptotically for vanishing sheet thickness $\u03f5$ where the skin depth is scaled like $\u03f5$. The first order condition is the perfect electric conductor boundary condition. The second order condition turns out to be a Poincaré-Steklov map between tangential components of the magnetic field and the electric field [49]. Numerical experiments have been performed to assess the accuracy of the second order model. Complementary simulations will be conducted to study the robustness with respect to the sheet conductivity and the convergence of the modelling error. These results have been obtained in collaboration with Mohammad Issa and Ronan Perrussel (LAPLACE, CNRS/IMPT/UPS, Univ. de Toulouse) and this work has been presented in the international conference ACOMEN 2017 [33].

#### Numerical robustness of single-layer method with Fourier basis for multiple obstacle acoustic scattering in homogeneous media

Participants : Hélène Barucq, Juliette Chabassier, Ha Pham, Sébastien Tordeux.

We investigate efficient methods to simulate the multiple scattering of obstacles in homogeneous media. With a large number of small obstacles on a large domain, optimized pieces of software based on spatial discretization such as Finite Element Method (FEM) or Finite Difference lose their robustness. As an alternative, we work with an integral equation method, which uses single-layer potentials and truncation of Fourier series to describe the approximate scattered field. In the theoretical part of the paper, we describe in detail the linear systems generated by the method for impenetrable obstacles, accompanied by a well-posedness study. For the numerical performance study, we limit ourselves to the case of circular obstacles. We first compare and validate our codes with the highly optimized FEM-based software Montjoie. Secondly, we investigate the efficiency of different solver types (direct and iterative of type GMRES) in solving the dense linear system generated by the method. We observe the robustness of direct solvers over iterative ones for closely-spaced obstacles, and that of GMRES with Lower–Upper Symmetric Gauss–Seidel and Symmetric Gauss–Seidel preconditioners for far-apart obstacles.

This work has been published in the journal Wave Motion, [15] and is also connected to the following conference presentations, [31], [41].

#### A study of the Numerical Dispersion for the Continuous Galerkin discretization of the one-dimensional Helmholtz equation

Participants : Hélène Barucq, Ha Pham, Sébastien Tordeux.

*This work is a collaboration with Henri Calandra (TOTAL).*

Although true solutions of Helmholtz equation are non-dispersive, their discretizations suffer from a phenomenon called numerical dispersion. While the true phase velocity is constant, the numerical one changes with the discretization scheme, order and mesh size. In our work, we study the dispersion associated with classical finite element. For arbitrary order of discretization, without using an Ansatz, we construct the numerical solution on the whole $\mathbb{R}$, and obtain an asymptotic expansion for the phase difference between the exact wavenumber and the numerical one. We follow an approach analogous to that employed in the construction of true solutions at positive wavenumbers, which involves Z-transform, contour deformation and limiting absorption principle. This perspective allows us to identify the numerical wavenumber with the angle of analytic poles. Such an identification is useful since the latter (analytic poles) can be numerically evaluated by an algorithm, which then yields the value of numerical wavenumber.

This work is detailed in the research report [44].