Section: New Results


Implementation of a non-reflecting boundary condition on ellipsoidal boundary

Participants : Hélène Barucq, Anne-Gaëlle Saint-Guirons, Sébastien Tordeux.

The modeling of wave propagation problems using finite element methods usually requires the truncation of the computational domain around the scatterer of interest. Absorbing boundary condition are classically considered in order to avoid spurious reflections. This year we have proposed and tested a formulation which allows to take into account with no extra-cost a quasi-exact radiation condition based on a non local Dirichlet to Neumann operator.

Modeling of small heterogeneities in the context of the time domain wave equation

Participants : Vanessa Mattesi, Sébastien Tordeux.

We have proposed an approximate model to take into account small heterogeneities for the three dimensional time dependent wave equation. One of the most important result of this work is the generalization of the multipole theory (classically written for the Helmholtz equation) to the wave equation. This work has been presented at Waves 2013 and at JSA 2013 [58] , [39]

A new modified equation approach for solving the wave equation

Participants : Hélène Barucq, Henri Calandra, Julien Diaz, Florent Ventimiglia.

In order to obtain high-order time-schemes, we are considering an alternative approach to the ADER schemes and to the modified equation technique described in section  3.2 . The two first steps of the construction of the schemes are similar to the previous schemes : we apply a Taylor expansion in time to the solution of the wave equation and we replace the high-order derivatives with respect to the time by high order space operators, using the wave equation. The difference is that we do not use auxiliary variables and we choose to discretize directly the high-order operators in space. These operators can not be discretized by classical finite elements. For the discretization of the biharmonic operator in an homogeneous acoustic medium, both C1 finite elements, such as the Hermite ones, and Discontinuous Galerkin Finite Elements (DGFE) can be used, while in a discontinuous medium, or for higher-order operators, DGFE should be preferred  [80] . We have applied this method to the second order wave equation [15] and the numerical results showed that this technique induced less computational burden than the modified equation scheme or the ADER scheme.

In the framework of the PhD thesis of Florent Ventimiglia, we have extended the new method involving p-harmonic operator to the first order formulation of the acoustic wave equation, which is the formulation discretized in the DIVA platform of TOTAL. In this case, the high order operators in space are not are not powers of the Laplace operator but powers of the gradient. Hence, we also had to adapt the space discretization, and we have extended the DG formulation with centered fluxes proposed in  [87] to higher order operators. A numerical analysis of performance in 2D indicates that, for a given accuracy, this method requires less computational costs and less storage than the High-Order ADER Scheme. These results have been presented to the SMAI conference [49] , Waves 2013 [61] , Numerico IV [61] and HF2013 [34] . A paper has been accepted in ESAIM Proceedings [49] .

Constructing and using Absorbing Boundary Conditions

Higher Order On-Surface Radiation Conditions for elastic scatterers

Participants : Hélène Barucq, Chokri Bekkey, Juliette Chabassier, Julien Diaz.

The numerical simulation of wave propagation is generally performed by truncating the propagation medium and the team works on new Absorbing Boundary Conditions (ABCs), trying to improve the performance of existing conditions. As we explained at Section  3.2 , item Boundary conditions, we are developing ABCs for curved boundaries, based on the full factorization of the wave equation. These ABCs should take propagating, grazing and evanescent waves into account. In [17] , we have considered the issue of constructing high-order ABCs taking into account both propagating and evanescent waves for the Helmholtz equation. In case of the simulation of acoustic waves diffracted by a solid immersed in a fluid, we investigate the performance of the new ABCs when used as On-Surface Radiation Conditions. The ABCs are set directly on the boundary of the solid. The unbounded problem is then replaced by a problem involving an acoustic pressure computed on the surface of the solid only. Preliminary results have been obtained by considering the toy problem where the scatterer is a disk. Analytic solutions are then available and we show that that taking into account evanescent waves in the ABC could improve the accuracy of classical ABCs by two orders of magnitude at mid-frequency range, for ka between 1 and 100, k being the frequency and a the typical size of the diffracting obstacle. These results have been presented to the Waves 2013 conference [47] .

Radiation boundary condition at high frequency

Participants : Hélène Barucq, Elodie Estécahandy, Juliette Chabassier, Julien Diaz.

Regarding the solution of the Helmholtz equation at high frequency with finite element methods, it is current to refine the mesh in order to limit the effect of numerical pollution. It is then interesting to dispose of radiation boundary conditions which do not require to set the artificial boundary far from the scatterer. In this work, we have investigated the possibility of bringing closer the artificial boundary when standard conditions are enhanced by the modeling of grazing waves. The preliminary results we obtained show that this new ABC outperforms classical ABCs at high-frequency, for ka>50, and that it is highly accurate, even when the boundary is very close to the scatterer. These results have been presented to the Waves 2013 conference [45] . Now, the next step is to consider an ABC taking the three types of waves into account.

Absorbing Boundary Conditions for Tilted Transverse Isotropic Elastic Media

Participants : Hélène Barucq, Lionel Boillot, Henri Calandra, Julien Diaz.

The simulation of wave propagation in geophysical media is often performed in domains which are huge compared to the wavelengths of the problem. It is then necessary to reduce the computational domain to a box. When considering acoustic or elastic isotropic media, this can be done by applying an Absorbing Boundary Condition (ABC) or by adding a Perfectly Matched Layer (PML). However, a realistic representation of the Earth subsurface must include anisotropy and, in particular, the so-called Tilted Transverse Isotropy. Perfectly Matched Layers are known to be unstable for this kind of media and, to the best of our knowledge, no ABC have been proposed yet. We have thus proposed a low-order ABC for TTI media. The construction is based on comparing and then connecting the slowness curves for isotropic and elliptic TTI waves. Numerical experiments illustrate the performance of the new ABC. They are performed by integrating the ABC in a DG formulation of Elastodynamics. When applied in a TTI medium, the new ABC performs well with the same level of accuracy than the standard isotropic ABC set in an isotropic medium. The condition demonstrates also a good robustness when applied for large times of simulation. These results have been presented to the Smai and to the Waves conferences [67] , [68] and a paper has been submitted.

Modeling of the damping factor of Multiperforated plates

Participants : Estelle Piot, Vincent Popie, Sébastien Tordeux.

Multiperforated plates are classically used as a damping material. Melling has proposed in [94] a model to estimate the energy dissipated by these devises. However its well-known result should be corrected by a factor two to fit with experimental data. We have proposed a correct way to compute the energy dissipated by the multiperforated plates. This work has been presented at the Fifth International Scientific Conference and Young Scientists School "Theory and Computational Methods for Inverse and Ill-posed Problems.

Performance Assessment of IPDG for the solution of an elasto-acoustic scattering problem

Participants : Hélène Barucq, Rabia Djellouli, Élodie Estécahandy.

We present a solution methodology for the direct elasto-acoustic scattering problem that falls in the category of Discontinuous Galerkin methods. The method distinguishes itself from the existing methods by combining high-order Discontinuous Galerkin approximations, local stabilizations for the coupled problem and the use of curved element edges on the boundaries. We present some numerical results that illustrate the salient features and highlight the performance of the proposed solution methodology on the resonance phenomenon existing in the elastic scatterer for simple geometries such as circles. Moreover, the designed method ensures a convergence order with a gain of two order of magnitude compared to polygonal boundaries, and a potential to address both mid- and high-frequency regimes. These results have been accepted for publication in International Journal for Numerical Methods in Engineering [19] .

On the influence of curvature on transmission conditions

Participants : Hélène Barucq, Martin Gander, Yingxiang Xu.

Domain decomposition methods are both highly successful parallel solvers and also important modeling tools, since problems in subdomains can be treated by adapted methods to the physics in each subdomain. Subdomain boundaries are therefore rarely straight lines. The focus of this paper is to study the influence of curvature on transmission conditions used in optimized Schwarz methods. For straight interfaces and simple geometries, optimized interface conditions are typically determined using Fourier analysis.Asymptotically, these optimized conditions are still valid for curved interfaces.Since however the curvature is the most important information for a smooth curve, we want to study in this paper if and how the interface curvature influences the constants in the optimized parameters.

We consider the model problem

( Δ - η ) u = f , on Ω = 2 , η > 0 , (2)

and we require the solution to decay at infinity.We decompose Ω into two overlapping subdomains Ω1=(-,a(y))× and Ω2=(b(y),)×, where Γ1 given by a(y) and Γ2 given by b(y) are smooth curves satisfying a(y)b(y). A general parallel Schwarz algorithm is then given by

( Δ - η ) u i n = f in Ω i , i ( u i n ) = i ( u j n - 1 ) on Γ i , 1 i j 2 , (3)

where i,i=1,2, are transmission conditions to be chosen. If i, i=1,2 are chosen as ni+DtNi, with DtNi the Dirichlet to Neumann operators, the iterates will converge in two steps. These operators are however non-local, and thus difficult to use in practice. Therefore, local approximations are used in optimized Schwarz methods. We presented two different approaches to take the curvature of interfaces into account in the transmission conditions of optimized Schwarz methods: micro-local analysis, and analysis using a circular model problem. In both cases, we obtained curvature dependent transmission conditions. A preliminary comparison shows that the transmission conditions based on optimization perform better on the model problem, and that it could be important to take the curvature into account in transmission conditions. In our opinion it is however essential to do a more thorough theoretical and numerical study on more general geometry, where micro-local analysis is still applicable, before we can definitely draw conclusions. This work has been published in the peer-reviewed proceedings of the conference Decomposition Methods in Science and Engineering XXI [51] .

Operator Based Upscaling for Discontinuous Galerkin Methods

Participants : Hélène Barucq, Théophile Chaumont, Christian Gout.

Realistic numerical simulations of seismic wave propagation are difficult to handle because they must be performed in strongly heterogeneous media. Two different scales must then be taken into account. Indeed, the medium heterogeneities are currently very small compared to the characteristic dimensions of the propagation medium. To get accurate numerical solutions, engineers are then forced to use meshes that match the finest scale representing the heterogeneities. Meshing the whole domain with a fine grid leads then to huge linear systems and the computational cost of the numerical method is then too high to consider 3D realistic simulations. To dispose of a numerical method allowing to represent the heterogeneity of the medium accurately while computing on a coarse grid is thus relevant. This is the challenge of multiscale approaches like homogenization or upscaling. The ultimate objective of this work is to develop a software package for the Helmholtz equation set in heterogeneous domains. We focus on the operator-based upscaling method. Operator-based upscaling methods were first developed for elliptic flow problems (see  [81] ) and then extended to hyperbolic problems (see  [90] , [104] , [103] ). Operator-based upscaling method consists in splitting the solution into a coarse and a fine part. The coarse part is defined on a coarse mesh while the fine part is computed on a fine mesh. In order to speed up calculations, artificial boundary conditions (ABC) are imposed. By enforcing suitable ABCs on the boundary of every cells of the coarse mesh, calculations on the fine grid can be carried out locally. The coarse part is next computed globally on the coarse mesh. Operator-based upscaling methods were so far developed in joint with standard finite element discretization strategy. In this work, we investigate the idea of combining an operator based upscaling method with discontinuous Galerkin finite element methods(DGFEM). During this year, we have performed the mathematical analysis of the Helmholtz equation set in a domain represented by a discontinuous velocity. The analysis has been achieved both for the continuous and the discretized problem which is based on a quadrature scheme which allows to take the discontinuities of the velocity into account. We get new stability results for stratified-like domains and numerical experiments show that in case of industrial benchmarks, the quadrature scheme leads to lower computational costs with a very good level of accuracy. A paper is in preparation and the results will be presented at the 2014 ECCOMAS conference in Barcelona.

Efficient solution methodology based on a local wave tracking strategy for high-frequency Helmholtz problems.

Participants : Mohamed Amara, Sharang Chaudhry, Julien Diaz, Rabia Djellouli, Steven Fiedler.

We have proposed a procedure for selecting basis function orientation to improve the efficiency of solution methodologies that employ local plane-wave approximations. The proposed adaptive approach consists of a local wave tracking strategy. Each plane-wave basis set, within considered elements of the mesh partition, is individually or collectively rotated to best align one function of the set with the local propagation direction of the field. Systematic determination of the direction of the field inside the computational domain is formulated as a minimization problem. As the resultant system is nonlinear with respect to the directions of propagation, the Newton method is employed with exact characterization of the Jacobian and Hessian. To illustrate the salient features and evaluate the performance of the proposed wave tracking approach, we present error estimates as well as numerical results obtained by incorporating the procedure into a prototypical plane-wave based approach, the least-squares method (LSM) developed by Monk and Wang  [95] . The numerical results obtained for the case of a two-dimensional rigid scattering problem indicate that (a) convergence was achievable to a prescribed level of accuracy, even upon initial application of the tracking wave strategy outside the pre-asymptotic convergence region, and (b) the proposed approach reduced the size of the resulting system by up to two orders of magnitude, depending on the frequency range, with respect to the size of the standard LSM system. These results has been presented to the Waves 2013 [43] and in a Research Report [71] . A paper has been submitted.

Mesh Free Frontier-based Formulation (MF3) for High Frequency Helmholtz Problems.

Participants : Mohamed Amara, Julien Diaz, Rabia Djellouli.

We have proposed a novel approach for solving efficiently Helmholtz problems. The proposed solution method employs a boundary-type formulation without however involving Green functions and/or incurring singular integrals. In addition, this approach does not necessitate the use of a mesh. For these reasons, the method is named Mesh Free Frontier-based Formulation (MF3). Furthermore, the sought-after field is locally approximated using a set of basis of functions that consists of a Bessel kind function computed at a prescribed finite set of points. The number of the functions determines mainly the size of the resulting system, the complexity, as well as the computational cost of the proposed method. Preliminary numerical results obtained in the case of 2D-Helmholtz problems in the high-frequency regime are presented to illustrate the computational efficiency of MF3 (the method delivers results with high accuracy level, about 108 on theL2 relative error, while requiring the solution of small linear systems). In addition, these results tend to suggest that MF3 is pollution free. These results has been presented to two conferences [32] , [44] . A paper is in preparation.

Energy based simulation of a Timoshenko beam in non-forced rotation. Application to the flexible piano hammer shank.

Participants : Juliette Chabassier, Marc Duruflé.

A nonlinear model for a vibrating Timoshenko beam in non-forced unknown rotation is derived from the virtual work principle applied to a system of beam with mass at the end. The system represents a piano hammer shank coupled to a hammer head. An energy-based numerical scheme is then provided, obtained by non classical approaches. A major difficulty for time discretization comes from the nonlinear behavior of the kinetic energy of the system. Numerical illustrations are obtained by coupling this new numerical scheme to a global energy-preserving numerical solution for the whole piano. These numerical results show that the pianistic touch clearly influences the spectrum of the piano sound of equally loud isolated notes. These differences do not come from a possible shock excitation on the structure, nor from a changing impact point, nor a “longitudinal rubbing motion” on the string, since neither of these features are modeled in our study.

This work has been submitted for publication in Journal of Sound and Vibration [79] . discretization

Simulating the propagation of ultra short laser pulses in a dispersive non linear medium.

Participants : Juliette Chabassier, Marc Duruflé, Nayla Herran.

This collaboration with CEA-CESTA aimed at evaluating the limits of validity of the MIRÓ software provided by CEA. The MIRÓ software implements the propagation of laser pulses with a non-linear Schrodinger-like equation obtained from Maxwell's equations in non-linear dispersive medium, assuming that the pulse spectrum is narrow and the non-linearity small enough. When considering intense ultra-short pulses, the spectrum bandwidth and the amplitude of the pulse may violate the constitutive assumptions of the model used by MIRÓ. This collaboration began with the internship of Nayla Herran (march 2013- august 2013). During this internship, different alternative models have been explored, and some of them were able to provide a solution much more accurate than MIRÓ's models. The comparisons between alternative models and MIRÓ's model have been performed in 1-D on small cases, for which the reference solution could be obtained from a direct numerical simulation of non-linear Maxwell's equations. An efficient and accurate solver for non-linear Maxwell's equations has been implemented in Montjoie software. This solver uses high order finite element in space, and high order Runge-Kutta-Nystrom scheme in time. The space grid moves with respect to the group velocity such that the computational domain stays relatively small and centered around the pulse. Thanks to this solver, we have been able to validate the different models and compare them. These 1-D promising results encourage us to continue this collaboration in order to obtain efficient and accurate numerical methods in 3-D. Our desire is also to be able to perform realistic cases involving physically relevant phenomena.

Asymptotic Modeling for Elasto-Acoustics

Participants : Julien Diaz, Victor Péron.

In the papers [31] , [74] , we derive equivalent conditions and asymptotic models for the diffraction problem of elastic and acoustic waves in a solid medium surrounded by a thin layer of fluid medium. Due to the thinness of the layer with respect to the wavelength, this problem is well suited for the notion of equivalent conditions and the effect of the fluid medium on the solid is as a first approximation local. We derive and validate equivalent conditions up to the fourth order for the elastic displacement. These conditions approximate the acoustic waves which propagate in the fluid region. This approach leads to solve only elastic equations. The construction of equivalent conditions is based on a multiscale expansion in power series of the thickness of the layer for the solution of the transmission problem.

Questions regarding the implementation of the conditions have been addressed carefully and the boundary conditions have been integrated without changing the structure of the code Hou10ni.

This work has been presented in two international conferences and Workshops : Workshop HPC-GA and WAVES'2013 [55] .

A paper with numerical results for the elasto-acoustic problem with several configurations (a thin layer of variable thickness; coupling with an exterior "acoustic" medium) is in preparation.

Thin layer models for electromagnetics

Participants : Marc Duruflé, Victor Péron, Clair Poignard.

We have considered the transmission of electromagnetic waves through a thin layer. This thin layer can be replaced by transmission conditions. Media with thin inclusions appear in many domains: geophysical applications, microwave imaging, biomedical applications, cell phone radiations, radar applications, non-destructive testing. Our application concerns also media for which the conductivity drops inside the layer, such that the low-frequency regime has a different behavior from the mid-frequency regime. Different models are compared for these two regimes, drawbacks and disadvantages of each model are detailed. This work has been accepted for publication [27] .

Corner Asymptotics of the Magnetic Potential in the Eddy-Current Model

Participants : Monique Dauge, Patrick Dular, Laurent Krähenbühl, Victor Péron, Ronan Perrussel, Clair Poignard.

The following results rely on a problematic developed in section  3.2 , item Asymptotic modeling.

In the paper [25] , we describe the magnetic potential in the vicinity of a corner of a conducting body embedded in a dielectric medium in a bidimensional setting. We make explicit the corner asymptotic expansion for this potential as the distance to the corner goes to zero. This expansion involves singular functions and singular coefficients. We introduce a method for the calculation of the singular functions near the corner and we provide two methods to compute the singular coefficients: the method of moments and the method of quasi-dual singular functions. Estimates for the convergence of both approximate methods are proven. We eventually illustrate the theoretical results with finite element computations. The specific non-standard feature of this problem lies in the structure of its singular functions: They have the form of series whose first terms are harmonic polynomials and further terms are genuine non-smooth functions generated by the piecewise constant zeroth order term of the operator. This work has been presented in the international conference JSA 2013 [38] .

Finite Element Subproblem Method

Participants : Patrick Dular, Christophe Geuzaine, Laurent Krähenbühl, Victor Péron, Ronan Perrussel.

In the paper [26] , we develop a finite element subproblem method to correct the inaccuracies proper to perfect conductor and impedance boundary condition models, in particular near edges and corners of conductors, for a large range of conductivities and frequencies. Successive local corrections, supported by fine local meshes, can be obtained from each model to a more accurate one, allowing efficient extensions of their domains of validity. This work has been presented in the international conference Compumag 2013 [57] .

We develop also a finite element subproblem method for progressive eddy current modeling. The modeling of eddy currents in conductors is split into a sequence of progressive finite element subproblems. The source fields generated by the inductors alone are calculated at first via either the Biot-Savart law or finite elements. The associated reaction fields for each added conductive region, and in return for the source regions themselves when massive, are then calculated with finite element models, possibly with initial perfect conductor and/or impedance boundary conditions to be further corrected. The resulting subproblem method allows efficient solving of parameterized analyses thanks to a proper mesh for each subproblem and the reuse of previous solutions to be locally corrected. This work has been presented in the international conference ISEF'2013 [56] .