Vahan Baronian

Anne-Sophie Bonnet

Ben Dhia

Patrick Joly

En collaboration avec Xavier Claeys (LJLL, Univ. Paris 6), nous avons donc revu le processus d'homogénéisation afin de proposer des modèles plus riches au niveau des interfaces. Ce travail a été réalisé, dans un premier temps, dans le cadre de la thèse de Valentin Vinoles, pour la propagation des ondes acoustiques en régime harmonique entre un demi-espace périodique classique (non métamatériau : les coefficients effectifs sont positifs) et un demi-espace homogène. En combinant les méthodes d'homogénéisation double échelle et celle des développements

Introduction

J'ai fait ma thèse de 2005 à 2009 avec Patrick Joly sur l'analyse de problèmes de diffraction dans des milieux périodiques. Arrivée à l'Ensta juste après, je me suis intéressée aux sujets de recherche de mes nouveaux collègues (théorie spectrale, guide d'ondes élastiques, problèmes inverses) tout en essayant de leur apporter mon expertise sur les milieux périodiques et sur les techniques et méthodes sous-jacentes. J'ai par exemple travaillé avec Laurent Bourgeois (POEMS) sur les méthodes d'échantillonnage dans des guides d'ondes périodiques et j'ai commencé ma collaboration avec Anne-Sophie Bonnet-Ben Dhia (POEMS) sur les problèmes de diffraction des plaques élastiques anisotropes, en collaboration avec le CEA. Il me semble que cette expertise acquise sur les milieux périodiques a intéressé aussi à l'exterieur. Si je ne devais citer que les travaux qui ont abouti à une publication, je parlerais de ma collaboration avec K. Ramdani (Inria Nancy), C. Besse (Univ. Toulouse) et I. Lacroix-Violet (Univ. Lille) et avec K. Schmidt et D. Klindworth (TU Berlin). Parallèlement, je me suis beaucoup interessée aussi aux travaux de l'équipe sur les métamatériaux, (sujet qui excite beaucoup les physiciens aujourd'hui du fait de leurs applications multiples) et en particulier les difficultés théoriques et numériques qui interviennent quand on utilise des modèles effectifs ou homogénéisés pour étudier la diffraction des ondes par ces milieux. Nous avons donc monté un projet ANR (le projet METAMATH 2011(le projet METAMATH -2016)), avec l'Univ. de Toulon, l'équipe DeFI de l'Inria et le LJLL, dont j'étais responsable et où il était question entre autres de revisiter les techniques d'homogénéisation pour proposer des modèles effectifs plus riches de ces métamatériaux, notamment au niveau des interfaces. J'ai commencé à travailler sur ces questions avec Xavier Claeys (LJLL) et beaucoup reste à faire.

Après cette courte présentation de mon parcours, je vais commencer ce document par présenter le contexte et certaines applications qui motivent l'étude de phénomènes de propagation d'ondes dans des milieux périodiques. Je résume ensuite mes contributions et le contenu de ce manuscrit.

Un milieu est dit "périodique" quand sa géométrie et/ou ses caractéristiques physiques (typiquement les coefficients élastiques quand il s'agit d' une structure mécanique ou la permittivité diélectrique quand il s'agit d'une structure optique) sont des fonctions périodiques d'une ou plusieurs variables d'espace. Les milieux périodiques apparaissent dans un grand nombre d'applications. Donnons deux exemples, un en mécanique et l'autre en optique.

Tout d'abord, le matériau à fibre ou composite est un assemblage souvent supposé périodique 6 Introduction d'au moins deux matériaux. (voir par exemple [Christensen, 1979[START_REF] Aboudi | Mechanics of composite materials[END_REF], Kaw, 2005] pour leur modélisation). Ces mélanges permettent d'améliorer la qualité de la matière pour des utilisations spécifiques (légèreté, résistance particulières à certains efforts, ...), ce qui explique l'utilisation croissante de ces matériaux dans différents secteurs, notamment l'industrie aéronautique et aérospatiale. On trouve souvent des défauts de structure, des ruptures géométriques ou des fissures, liés à une fabrication défectueuse de la pièce ou à la fatigue durant son utilisation. Pour détecter ces défauts et en déterminer certaines propriétés, les techniques de contrôle non destructif (CND), parmi lesquelles les méthodes ultrasonores, sont très utilisées car elles permettent de vérifier l'intégrité des structures sans les altérer. Mais il est souvent difficile d'interpréter les signaux du fait de l'anisotropie des milieux, de la présence de raidisseurs. Des outils de simulation peuvent aider à maitriser cette complexité et peuvent donc être utilisés pour aider au positionnement des capteurs, à l'amélioration des techniques de contrôle et au diagnostic des défauts.

Les milieux périodiques présentent également des propriétés très intéressantes en optique, en particulier en micro-technologie et nano-technologie. En effet, un des récents sujets porteurs en optique concerne les matériaux à bandes interdites de photons connus sous le nom de cristaux photoniques. [Kuchment, 2001, Joannopoulos et al., 1995, Johnson and Joannopoulos, 2002, Sakoda, 2001]. Ces derniers sont des structures périodiques composées de matériaux diélectriques qui présentent souvent un fort contraste d'indice. En particulier, de telles structures permettent de sélectionner les bandes de fréquences pour lesquelles les ondes peuvent ou non se propager dans le milieu en question. Plusieurs travaux mathématiques et numériques indiquent qu'un choix adéquat de la structure du cristal et des matériaux diélectriques le composant permettraient de créer des bandes interdites particulières et donc d'un point de vue pratique, bannir du cristal certaines ondes électromagnétiques [Kuchment, 2004, Kuchment, 2001]. De ce fait, ces milieux pourraient être utilisés dans la réalisation de filtres, d'antennes et de différents composants utilisés en télécommunications. En outre, des défauts peuvent être introduits volontairement dans le milieu pour en changer les propriétés. Ainsi, en optique, dans le but de réaliser des lasers, des filtres, des fibres ou des guides d'ondes, des défauts linéiques sont introduits pour permettre la propagation de modes guidés. Là encore, des travaux mathématiques et numériques indiquent comment optimiser la structure pour créer un mode guidé.

Positionnement et présentation de mes travaux de recherche

Mes travaux portent sur l'analyse mathématique et numérique des problèmes de diffraction (principalement en régime harmonique ou régime periodique établi) dans des milieux périodiques présentant des perturbations ou non (voir les chapitres 1 et 2). Je me suis également intéressée à des problèmes spectraux (voir le chapitre 3) dans ces milieux qui permettent de caractériser et calculer les modes (généralisés, guidés ou localisés suivant le milieu). Dans chacun de ces problèmes, ma démarche est toujours la même : après une analyse mathématique des problèmes (étude théorique des modèles mathématiques relativement à l'existence, l'unicité et aux propriétés qualitatives des modèles), je conçois, analyse et met en oeuvre une méthode numérique pour résoudre ces problèmes.

Ces problèmes font intervenir plusieurs échelles (la dimension du milieu, la longueur d'onde, la période,...) et des questions de modélisation se posent également. Dans toutes les applications que j'ai citées, la taille du milieu est beaucoup plus grande que la longueur d'onde et la période.

Il parait donc légitime de supposer que le milieu est infini. En outre, on sait que quand la longueur d'onde est de l'ordre de la période, il faut considérer le milieu périodique en tant que tel (c'est souvent le cas dans les applications en photonique). Mais quand elle est grande devant la période (comme souvent dans les techniques de CND par ondes ultrasonores des matériaux composites), le milieu peut être homogénéisé c'est-à-dire remplacé par un milieu dit effectif ou homogène (souvent anisotrope) équivalent. J'ai commencé à m'intéresser assez récemment à la construction et la justification de modèles effectifs, obtenu à l'aide de méthodes asymptotiques (voir le chapitre 4).

Chapitres 1 et 2 -Problèmes de diffraction des ondes acoustiques et élastiques en régime harmonique dans des milieux périodiques ou anisotropes : depuis ma thèse, je me suis intéressée à l'analyse mathématique et numérique d'équations d'ondes acoustiques, électromagnétiques ou élastiques en régime harmonique dans des milieux périodiques infinis dans une direction (dans ce cas c'est un guide d'onde) ou plusieurs directions. Ces équations sont obtenues à partir de l'équation des ondes en régime temporel en supposant que le terme source est harmonique en temps et en cherchant la solution elle aussi harmonique en temps. L'intérêt est que la solution recherchée ne dépend plus que des variables d'espace puisque la dépendance en temps peut être factorisée. Bien que d'apparence plus simple, les équations harmoniques, dont l'exemple type est l'équation de Helmholtz, posées dans un domaine non borné, soulèvent de nombreuses difficultés. La première difficulté est que ces équations posées dans un milieu infini sont en général mal posées dans le cadre classique des solutions d'énergie finie : il est assez facile de vérifier que les solutions physiques ne sont en général pas d'énergie finie. Si on cherche des solutions ayant une régularité seulement locale, le problème est également mal posé : il est assez facile dans ce cas de construire une infinité de solutions. Il faut en général imposer en plus, une condition sur le comportement de la solution à l'infini, appelée condition de rayonnement ou radiation à l'infini. La deuxième difficulté est d'ordre numérique : comment calculer la solution alors que le milieu est infini ?

Dans le chapitre 1, je donne quelques éléments sur mes travaux qui traitent des problèmes de diffraction dans des guides d'ondes. La spécificité des guides d'ondes est que le milieu est infini dans seulement une direction (ce qui peut paraitre plus simple !) mais la solution physique se propage, en général, jusqu'à l'infini sans décroitre. Dans ce chapitre, je commence par traiter le cas des guides d'ondes acoustiques homogènes pour expliquer les spécificités du guide d'onde et pour rappeler les résultats classiques. En effet, dans ce cas, des techniques de type séparation de variables permettent de résoudre les problèmes d'ordre théorique et numérique évoqués plus haut. Cependant, ces techniques ne peuvent pas s'appliquer ou s'étendre à des guides d'ondes plus généraux par exemple si le milieu est périodique ou si il s'agit d'ondes élastiques. Depuis ma thèse, je travaille sur les guides d'ondes périodiques avec Patrick Joly. Tout d'abord, en utilisant des outils adaptés aux milieux périodiques comme la transformée de Floquet Bloch, nous avons proposé des conditions de radiation pour les guides d'ondes périodiques [START_REF] Fliss | Solutions of the time-harmonic wave equation in periodic waveguides: asymptotic behaviour and radiation condition[END_REF]. Ces conditions permettent de rendre le problème bien posé en général. Néanmoins, elles ne peuvent pas être utilisées directement numériquement. Nous avons répondu aussi aux Introduction difficultés numériques en construisant des conditions transparentes imposées sur une frontière artificielle, qui sont basées sur des opérateurs de type Dirichlet-to-Neumann (DtN). Même s'il reste quelques questions ouvertes, l'utilisation des opérateurs DtN répond à la fois aux difficultés théoriques et numériques [1] [24]. Toutes ces techniques ont été utilisées avec Laurent Bourgeois (POEMS), pour résoudre un problème inverse dans un guide d'onde périodique (quand la cellule de périodicité est connue) via une méthode d'échantillonnage de type Linear Sampling Method [START_REF] Bourgeois | On the identification of defects in a periodic waveguide from far field data[END_REF]. Ces méthodes ont également été utilisées pour résoudre des problèmes de diffraction entre un demi-espace homogène et un demi-espace périodique [START_REF] Fliss | Computation of light refraction at the surface of a photonic crystal using dtn approach[END_REF]. Plus récemment, dans le cadre de la thèse d'Antoine Tonnoir, avec Anne-Sophie Bonnet Ben Dhia (POEMS) et Vahan Baronian (CEA-LIST), nous avons considéré la diffraction des ondes élastiques dans un guide d'ondes isotrope ou anisotrope (correspondant à un matériau composite homogénéisé). Une des applications visées est le contrôle non destructif de tubes ou de câbles. Les difficultés ici sont induites par le caractère vectoriel des équations qui ne permet pas d'étendre les techniques du cas scalaire. Nous avons proposé de nouvelles conditions transparentes dans [START_REF] Vahan Baronian | Iterative methods for scattering problems in isotropic or anisotropic elastic waveguides[END_REF] ce qui n'avait jamais été fait pour des guides élastiques anisotropes. Ces conditions utilisent une décomposition de la solution en fonction des modes du guide (appelés modes de Lamb) pour lesquels la propagation le long du guide est explicite. Ces conditions nécessitent l'introduction de 2 frontières artificielles (comme dans certaines méthodes de couplage Eléments Finis-Equations intégrales ou dans les méthodes de décomposition de domaines avec recouvrement). Nous nous sommes ensuite rendus compte que cette technique de construction de conditions transparentes avaient un portée plus générale. Par exemple, pour les guides périodiques, cela permet de construire d'autres conditions transparentes basées sur une décomposition sur les modes de Floquet. De manière peut être surprenante, la diffraction dans des guides d'ondes est un sujet ancien, bien connu mais loin d'être clos. En effet, dès qu'il s'agit de traiter des guides d'ondes autres que scalaires et homogènes, des méthodes un peu ad-hoc voire bricolées sont développées. Nous pensons que l'approche présentée dans la dernière section de ce chapitre est très générale et il est très tentant de la tester sur toutes les situations (guides d'ondes électromagnétiques, de type plaque (voir [20]), ou guides d'ondes ouverts par exemple) pour lesquelles il n'existe peu ou rien.

Dans le chapitre 2, j'évoque mes travaux sur les problèmes de diffraction dans des milieux infinis dans au moins 2 directions. Les questions théoriques concernant les conditions de radiation pour les milieux périodiques restent encore ouvertes. Pour ces problèmes, je me suis concentrée pour l'instant sur les difficultés numériques en me plaçant dans un cadre coercif (en rajoutant de la dissipation au modèle). Si le milieu est périodique ou si le milieu est homogénéisé mais qu'il s'agit d'ondes élastiques, les méthodes classiques pour restreindre les calculs autour de zones d'intérêt (sources, obstacles, perturbations) ne fonctionnent pas (même avec de la dissipation). Nous avons développé la méthode dite des demis-espaces raccordés ou la Half-Space Matching Method. La méthode part d'une idée assez simple, que le milieu soit homogène ou périodique, on peut exprimer la solution dans un demi-plan, de manière explicite ou semi-explicite, à partir de sa trace sur le bord en utilisant une transformation adaptée (Fourier ou Floquet-Bloch). L'approche consiste ensuite à coupler les représentations de la solution dans plusieurs demi-plans (au moins 3) avec une représentation de la solution autour du défaut. En assurant la compatibilité des représentations dans les zones de recouvrement, on aboutit à une formulation couplant, via des opérateurs intégraux, la solution dans un domaine borné contenant le défaut et ses traces sur les bords des demi-plans. Pendant ma thèse avec Patrick Joly, nous avons introduit la méthode dans le cas d'un milieu périodique avec des directions de périodicité orthogonales [START_REF] Fliss | Exact boundary conditions for time-harmonic wave propagation in locally perturbed periodic media[END_REF]3,[START_REF] Fliss | Wave propagation in locally perturbed periodic media (case with absorption): Numerical aspects[END_REF], le cas d'un réseau hexagonal a été traité dans [6]. Cette méthode a ensuite été utilisée dans le cadre d'une collaboration avec Vahan Baronian du CEA, de la thèse d'Antoine Tonnoir (soutenue en 2015) et de celle de Yohanes Tjandrawidjaja (débutée en 2016), toutes deux co-encadrées avec Anne-Sophie Bonnet-Ben Dhia (POEMS), pour étudier la propagation des ondes élastiques dans des plaques dites épaisses (bornées dans une direction et infinies dans les 2 autres) de matériaux composites. Pendant la thèse d'Antoine, nous avons traité des milieux 2D et le cas d'ondes acoustiques et élastiques. Le cas de la diffraction des ondes acoustiques dans des milieux homogènes même anisotropes, que nous traitons dans [START_REF] Bonnet-Ben Dhia | The halfspace matching method: a new method to solve scattering problem in infinite media[END_REF] et qui fait l'objet du début du chapitre 2, est évidemment un problème beaucoup plus simple que tous ceux que je viens d'évoquer. Et il semble sans intérêt d'utiliser cette méthode d'une telle complexité (avouons le !) pour un cas où toutes les méthodes classiques marchent très bien. Mais sur ce cas simple, il est plus commode, je pense, de comprendre le principe de la méthode. Et surtout, sur ce cas nous pouvons faire l'analyse théorique [START_REF] Bonnet-Ben Dhia | The halfspace matching method: a new method to solve scattering problem in infinite media[END_REF], l'analyse numérique [START_REF] Bonnet-Ben Dhia | Numerical analysis of the halfspace matching method with robin traces on a convex polygonal scatterer[END_REF], déduire assez simplement les analyses dans le cas périodique et avoir des pistes pour traiter le cas des ondes vectorielles. Comme je le mentionne à la fin du chapitre, il reste des questions ouvertes, notamment pour analyser le cas sans dissipation. Nous avons posé une première brique avec un résultat théorique qui généralise le théorème de Rellich [START_REF] Bonnet-Ben Dhia | A rellich type theorem for the helmholtz equation in a conical domain[END_REF]. Et la méthode marche très bien même en l'absence de dissipation, il suffit d'imposer une condition de radiation dans chaque demi-espace. Enfin avec cette méthode, nous pouvons envisager de traiter des situations encore plus délicates pour lesquelles il n'existe aucune analyse ou méthode numérique : citons par exemple les plaques élastiques (qui fait l'objet de la thèse de Yohanes Tjandrawidjaja, débutée en octobre 2016), les jonctions de fibres optiques périodiques et l'étude de toutes ces équations en régime temporel (traitée en partie dans la thèse de Hajer Methenni, débutée en novembre 2017.).

Chapitre 3 -Ondes piégées et guidées dans les milieux périodiques : ce deuxième sujet a été financé en partie par un projet DGA (via le financement de 18 mois de post-doc, 6 mois pour Bérangère Delourme (aujourd'hui au LAGA, Univ. Paris 13) Alors que dans les résultats théoriques existants, les conditions d'existence de gaps pour des milieux périodiques et de modes guidés ont été obtenus pour des variations de coefficients, nous avons proposé de telles conditions en gardant des coefficients constants mais en faisant varier la géométrie [15] [19, 23] (dans le premier cas, l'opérateur perturbé a le même domaine que l'opérateur sans perturbation, ce qui n'est plus vrai dans le second cas). Pour cela, nous avons considéré une géométrie qui est asymptotiquement proche d'un graphe pour lequel des calculs analytiques sont possibles. Des développements asymptotiques nous ont permis de conclure. J'ai également proposé une méthode numérique exacte pour le calcul des modes guidés qui soit indépendante du confinement. Le principe de la méthode consiste à restreindre les calculs autour de la perturbation en imposant des conditions aux limites transparentes dont la construction

Introduction

est la même que pour les problèmes de diffraction. La méthode est exacte (comparée à des méthodes dites de Supercell) mais le prix à payer est que le problème aux valeurs propres posé en domaine borné devient non linéaire. L'étude mathématique et numérique de la méthode a fait l'objet de plusieurs articles [START_REF] Fliss | A dirichlet-to-neumann approach for the exact computation of guided modes in photonic crystal waveguides[END_REF]9,[START_REF] Fliss | Robin-to-robin transparent boundary conditions for the computation of guided modes in photonic crystal wave-guides[END_REF], en collaboration notamment avec Kirsten Schmidt et Dirk Klindworth (TU Berlin) et un article sur des configurations 3D plus réalistes [25]

Introduction

In this chapter, I consider time harmonic scattering problems in unbounded closed waveguides.

Let me begin by a general introduction on waveguides. A waveguide is a structure that guides waves in one direction with minimal loss of energy. A waveguide is called closed if it is bounded in the directions which are orthogonal to the propagation direction by an impenetrable wall (such as pipelines). A waveguide is called open if its cross section can be considered unbounded, as for instance optical fibers or immersed pipes. In such devices, guided waves are time-harmonic waves that propagate without attenuation in the longitudinal direction and remain confined in the other transverse directions. For a closed waveguide, such a confinement is simply due to the boundedness of the cross section. And for an open waveguide, this confinement results from a particular layout of the various materials which compose the waveguide. In this chapter, I detail my work on closed waveguides. More precisely, I have considered two types of waveguide problems (periodic and elastic waveguides). Some specific applications are the study and the computation of ( 1) the light propagation in optical fibers for the communication industry and

(2) the ultrasonic wave propagation in cables and pipelines for non destructive testing purposes. Concerning open waveguides, I haven't addressed yet scattering problems but I have some perspectives in this direction (see Section 1.5) and in Chapter 3, I detail my work devoted to the existence, the properties and the computation of guided modes in particular open waveguides.

Time harmonic scattering problems in unbounded waveguides raise several difficulties of theoretical and numerical nature which are intricately linked. From a theoretical point of view, the difficulty concerns the definition of the physical solution. It is often defined as the unique solution of a well-posed problem. However, the time harmonic scattering problems in waveguides are in general not well posed in the classical L 2 framework. This is linked to the fact that the physical solution is in general not of finite energy since a propagation without attenuation is possible in the direction of the waveguide. In the L 2 loc framework, an infinity of solutions can be 14

Chapter 1. Scattering problem in waveguides found. Usually radiation conditions which characterize the behaviour at infinity of the physical solution have to be determined and added to the problem in order to recover well-posedness. From a numerical point of view, the domain being unbounded, the difficulty is to compute the physical solution.

These difficulties are well known and solved for homogeneous acoustic waveguides. Indeed, from a theoretical point of view, the classical approach is to remark that the equivalent problem with a dissipation term is well posed in a classical setting and then define the physical solution as the limit when the dissipation tends to 0 of the family of solutions of the problem with dissipation. This is the limiting absorption principle. When the waveguide is homogeneous and the time harmonic scalar wave equation is considered, an explicit expression of the solutions enables to show that the family of solutions has a limit when the dissipation tends to 0 and it is a solution of the problem without dissipation. This is by definition the physical solution. At infinity, the physical solution can be decomposed as the sum of the so-called outgoing propagative modes. Moreover, from the behaviour at infinity of the physical solution, one can derive some radiation conditions. Adding these radiation conditions to the original problem makes it well-posed (except for a countable set of frequencies which corresponds to the so-called resonances of the problem). This approach cannot be extended directly to periodic or elastic waveguides since the solution, even in the dissipative case, cannot be expressed as explicitly as for the homogeneous acoustic case. I have worked on the limiting absorption principle and the derivation of radiation conditions for periodic waveguides (see Section 1.3.2), for another approach, see [Nazarov andPlamenevsky, 1990,Nazarov andPlamenevsky, 1991], which was already applied to periodic and/or elastic waveguides [Nazarov, 2013, Nazarov, 2014b].

From a numerical point of view, there exist two classes of method to compute the solution of the time harmonic scalar wave equation in a locally perturbed homogeneous waveguide. A first class of methods consists in putting on each side of the computational domain an absorbing layer in which the Perfectly Matched Layer (PML) technique [Berenger, 1994] is applied. Physically the method can be interpreted as letting a wave coming from the computational domain enters the layer without reflexion and absorbs the wave inside the layer preventing it to come back in the computational domain. This technique is easy to integrate in any Finite Element codes (no specific implementation is needed, adding a PML layer corresponds to making a (complex) change of variable) and leads to solve a classical sparse linear system. Unfortunately, it is well-known that PMLs do not work in elastic and periodic waveguides (see for instance [START_REF] Skelton | Guided elastic waves and perfectly matched layers[END_REF] for elastic waveguides). Indeed, PMLs absorb the waves having positive phase velocities. However, the physical or outgoing solution can be characterized through its positive group velocity. Even if in scalar homogeneous waveguide, phase and group velocities are of same sign, in general (as in elastic or periodic waveguides for instance) they are not. A remedy has been proposed and analyzed for elastic waveguides in [Bonnet-Ben Dhia et al., 2014] where the physical solution is reconstructed a posteriori by combining several "wrong fields" computed with PMLs. An alternative consists in using adiabatic viscoelastic absorbing layers [START_REF] Drozdz | Efficient numerical modelling of absorbing regions for boundaries of guided waves problems[END_REF] which are not perfectly matched and need to be sufficiently large to avoid spurious reflections. The main drawback of this approach is then its computational cost. Moreover, let us point out that absorbing layer techniques (perfectly matched or not) require a fine adjustment of some parameters, which may limit their systematic use. Let us finally mention the Hardy space infinite elements method known also as the pole condition [Hohage et al., 2003a[START_REF] Hohage | Solving timeharmonic scattering problems based on the pole condition. II. Convergence of the PML method[END_REF] which was recently extended to elastic waveguides [Halla andNannen, 2015, Halla et al., 2016].

A second possibility consists in using the modal decomposition of the field outside the perturbed area to derive transparent boundary conditions which can be set on the artificial boundaries of the finite element domain. The advantage is that such conditions are exact (and with an exponentially small error at the discrete level if enough modes are kept in the modal expansion). However, it requires specific properties of the modes of the isotropic homogeneous acoustic waveguide (in particular the orthogonality of the modes in the transverse section of the waveguide), which does not hold for instance in general waveguides (stratified, anisotropic or periodic waveguides or in (even isotropic) elastic waveguides). Let me finally mention that for isotropic elastic waveguides, an alternative has been proposed in [START_REF] Baronian | Couplage des méthodes modale et éléments finis pour la diffraction des ondes élastiques guidées[END_REF], Baronian et al., 2010] based on bi-orthogonality relations but it cannot be extended to general anisotropic elastic waveguides.

In the first section of this chapter, I will recall quickly the principal results that I have mentioned, for the homogeneous acoustic waveguide. I present my contributions on periodic waveguide in Section 1.3. As you will see, there are theoretical and numerical results. One of the most original result is maybe the construction of the DtN operators. Let me add that this method can be extended directly to scalar or vectorial equations and of course isotropic/anisotropic homogeneous or periodic media. Moreover, this construction is not based on a modal decomposition of the solution far from the defect. Of course there exists a link but from a numerical point of view, this does no require the computation of the modes.

For some applications (non destructive testing of pipes for instance), it could be important to keep the decomposition in terms of the modes of the guides. However, the modes in general do not have the magic properties (orthogonality for instance) of the modes of the isotropic homogeneous acoustic waveguide. The construction of transparent boundary conditions can then be really intricate (for isotropic elastic waveguides) or impossible (anisotropic elastic waveguides). By working on elastic waveguides, we have constructed new transparent boundary conditions whose construction has, I think, a general range. For pedagogical purposes, I will describe these new transparent boundary conditions in Section 1.4 in the simple case of the isotropic homogeneous acoustic waveguide and I will extend them for anisotropic elastic waveguides and periodic waveguides. This transparent boundary conditions require a prior computation of the modes. Finally, my ongoing works and perspective for this subject are described in Section 1.5.

Classical results for homogeneous acoustic waveguides

Let me here recall some well-known results for the acoustic isotropic homogeneous waveguides that can be found for instance in [START_REF] Harari | Dirichlet-to-Neumann maps for unbounded wave guides[END_REF], Hagstrom, 1999, Lenoir and Tounsi, 1988].

Let us consider a diffraction problem in an acoustic isotropic straight guide Ω = S × R where S ⊂ R 2 (3D waveguides) or S ⊂ R (2D waveguides) denotes the bounded cross-section of the Chapter 1. Scattering problem in waveguides guide (see Figure 1.1(left)). We want to study in this section the physical solution u (also called

∂ ν u = 0 on ∂Ω, (1.1)
where ω is the frequency, ν is the exterior normal to ∂Ω and the source term f is supposed to be compactly supported in, let us say, {|z| < a}.

We want to define the physical solution of (1.1) and derive an equivalent problem which is suitable for numerical computations.

A general approach is to consider a similar problem adding some absorption (the presence of the damping term ε > 0 guaranteeing the well-posedness of this problem in H 1 ). More precisely, let us consider the unique solution

u ε in H 1 (Ω) of -∆u ε -(ω 2 + ıε)u ε = f in Ω, ∂ ν u ε = 0 on ∂Ω. (1.2)
If u ε has a limit in a certain sense (not necessarily in H 1 (Ω)), then the physical solution is defined as this limit. This definition is not explicit so the difficulty now is to derive a more explicit characterization of this physical solution.

Let us concentrate on the problem with dissipation. We know that u ε is solution of (1.1) if and only if

(u a ε = u ε | Ω a , u + ε = u ε | Ω + , u - ε = u ε | Ω -) is solution of -∆u a ε -(ω 2 + ıε)u a ε = f in Ω a , ∂ ν u a ε = 0 on ∂Ω a ∩ ∂Ω, -∆u ± ε -(ω 2 + ıε)u ± ε = 0 in Ω ± , ∂ ν u ± ε = 0 on ∂Ω ± ∩ ∂Ω, u a ε | Γ a,± = u ± ε | Γ a,± , ∂ z u a ε | Γ a,± = ∂ z u ± ε | Γ a,± ,
where

Ω a = S × (-a, a), Ω + = S × (a, +∞), Ω -= S × (-∞, -a) and Γ a,± = {(x s , z) ∈ Ω, z =
u ± ε (x s , z) = k∈N a ± k,ε ϕ k (x s )e ∓β k,ε (z-a) ,
where the {ϕ k , k ∈ N} are an orthonormal basis of eigenvectors (also called transverse modes)

of the transverse laplacian operator -S defined on S with Neumann boundary conditions, that are associated to the eigenvalues, denoted λ k which are positive and tends to +∞. We order the eigenvalues (and then the eigenvectors) increasingly. For instance for S = (0, 1), λ k = (kπ) 2 and ϕ k (x s ) = c k cos(kx s ). The coefficients β k,ε are given by

β k,ε = λ k -(ω 2 + ıε),
where by convention Re √ ≥ 0.

From the properties of the transverse modes, we have

a ± k,ε = (u ± ε | Γ a,± , ϕ k ) L 2 (Γ a,± ) .
This explicit expression for u ± ε with respect to their traces on Γ a,± enables to introduce the problem satisfied by u a ε (by eliminating u ± ε ):

-∆u a ε -(ω 2 + ıε)u a ε = f in Ω a , ∂ ν u a ε = 0 on ∂Ω a ∩ ∂Ω, ±∂ z u a ε + Λ ± ε u a ε = 0 on Γ a,± , (1.3)
where Λ ± ε are the so-called Dirichlet-to-Neumann operators given here by

∀ϕ ∈ H 1/2 (Γ a,± ), Λ ± ε ϕ = k∈N β k,ε (ϕ, ϕ k ) L 2 (Γ a,± ) ϕ k ∈ H -1/2 (Γ a,± ).
The last relations in (1.3) satisfied by u a ε are called transparent boundary conditions since they represent exactly the behaviour of u ± ε in Ω ± . It is easy to show that Problem (1.3) is coercive in H 1 (Ω a ) and well-posed.

Moreover, Problem (1.3) is equivalent to Problem (1.2) since the restriction to Ω a of the solution of (1.2) is solution of (1.3) and from a solution u a ε of (1.3), we can construct the solution -a) .

u ε of (1.2) by u ε | Ω a = u a ε u ε | Ω ± = k∈N (u a ε | Γ a,± , ϕ k ) L 2 (Γ a,± ) ϕ k (x s )e ∓β k,ε (z
(1.4)

We have then derived an equivalent problem which is suitable for numerical computations (and the numerical analysis of the associated numerical method can be done since the problem is coercive).

Let us now study the case without dissipation by applying the limiting absorption principle. We can show that the family {u ε , ε > 0} has a limit in H 1 loc (Ω) (except for a countable set of frequencies). We define this limit as the physical solution of (1.1) and we can derive a well-posed problem satisfied by the physical solution which is suitable for numerical computations.

1. First, let us remark that the Dirichlet-to-Neumann operator Λ ± ε ∈ L(H 1/2 (Γ a,± ), H -1/2 (Γ a,± )) has a limit in norm towards the operator Λ ± defined by

∀ϕ ∈ H 1/2 (Γ a,± ), Λ ± ϕ = k∈N β k (ϕ, ϕ k ) L 2 (Γ a,± ) ϕ k ∈ H -1/2 (Γ a,± ),
where

β k = -ı ω 2 -λ k , if λ k < ω 2 , λ k -ω 2 , if λ k ≥ ω 2 .
(1.5)

2. Let us introduce the following problem set in

Ω a -∆u a -ω 2 u a = f in Ω a , ∂ ν u a = 0 on ∂Ω a ∩ ∂Ω, ±∂ z u a + Λ ± u a = 0 on Γ a,± .
(1.6)

Noticing that Re(Λϕ, ϕ) ≥ 0 for any ϕ, we can show that the bilinear form associated to (1.6) in H 1 (Ω a ) is the sum of a coercive bilinear form and a compact one. Then Fredholm alternative holds : uniqueness implies well-posedness of the problem. To show uniqueness, let us consider the solution of (1.6) with f = 0. Multiplying by u a , integrating by parts and taking the imaginary part of the obtained expression, we find

(u a , ϕ k ) L 2 (Γ a,± ) = 0 if λ k < ω 2 .
From this u a , the function u defined by

u| Ω a = u a u| Ω ± = λ k >ω 2 (u a | Γ a,± , ϕ k ) L 2 (Γ a,± ) ϕ k (x s ) e ∓β k (z-a)
solves (1.1) with f = 0. If ω 2 is not a cutoff frequency (i.e. ω 2 is not an eigenvalue of the transverse laplacian or equivalently ∀, β k = 0), u is in L 2 (Ω) and if ω 2 is a cut-off frequency, u tends to a constant when x tends to infinity. The question of uniqueness is then linked to the discrete spectrum of the laplacian in the whole waveguide defined in a (weighted) Sobolev spaces. Consequently, uniqueness can be established except for at most a countable set of frequencies.

3. Finally, except for this countable set of frequencies, we conclude by showing that u a ε tends to u a in H 1 (Ω a ) when ε tends to 0. And taking the limit in (1.4), we have that u ε tends to the function u in H 1 loc defined by

u| Ω a = u a , u| Ω ± = k∈N a ± k ϕ k (x s ) e ∓β k (z-a) , with a ± k = (u a | Γ a,± , ϕ k ) L 2 (Γ a,± ) (1.7)
where β k is given in (1.5). This is by definition the physical or outgoing solution of Problem (1.1). Let us remark that far from the perturbations (i.e. in Ω ± ), the physical solution is the combination of the so-called modes of the guide, i.e. particular solutions in the distributional sense of the homogeneous problem (f = 0) : z-a) .

∀k ∈ N, w k (x s , z) = ϕ k (x s ) e ∓β k (
A finite number of the modes are propagating along the guide (if λ k < ω 2 ), they are called guided modes, and the others are exponentially decaying at infinity.

Except for a countable set of frequencies, we have then defined the physical solution and derived the well-posed problem (1.6) which is suitable for numerical computations (and the numerical analysis can be done since the problem is "coercive + compact"). Moreover, after computing u in Ω a (using the FE Method for instance), by (1.7) u can be computed in the whole waveguide.

Remark 1.2.1

This can be easily extended to some more general problems : to domains which are local perturbations of this straight domain Ω, to the presence of local penetrable or non penetrable obstacles, to junctions of different waveguides (since the Dirichlet-to-Neumann operators are defined independently the one from the other), to Dirichlet, Robin or periodic boundary conditions. See Figure 1.1(right).

Let us finally make some remarks about the countable set of frequencies for which the problem in the general case is not well-posed. Indeed, there may exist frequencies for which there exists L 2 -solutions, the so-called trapped modes, of Problem (1.6) with f = 0. These frequencies corresponds to the well-known resonances, i.e. eigenvalues of the laplacian and the trapped modes are associated eigenvectors. We know that in the simple problem (1.6) Let us emphasize that this approach allows to justify the use of PML techniques but since these techniques cannot be extended, to our knowledge, for periodic media, we do not give more details.

The key property of the Helmholtz equation set in homogeneous waveguides that makes the separation of variables possible and enables an explicit expression of the solution, is that the transverse operator is self-adjoint. In the periodic case and/or the elastic equation, this is not true anymore and alternatives have to be found.

Periodic waveguides

Model problem

The model problem that we consider in this section is that of a periodic waveguide Ω ⊂ S×R, and a simple scalar model. The treatment of other equations (Maxwell's equations or elastodynamic equations) would be similar in principle. As in the previous section, we will consider a perfectly periodic problem and a more general problem. In the first case, the geometry Ω as well as the material properties, typically the refraction index n, are periodic in the z-direction

(z, x s ) ∈ Ω ⇒ (z ± L, x s ) ∈ Ω and n(z ± L, x s ) = n(z, x s ) (1.8)
and in the second case they are z-periodic except in a bounded region and the periodicity could be different on each side of this bounded region where Ω ± p and n ± p are periodic in the z-direction with period L ± (see Figure 1.3). We suppose also that the boundary of Ω is smooth enough.

Ω ∩ {±z ≥ a} = Ω ± p ∩ {±z ≥ a} ± z ≥ a, n(z, x s ) = n ± p (z, x s ),
Let us consider then the propagation of a time harmonic scalar wave in Ω due to a given compactly supported source f ∈ L 2 (Ω)

-u -n 2 ω 2 u = f in Ω ∂ ν u = 0 on ∂Ω (1.9) ν being the exterior normal of Ω, n ∈ L ∞ (Ω) with n ≥ c > 0.
In general (more precisely if the frequency ω is such that ω 2 is in the spectrum of the corresponding periodic differential operator

A = -n -2 , D(A) = {u ∈ H 1 (Ω), u ∈ L 2 (Ω), ∂ ν u ∂Ω = 0})
, waves can propagate in the waveguide without attenuation in the longitudinal direction. In this case, as we have seen in the homogeneous case, the H 1 classical framework is not appropriate anymore : solutions cannot exist in this space. On the other hand, in H 1 loc = {u, ∀ϕ ∈ D(R) ϕ(x 1 )u(x 1 , x s ) ∈ H 1 (Ω)}, one can define several solutions to the original problem. Additional conditions are required to define the "good" physical solution of (1.9). Let us remark that for the homogeneous case, the spectrum of the differential operator is R + and then the problem is never well posed in H 1 . For the periodic case, the spectrum may contain gaps and when the frequency lies is the gap, the problem is well posed in H 1 (except for the resonances of the problem).

Again, we can characterize the physical solution by the limiting absorption principle. It is the limit, when there exists, of the unique solution in H 1 (Ω) of the corresponding mathematical model with absorption, when this absorption tends to 0.

The limiting absorption principle has been extensively developed in the literature in various situations such as locally perturbed homogeneous media (see for instance [START_REF] Èidus | On the principle of limiting absorption[END_REF],Wilcox, 1967,Agmon, 1982]) or locally perturbed stratified media (see for instance [Wilcox, 1984[START_REF] Eidus | The limiting absorption and amplitude principles for the diffraction problem with two unbounded media[END_REF], Weder, 1990, Bonnet-Ben Dhia and Tillequin, 2001b]). There are mush less results in periodic media. In [START_REF] Gérard | The Mourre theory for analytically fibered operators[END_REF]] , the authors use Mourre's theory [Mourre, 1981, Jensen et al., 1984] to prove the limiting absorption principle for analytically fibered operators (differential operators with periodic coefficients are a particular case). In [Iftimie, 2003], the author uses the same theory to deal with the Laplace-Dirichlet operator in an (n +1)dimensional homogeneous layer with periodically shaped boundary. In [Levendorskiȋ, 1998], a limiting absorption principle for a two-dimensional periodic layer with perturbation is proven. In all these works, the authors are interested mainly in the existence of the limit of the solution with absorption which amounts to study the limit of the resolvent of the corresponding periodic differential operator, and in which sense this limit holds. Their approach is not constructive and therefore cannot be used to characterize more explicitly the physical solution.

For the perfectly periodic problem, we consider a more constructive approach, using the Floquet-Bloch theory, that leads to a semi-analytical expression of the limit, the physical solution, as in [Hoang, 2011] for the periodic half-waveguide and [Radosz, 2009] for more general infinite periodic media. Actually, it is easy to construct other solutions of (1.9). How can we distinguish this physical solution to other solutions? This is from its behaviour at infinity that we can deduce a radiation condition that characterizes uniquely the physical solution. This is explained in Section 1.3.2.

It is not obvious to derive a numerical method to compute the physical solution from this approach. Similarly to what it is done for homogeneous acoustic waveguides (explained in Section 1.2), we have constructed transparent boundary conditions involving DtN operators. This answers both to the theoretical and the numerical difficulties. Of course, the construction of the DtN operators, comparing to Section 1.2, is more involved since separation of variables can no longer be used. As you will see in Section 1.3.3, the expressions of the solution and the Dirichletto-Neumann operator even in the case with dissipation are indeed less explicit. The case without absorption will be the subject of Section 1.3.4.

Remark 1.3.1 (Extensions)

Most of the results can be extended to

• more general symmetric, second order elliptic differential operators with real periodic coefficients (with the same period);

• other boundary conditions as soon as they satisfy the same periodicity properties than the geometry and the coefficients.

Limiting absorption principle and radiation condition for perfectly periodic waveguides

This is a summary of the results that you can find in [START_REF] Fliss | Solutions of the time-harmonic wave equation in periodic waveguides: asymptotic behaviour and radiation condition[END_REF]. We want to characterize the physical solution of (1.9) in the perfectly periodic case. Let us define the periodicity cell of the medium by

C = {(z, x s ) ∈ Ω, -L/2 < z < L/2}
and let us study the unique solution in

H 1 (Ω) of -u ε -n 2 (ω 2 + ıε) u ε = f in Ω ∂ ν u ε = 0 on ∂Ω, (1.10)
where n and Ω are perfectly periodic, i.e. they satisfy (1.8).

Using the Floquet Bloch Transform in the z-direction and the well posedness in L 2 (Ω) of Problem (1.10), it is easy to show that the solution u ε of (1.10) is given by

∀(z, x s ) ∈ C, ∀p ∈ Z, u ε (z + pL, x s ) = n∈N L √ 2π π/L -π/L P n (f )(z, x s ; k) λ n (k) -(ω 2 + ıε) e ıpk dk in H 1 (Ω)
where .11) where, f is the Floquet-Bloch transformation of f (see [START_REF] Fliss | Étude mathématique et numérique de la propagation des ondes dans des milieux périodiques localement perturbés[END_REF] and [START_REF] Fliss | Solutions of the time-harmonic wave equation in periodic waveguides: asymptotic behaviour and radiation condition[END_REF] for the definition and the properties), for all k ∈ (-π/L, π/L), (•, •) C is the scalar product in L 2 (C), λ n (k) is an eigenvalue and ϕ n (•; k) an associated eigenvector of the self-adjoint and positive operator

P n (f )(• ; k) = f (• ; k), ϕ n (• ; k) C ϕ n (• ; k). ( 1
A(k) = - 1 n 2 D(A(k)) = {u ∈ H 1 (C), u ∈ L 2 (C), ∂ ν u ∂C∩∂Ω = 0, u z=L/2 = e ıkL u z=-L/2 ∂ z u z=L/2 = e ıkL ∂ z u z=-L/2 }.
(1.12) such that (using [Kato, 1995]), each dispersive curve k → λ n (k) is analytic. Let me emphasize that the eigenvalues λ n (k) are not the eigenvalues of A(k) which are ordered increasingly (these ones are only piecewise analytic). The λ n (k)'s are a rearrangement of the sequence of the eigenvalues which are ordered increasingly such that they are analytic with respect to k. See [START_REF] Fliss | Solutions of the time-harmonic wave equation in periodic waveguides: asymptotic behaviour and radiation condition[END_REF] where we show that we can apply the perturbation theory of operators of [Kato, 1995].

To pass to the limit when ε goes to 0, let us define the finite sets .14) We can show, using the properties of the dispersion curves, that

I(ω) = {n ∈ N, ∃ ξ ∈ (-π/L, π/L), λ n (ξ) = ω 2 } (1.13) and for n ∈ I(ω) Ξ n (ω) = {ξ ∈ (-π/L, π/L), λ n (ξ) = ω 2 }. ( 1
ξ ∈ {Ξ n (ω), n ∈ I(ω)} ⇔ -ξ ∈ {Ξ n (ω), n ∈ I(ω)} Let A = -n -2 , D(A) = {u ∈ H 1 (Ω), u ∈ L 2 (Ω), ∂ ν u ∂Ω = 0}.
Using the Floquet-Bloch Theory [Kuchment, 2004], one can show that σ

(A) = ∪{σ(A(k)), k ∈ (-π/L, π/L]}. If ω 2 / ∈ σ(A) we have I(ω) = ∅. If ω 2 ∈ σ(A)
, for all n ∈ I(ω) and all ξ ∈ Ξ n (ω), we can introduce the corresponding propagating Floquet modes namely the functions in H 1 loc (Ω) which are the ξ-quasiperiodic extensions of the eigenvectors ϕ n (• ; ξ), that we still denote ϕ n for simplicity. Introducing for all n ∈ I(ω) and for all ξ ∈ Ξ n (ω), ω n (ξ) = λ n (ξ), one can reinterpret the propagative Floquet mode as the time harmonic wave obeying the dispersion relation

ω = ω n (ξ).
We can then define the group velocity by

V n (ξ) := dω n dξ (ξ) = 1 2 λ n (ξ) -1/2 λ n (ξ).
The sign of the group velocity indicates how the energy propagates in the waveguide.

Using the abstract result of [Levendorskiȋ, 1998], or the more explicit result of [Hoang, 2011] or [START_REF] Fliss | Solutions of the time-harmonic wave equation in periodic waveguides: asymptotic behaviour and radiation condition[END_REF], the limiting absorption principle can be shown except for the countable set of frequencies

σ 0 = ω ∈ R + , ∃ n ∈ I(ω), ∃ k ∈ Ξ n (ω), λ n (k) = 0 .

Periodic waveguides

Theorem 1.3.1 (Limiting absorption principle [START_REF] Fliss | Solutions of the time-harmonic wave equation in periodic waveguides: asymptotic behaviour and radiation condition[END_REF])

For all ω / ∈ σ 0 then ∀ p ∈ Z, lim ε→0 u ε -u H 1 (C+pL) = 0,
where u is solution of the Helmholtz equation (1.9) and is given by 1.15) where P n (f ) is defined in (1.11).

∀(z, x s ) ∈ C, ∀p ∈ Z, u(z + pL, x s ) = n / ∈I(ω) L 2π π/L -π/L P n (f )(z, x s ; k) λ n (k) -ω 2 e ıpkL dk + n∈I(ω) L 2π   p.v. π/L -π/L P n (f )(z, x s ; k) λ n (k) -ω 2 e ıpkL dk + ıπ ξ∈Ξn(ω) P n (f )(z, x s ; ξ) |λ n (ξ)| e ıpξL   (
Remark 1.3.2

When ω 2 / ∈ σ(A) (i.e. when I(ω) = ∅), we can show that u ∈ H 1 (Ω).

The set σ 0 corresponds to the frequencies for which there exists a Floquet mode whose group velocity can vanish. This corresponds for instance to the cut-off frequencies for the homogeneous waveguide for which we know that the limiting absorption principle does not hold in the perfectly straight waveguide.

We know now that for ω 2 / ∈ σ 0 , the physical solution exists. It can be shown easily that it is a solution of (1.9) in H 1 loc . When ω 2 ∈ σ(A), another solution can be constructed taking the limit, when ε goes to 0, of the solution of the dissipative Helmholtz equation (1.10) replacing ε by -ε. The limit has the same expression than (1.15) replacing in (1.15), ıπ by -ıπ. And we could even take an linear combination of these two solutions to obtain other solution of (1.9) in H 1 loc . A condition has to be added to Problem (1.9) to recover the uniqueness and then the well-posedness. This is the asymptotic behaviour at infinity of the physical solution which allows to distinguish it from the other solutions of (1.9).

To find the asymptotic behaviour of the physical solution at infinity, we use the analyticity of the eigenvalues k → λ n (k) and the corresponding eigenvectors k → ϕ n (•; k) with respect to k for n ∈ I(ω) (see [Kato, 1995]). We can show that Theorem 1.3.2 (Asymptotic behaviour of the physical solution [START_REF] Fliss | Solutions of the time-harmonic wave equation in periodic waveguides: asymptotic behaviour and radiation condition[END_REF])

u ∼ z→+∞ ı √ 2πL n∈I(ω) ξ∈Ξn(ω) λ n (ξ)>0 f (• ; ξ), ϕ n (• ; ξ) L 2 (C) λ n (ξ) ϕ n (•; ξ) u ∼ z→-∞ ı √ 2πL n∈I(ω) ξ∈Ξn(ω) λ n (ξ)<0 f (• ; ξ), ϕ n (• ; ξ) L 2 (C) |λ n (ξ)| ϕ n (•; ξ)
where the notation u ∼ z→±∞ v ± means that uv ± tends exponentially to 0 when z tends to ±∞.

Remark 1.3.3

For the homogeneous case, a similar radiation condition can be obtained thanks to the decomposition Chapter 1. Scattering problem in waveguides (1.7) of the solution as a infinite sum of evanescent modes and a finite sum of the outgoing propagative modes. The situation for a periodic closed waveguide is similar in nature to the case of the homogeneous waveguide. However, as you have noticed, the notion of outgoing modes is much more delicate and the analysis relies on quite different mathematical tools (Floquet-Bloch transform, spectral theory of operators depending analytically on a parameter, complex contour integral techniques,...).

Finally, we show that this behaviour characterizes the physical solution. In other words, there exists a unique solution of (1.9) having this behaviour at infinity. In order to simplify the expression of the result, let us introduce a relabelling of the set

n∈I(ω) {ξ ∈ Ξ n (ω), λ n (ξ) > 0} = {ξ + 1 ≤ . . . ≤ ξ + ≤ . . . ≤ ξ + N (ω) }
where N (ω) is the number of Floquet modes propagating to the right. Using the properties of the dispersion curves, we have

n∈I(ω) {ξ ∈ Ξ n (ω), λ n (ξ) < 0} = {-ξ + 1 ≥ . . . ≥ -ξ + ≥ . . . ≥ -ξ + N (ω) }
We relabel accordingly the set {ϕ n (•; ξ), ξ ∈ Ξ n (ω), n ∈ I(ω) } and renormalize them for convenience of notation. More precisely we introduce ϕ 1 , . . . , ϕ N (ω) respectively ϕ -1 , . . . , ϕ -N (ω) defined as

∀n ∈ I(ω), ∀ξ ∈ Ξ n (ω), ξ = ±ξ + m ⇒ ϕ ±m = ϕ n (• ; ξ) |λ n (ξ)|
Then the result of Theorem 1.3.2 becomes

u ∼ z→±∞ ı √ 2πL N (ω) m=1 f (• ; ±ξ + m ), ϕ ±m C ϕ ±m .
The propagative Floquet modes {ϕ m , 1 ≤ m ≤ N (ω)} are outgoing regarding the propagation towards +∞ and ingoing regarding the propagation towards -∞. Conversely, the propagative Floquet modes {ϕ -m , 1 ≤ m ≤ N (ω)} are ingoing regarding the propagation towards +∞ and outgoing regarding the propagation towards -∞.

We can now define a radiation condition and establish the well-posedness of Problem (1.9) set in the perfectly periodic waveguide.

Definition 1.3.3 (The outgoing radiation condition)

We say that u satisfies the outgoing radiation condition if and only if there exist

2N (ω) complex numbers (u ± n ) n∈ 1,N (ω) such that u ∼ z→±∞ N (ω) m=1 u ± m ϕ ±m . (1.16)
Theorem 1.3.4 (Well-posedness of the problem [START_REF] Fliss | Solutions of the time-harmonic wave equation in periodic waveguides: asymptotic behaviour and radiation condition[END_REF]) Suppose ω 2 / ∈ σ 0 . There exists a unique solution of problem (1.9) which satisfies the outgoing radiation condition.

The proof is done in two steps. Firstly, we show that the outgoing Floquet modes have a positive energy flux. By energy conservation, we can deduce that the outgoing solution of the problem with f = 0, has no contribution on the propagative modes and is then necessarily exponentially decaying at infinity. Secondly, it suffices to use the spectral property of the operator A to conclude that except for the cut-off frequencies, the physical solution necessarily vanishes.

By using truncature functions as in [Nazarov, 2014b], this result can be extended to more general periodic waveguides (See Figure 1.3 (right)). Let me mention the work of [Nazarov, 2014b, Kirsch andLechleiter, 2018] where another form of these radiation conditions were proposed. We have then a problem which characterizes uniquely the physical solution. However it is not obvious to derive a numerical method from this non local radiation condition to compute the physical solution (in the radiation conditions, the values (u ± n ) n are not known and the dependence with respect to u is not straighforward). The Dirichlet-to-Neumann operators or more generally transparent boundary conditions, are a natural alternative.

Construction of transparent boundary conditions in presence of dissipation

In [START_REF] Joly | Exact boundary conditions for periodic waveguides containing a local perturbation[END_REF], we propose a method for constructing DtN operators by solving local problems on a single periodicity cell. This is closely connected to operator-valued Riccati equations (here, of stationary nature), a topic which is already present in many problems concerning artificial boundary conditions (see, for instance, [START_REF] Lu | The riccati method for the helmholtz equation[END_REF], Henry and Ramos, 2004[START_REF] Champagne | Riccati equation for the impedance of waveguides[END_REF]). It appears also that our method is similar to the matrix transfer approach developed for ordinary equations with periodic coefficients [START_REF] Magnus | Hill's equation[END_REF]. However, except in the 1D case ( [Potel et al., 2001, Figotin andGorentsveig, 1998]), this theory cannot be applied directly to our problem due to the fact that the Cauchy problem for the Helmholtz equation is ill-posed in higher dimensions.

To treat the general periodic waveguide problem (1.9), let us consider first the problem with dissipation

-u ε -n 2 (ω 2 + ıε) u ε = f in Ω ∂ ν u ε = 0 on ∂Ω (1.17)
The support of the source term f is supposed to be compactly supported in Ω a = Ω ∩ {-a < z < a} (see Figure 1.4). The two infinite periodic sub-domains Ω ± = Ω ∩ {±z > a} are of the Chapter 1. Scattering problem in waveguides form:

Ω ± = ∞ j=0 C ± ± (j L ± , 0) ,
where the unit periodicity cells are

C ± = Ω ∩ { ± a ≤ ±z ≤ a + L ± }.
The function n(x) is "z-periodic" as well:

n(z, x s ) = n(z ± L ± , x s ), (z, x s ) ∈ Ω ± .
As you have noticed, the periodic media at each side of Ω a are not necessarily the same.

Remark 1.3.4

The boundary condition on ∂Ω can be Dirichlet, Neumann, or any combination, but they need to be compatible with the periodicity of Ω -and Ω + .

We want to find transparent boundary conditions on the two "vertical" boundaries .18) They are transparent in the sense that they are satisfied by the solution u ε of (1.17). Thus, the restriction of u ε to Ω a is solution of (1.17) in Ω a and satisfies the transparent boundary conditions on Γ ± .

Γ ± = Ω∩{z = a ± } of the form ± ∂u ε ∂z + Λ ± ε u ε = 0. ( 1
It is easy to show that the Dirichlet-to-Neumann operators involved in (1.18) are defined by

∀ϕ ∈ H 1/2 (Γ ± ), Λ ± ε ϕ = ∓ ∂ ∂z u ± ε (ϕ) Γ ± where u ± ε (ϕ) is the unique solution in H 1 (Ω ± ) of -u ± ε -n 2 p (ω 2 + ıε) u ± ε = 0 in Ω ± ∂ ν u ± ε = 0 on ∂Ω ± ∩ ∂Ω u ± ε = ϕ on Γ ± (1.19)
The constructions of Λ + ε and Λ - ε are done independently but similarly.

Note that, by periodicity, all the "vertical" interfaces

Γ ± j = Γ ± ± (jL ± , 0) can be identified to Γ ± (= Γ ± 0 ) and all the cells C ± j = C ± ± (jL ± , 0) to C ± (= C ± 0 ).
First let us give a basic result on the structure of u ± ε (ϕ) for all ϕ which use the periodic structure of the problem in Ω ± . Let P ± ε be the operator defined by

∀ϕ ∈ H 1 2 (Γ ± ), P ± ε ϕ := u ± ε (ϕ)| Γ ± 1 By identifying Γ ± 1 with Γ ± , we can consider P ± ε ∈ L(H 1 2 (Γ ± )).
The operator P ± ε is compact, injective, and its spectral radius ρ(P ± ε ) < 1. Moreover, we have, by well-posedness of the halfguide problem with periodic coefficients (and with dissipation) (1.19) in .20) This property explains the name that we often give to this operator : the propagation operator.

Ω ± that ∀ϕ ∈ H 1 2 (Γ ± ), ∀j ∈ N, u ± ε (ϕ) C ± j = u ± ε (P ± ε ) j ϕ C ± . ( 1
As a consequence, once we know the operator P ± ε , we can reconstruct the solution u ± ε (ϕ) in the whole domain Ω ± . More precisely, let us define for all ϕ ∈ H 1 2 (Γ ± ), the solutions e 0,± ε (ϕ) and e 1,± ε (ϕ) of the unit cell problems

∈ {0, 1}, -e ,± ε -n 2 p (ω 2 + ıε) e ,± ε = 0 in C ± ∂ ν e ,± ε = 0 on ∂C ± ∩ ∂Ω (1.21)
with the boundary conditions on

Γ ± 0 and Γ ± 1 e 0,± ε (ϕ) Γ ± 0 = ϕ and e 0,± ε (ϕ) Γ ± 1 = 0 e 1,± ε (ϕ) Γ ± 0 = 0 and e 1,± ε (ϕ) Γ ± 1 = ϕ. (1.22)
By identifying all the cells C ± j to C ± , using (1.20), we have that .23) In order to determine P ± ε , it suffices to note that the function defined on each cell C ± j by (1. 23) is solution of (1.19) in Ω ± only if its normal derivative across each Γ ± j is continuous. By injectivity of P ± ε , this will be satisfied if and only if it is continuous across Γ ± 1 . This corresponds to

u ± ε (ϕ) C ± j = e 0,± ε ((P ± ε ) j-1 ϕ) + e 1,± ε ((P ± ε ) j ϕ). ( 1
∂ z e 0,± ε (ϕ) Γ + 1 + ∂ z e 1,± ε (P + ε ϕ) Γ + 1 = ∂ z e 0,± ε (P + ε ϕ) Γ + 0 + ∂ z e 1,± ε ((P + ε ) 2 ϕ) Γ + 0 .
By defining the local DtN operators for , k ∈ {0, 1},

T k,± ε ∈ L(H 1 2 (Γ ± ), H -1 2 (Γ ± ))
T k,± ε ϕ = ∓(-1) k ∂ z e ,± ε (ϕ) Γ ± k (1.24)
we obtain that P ± ε satisfies the stationary Riccati equation

T 10,± ε (P ± ε ) 2 + T 00,± ε + T 11,± ε P ± ε + T 01,± ε = 0 (1.25)
This equation has an infinity of solutions but there is only one solution whose spectral radius is strictly less than 1. This can be shown by well-posedness of the half-guide problem (1.19). This characterizes uniquely the propagation operator P ± ε . The solution u ± ε (ϕ) of the half-guide problem (1.19) can then be constructed cell by cell and in particular we have

Λ ± ε = T 00,± ε + T 10,± ε P ± ε . (1.26)
We can solve (1.17) by using the following algorithm.

Chapter 1. Scattering problem in waveguides

• Compute the DtN operators Λ ± :

solve the two cell problems (1.21)- (1.22);

-compute the local DtN operators T k,± ε , , k ∈ {0, 1} defined in (1.24);
compute the unique solution P ± ε of spectral radius strictly less than 1 of the stationary Riccati equation (1.25);

compute the DtN operators by (1.26).

• Solve the coercive problem

-u a ε -n 2 (ω 2 + ıε) u ε = f in Ω a ∂ ν u a ε = 0 on ∂Ω a ∩ ∂Ω ± ∂u a ε ∂z + Λ ± ε u a ε = 0 on Γ ± .
• The solution of (1.17) is given by

u ε Ω a = u a ε u ε Ω ± = u ± ε (ϕ ± ), where ϕ ± = u a ε Γ ±
where u ± ε (ϕ ± ) is computed cell by cell as in (1.23) by using the associated solutions of cell problems and the propagation operators P ± ε .

From a numerical point of view, all the steps can be handled classically except the solution of the Riccati equation. To solve it, we use a modified Newton algorithm (modified in order to take into account the condition on the spectral radius) or a spectral method. This last method is based on the following result

ϕ = 0, P ± ε ϕ = pϕ ⇔ ϕ ∈ Ker(p 2 T 10,± ε + p T 00,± ε + T 11,± ε + T 01,± ε ) and |p| < 1
This characterizes the eigenvalues and the eigenvectors of P ± ε and the Jordan blocks can be determined in a similar way (in [START_REF] Hohage | Riesz bases and jordan form of the translation operator in semi-infinite periodic waveguides[END_REF], the authors have shown that P ± ε has a Jordan form.). Moreover, we can show that if p is solution of the quadratic eigenvalue problem then 1/p too. It suffices then to solve the quadratic eigenvalue problem, couple the solutions by pair (p, 1/p) and select the eigenvalue of modulus strictly less than one (see Figure 1.6). We represent a solution which was computed using this algorithm in Figure 1.7.

Transparent boundary conditions in absence of dissipation

We are interested now in defining and computing the physical solution of

-u -n 2 ω 2 u = f in Ω ∂ ν u = 0 on ∂Ω.
(1.27)

A natural extension of the homogeneous acoustic waveguide (described in Section 1.2) which is based on the limiting absorption principle would be 1. Limiting absorption principle for the periodic half guide problems : prove that for any ϕ, u ± ε (ϕ) has a limit in H 1 loc (Ω ± ) that we denote u ± (ϕ) and deduce that the DtN operators Λ ± ε has a limit in operator norm when ε tends to 0 that we denote in the following Λ ± ;

2. Characterize for any ϕ, u ± (ϕ) and Λ ± ϕ in order to compute them;

3. Show that the problem

-u a -n 2 ω 2 u a = f in Ω a , ∂ ν u a = 0 on ∂Ω a ∩ ∂Ω, ± ∂u a ∂z + Λ ± u a = 0 on Γ ± ,
enters in the framework of the Fredholm alternative;

4. Show that the problem (1.27) has at most one solution ;

5. Show that u a ε tends to u a in H 1 (Ω a ) when ε tends to 0;

Chapter 1. Scattering problem in waveguides 6. The outgoing solution is then given by

u Ω a = u a u Ω ± = u ± (ϕ ± ), where ϕ ± = u a Γ ± .
Of course, as in the homogeneous case, this would be valid except for at most a countable set of frequencies.

However, each step raises several difficulties that we describe in the sequel. We mention also how to overcome these difficulties from a theoretical and a numerical point of view. This is the subject of a paper that I am finishing with Vincent Lescarret (Supelec).

Limiting absorption principle for the periodic half guide problems . First, let me mention that it may exist frequencies ω for which Problem

-u ± -n 2 ω 2 u ± = 0 in Ω ± , ∂ ν u ± = 0 on ∂Ω ± ∩ ∂Ω, u ± = 0 on Γ ± ,
admits a non trivial solution in H 1 (Ω). They correspond to edge resonances (with the associated edge mode). The existence and the value of these edge resonances depend on the periodic medium and the position of the boundary. Moreover it is possible to show that there exists only at most a countable set of such frequencies. When the periodicity cell is symmetric with respect to the axis z = L ± /2, there is no edge resonances. But if such resonance exist, we cannot hope defining the corresponding DtN operators at these frequencies so we have to exclude them or at these fixed frequencies, we could move the boundary Γ ± whose position is artificial. Finally, another alternative to avoid these artificial forbidden frequencies is to construct Robin-to-Robin (instead of Dirichlet-to-Neumann) transparent boundary conditions that have the form

(± ∂u ∂z + ıαu) + Λ ± (∓ ∂u ∂z + ıαu) = 0.
with typically α = ω. The corresponding half-guide problem is of the same type than before with the Dirichlet boundary condition on Γ ± replaced by a Robin boundary conditions. In that case, of course, edge resonances associated to real frequencies cannot exist. Let us just emphasize that the construction of these operators is a little bit more technical (see [3,[START_REF] Fliss | Robin-to-robin transparent boundary conditions for the computation of guided modes in photonic crystal wave-guides[END_REF]) than the construction of the DtN operators explained in the previous section.

When the periodicity cell is symmetric with respect to the axis z = L ± /2, the limiting absorption principle for the halfguide problem can be deduced from the one of a guide problem (the symmetrized one) which is perfectly periodic. In the case on a non symmetric periodicity cell, one can use the result of [Hoang, 2011]. In both cases, one has to exclude a countable set of frequencies. For any ϕ ∈ H 1/2 (Γ ± ), the limit function is solution in

H 1 loc (Ω ± , ) of -u ± -n 2 ω 2 u ± = 0 in Ω ± ∂ ν u ± = 0 on ∂Ω ± ∩ ∂Ω u ± = ϕ on Γ ± .
(1.28)

Periodic waveguides

It is by definition the outgoing solution of this half-guide problem. The associated DtN operator Λ ± is the limit in operator norm of the DtN operators Λ ± ε when ε tends to 0.

Characterization of the limit solution of the half-guide Let us suppose now that the limiting absorption principle holds. The matter now is to characterize the outgoing solution of the half-guide problem and the associated DtN operator in order to compute them. Of course, the propagation operators P ± ε has a limit in operator norm which is defined by

∀ϕ ∈ H 1/2 (Γ ± ), P ± ϕ = u ± (ϕ) Γ ± 1 .
It is then easy to show that the spectral radius of the propagation operator P ± is less or equal to 1. How to characterize uniquely this operator?

For any data ϕ in H 1/2 (Γ ± ), the cell solution e ,± ε (ϕ) has a limit in H 1 (C ± ) if and only if ω 2 is not an eigenvalue of the Dirichlet cell problem. More precisely, except for these resonances, e ,± ε (ϕ) are differentiable in the neighborhood of ε:

e ,± ε (ϕ) = e ,± (ϕ) + ε e ,± (1)
(ϕ) + O H 1 (ε 2 )
where e ,± (ϕ) is solution of the Dirichlet cell problem (1.21-1.22) with ε = 0 and e ,± (1) (ϕ) is solution of a similar problem with homogeneous Dirichlet boundary conditions on Γ ± 0 and Γ ± 1 and with a term source which is ın 2 p e ,± (ϕ).

Except for the resonances of the Dirichlet cell problem, the local DtN operators T k,± ε are differentiable in operator norm in the neighborhood of ε:

T k,± ε = T k,± + ε T k,± (1) + O(ε 2 ) where T k,± ϕ = ∓(-1) k ∂ z e ,± (ϕ) Γ ± k and T k,± (1) ϕ = ∓(-1) k ∂ z e ,± (1)
(ϕ) Γ ± k .
Let us remark that we have to exclude the countable set of resonances of the Dirichlet cell problem to define these limits. This is as artificial as in the first step and this can be avoided by changing the periodicity cell -which corresponds to move the artificial boundaries on which the DtN conditions is imposed-or by introducing cell problem with Robin boundary condition instead of Dirichlet boundary conditions. It is a little bit more technical so we have privileged the simplicity for the presentation of the method. In the sequel, we suppose that the Dirichlet cell problem is well posed.

In consequence, the limit propagation operator P ± is solution of the stationary Riccati equation

T 10,± (P ± ) 2 + T 00,± + T 11,± P ± + T 01,± = 0 (1.29)
The difficulty now is that there does not exist necessarily only one operator solution of (1.29) whose spectral radius is less or equal to 1. On one hand, if ω 2 is not in the spectrum of the operator with periodic coefficients associated to the half periodic waveguide Ω ± , there exists a unique solution of spectral radius less or equal to 1 and its spectral radius is strictly less than one. Moreover, the solution of the half-guide problem u ± (ϕ) is in H 1 . On the other hand, if ω 2 is in the spectrum, there is no uniqueness. This is linked to the non uniqueness of the solutions in H 1 loc (Ω) of (1.27). This is illustrated in Figure 1.6 : some solutions (which can be also coupled by pairs (p, 1/p)) of the corresponding quadratic eigenvalue problem tend to the unit circle when ε tends to 0. At the limit ε = 0, it becomes impossible to select the eigenvalue of the outgoing propagation operator without knowing the behaviour of these eigenvalues in the neighborhood of ε. We have then to add some conditions to uniquely characterize the "outgoing" propagation operator. Let me give here the principle. Let p be an eigenvalue of the outgoing propagation operator P ± of modulus 1 then it is the limit of p ε , where p ε is an eigenvalue of

P ± ε . As ε → |p ε | is a decaying function at ε = 0, we obtain that Re dp ε dε ε=0 p < 0 (1.30)
where the derivative of p ε at ε = 0 can be computed thanks to p, the local DtN operators T k,± and T k,±

(it suffices to differentiate the quadratic eigenvalue problem with respect to ε). This condition is equivalent to the fact that the derivative of p ε at ε = 0, which is a complex, points inside the unit circle. The property (1.30) is not satisfied by 1/p (see Figure 1.8). Moreover, it can be shown that the associated eigenvector is linked to the Floquet mode of the periodic waveguide with positive group velocity for Ω + and negative group velocity for Ω -. When these conditions for the eigenvalues of modulus 1 are added to the Ricatti equation, we can show that this uniquely defines the outgoing propagation operator.

When the outgoing propagation operator is computed, then the outgoing solution of the halfguide problem and the associated DtN operator can be constructed.

The problem in the bounded domain Ω a For Problem (1.28), does the Fredholm alternative hold? It is less obvious than for the homogeneous case for which it suffices to use the explicit expression of the DtN operators and show that it is a positive operator. In the case of a periodic media that is homogeneous in a small neighborhood of Γ ± , we are able to show that the operators can be decomposed as the sum of a positive operator (the DtN operators of an homogeneous half-waveguide) and a compact operator. This is enough to show that Fredholm alternative holds.

Remark 1.3.5

There is another way to propose a problem in a bounded domain for which stability can be proven. It suffices to introduce DtN operators with overlap, i.e. defined thanks to two boundaries. For instance,

∀ϕ ∈ H 1/2 (Γ ± ), Λ± ϕ = ∓ ∂ ∂z u ± (ϕ) Γ ± 1 = Λ ± P ±
This overlap makes the associated operator always compact. The associated problem defined in

Ω a+1 = Ω ∩ {|z| < a + 1} given by            -u a+1 -n 2 ω 2 u a+1 = f in Ω a+1 , ∂ ν u a+1 = 0 on ∂Ω a+1 ∩ ∂Ω, ± ∂u a+1 ∂z + Λ ± u a+1 Γ ± = 0 on Γ ± 1 is naturally of Fredholm type.
For the uniqueness, it is possible to show that the DtN conditions are equivalent to the radiation conditions (1.16), for which uniqueness can be proven, except for a countable set of frequencies corresponding to the resonances of the problem. Finally, using classical arguments, we can show that u a ε tends to u a in H 1 (Ω a ).

Except for a countable set of frequencies containing the set of the cut-off frequencies σ 0 and the resonances of the problem, we were able to characterize uniquely the outgoing solution of the Problem. This characterization leads naturally to a numerical method to compute this solution.

This algorithm can be also used to solve scattering problems which are defined as follows.

Let u be the total field satisfying the problem

-u -n 2 ω 2 u = 0 in Ω ∂ ν u = 0 on ∂Ω u -u i is outgoing,
where Ω is as described at the beginning of Section 1.3.1 (see for an example Figure 1.3 (right)) and u i is the incident field, i.e. an ingoing floquet mode of one of the half guide (see the end of Section 1.3.2 for a precise definition). For the incident field represented in Figure 1.9 (middle figure), we have computed and represented the total field (bottom figure).

Application to the transmission problem between periodic half-spaces

Here we consider, the transmission problem between two periodic halfspaces of commensurate period. This is an extension of the more classical diffraction problem by a periodic grating or layer which has been the subject of a huge amount of papers since the beginning of the 80's (see for instance [Petit, 1980]). In that case, the complete mathematical analysis of this problem was done using variational techniques (see for instance From a numerical point of view, for incident plane waves, the problem can be restricted to a band problem with quasi periodicity conditions and can then be solved using transparent boundary conditions in homogeneous waveguides [Abboud, 1993,Bao, 1997] or integral equation techniques [START_REF] Meier | A Nyström method for a class of integral equations on the real line with applications to scattering by diffraction gratings and rough surfaces[END_REF],Arens et al., 2006]. The case of more general sources and the presence of perturbations have been handled more recently by using the Floquet-Bloch Transform in [Lechleiter, 2017]. In presence of guided waves, the study is much more intricate since appropriate radiation conditions (which are naturally not the same in the longitudinal and transverse directions) have to be added in order to make the problem well-posed. See [Bonnet-Ben Dhia and Tillequin, 2001a, Bonnet-Ben Dhia and Tillequin, 2001b, Bonnet-Ben Dhia et al., 2000, Bonnet-Ben Dhia and Ramdani, 2002] for stratified media and more recently [START_REF] Kirsch | The limiting absorption principle and a radiation condition for the scattering by a periodic layer[END_REF] for periodic layer.

This work is an extension in the sense that the media in each side of the interface can be periodic in the two directions. The restriction is that the periodicity in the direction of the interface has to be the same (or commensurable of course), see Figure 1.10. We consider the case of absence of guided waves by the interface. We can then apply the Floquet Bloch Transform in the direction of the interface in a classical framework. The problem is then equivalent to a family of independent waveguide problems with quasi-periodic boundary conditions, the quasi-period being the dual Floquet variable. It suffices then to use the method to solve waveguide problem described in the previous sections. The Floquet-Bloch inverse allows to recover the solution in the whole domain. The discretization of the dual Floquet variable is studied in [Coatléven, 2012, Lechleiter, 2017].

In collaboration with physicists from l'Institut d'Electronique Fondamental of Orsay University, we have made some computations in realistic periodic media to illustrate negative refraction and superlens phenomena of some photonic crystals. The description of the numerical method and the numerical results were the topic of the paper [START_REF] Fliss | Computation of light refraction at the surface of a photonic crystal using dtn approach[END_REF].

A modal transparent boundary conditions for general waveguide problems

The method of construction of transparent boundary conditions which I have explained in the previous sections can be extended easily to vectorial equations (Maxwell's or elasticity for instance) for locally perturbed, homogeneous or stratified, isotropic or anisotropic media. This construction is not based on a modal decomposition and in some applications, it is important to keep the decomposition in terms of the modes of the guides, as it is done in Section 1.2. The extension of the results of Section 1.2 for isotropic acoustic waveguide is really intrincate for isotropic elastic waveguide (see for instance [START_REF] Baronian | Couplage des méthodes modale et éléments finis pour la diffraction des ondes élastiques guidées[END_REF], Baronian et al., 2010]) and impossible for anisotropic elastic waveguide. The difficulty comes from the fact that the modal amplitudes (corresponding to the coefficients involving in the decomposition of the outgoing solution in terms of the modes) can be directly linked to the solution in the isotropic acoustic case (see (1.7)) because the modes are orthogonal in any transverse sections. This magical property is not true for general waveguide problem. In other words, constructing Dirichlet-to-Neumann operators fail in general.

By working on elastic waveguides, we have constructed new transparent boundary conditions linked to what we have called the Poynting-to-Neumann operator because it is linked to the energy flux and the Poynting vector. For pedagogical purposes, I begin by explaining the conditions on the simple case of the isotropic acoustic waveguide. Then I explain how they can be extended to elastic waveguides and also to periodic waveguides, in order to show the general nature of these conditions.

The Poynting-to-Neumann map for isotropic homogeneous acoustic waveguides

Let me here consider a diffraction problem in an acoustic isotropic half-guide Ω = S×] -a, +∞[ where S ⊂ R 2 denotes the bounded cross-section of the guide. We look for the outgoing solution u of

-∆u -ω 2 u = f in Ω, ∂ ν u = 0 on ∂Ω,
where ω is the frequency, ν is the exterior normal to ∂Ω and the source term f is supposed to be compactly supported in {z < 0}. We denote by u a (resp. u + ) the restriction of u to the subdomain Ω a = Ω ∩ {z < a} (resp. Ω + = Ω ∩ {z > 0}) and we want to derive transparent boundary conditions for u a on Γ a = {(x s , z) ∈ Ω, z = a}. In Section 1.3.1, we have considered the case where Ω a and Ω + do not overlap. Here we consider the general case in the sense that a ≥ 0 and for a > 0, the two domains overlap. Looking for the outgoing solution, u + admits the following expression

u + (x s , z) = k≥0 a + k w k (x s , z) (1.31)
involving the right-going modes w k (x s , z) = ϕ k (x s )e ıβ k z . We recall the reader (see Section 1.3.1) that a finite number (N) of them are propagative (Im(β k ) = 0 and Re(β k ) > 0) and the rest are evanescent (Im(β k ) > 0 and Re(β k ) = 0). Here, the sequence {ϕ k , k ∈ N} forms an orthonormal basis of L 2 (S). The a + k are the unknown modal amplitudes.

Imposing the following matching conditions on the boundaries Γ 0 and Γ a :

u a Γ 0 = u + | Γ 0 and ∂ z u a | Γa = ∂ z u + | Γa ;
using the formula (1.31) and the orthogonality of the ϕ k , one can derive another transparent condition for u a -except for a countable set of frequencies-involving a Dirichlet-to-Neumann operator with overlap Λ DtN,a :

Λ DtN,a u a = ∂ z u a Γa = k≥0 iβ k (u a , ϕ k ) Γ 0 ϕ k e ıβ k a .
This last operator is compact thanks to the exponentially decaying factors e ıβ k a . And it is easy to show that the problem satisfied by u a is of Fredholm type. The countable set of frequencies that we have to exclude are linked to the resonances of the overlap box Ω a ∩ Ω + . To avoid them, it suffices to replace the matching conditions on Γ a by a Robin type matching condition.

This cannot be extended to stratified, anisotropic or periodic waveguides or in (even isotropic) elastic or electromagnetic waveguides principally for two reasons :

• the modal decomposition (1.31) is obtained by using a separation of variables technique which is possible because the transverse operator is self-adjoint. In general, the transverse operator has no meaning (in the periodic case) or is not self-adjoint (in the elastic case).

• the way to obtain the modal amplitudes thanks to the trace of u + requires the orthogonality of the modes in L 2 (S), which does not hold for instance in stratified, anisotropic or periodic waveguides or in (even isotropic) elastic or electromagnetic waveguides.

Let us explain now how to derive the modal amplitudes using a more general framework.

Using the expression of w k , it is easy to see that q(w j , w k ) = 0 if j = k q(w j , w j ) = 0 for evanescent waves 2ıβ j for propagative waves (1.32) where q is a (energy flux) sesquilinear form defined by

∀u, v ∈ H 2 loc , q(u, v) = Γ ∂ z u v -u ∂ z v.
An important property of q is that if u and v are two different modes then, by Green's formulas, q(u, v) = 0. This gives a bi-orthogonality property for the modes which is extendable to more Chapter 1. Scattering problem in waveguides general waveguide since it can be proven using only Green's formula.

We deduce that, except for the cut-off frequencies (the frequencies for which one β k vanishes), u + is given by

u + = k≤N q(u + , w k ) 2ıβ k w k + u evan .
where u evan is exponentially decaying at +∞.

Imposing that only the propagative modal amplitudes of u + and u a match as well as their normal derivative on Γ a leads to introduce the so-called Poynting-to-Neumann (PtN) operator T P tN defined by

T P tN u a = k≤N q(u a , w k ) 2ıβ k ∂ z w k Γ a .
Using the properties of q, we can show that Theorem 1.4.1

The operator T P tN is of finite rank from V to H -1/2 (Σ ) where

V = {u ∈ H 2 (Ω a \ Ω 0 ), ∆u + ω 2 u = 0}. and ∀u ∈ V, Σ u T P tN u -u T P tN u ∈ ıR + .
The first property gives that the problem

∆ũ a + ω 2 ũa = f in Ω a , ∂ ν ũa = 0 on ∂Ω a ∩ Ω, ∂ z ũa = T P tN ũa on Γ a (1.33)
is coercive + compact. From the second property, we can prove that this problem has at most one solution except for a countable set of frequencies corresponding to the problem in Ω a with Neumann boundary conditions (to avoid them it suffices to construct a Poynting-to-Robin operator instead of the Poynting-to-Neumann one). This problem is however not equivalent to the initial one since we have neglected the evanescent part in u a but by stability of (1.33) (except for a countable set of frequencies), we can show that .34) In Figure 1.13, I compare a reference solution computed using DtN operator with the solution computed using PtN transparent boundary conditions at a frequency for which there are 5 propagative modes. The difference between the two solutions is linked to the first evanescent modes which is neglected in the construction of the PtN operator. Of course, as indicated in (1.34) if ω is close to a cut-off frequency (i.e. ω < β N +1 ) or if the straight part of the waveguide is small (i.e. a is small), then the approximation is poor. A natural idea is to consider in the construction of the PtN operator the first evanescent modes. The difficulty is that as indicated in (1.32) the energy flux of the evanescent modes vanishes q(w k , w k ) = 0 for k ≥ N + 1. But if we consider the associated left going modes w - k (x s , z) = ϕ k (x s )e -ıβ k z which exponentially grows at +∞, we find q(w k , w -

∀ã < a, ∃C > 0, u a -ũa ≤ C e -Im(β N +1 )ã . ( 1
k ) = 2ıβ k ∈ R
which vanishes at the associated cut-off frequency. So one can consider as many evanescent modes as for the DtN operator and define the complete PtN operator TP tN defined by

TP tN u a = k≤N q(u a , w k ) 2ıβ k ∂ z w k Γ a + k>N q(u a , w - k ) 2ıβ k ∂ z w k Γ a .
The function u a is then solution of

∆u a + ω 2 u a = f in Ω a , ∂ ν u a = 0 on ∂Ω a ∩ Ω, ∂ z u a = TP tN u a on Γ a ,
this problem being of Fredholm type and well-posed except for a countable set of frequencies corresponding to the resonances of the problem.

As we are going to illustrate now, this PtN operator can be constructed for more general waveguide problems, as it is based on general properties of the modes (linked to their energy flux). Indeed, by definition the rightgoing modes w k 's are solutions of the homogeneous equations far from the perturbations whose group velocity is strictly positive. It can be proven that they verify

q(w j , w k ) = 0 if j = k, q(w j , w j ) = ıλ j with λ j ≥ 0
where q is a (energy flux) sesquilinear form derived from the Green's formula associated to the problem and λ j corresponds to the group velocity of the mode w j . The most delicate work is to show that the outgoing solution can be decomposed as a linear combination of a finite set of propagative rightgoing modes up to an exponentially decaying function at +∞, as in (1.31). This has to be investigated case by case. Let us however mention that the tools introduced in [START_REF] Nazarov | Elliptic problems in domains with piecewise smooth boundaries[END_REF], Nazarov, 2013, Nazarov, 2014b] and based on the use of an appropriate transformation in the direction of the guide (Fourier for homogeneous or stratified and Floquet for periodic) in weighted Sobolev spaces seem really relevant and general.

As soon as one can show that the associated u + (which can be vectorial) can be decomposed as

u + = k≤N a + k w k + u evan .
the PtN operator can be derived in the same way than for the isotropic acoustic problem and the properties of the operator remain the same.

Remark 1.4.1

Let us emphasize that this method requires the a priori knowledge or computation of the modes.

Let us explain more precisely the extension to periodic waveguides and anisotropic elastic waveguides.

Construction of the PtN operator for periodic waveguides

We have already introduced all the tools in Section 1.3.2. Let us consider to simplify

Ω ⊂ S×] - A, +∞[ where S ⊂ R 2 is bounded, Ω being periodic for z > 0, i.e. Ω ∩ {z > 0} = Ω p ∩ {z > 0}
where Ω p is a L-periodic domain. We look for the outgoing solution u of

-∆u -ω 2 n 2 u = f in Ω, ∂ ν u = 0 on ∂Ω, (1.35)
where n(z, x s ) = n p (z, x s ) is a L-periodic function for z > 0.

The modes are particular solution of

-∆u -ω 2 n 2 p u = 0 in Ω p , ∂ ν u = 0 on ∂Ω p , of the form w(z, x s ) = ϕ(z, x s ) e ıβz .
where ϕ is L-periodic in the z-direction. This leads to solve a quadratic eigenvalue problem for (β, ϕ):

-(∇ -ıβe z ) 2 ϕ -ω 2 n 2 p ϕ = 0 in C ∂ ν ϕ = 0 on ∂C ∩ ∂Ω p ,
where C is the periodicity cell of Ω p . From the Steinberg'a analytic Fredholm theorem [Steinberg, 1968] (see also [START_REF] Reed | [END_REF]), we conclude that there exists a countable set of modes which can be classified as follows: (a finite number) of the propagative modes (Im (β)=0), the evanescent ones (Im (β)>0) and the exponentially growing ones. Moreover, one can show easily that if (β, ϕ) defines a propagative modes, (-β, ϕ) as well. The propagative modes are given by β ∈ {Ξ n (ω), n ∈ I(ω)} where Ξ n (ω) and I(ω) are defined in (1.13)- (1.14) and ϕ is the eigenvector of A p (β) associated to the eigenvalue λ n (β) (see (1.12) for the definition). Finally, let us introduce q the (energy flux) sesquilinear form defined by

∀u, v ∈ H 2 loc , q(u, v) = Γ ∂ z uv -u∂ z v.
We can show using Green's formula that if w 1 = e ıβ 1 z ϕ 1 (x s , z) and w 2 = e ıβ 2 z ϕ 2 (x s , z) are two modes such that β 1 -β 2 = 0, they satisfy the bi-orthogonality condition

q(w 1 , w 2 ) = 0,
(the proof is straightforward) and if w = e ıβz ϕ(z, x s ) is a propagative mode then q(w, w) = ıλ n (β)

(here this is not straightforward). We have then a simple relation between the energy flux and the group velocity of the modes. The rightgoing propagative modes w + k are the ones for which the imaginary part of the energy flux (equivalently the group velocity) is positive. The set of cut-off frequencies (the frequencies for which the energy flux of one propagative mode vanishes) which has to be excluded, corresponds to the set σ 0 .

Finally, in [Nazarov, 2014b, Hoang, 2011] and by adapting Section 1.3.4, it is shown that there exists a unique solution u of (1.35), called outgoing solution, such that u

+ = u| Ω + (Ω + = Ω ∩ {z > 0}
) is a linear combination of the finite set of the rightgoing modes up to an exponentially decaying function at +∞

u + = k≤N a + k w + k + u evan .
From all these results, we can construct the associated PtN operator

T P tN u a = k≤N q(u a , w + k ) q(w + k , w + k ) ∂ z w + k Γ a .
It is possible, as in the homogeneous case, to add the first evanescent modes.

The properties of the PtN operator are the same than the one given in Theorem 1.4.1 and the corresponding problem with PtN transparent boundary conditions are of the same type than the one of the homogeneous case (1.33).

The PtN operator for anisotropic elastic waveguide

These results are part of the paper [START_REF] Vahan Baronian | Iterative methods for scattering problems in isotropic or anisotropic elastic waveguides[END_REF]. Let us now consider a diffraction problem in an elastic half-guide Ω = S×] -a, +∞[ where S ⊂ R 2 denotes the bounded cross-section of the guide. The density of the material denoted by ρ is supposed to be independent of z for z > 0, but it can depend on the transverse variable x s . We suppose that the waveguide has a stress-free boundary ∂S × R. In time harmonic regime (of pulsation ω > 0), the propagation in the waveguide is modeled by the following equations

-div σ(u) -ω 2 ρ u = 0 in Ω, σ(u) • ν = 0 on ∂Ω, (1.36)
where u represents the displacement field (u = (u x , u z ) in the 2D case and u = (u x , u y , u z ) in the 3D case) and σ(u) the stress tensor which is related to the strain tensor ε(u) = 1/2 (∇u + ∇ T u)(in the small deformation assumption) through the stiffness tensor C by the Hooke's law

σ(u) ij = C ijkl ε kl (u) with i, j, k, l ∈ {x, y, z} in 3D
{x, z} in 2D.

(1.37)

The stiffness tensor is supposed to be independent of z for z > 0, but it can depend on x s . The notion of outgoing solution is defined below.

The modes are particular solutions of

-div σ(w) -ω 2 ρ w = 0 in S × R, σ(w) • ν = 0 on ∂S × R, of the form w(z, x s ) = ϕ(x s ) e ıβz .
This leads to solve a quadratic eigenvalue problem for (β, ϕ):

-div β σ β (ϕ) -ω 2 ρ ϕ = 0 in S, σ β (ϕ) • ν = 0 on ∂S,
where div β and σ β correspond to the operators div and σ replacing ∂ z by ıβ. From the Steinberg'a analytic Fredholm theorem [Steinberg, 1968] (see also [START_REF] Reed | [END_REF]), we conclude that there exists a countable set of modes which can be classified as follows: (a finite number) of the propagative modes (Im(β)=0), the evanescent ones (Im(β)>0) and the exponentially growing ones.

For orthotropic media, the modes satisfy a bi-orthogonality relation, known as Fraser biorthogonality relation (see [Fraser, 1976, Gregory, 1983]), which was used in [START_REF] Pagneux | Lamb wave propagation in inhomogeneous elastic waveguides[END_REF], Pagneux and Maurel, 2004, Pagneux and Maurel, 2006] for 2D waveguides and extended in [START_REF] Baronian | Couplage des méthodes modale et éléments finis pour la diffraction des ondes élastiques guidées[END_REF], Baronian et al., 2010] for 3D waveguides to construct transparent boundary conditions. In the case of a general anisotropy, Fraser's relation does not hold anymore. However, the general bi-orthogonality relations, which rely on simple Green's formulas, hold (see [Auld, 1973]). Let us introduce the sesquilinear form q defined by

∀ u, v, q(u, v) = Γ u S v -S u v 1.

A modal transparent boundary conditions for general waveguide problems

where S u denotes the normal stress associated to u:

S u = σ (u) • e z .
And we have that if w 1 (z, x s ) = e ıβ 1 z ϕ 1 (z) and w 2 (z, x s ) = e ıβ 2 z ϕ 2 (z) are two modes such that β 1 -β 2 = 0, they satisfy the biorthogonality condition q(w 1 , w 2 ) = 0;

and if w(z, x s ) = e ıβz ϕ(z) is a propagative mode then

q(w, w) = ı J β .
It is easy to show that if J β ≥ 0 then J -β ≤ 0. The cut-off frequencies correspond to the frequencies for which one J β vanishes. The rightgoing (resp. leftgoing) modes w + k are the ones for which the energy flux is strictly positive (resp. negative). Finally, in [Nazarov, 2013], it is shown that there exists a unique solution of (1.36), called the outgoing solution, such that

u + = u| Ω + (Ω + = S × R +
) is a linear combination of the finite set of the rightgoing modes up to an exponentially decaying function at +∞

u + = k≤N a + k w + k + u evan
From all these results, we can construct the associated PtN operator

T P tN u a = k≤N q(u a , w + k ) q(w + k , w + k ) ∂ z w + k Γ a .
To illustrate these new conditions, let me present numerical simulations done by Antoine Tonnoir (2D) and Vahan Baronian (3D). We have solved problem (1.36) for artificial 2D and 3D unperturbed anisotropic waveguides.

1. 2D example : S = (0, 1) and Ω = S × R, the frequency is such that there exists 4

propagative modes, we have computed the solution for two sizes of the computational area Ω a . Figure 1.14 represents the modulus of the displacement field obtained in these two cases. We also represent the whole solution reconstructed in the half-guides Ω ± thanks to the modal expansion.

2. 3D example : S is a rectangle and Ω = S × R, the frequency is such that there exists 24 propagative modes. Figure 1.15 represents the modulus of the displacement field three different sizes of the computational domain.

As we can see in the 2D and 3D configurations, the restrictions of the solutions to the smallest domain match very well. It is the first time that a transparent boundary condition for anisotropic waveguides is proposed.

Let me finally emphasize that these new conditions have also some numerical advantages for isotropic elastic waveguide for which other methods can be used. Let me first recall that the use of transparent boundary conditions leads to a partially dense linear system since the conditions are non local. Such system can be difficult to invert directly, in particular for vectorial equations and 3D domains. A natural alternative for the inversion of such system is the use of iterative algorithms such as preconditioned Krylov methods like GMRES [START_REF] Saad | Gmres: A generalized minimal residual algorithm for solving nonsymmetric linear systems[END_REF]. Roughly speaking, the preconditioner is chosen as a sparse part of the complete system and the dense part of the matrix coming from the non local transparent boundary condition is involved only in the matrix vector product step. The choice of the transparent boundary conditions has an influence on the convergence of the associated iterative algorithm. Comparing to other transparent boundary conditions constructed even with an overlap, this transparent boundary conditions based on the Poynting-to-Neumann operator is much more efficient.

Ongoing works and perspectives

Many theoretical questions remain open today. But it is worth continuing working on the transparent boundary conditions (based on DtN or PtN operators) that I have described in this chapter since they seem to work, at least from a numerical point of view where all existing methods fail in environments of such complexity. Moreover, they allow us to approach from another angle the theoretical questions related to the definition of the physical solution.

For the periodic waveguides, apart of finishing the paper whose content is summarized in Section 1.3.4, the precise analysis of the problem at the cutoff frequencies, will allow us to justify theoretically the numerical method described in Section 1.3.5 for the transmission problem between two periodic media which have the same period in the direction of the interface. Moreover, I want to study such transmission problem in presence of local perturbations and with possible presence of guided waves by the interface. The difficulties are theoretical (derive a radiation condition which will be different in the direction of the interface and the other directions and show well-posedness of the problem) and numerical (compute the solution by restricting the problem around the perturbation). One (long-term) application of this work is the computation of topologically protected edge states (See Chapter 3, Section 3.4). Finally, even in presence of dissipation, we do not know how to study the transmission problem when the periods in both side of the interface are not commensurate. A particular case is the transmission between an homogeneous medium and a periodic medium in the case where the interface cuts the periodic medium in a direction for which the medium is not periodic. I want to extend the tools which were developed for the homogenization of this kind of problem [Gérard-Varet andMasmoudi, 2011, Gérard-Varet andMasmoudi, 2012] and derive an associated numerical method.

The method explained in Section 1.4 can be applied, I think, to construct transparent boundary conditions for general waveguides for which classical methods fail. I want to apply this method to (1) Maxwell's equations for anisotropic or periodic waveguides (with M. Kachanovska and E. Becache (POEMS)), [START_REF] Fliss | Exact boundary conditions for time-harmonic wave propagation in locally perturbed periodic media[END_REF] to Kirchoff Love equations for isotropic or anisotropic waveguides (with L. Chesnel (Inria Saclay, CMAP) and L. Bourgeois (POEMS), see [20] for the derivation of the radiation condition and the well-posedness of the Kirchoff Love equations in the time harmonic regime in 2d waveguides), (3) to dispersive waveguides. The most intrincate part is to show that the physical solution can be decomposed in terms of the outgoing propagative modes. This can be done, we think, by extending the framework of [START_REF] Nazarov | Elliptic problems in domains with piecewise smooth boundaries[END_REF]. [Goursaud, 2010] for instance) but again the method described in Section 1.4 is promising.

1.5 Ongoing works and perspectives 47

Introduction

In this chapter, we are interested in the diffraction of time-harmonic waves in a homogeneous 2D or 3D infinite anisotropic elastic or periodic medium. The difficulties are similar to the ones describes in the previous chapter about the waveguide problems : since the medium is infinite, there are theoretical difficulties -how to define the so called outgoing solution of such problem? -and numerical difficulties -can we introduce an equivalent formulation which is suitable for numerical purposes ? In contrast with waveguide configurations, the solution can propagate in all the directions of the medium and have in consequence a certain decay at infinity. In this chapter, we will focus on the numerical difficulties and we will avoid the theoretical ones by adding a small dissipation term to our model. The problem is then well-posed in the classical L 2 setting. I will explain in a dedicated section (see Section 2.3) the theoretical difficulties raised by the case without dissipation.

Solving time harmonic scalar waves equations in infinite homogeneous media is an old topic [Givoli, 1992] and there exist several methods. They are all based on the natural idea of reducing the pure numerical computations to a bounded domain containing the perturbations (achieved using for instance Finite Element methods). A first class of methods consists in applying an artificial boundary condition, around the bounded domain, which is transparent or approximately transparent as in: (1) integral equation techniques, (2) Dirichlet-to-Neumann approaches providing that the boundary is properly chosen to allow separation of variables and (3) local radiation conditions at finite distance constructed as local approximations at various order of the exact non local condition. These methods were first introduced for the time harmonic scalar wave equation -the Helmholtz equation -and then extended to 2D isotropic elasticity problems by simply using the Helmholtz decomposition of the displacement field in terms of potentials (see for instance [START_REF] Givoli | Non-reflecting boundary conditions for elastic waves[END_REF]). However it seems that all these methods 50 Chapter 2. Scattering problem in infinite complex media either do not extend to anisotropic elastic or periodic media -the separation of variables is not possible anymore to determine the Dirichlet-to-Neumann (DtN) operator and the Green's function is not known for periodic media -or do extend but with a tremendous computational cost -for the integral equation techniques, the Green tensor for anisotropic elastic media depends not only on the distance between two points but also on the orientation [START_REF] Wang | Three-dimensional timeharmonic elastodynamic green's functions for anisotropic solids[END_REF]. A second class of methods consists in surrounding the computational domain by a Perfectly Matched absorbing Layer (PML). PML techniques are very popular because they are efficient and easy to implement in a large class of problems. But they may be inoperant. Roughly speaking, the PML absorbs the wave with an outgoing phase velocity, preventing them to come back in the computational domain, while in order to catch the physical solution, it should absorb the waves with outgoing group velocities. That is why to our knowledge the standard PML technique works for isotropic elastic media (in which the waves with outgoing phase velocities have outgoing group velocities and vice versa) but may fail for general anisotropic elastic or periodic media where the two velocities differ [START_REF] Bécache | Stability of perfectly matched layers, group velocities and anisotropic waves[END_REF].

By contrast, our method is based on a simple and quite general idea: the solution of homogeneous -isotropic or anisotropic, acoustic or elastic -halfspace problems can be expressed thanks to its trace on the halfspace boundary. As several halfspaces surrounding the perturbations are needed to recover the whole domain, they will necessarily overlap. The second step is then to find conditions to ensure the compatibility of the representations in the overlapping zones. This method has links with domain decomposition methods with overlap [START_REF] Lions | On the schwarz alternating method. i[END_REF], Dryja and Widlund, 1994, Toselli and Widlund, 2005], with the specific difficulty that the overlapping zones are unbounded. More precisely, the idea in 2D is to split the whole domain into five parts (see Figure 2.1):

• a square that includes the defect (and all the inhomogeneities) in which we will use a Finite Elements representation of the solution,

• and 4 half-planes, parallel to the four edges of the square in which the medium is homogeneous.

Taking advantage of the medium properties in a half-plane, we can give a semi-explicit (integral) expression of the solution given (for instance) its trace on the edge of the half-plane, via the Fourier transform in the transverse direction in the homogeneous case or via the Floquet-Bloch Transform in the periodic case. With these integral representations and the Finite Element representation of the solution in the square, we can formulate a coupled problem. To ensure the compatibility of the different representations, as in domain decomposition methods, we impose transmission conditions on the edges of the subdomains. This leads us to a system of coupled equations where the unknowns are the solution in the bounded square and the traces of the solution on the edges of the half-planes.

When compared to absorbing layers methods, this approach is obviously more costly due to the additional unknowns (the traces) linked by non-local integral equations. However, one benefit is that this additional computation of the traces enables to reconstruct a posteriori the solution in the half-planes (and therefore in the whole domain), which is impossible for instance when using non exact absorbing boundary conditions or PML.

I first present in Section 2.2 the principles of the method and the main results for a simple model problem : a dissipative anisotropic Helmholtz equation. Let me underline that thought the mathematical analysis holds only for dissipative media, the method gives good numerical results also in the non dissipative case, as I will show in Section 2.3. I will explain in this section the theoretical difficulties raised by the case without dissipation. This method applies also to anisotropic elastic media and periodic media as I will show respectively in Sections 2.4 and 2.5. I will finally describe my ongoing work and my perspectives for this subject.

2.2

The Halfspace Matching Method (The HsMM) on a toy problem

The model problem

The model problem that I consider in this section is

-div(A ∇u) -ω 2 ε ρ u = f in Ω, (2.1)
in the time harmonic regime at the frequency Re(ω ε ) = ω with a small absorption Im(ω ε ) = ε > 0, where A is a symmetric positive definite matrix of (L ∞ (Ω)) 2×2 modeling the anisotropy and ρ is a strictly positive function of L ∞ (Ω).

The propagation domain Ω is typically R 2 , or R 2 minus a set of obstacles which are included in a bounded region

∃a > 0, ∂Ω ⊂ Ω a ≡ (-a, a) 2 .
In presence of obstacles, some boundary conditions have to be added to the model.

Remark 2.2.1

The principle of the method extends easily to 3D propagation domain that is infinite in the 3 directions.

The source term f is supposed to be a function of L 2 (Ω) with a compact support included in Ω a . Finally, the matrix A is a local perturbation of a constant matrix A 0

supp(A -A 0 ) ⊂ Ω a , where A 0 = c 1 c 3 c 3 c 2 with c 1 , c 2 > 0, c 1 c 2 -(c 3 ) 2 > 0,
and the function ρ is a local perturbation of a constant function, which is taken, without loss of generalities, equal to 1 supp(ρ -1) ⊂ Ω a .

For variational boundary conditions on ∂Ω -for instance Neumann or Dirichlet conditions-it is well known that thanks to the dissipation, this problem admits a unique solution in H 1 (Ω).

To clarify the presentation of the method, I will consider three situations of increasing difficulty.

(Case 1) The HsMM formulation. The propagation medium is Ω = R 2 \ Ω a and non homogeneneous Dirichlet boundary conditions are imposed on its boundary. The source term f = 0, A = A 0 = Id and ρ = 1 in Ω. Coupling analytical representations of the solution in the 4 halfspaces surrounding the obstacle O = Ω a and ensuring that all representations match, we end up with a system of integral equations whose unknowns are the 4 traces of the solution on the edges of the halfspaces. We show stability and well-posedness of this formulation in a framework which is suitable for numerical simulations.

(Case 2) Coupling the HsMM with FE method. We consider here the general isotropic case A 0 = Id, in presence of source terms and possible perturbations of the geometry, the matrix A and the coefficient ρ. The idea here is that any perturbation can be taken into account using a Finite Element (FE) method, as soon as it is contained in a bounded region. Here, we examine the coupling between the FE representation of the solution and the system of integral equations obtained in the previous case. We point out the importance of the presence of an overlap between the FE box and each halfspace. Thanks to this overlap, we show again stability properties and well-posedness for this problem.

(Case 3) Case of a general anisotropic case. We consider then the general anisotropic case.

I explain formally how to derive the new formulations for the 3 cases in the next subsections and I give without proof the associated stability results. For the proof, see [START_REF] Bonnet-Ben Dhia | The halfspace matching method: a new method to solve scattering problem in infinite media[END_REF].

Case 1: the HsMM for an exterior isotropic problem

We consider the following problem: let Ω = R 2 \ Ω a , g ∈ H 1/2 (∂Ω) and find the unique solution

u ∈ H 1 (Ω) of -∆u -ω 2 ε u = 0 in Ω, u = g on ∂Ω.
(2.

2)

The domain Ω is the union of 4 half-planes Ω j a that lie on the 4 edges of the square Ω a . Using the following local coordinates for all j ∈ 0, 3

x j y j = cos(θ j ) -sin(θ j ) sin(θ j ) cos(θ j ) x y
, where θ j = jπ 2 .

(2.3) the half-planes are defined as follows for all j ∈ 0, 3 Ω j a = {x j ≥ a} × {y j ∈ R} Σ j a ≡ ∂Ω j a := {x j = a} × {y j ∈ R}.

(2.4)

Finally, we denote Σ aa = ∂Ω a and Σ j aa = Σ aa ∩ Σ j a .

(2.5)

These notations are summarized on Figure 2.1.

As explained previously, the formulation uses the representation of the solution in each halfspace Ω j a surrounding the obstacle. Let us then introduce the halfspace problem : for any ψ ∈ H 1/2 (Σ j a ), find the unique solution

U j (ψ) ∈ H 1 (Ω j a ) -∆U j -ω 2 ε U j = 0 in Ω j a , U j = ψ on Σ j a (2.6)
Using the Fourier transform in the y j -direction, it is easy to see that, the solution of (2.6) is given by

∀(x j , y j ) ∈ [a, +∞[×R, U j (ψ)(x j , y j ) = 1 √ 2π R ψ(ξ)e ı √ ω 2
ε -ξ 2 (x j -a) e ıξy j dξ, (2.7)

where the square root is defined with the convention Im √ ≥ 0 and ψ is the Fourier transformation of ψ using the convention,

∀ ξ ∈ R, ψ(ξ) = 1 √ 2π R
ψ(y j ) e -ıξy j dy j .

Let us now derive the halfspace matching formulation which involves only the traces of the solution u on the edges Σ j a of Ω j a that we denote

ϕ j = u Σ j a .
(2.8)

We have thanks to the boundary conditions satisfied by u in (2.2)

ϕ j Σ j aa = g Σ j aa .
(2.9)

Moreover, because u is solution of (2.2), its restriction to Ω j a is solution of (2.6) with ψ = ϕ j so u Ω j a = U j (ϕ j ).

(2.10)

The trace of u on Ω 0 a ∩ Σ 1 a is given by ϕ 1 by (2.8) and by U 0 (ϕ 0 ) by (2.10) so we have necessarily

ϕ 1 Σ 1 a ∩Ω 0 a = U 0 (ϕ 0 ) Σ 1 a ∩Ω 0 a .

Chapter 2. Scattering problem in infinite complex media

Applying this reasoning to the other half lines Σ j±1 a ∩ Ω j a , we get

ϕ j±1 Σ j±1 a ∩Ω j a = U j (ϕ j ) Σ j±1 a ∩Ω j a , ∀j ∈ Z/4Z.
(2.11)

The system of coupled equations (2.11) and (2.9) constitutes the halfspace matching formulation for the problem (2.2). To understand the nature of this system of equations and to analyze it, let us introduce the operators

∀ψ ∈ H 1/2 (Σ j a ), D j j±1 ψ := U j (ψ) Σ j±1 a ∩Ω j a .
(2.12)

By classical trace theorems, the operators D j j±1 are continuous operators from

H 1/2 (Σ j a ) to H 1/2 (Σ j±1 a ∩ Ω j a ).
In the isotropic case, the expressions of the operators D j j±1 derive directly from the expressions of two operators

D ± ∈ L(H 1/2 (R), H 1/2 (a, +∞)) as follows ∀ψ ∈ H 1/2 (Σ j a ), D j j±1 ψ j (x j , y j = a) = D ± ψ(x j ) for x j > a (identifying H 1/2 (Σ j a ) to H 1/2 (R))
where D ± are defined by

D ± ψ(x) = 1 √ 2π R ψ(ξ)e ı √ ω 2
ε -ξ 2 (x-a) e ±ıξa dξ, for x ≥ a.

Gathering (2.9), (2.11) and (2.12), we have shown that the set of traces

(ϕ j ) j∈ 0,3 ∈ 3 j=0 H 1/2 (Σ j a )
is solution of ∀j ∈ Z/4Z,

ϕ j = D j-1 j ϕ j-1 on Σ j a ∩ Ω j-1 a , ϕ j = D j+1 j ϕ j+1 on Σ j a ∩ Ω j+1 a , ϕ j Σ j aa = g Σ j aa . (2.13)
We have then derived a system of coupled integral equations satisfied by the traces ϕ j of the solution u of (2.2). Moreover, there is an equivalence result between problem (2.2) and system (2.13). Indeed, from a solution (ϕ j ) j∈ 0,3 in H 1/2 (Σ j a ) of (2.13), we can use (2.7) to construct a solution of the Helmholtz equation in each halfspace. To deduce a solution u of (2.2) unequivocally (imposing that its restriction in each halfspace is equal to the halfspace representation), it suffices to show that the halfspace representations match in their intersections -which are quarter planes. The system of coupled integral equations implies that for each quarter planes, the corresponding two halfspace representations match on the boundary. By the well-posedness of the Helmholtz equation with dissipation in the quarter plane, they match also in the quarter plane.

We can state then the following proposition whose proof is detailed in [START_REF] Bonnet-Ben Dhia | The halfspace matching method: a new method to solve scattering problem in infinite media[END_REF].

Proposition 2.2.1 ( [17]) Let g ∈ H 1/2 (Σ aa ). If u ∈ H 1 (Ω) is solution of (2.2) then (ϕ j = u Σ j a ) j∈ 0,3 is solution in H 1/2 (Σ j a ) of (2.13).
Conversely, if (ϕ j ) j∈ 0,3 in H1/2 (Σ j a ) is solution of (2.13) then u defined by

∀j ∈ 0, 3 , u Ω j a = U j (ϕ j ),
where U j (•) is solution of the halfspace problem (2.6) (see also the expression (2.7)), is a function defined "unequivocally", is in H 1 (Ω) and is solution of (2.2).

This system (2.13) is the one that we want to discretize. If we want to use a FE method, a variational formulation has to be derived and using the functional framework H 1/2 (Σ j a ) could be intrincate 1 . That is why we consider (2.13) when looking to the traces in L 2 (Σ j a ). The operators D j j±1 that intervene in the formulation are well-defined and even continuous from

L 2 (Σ j a ) to L 2 (Σ j±1 a
∩ Ω j a ) (see [START_REF] Bonnet-Ben Dhia | The halfspace matching method: a new method to solve scattering problem in infinite media[END_REF] for more details). Writing for all j ∈ 0, 3 , ϕ j = ϕ j 0 + g Σ j aa , we can easily show that Problem (2.13) in the L 2 -framework is equivalent to

Φ 0 = (ϕ j 0 ) j∈ 0,3 ∈ V 0 = {Ψ 0 = (ψ j 0 ) ∈ 3 j=0 L 2 (Σ j a ), ∀j, ψ j 0 Σ j aa = 0} (2.14) is solution of A Φ 0 = B g, where A ∈ L(V 0 ), B ∈ L(L 2 (Σ aa ), V 0 ) and ∀j ∈ Z/4Z, (AΦ 0 ) j = ϕ j 0 -D j-1 j ϕ j-1 0 -D j+1 j ϕ j+1 0 (Bg) j = D j-1 j g Σ j-1 aa + D j+1 j g Σ j+1 aa (2.15)
where for all j, we define the functions D j±1 j ϕ j±1 0 , D j±1 j g Σ j±1 aa ∈ L 2 (Σ j a ∩ Ω j±1 a ) as functions of L 2 (Σ j a ) by extending them by 0. We are able to show stability property and well posedness of problem (2.15) in V 0 . Theorem 2.2.2 (Stability result for the HsMM formulation [START_REF] Bonnet-Ben Dhia | The halfspace matching method: a new method to solve scattering problem in infinite media[END_REF]) (1) The operator A ∈ L(V 0 ) is the sum of a coercive operator and a compact one. Thus, for Problem (2.15), Fredholm alternative holds.

(2) Problem (2.15) is well posed in V 0 .

The operator can be rewritten in the form I -D where the operator D is defined thanks to the operator D j±1 j . The difficulty comes from the fact that the operators D j±1 j are not compact. This is due to the intersection points between Σ j a and Σ j±1 a ∩ Ω j a as it can be shown that χ x j ≥b D j±1 j with b > a is a compact operator. Inspired by the singularity theory [START_REF] Kozlov | Elliptic boundary value problems in domains with point singularities[END_REF], Dauge, 2006], the idea of the proof is to introduce the static equivalent of D j±1 j for ω ε = 0, denoted L j±1 j . The difference of the two operators is compact so it suffices to study the property of L j±1 j , which is done by using the Mellin transform.

To show the second point of Theorem 2.2.2, we cannot use the equivalence result since the Chapter 2. Scattering problem in infinite complex media functional framework has changed (L 2 instead of H 1/2 ). However, uniqueness (and then wellposedness) can be proven in this framework.

Remark 2.2.2 1. Let us point out that what is done here with Ω a being a square can be extended naturally to the case of any convex polygon Ω a . The unknowns, which are involved in the corresponding system of equations, correspond to the traces of the solution on the boundary of the halfspaces supported by each edge of the polygon. See [START_REF] Bonnet-Ben Dhia | Numerical analysis of the halfspace matching method with robin traces on a convex polygonal scatterer[END_REF][START_REF] Bonnet-Ben Dhia | The halfspace matching method: a new method to solve scattering problem in infinite media[END_REF] for more details.

2. The method can be extended to any other boundary conditions on ∂Ω a by taking the unknowns on the lines Σ j a of the same nature than the boundary conditions. Stability results can be shown as well. See [START_REF] Bonnet-Ben Dhia | Numerical analysis of the halfspace matching method with robin traces on a convex polygonal scatterer[END_REF] for more details.

3. If g ∈ H 1/2 (Σ aa ), the last theorem gives that the problems (2.2) and (2.15) are equivalent. If the data g is only in L 2 (Σ aa ), we cannot expect that the solution u is in H 1 (Ω). The functional framework of such problem and its discretization is not standard. By using the transposition method, in [START_REF] Lions | Problèmes aux limites non homogènes et applications[END_REF], the authors analyze how to understand the solution of such boundary value problem in bounded convex domain. They introduce a very weak formulation (the unknown is only in L 2 but the test functions is much more regular). In [START_REF] Apel | Discretization of the poisson equation with non-smooth data and emphasis on non-convex domains[END_REF], an extension to non convex polygonal bounded domain which requires an involved analysis is proposed. For the particular case of Problem (2.15) (and more generally an elliptic PDE set outside a convex polygon), we have introduced a simple formulation which is stable and whose unique solution is such that its trace on each line Σ j a is L 2 (or equivalently its restriction to each halfspace Ω j a is H 1/2 ).

Case 2 : coupling the HsMM with the FE method in the isotropic case

We consider now Problem (2.2), still in the isotropic case : A 0 = Id. Without loss of generality, for conciseness of the statement, we suppose A = A 0 and Ω = R 2 (we only need to adapt the variational formulation and/or the functional framework to take into account the presence of perturbations, obstacles and the associated boundary conditions). Thus the problem simplifies into

-u -ω 2 ε ρ u = f in R 2 (2.16)
where Supp(ρ -1) and Supp(f ) are included in Ω a .

Let us first introduce some new notations. For b ≥ a, we denote (see Figure 2.2)

Ω b = (-b, b) 2 , Σ bb = ∂Ω b and Σ j bb = {(x, y), x j = b and y j ∈ (-b, b)} (2.17)
where the (x j , y j ) are defined in (2.3). In Figure 2.2, five domains appear: a bounded one Ω b and the four halfspaces Ω j a . If b = a, there is no overlap between the bounded domain and the halfspaces whereas if b > a, there exists an overlap. This overlap has an important role from a theoretical point of view.

Let us now derive the new formulation. Thanks to the previous section, it is easy to see that the traces ϕ j of u on the lines Σ j a have to satisfy (2.11). On the other hand, if u b denotes the restriction of u to Ω b , then it is solution of

-u b -ω 2 ε ρ u b = f in Ω b .
(2.18) For all j, in Ω b ∩ Ω j a , we have introduced two representations of the same function: u b solution of (2.18) and U j (ϕ j ) solution of (2.6). We have then in particular that

u b Σ j aa = ϕ j Σ j aa (2.19) and ∇u b • n j Σ j bb = ∇U j (ϕ j ) • n j Σ j bb , (2.20)
where n j = (x j = 1, y j = 0) is the normal to Σ j bb . Let us introduce the operators

∀j ∈ 0, 3 , ∀ψ ∈ H 1/2 (Σ j a ), Λ j ψ = ∇U j (ψ) • n j Σ j bb (2.21)
By classical trace theorem, the operators Λ j are continuous operators from H 1/2 (Σ j a ) to H -1/2 (Σ j bb ) where H -1/2 (Σ j bb ) is defined as the dual of H1/2 (Σ j bb ) which contains the functions of H 1/2 (Σ j bb ) which, when extending by 0, are H 1/2 (Σ bb ). In the isotropic case, we have that for all j ∈ 0, 3 , Λ j can be expressed directly from the operator

Λ ∈ L(H 1/2 (R), H -1/2 (-b, b)) by ∀ψ ∈ H 1/2 (Σ j a ), Λ j ψ (x j = b, y j ) = Λψ (y j ) for y j ∈ (-b, b) (identifying H 1/2 (Σ j a ) and H 1/2 (R))
where Λ is defined by

Λψ(y) = 1 √ 2π R ı ω 2 ε -ξ 2 ψ(ξ) e ı √ ω 2 ε -ξ 2 (b-a) e ıξy dξ. for y ∈ (-b, b).
Gathering (2.11)

u b ∈ H 1 (Ω b ) and (ϕ j ) j∈ 0,3 ∈ 3 j=0 H 1/2 (Σ j a ) is solution of -u b -ω 2 ε ρ u b = f in Ω b ∇u b • n j Σ j bb = Λ j ϕ j , ∀j ∈ 0, 3 ϕ j Σ j aa = u b Σ j aa , ∀j ∈ 0, 3 ϕ j = D j-1 j ϕ j-1 a on Σ j a ∩ Ω j-1 a , ∀j ∈ Z/4Z ϕ j = D j+1 j ϕ j+1 a on Σ j a ∩ Ω j+1 a , ∀j ∈ Z/4Z
(2.22)

Chapter 2. Scattering problem in infinite complex media

There is here also an equivalence result between problem (2.2) in the isotropic case and the system (2.22). Indeed, from the trace unknowns, we can construct a solution u ext of the homogeneous isotropic Helmholtz equation outside Ω a , as in Case [START_REF] Joly | Exact boundary conditions for periodic waveguides containing a local perturbation[END_REF]. Now the difficulty is to show that this function u ext coincides with u b in the overlapping zone Ω b \ Ω a . The well-posedness of the Helmholtz equation with dissipation in the ring Ω b \ Ω a is used.

Proposition 2.2.3 If u ∈ H 1 (Ω) is solution of (2.16) then u b , (ϕ j ) j∈ 0,3 where u b = u Ω b and for all j, ϕ j = u Σ j a is solution in H 1 (Ω b ) × H 1/2 (Σ j a ) of (2.22). Conversely, if u b , (ϕ j ) j∈ 0,3 in H 1 (Ω b ) × H 1/2 (Σ j a ) is solution of (2.22) then u defined by u Ω b = u b , and ∀j ∈ 0, 3 , u Ω j a = U j (ϕ j ), (2.23)
where U j (•) is solution of the halfspace problem (2.6) (see also the expression (2.7)), is a function defined "unequivocally", is in H 1 (Ω) and is solution of (2.16).

We can derive a variational formulation only for b > a, the problem, being as in Domain Decomposition methods [START_REF] Gander | Cross-points in domain decomposition methods with a finite element discretization[END_REF], the cross points. Moreover, as for the case (1), we want to consider the trace unknowns in L 2 for numerical purposes. For b > a, the operators Λ j which intervene in the formulation are well-defined and continuous from L 2 (Σ j a ) to L 2 (Σ j bb ). Thus, for b > a, and by using, as for Case [START_REF] Joly | Exact boundary conditions for periodic waveguides containing a local perturbation[END_REF], for all j, ϕ j 0 = ϕ ju b Σ j aa , we introduce a variational formulation associated to (2.22) :

Find u b , Φ = (ϕ j 0 ) j∈ 0,3 ∈ H 1 (Ω b ) × V 0 such that ∀ v b , Ψ ∈ H 1 (Ω b ) × V 0 Ω b ∇u b • ∇v b -ω 2 u b v b + (AΦ, Ψ) V 0 -(B γ u b , Ψ) V 0 - 3 j=0 (Λ j γ j u b , v b ) Σ j bb -(Λ j ϕ j 0 , v b ) Σ j bb = Ω b f v b (2.24)
where V 0 is defined in (2.14), the operators A and B are defined in (2.15), the operator γ (resp.

γ j ) is the trace operator from H 1 (Ω b ) to L 2 (Σ aa ) (resp. from H 1 (Ω b ) to L 2 (Σ j aa )), (•, •) Σ j bb is the scalar product in L 2 (Σ j bb ), (•, •) V 0 is the scalar product in L 2 (Σ j a
) and where we consider the functions of L 2 (Σ j aa ) and L 2 (Σ j a ∩ Ω j±1 a ) as functions of L 2 (Σ j a ) (by extending them by 0).

Extending the proof of Theorem 2.2.2, we can show a stability result as well.

Theorem 2.2.4 (Stability result for the coupling HsMM-FE -Isotropic case [START_REF] Bonnet-Ben Dhia | The halfspace matching method: a new method to solve scattering problem in infinite media[END_REF]) [START_REF] Joly | Exact boundary conditions for periodic waveguides containing a local perturbation[END_REF] The bilinear form associated to (2.24) in (H 1 (Ω b ) × V 0 ) 2 is the sum of a coercive bilinear form and a compact one. Thus, for Problem (2.24), Fredholm alternative holds.

(2) Problem (2.24) is well posed and if u b , Φ 0 = (ϕ j 0 ) j∈ 0,3 in H 1 (Ω b ) × V 0 is solution of (2.24) then u defined as in (2.23) with for all j, ϕ j = ϕ j 0 + u b Σ j aa is a function defined "unequivocally", is in H 1 (Ω) and is solution of (2.16).

2.2

The Halfspace Matching Method (The HsMM) on a toy problem 59 2.2.4 Case 3 : coupling the HsMM and the FE method in the anisotropic case

Now, we consider Problem (2.1) in the general anisotropic case. Again without loss of generality, we suppose Ω = R 2 . The method to derive the system of coupled equations is exactly the same than in Case ( 2), the only difference being the expression of the solution of the corresponding halfspace problems and then the expression of the integral operators defined in (2.12) and (2.21).

The solution U j (ψ) ∈ H 1 (Ω j a ), for any ψ ∈ H 1/2 (Σ j a ) of the halfspace problems

-∇ • (A 0 ∇U j ) -ω 2 ε U j = 0 in Ω j a , U j = ψ on Σ j a is now given by U j (ψ)(x j , y j ) = 1 √ 2π R ψ(ξ)e r j (ξ)(x j -a) e iξy j dξ, ∀(x j , y j ) ∈ [a, +∞[×R.
where the coefficients r j (ξ) are defined by :

r 0 (ξ) = r 2 (ξ) = -iξc 3 c 1 + i d 1 (ξ) with d 1 (ξ) = ω 2 ε c 1 -d ξ 2 (c 1 ) 2 , r 1 (ξ) = r 3 (ξ) = iξc 3 c 2 + i d 2 (ξ) with d 2 (ξ) = ω 2 ε c 2 -d ξ 2 (c 2 ) 2 .
where d := c 1 c 2 -(c 3 ) 2 . We deduce then the expressions of the operators D j j±1 and Λ j for all j ∈ Z/4Z defined respectively in (2.12) and (2.21). Contrary to the isotropic case, the operators are a priori different from each other. Using the same ideas than in Case (2), we can show that for b > a

u b = u Ω b ∈ H 1 (Ω b ) and (ϕ j = u Σ j a ) j∈ 0,3 ∈ 3 j=0 L 2 (Σ j a )
is solution of (2.24) with for all j, ϕ j 0 = ϕ ju b Σ j aa where, here, the operators are given by for

all ψ ∈ H 1/2 (Σ j a ) D j j±1 ψ(x j , y j = ±a) = 1 √ 2π R ψ(ξ)e ı √ r j (ξ)(x j -a) e ±ıξa dξ, for x j ≥ a and Λ j ψ(x j = b, y j ) = 1 √ 2π R
r j (ξ) ψ(ξ)e r j (ξ)(b-a) e ıξy dξ, for y j ∈ (-b, b).

Proposition 2.2.3 extends easily to the anisotropic case. Indeed, for the proof, we do not use the expression of the solution of the halfspace problems or the expression of the operators but only that the halfspace representations, two by two, are solution of the same equations in the quarter planes where they coexist and that their traces coincide on the boundary of the quarter plane.

We can also introduce a variational formulation that is similar to (2.24). The main originality in the anisotropic case concerns the stability results of the variational formulation which are linked to the properties of the operators D j j±1 and Λ j .

Theorem 2.2.5 (Stability result for the anisotropic case [START_REF] Bonnet-Ben Dhia | The halfspace matching method: a new method to solve scattering problem in infinite media[END_REF]) Let b > a and θ defined by tan

(θ) = c 1 c 2 -c 2 3 /c 3 if c 3 = 0 and θ = π/2 if c 3 = 0. (i) If max(sin(θ/2), cos(θ/2)) max( √ c 1 , √ c 2 ) < (c 1 c 2 -c 2 3 ) 1/4 , (2.25)
the bilinear form associated to the problem is the sum of a coercive bilinear form and a compact one in (H 1 (Ω b ) × V 0 ) 2 ;

(ii) For any values of c i , Fredholm alternative holds for the problem and it is well-posed.

As you can notice, the first result was proved only for a certain class of moderate anisotropy.

In particular if the x 0 and y 0 directions are chosen regarding the directions of anisotropy of A 0 then c 3 = 0 and the condition (2.25) reduces to

max c 1 c 2 , c 2 c 1 < 4.

Remark 2.2.3

Let us emphasize than in the general anisotropic case, the operator associated to the bilinear form of the problem cannot be decomposed as the sum of a coercive operator and a compact one but only as the sum of an invertible one and a compact one. The well posedness is then ensured by the uniqueness. The consequence of this weaker result is that a priori, the numerical analysis cannot be done in the general case. However, the problem can be adapted to recover the stability property. Indeed, we can show that a solution of (2.24) is also solution of a similar problem where the A is replaced by an operator à = I -(I -A) 2 and the operator B by B = (2I -A)B. In the general anisotropic case, the operator à is always the sum of a coercive operator and a compact one.

Numerical results

To end this section, let us briefly describe how to discretize this new formulation and show some numerical results to illustrate and validate the method.

For the approximation of the formulation (2.15) (resp. (2.24)), we use 1D Lagrange FE method for the trace unknowns (resp. 2D Lagrange FE method for the volume unknown). More precisely, the finite dimensional space V h,T ⊂ V 0 contains piecewise polynomial functions of L 2 (-T, T ) 4 .

The difficulty in practice is to handle the integral operator terms. For instance, to compute

D + ϕ, ψ (a,+∞) = 1 √ 2π +∞ a ψ(x) ∞ -∞ ϕ(ξ)e ı √ ω 2 -ξ 2 (x-a) e ıξa dξ dx
where ϕ and ψ are functions of L 2 (a, +∞),

• we reduce the inner integral in ξ to the interval [-T , T ] where T > 0,

• we use a quadrature formula to discretize the integral in ξ.

Doing so, the previous integral is approximated by

N ξ n=0 1 √ 2π ϕ(ξ n )e ıξna +∞ 0 e ı √ ω -ξ 2 n (x-a) ψ(x)dx
where the {ξ n } n∈ 1,N ξ correspond to the quadrature nodes. Since ϕ and ψ are piecewise polynomial (and compactly supported), we can compute exactly ϕ(ξ n ) and the integral term in x.

We can treat similarly the integral terms related to the operators D j j±1 and Λ j .

The numerical analysis of the method, before using the quadrature formulas in ξ, was done in [START_REF] Bonnet-Ben Dhia | Numerical analysis of the halfspace matching method with robin traces on a convex polygonal scatterer[END_REF]. The difficulties are [START_REF] Joly | Exact boundary conditions for periodic waveguides containing a local perturbation[END_REF] to prove the problem is coercive+compact and the bilinear form is approximated because of the truncature of the Fourier variable and ( 2) to show the stability properties of the approximated bilinear form, the equivalent of Theorem 2.2.2 has to be shown and the properties of the operators L j±1 j,

T (with the truncature of the fourier variable) cannot be deduced from the properties of L j±1 j . The error estimates is given by the following result.

Theorem 2.2.6 (Error estimates) Let the ϕ j 's (resp. the ϕ j h,T, T ) be the solution of the continuous HsMM formulation (resp. the discrete approximated formulation after quadrature in the Fourier variable). We suppose that the ϕ s are in H s (Σ j ) then ∀j,

ϕ j -ϕ j h,T, T L 2 (Σ j ) ≤ C T s + Ce -εT + C h min(s,l+1)
You can find the error plots with respect to each discretization parameters in [START_REF] Bonnet-Ben Dhia | Numerical analysis of the halfspace matching method with robin traces on a convex polygonal scatterer[END_REF].

To validate the method, we have considered in Case ( 1) the particular data of the Hankel function on the boundary of the square:

u(x, y) Σaa = 1 4ı H(ω x 2 + y 2 ) Σaa .
In that situation, we can validate the results since we know the exact solution is given by u(x, y) = 1 4ı H(ω x 2 + y 2 ). On Figure 2.3 (left), we have represented on (-T, T ) the real part of ϕ 0 and its approximation taking ω = 10+0.01ı and using 1D P2 finite elements with h = 0.05, T = 12, T = 20 and a Gauss quadrature formula of order 4 in a regular mesh of size 0.025. We get a relative error in the L 2 (-T, T ) of 0.06%. Surprisingly, as we can see, even though the trace ϕ j are not close to zero at y j = ±T , the results are quite good.

Let us emphasize that with this method, we can reconstruct the solution in the exterior domain. Indeed, we can compute a posteriori the solution in each halfspace Ω j a using the formulae (2.6) as represented on Figure 2.4. In the left figure, we have represented the reconstructed halfspace solution in the right halfspace, in the right figure, we have added the halfspace representation in the top halfspace. As you can notice, the solutions coincide in the quarter plane where they coexist. We have finally represented the solution in the whole exterior domain. This representation is an a posteriori validation of the discretization since when the discretization is not fine enough, the halfspace solutions do not coincide (see Figure 2.4 Bottom right) Let us now show a numerical simulation in the anisotropic case. We have solved (2.1) in (2.26)

Ω = R 2 \ O (O

The case without dissipation 63

It is a case of strong anisotropy for which the condition (2.25) is not satisfied. We can see on Figure 2.5 that the method gives good numerical results. One simple way to validate these results is to see if the reconstructed solution in the halfspaces match in the quarters of plane.

The case without dissipation

Up to now, we have considered the dissipative problem (4) (with Im(ω ε ) > 0). But we claim that in practice, our method also works in the non dissipative case Im(ω) = 0 (See Figure 2.6).

It is quite easy to extend the discrete formulation to this case, selecting carefully the outgoing solution of the halfspace problem. The price to pay in order to get accurate results is then to use a refined discretization in the Fourier variable ξ. Indeed, the solution decaying much more slowly than in the dissipative case (the solution in the non-dissipative isotropic case behaves like e iωr √ r at infinity, while the solution in the dissipative case decays exponentially), its Fourier transform is less regular. The main difficulty lies at the theoretical level: the extension of the theoretical results to the non-dissipative case still raises many open questions, that we briefly discuss below.

Let us consider to fix ideas the simplest case of the (isotropic non-dissipative) Helmholtz equation. It is well-known that the well-posedness of the associated boundary value problem is ensured if one imposes to the solution to be "outgoing", which means for instance that the solution satisfies the Sommerfeld radiation condition at infinity. Our objective is to reformulate this problem using our halfspace matching formulation. The first difficulty is that the traces ϕ j of the solution are no longer in L 2 (again, they behave like e iωr √ r ). An appropriate functional framework has been introduced in [Bonnet- Ben Dhia andTillequin, 2001a,Bonnet-Ben Dhia andTillequin, 2001b], where the space H 1/2 (R), convenient for the dissipative case, is replaced by the space V (ω; R) of the functions ϕ such that |ξ 2 -ω 2 | 1/4 φ ∈ L 2 (R). The functions of V (ω; R) are locally in H 1/2 , but they may not belong to L 2 (R), due to their behavior at infinity. Using this framework we can state the following result:

Let g ∈ H 1/2 (Σ aa ). If u is the outgoing solution of -∆u -ω 2 u = 0 in Ω = R 2 \Ω a u = g on Σ aa
then {ϕ j a , j ∈ {0, 1, 2, 3}} is solution of (2.13) in the space 3 j=0 V (ω; Σ j a ) . But we are not able to prove the converse statement corresponding to the second part of Proposition 2.1 (even if we conjecture that it is true). Indeed, we did not succeed in proving that the different halfspace representations of the solution match in the quarters of planes. For instance, we are not able to prove like in paragraph 3.

1 that V = U 0 (ϕ 0 ) -U 1 (ϕ 1 ) vanishes in Ω 0 a ∩ Ω 1 a .
The reason is that we cannot prove that V is outgoing (in the sense of Sommerfeld radiation condition at infinity), since we just know that U 0 (ϕ 0 ) is propagating in the direction of positive x, while U 1 (ϕ 1 ) is propagating in the direction of positive y.

Another difficulty arises when we try to extend stability results like Theorem 2.2 to the nondissipative case. All the theory has been done with traces in L 2 , which is no longer the appropriate space. An idea could be to introduce the space L(ω; R) of the functions ϕ such that

|ξ 2 -ω 2 | 1/4 |ξ 2 +ω 2 | 1/4 φ ∈ L 2 (R)
, because the functions of L(ω; R) are locally in L 2 . But it is far from obvious to prove the different properties of the integral operators in this space. Finally, proving uniqueness will be of course much more intricate in the non-dissipative case than in the dissipative case, where we mainly used the coerciveness in H 1 . A first step has been done by proving a related uniqueness result: in [START_REF] Bonnet-Ben Dhia | A rellich type theorem for the helmholtz equation in a conical domain[END_REF], overlapping halfspaces representations are used to prove the absence of trapped modes (i.e. L 2 solutions of the homogeneous non-dissipative equation) under very weak hypotheses.

Les us emphasize that our approach seems to be well-suited to formulate a large class of prob-2.4 The case of an anisotropic elastic media lems, for which no equivalent of the Sommerfeld radiation condition is available. This is typically the case of anisotropic elastic media or periodic media. For these problems, we aim at using our halfspace formulation, not only to solve numerically the problem, but also to define the notion of outgoing solution and to prove existence and uniqueness of this outgoing solution. A similar idea has been used for instance in [Bonnet-Ben Dhia et al., 2011] using a formulation with non-overlapping halfspace representations for the junction of open waveguides.

The case of an anisotropic elastic media

The interest of the HsMM is that it exploits only the properties of the medium in each halfspace. This method extends then naturally to the time-harmonic elasto-dynamic equations even for anisotropic media. We begin this section by explaining the extension of the HsMM to that case. Numerical results will be shown at the end of this section to validate and illustrate the efficiency of the method. Moreover, let us recall that the PMLs fail for general anisotropic media and this will be also shown.

Let us consider the elastic diffraction problem in 2D

-∇ • σ(u) -ω 2 ε ρ u = 0 in Ω = R \ O, σ(u) • ν = 0 on ∂O (2.27)
where u represents the displacement field (u = (u x , u y ), σ(u) the stress tensor which is related to the strain tensor ε(u) = 1/2 (∇u + ∇ T u) (in the small deformation assumption) through the stiffness tensor C by the Hooke's law (as in (1.37)), C is supposed to be constant outside Ω a :=] -a, a[ 2 , ρ is the density of the material which si supposed to be constant outside Ω a , O is a set of bounded obstacles contained in Ω a , the source term f is of compact support and finally we consider the case with dissipation Im ω ε > 0.

Using the same approach and the same notations than previously, we can derive a system of coupled equations whose unknowns are vectorial and are the restriction of the displacement field in

Ω b :=] -b, b[ 2 (b > a)
, the traces of the displacement field on each lines Σ j a . The only difference is the expression of the solution of the halfspace problems and then the expression of the integral operators defined in (2.12) and (2.21).

We then explain in the rest of this section how the solution of the halfspace problems can be expressed thanks to its trace on the boundary. To simplify the notations and the explanations, we focus on the halfspace Ω 0 and we omit the subscript 0 for the halfspace solution

U := U 0 ∈ H 1 (Ω 0 a ) × H 1 (Ω 0 a ) of divσ(U) + ω 2 ε U = 0 in Ω 0 a , U = Ψ on Σ 0 a , where Ψ = (Ψ x , Ψ y ) ∈ H 1 2 (Σ 0 a ) × H 1 2 (Σ 0 a
). This can be rewritten using the Voigt-Mandel notations 2.28) where the C i are symmetric 2x2 matrices whose coefficients are related to the coefficients of C (see for instance [START_REF] Royer | Elastic waves in solids II: generation, acousto-optic interaction, applications[END_REF]).

C 1 ∂ 2 ∂x 2 U + C 2 ∂ 2 ∂x∂y U + C 3 ∂ 2 ∂y 2 U + ω 2 ε U = 0, in Ω 0 a , U = Ψ on Σ 0 a , (
By applying the Fourier Transform in the y-direction to the previous equations, we obtain for almost every ξ ∈ R,

C 1 ∂ 2 ∂x 2 U -ıξ C 2 ∂ ∂x U -ξ 2 C 3 U + ω 2 ε U = 0, for x ≥ a, U(a; ξ) = Ψ(ξ) (2.29)
which is a system of two coupled ODEs of order 2 with constant coefficients and where the Fourier variable ξ plays the role of a parameter. The solution of (2.29) can be written for a.e. ξ, as

U(x; ξ) = α 1 (ξ)U 1 (ξ)e ıq 1 (ξ)(x-a) + α 2 (ξ)U 2 (ξ)e ıq 2 (ξ)(x-a)
where 1. the q j (ξ)'s are the two only solutions having a positive imaginary part, of the quadratic eigenvalue problem that is obtained by replacing in (2.29), ∂ ∂x by ıq -the quadratic eigenvalue problem has 4 solutions : two of them have a positive imaginary part (corresponding to exponentially decaying functions) and the other two have a negative imaginary part (corresponding to exponentially growing functions);

2. U j (ξ)'s are the associated eigenvectors;

3. the α j (ξ)'s are such that U(x; ξ) satisfies the boundary condition at x = a.

The solution of (2.28) is finally obtained by applying the inverse Fourier transform and the expressions (2.12) and (2.21) of the operators appearing in the HsMM formulation (the equivalent of (2.24)).

Using the same arguments than for the scalar case, the HsMM formulation is equivalent to the original problem. The theoretical analysis of the formulation and its numerical analysis are in progress.

Even if the case without dissipation still raises open questions, the method works well. One has only to construct the "outgoing" solution of each halfspace problem. However, here this is less obvious than in the scalar case since this is linked to the choice of the "physical" solutions of the quadratic eigenvalue problem. In the non dissipative case, it can be shown that, for small values of ξ, the solutions of the quadratic eigenvalue problem may be real (and the associated eigenvector propagative) so it is not possible to select the good solutions directly. An additional criterion has to be considered. As explained in Chapter 1, the physical criterion is the energy flux so it suffices to select the solutions q's such that

Im(σ(U q (ξ, •))e x • U q (ξ, •)) > 0.
where U q is the associated eigenvector. It can be shown that at most two solutions q satisfy this criterion.

We validate the method by applying it for an isotropic media and by choosing the source term with a null rotational (resp. a null divergence) so that the solution of (2.27) is a P-wave (resp. a S-wave), see Figure 2.7. We can then apply the method for a general anisotropic media still by choosing the source term with a null rotational (resp. a null divergence): the solution is called Quasi-P (Quasi-S) wave even if in general the waves cannot be decoupled as for the isotropic case, see Figure 2.8.

Finally, as we have explained in the introduction of this chapter, the PML does not work for general anisotropic media (see Figure 2.10) even if they do for isotropic ones (see Figure 2.9). As you have understood the HsMM works for any general anisotropic elastic media.

The case of a periodic media

We consider a 2D anisotropic scalar acoustic model

-div(A∇u) -ω 2 ε ρ u = f in Ω (2.30)
in time-harmonic regime with a small absorption, Im(ω 2 ε ) = ε > 0, where

• the domain of propagation Ω is R 2 (geometrical defects included in Ω a can also be considered, boundary conditions have to be added to the model);

• the support of the source term is supposed to be compact and strictly included in Ω a ;

• A is a symmetric, positive definite matrix of (L ∞ (Ω)) 2×2 and ρ is a strictly positive function of L ∞ (Ω) and outside Ω a , the matrix A and the coefficient ρ are (L x , L y )-periodic in the two directions. Without loss of generality, we will suppose that L x = L y = 1. Our objective is to propose a numerical solution for the computation of this solution, the specific difficulty being the periodic nature of the infinite media far from the perturbation.

Remark 2.5.1

Here the periodic media have a square structure (the directions of periodicity are orthogonal) but a hexagonal structure can also be considered (see [6] for more details)

It seems that there are very few works for the simulation of wave propagation in infinite periodic media. A first class of methods covers problems where the periodicity can be treated by homogenization techniques [START_REF] Bensoussan | Asymptotic analysis for periodic structures[END_REF], Allaire, 1992], typically when the wavelength is much larger than the period. The unboundedness of the homogenized and often anisotropic media can then be handled using classical methods mentioned in the Introduction of this chapter.

A second class of methods considers the periodicity as such but only for finite media [START_REF] Yuan | An efficient bidirectional propagation method based on dirichlet-to-neumann maps[END_REF], Ehrhardt et al., 2008, Ehrhardt and Zheng, 2008] or periodic media which can be reduced to finite domain when the frequency is in a band gap of the underlying operator (by using the supercell method [Soussi, 2005]). In our work, we consider infinite periodic media where the periodicity has to be considered as such.

The Halfspace Matching Method was first conceived and introduced for periodic media during my phD with Patrick Joly. Now I often present the periodic case after the homogeneous case because the notations and the computations are much more involved. The idea of the method is the same : exploit the -namely here periodic -properties of the media in some halfspaces surrounding the perturbation.

To simplify the presentation, we assume that a = 1/2 and b = 3/2 but any a, b ∈ 1/2Z * can be considered. Using the same approach and the same notations than previously, we can derive a system of coupled equations whose unknowns are the restriction of the solution in

Ω b :=] -b, b[ 2 ,
the traces of the solution on each lines Σ j a . The only difference is the expression of the solution of the halfspace problems and then the expression of the integral operators defined in (2.12) and (2.21).

We explain then in the rest of this section how the solution of the halfspace problems can be expressed thanks to its trace on the boundary. To simplify the notations and the explanations, we focus on the halfspace Ω 0 and we omit the subscript 0 for the halfspace solution

U := U 0 ∈ H 1 (Ω 0 a ) of -div(A∇U ) -ω 2 ε ρ U = 0 in Ω 0 a , U = ψ on Σ 0 a .
(2.31)

where A and ρ are 1-periodic in the two directions x and y.

To solve this halfspace problem, we replace the Fourier Transform used for homogeneous media with the privileged tool for the study of PDEs with periodic coefficients: The Floquet Bloch Transform (FBT) (see also Section 1.3). The FBT of period 1 is defined by (see [Kuchment, 1993]):

F : ϕ ∈ L 2 (R) → Fϕ ∈ L 2 (K)
where for a.e. (y, ξ)

∈ K = (-1/2, 1/2) × (-π, π) Fϕ(y, ξ) = 1 √ 2π n∈Z ϕ(y + n)e -ınξ .
Moreover we have the inversion formula: for a.e.y ∈ [0, 1] ∀n ∈ Z,

ϕ(y + n) = 1 2π π -π Fϕ(y, ξ)e ınξ dξ
We denote F y u the FBT of a function u of R 2 applied in the y-direction.

By using the properties of the FBT (see [Kuchment, 1993] and [START_REF] Fliss | Étude mathématique et numérique de la propagation des ondes dans des milieux périodiques localement perturbés[END_REF]), we apply the FBT in the y-direction to the halfspace problem (2.31) and find that the FBT in the y-direction

Û (ξ) = F y U (•; ξ) of U is solution, for any dual variable ξ ∈ (-π, π), of the strip problem -div(A p ∇ Û (ξ) -ω 2 ε ρ p Û (ξ) = 0, in B 0 Û (ξ)| x=a = F y ψ(ξ)| x=a , Û (ξ)| y=1/2 = e ıξ Û (ξ)| y=-1/2 A p ∇ Û (ξ) • e y | y=1/2 = e ıξ A p ∇ Û (ξ) • e y | y=-1/2
where

B 0 = Ω 0 a ∩ {-1/2 ≤ y ≤ 1/2}.
This is a semi-infinite periodic waveguide problem with quasi-periodic boundary conditions parametrized by the Floquet dual variable ξ. As explained in Chapter 1 and Section 1.3.3, we can characterize and compute periodicity cell by periodicity cell the solution of this semi-infinite periodic waveguide thanks to the solutions e 0 ξ and e 1 ξ of two cell problems and the propagation operator P (ξ). This operator is the unique operator of spectral radius strictly less than 1 to the stationary Riccati equation

T 10 (ξ)P (ξ) 2 + (T 00 (ξ) + T 11 (ξ)) P (ξ) + T 01 (ξ) = 0,
where the operators T ij (ξ) are the local DtN operators. Applying the FBT inverse, U and the operators Λ 0 (we recall the reader that ∀ψ, Λ 0 ψ = A p ∇U •e x Σ 0 bb where here Σ 0 bb = {(3/2, y), y ∈ (-3/2, 3/2)), D 0 1 and D 0 2 (we recall the reader that ∀ψ, D 0 1 ψ = U y=1/2 and D 0 2 ψ = U y=-1/2 ) can then be expressed semi-analytically. In particular, we have for all n ∈ {-1, 0, 1}, for a.e.

1/2 + n < y < 1/2 + n + 1, (Λ 0 ψ) (y) = - 1 √ 2π π -π
(T 00 (ξ) + T 10 (ξ)P (ξ)) P (ξ) F y ψ(•; ξ) e ınξ dξ, and for all n ∈ N and for a.e. 1/2 + n < x < 1/2 + n + 1

(D 0 1 ψ) (x) = 1 √ 2π π -π D(ξ)P (ξ) n F 2 ψ(•; ξ) dξ,
where D(ξ) = D 0 (ξ) + D 1 (ξ)P (ξ) and the operators D 0 (ξ) ψ and D 1 ψ are respectively the traces on y = 1/2 of the cell solutions e 0 ξ ( ψ) and e 1 ξ ( ψ).

This gives an idea of the expressions of all the operators appearing in the HsMM formulation (2.24) applied to the periodic case. It is easy to show that this formulation is equivalent to the original problem (using the same arguments than for the homogeneous case). Using Theorems 2.2.4 and 2.2.5, we are able to show equivalent results.

Theorem 2.5.1 (Stability result for the periodic case) If A = I outside Ω a , the bilinear form associated to (2.24) written for the periodic case in

(H 1 (Ω b ) × V 0 ) 2
is the sum of a coercive bilinear form and a compact one.

In general, Problem (2.24) written for the periodic case, Fredholm alternative holds.

From the numerical point of view (see [START_REF] Fliss | Wave propagation in locally perturbed periodic media (case with absorption): Numerical aspects[END_REF] for more details), • we use 1D Lagrange FE method for the trace unknowns and a 2D Lagrange FE method for the volume unknown and for the periodicity cell;

• to handle the integral operators we use a quadrature formula to discretize the dual variable ξ ∈ (-π, π);

• to compute the integral operators, for each ξ, we have to (1) solve the cell problems, (2) compute the local DtN operators T ij (ξ) and the local DtD D j (ξ) and (3) solve the stationary Riccati equation to compute P (ξ).

Let me mention that the integral equations can be considered in the space variable (the unknowns are the traces) or in the Floquet variable (the unknowns are the Floquet Bloch Transform of the traces) -same remarks apply for the homogeneous case. I have implemented the second case but the first case could be done. For instance Figure 2.12 represents the solution (right figure) for a particular periodic media, whose ρ is represented in the left figure. In this computation, Ω a is a small region containing the perturbation and using the traces, I have then reconstructed the solution everywhere by using the solution of the halfspace problems.

Ongoing works and perspectives

Although the method works very well in practice in the non-dissipative case, we were not able to extend the theoretical analysis to that case. Making progress in that direction is a challenging objective.

In collaboration with CEA LIST, we want to propose numerical methods in the context of Non Destructive Testing or Structural Health Monitoring of elastic plates made of composite materials (in aeronautics or in civil engineering structures for instance). In these applications, the medium is infinite in two directions and bounded in the other direction. The extension to 3D elastic plates -elasticity equations in such 3D structures are considered and not asymptotic plate equations in a limit 2D structure-of the Half-Space Matching method requires some generalizations of the formulation which are not obvious. In 3D, the representation of the solution in a half-plate combines a Fourier transform in the parallel direction to the boundary of the half-plate, and representations in series expansion of Lamb modes in the thickness of the plate.

In order to exploit bi-orthogonality relations of Lamb modes, one has to use both Dirichlet and Neumann traces of the solution (displacement and normal stresses) (see Section 1.4.3), which raises difficult theoretical questions, even in the dissipative case. This is done in the framework of Yohanes Tjandrawidjaja PhD thesis. From the computational point of view, the efficiency of the method will probably be compromised by the presence of large full matrices, associated with the integral operators. One of the objectives is to adapt to the discretized Half-Space Matching method some compressing techniques.

The HsMM has been already used by Julian Ott (KIT, Karlsruhe) to deal with the junction of optical fibers (the algorithm was then integrated in an optimization process in order to optimize the efficiency of the junction), see Section 3.4 and for instance the configuration of Figure 1.17 top left. I would like to consider more general configurations of junctions of fibers, for instance with the presence of periodic structures (see Figure 1.17 top right). For this last configuration, the determination of the halfspace solution has to be done.

One limitation of the Half-Space Matching method is precisely that the use of a partial Fourier transform imposes the geometry of Half-Spaces for the exterior sub-domains. We plan to investigate the possibility of using more general unbounded overlapping sub-domains, replacing the Fourier representations by classical integral representations. The method would be therefore applicable in configurations where a Green function is available for each subdomain albeit not globally for the total exterior domain. This is the case for instance in 2D if each subdomain has its own stratification, which does not correspond to a global stratification. This should also allow to treat the junction of several embedded pipes or optical fibers (see Section 1.5 for a presentation of the problem and for an alternative method).

Another perspective consists in extending the Half-Space Matching method to the time domain. As for convolution quadrature methods, the idea is to use an implicit scheme. Then, the problem at each time step is solved by using the Half-Space Matching method. The difficulty is the presence of a source term which depends on the previous time steps and which is no longer compactly supported. This is in progress in the framework of Hajer Methenni's PhD thesis.

In time domain or in frequency domain, the extension of the method to multiple defects and/or lineic defects seems to be important for the applications in ultrasonic Non Destructive Testing.

Let us finally mention that most of the subjects described above for ultrasonic waves could Chapter 2. Scattering problem in infinite complex media be also relevant for electromagnetic waves. For instance, the Half-Space Matching method could be a way to solve Maxwell equations in an unbounded anisotropic medium, where PMLs can be unstable.

CHAPTER 3

Trapped and guided modes in periodic structures

Collaborations : Bérangère Delourme (LAGA, Paris 13), Patrick Joly (POEMS), Elizaveta Vasilevskaya (LAGA, Paris 13) Supervising : Berangere Delourme's Post-doc (2012), Khac-Long Nguyen 'Post-doc (2014'Post-doc ( -2015) ) 3

.1 Introduction

As I said in the introduction, periodic media play a major role in applications, in particular in optics for micro and nano-technology [START_REF] Joannopoulos | Photonic Crystal -Molding the Flow of Light[END_REF], Johnson and Joannopoulos, 2002, Kuchment, 2001, Sakoda, 2001]. From the point of view of applications, one of the main interesting features is that it can exist intervals of frequencies, called band gaps, for which the propagative waves cannot exist in the media. This phenomenon is due to the fact that a wave on the media is multiply scattered by the periodic structure, which can lead, depending on the characteristics of the media and the frequency, to possibly destructive interferences. It seems then that an appropriate choice of the structure and the dielectric materials of the photonic crystal can create particular band gap and then, from a practical point of view, banish some monochromatic electromagnetic waves. Thus, the periodic media could be used to several potential applications such as in the realization of filters, antennas and more generally, components used in telecommunications.

Mathematically, this property is linked to the gap structure of the spectrum of the underlying differential operator appearing in the model. For a complete, mathematically oriented presentation, we refer the reader to [Kuchment, 2001,Kuchment, 2004]. Without being exhaustive, let us review a few important results on the topic. In the one dimensional case, it is well known [Borg, 1946] that a periodic material can have (often infinitely) many gaps unless it is constant. By contrast, in 2D or 3D, a periodic medium might or might not have gaps and necessary conditions for the existence of band gaps are not known. Nevertheless, sufficient conditions where at least one gap exists can be found. For instance, Figotin and Kuchment have given examples of high contrast medium for which as many band gaps as one wants can be created [Figotin andKuchment, 1996a, Figotin andKuchment, 1996b]. Using asymptotic arguments, Nazarov and co-workers have established that a small perturbation of an homogeneous waveguide can open a gap in the continuous spectrum of the operator [Nazarov, 2010, Cardone et al., 2010]. Moreover, other band gap structures have been characterized through numerical approaches in [Figotin

76

Chapter 3. Trapped and guided modes in periodic structures and Godin, 1997]. In any case, except in dimension one, the number of gaps is expected to be finite. This statement, known as the Bethe-Sommerfeld conjecture, is fully demonstrated for periodic Schrödinger operators in 2D [Skriganov, 1979, Popov andSkriganov, 1981] and in 3D [Skriganov, 1983, Skriganov, 1985, Karpeshina, 1997] and for particular 2D periodic Maxwell operators [Vorobets, 2011].

Besides, a local perturbation of the crystal can produce defect mid-gaps modes, that is to say solutions to the homogeneous time-harmonic wave equation, at a fixed frequency located inside one gap, that remains localized in the vicinity of the perturbation. This localization phenomenon is of particular interest for a variety of promising applications in optics to design lasers, fibers or efficient waveguides in general. This localization effect is directly linked, from a mathematical point of view, to the possible presence of isolated eigenvalues of finite multiplicity inside the gaps appearing when perturbing the underlying periodic operator. Several papers exhibit situations where a compact (resp. lineic) perturbation of a periodic medium give rise to localized (resp. guided) modes. It seems that the first results concern strong material perturbations : for local perturbations [START_REF] Figotin | Localized classical waves created by defects[END_REF], Figotin and Klein, 1998a, Figotin and Klein, 1998b] and for lineic perturbation [START_REF] Kuchment | On guided waves in photonic crystal waveguides[END_REF]. There exist fewer results about weak material perturbations: [START_REF] Brown | Spectrum created by line defects in periodic structures[END_REF],Brown et al., 2015] deal with 2D lineic perturbations. Finally geometric perturbations are considered in [START_REF] Nazarov | On the spectrum of the laplace operator on the infinite dirichlet ladder[END_REF], Nazarov, 2014a].

A first part of our contributions to this subject, is to complement the references mentioned above by proving the existence of guided modes created by a lineic geometrical perturbation of a particular periodic medium. A summary of the results and the method of study is presented in Section 3.3.

From a numerical point of view, there exist only few methods to compute trapped or guided modes. The most well known is the Supercell method. It consists in making computations in a bounded domain of large size with periodic boundary conditions, the resulting solution converging to the true solution when the size tends to infinity. The convergence for the computation of defect modes has been shown for 2D problems and compact perturbations in [Soussi, 2005] and generalized to 3D problems and to exponentially decaying perturbations in [START_REF] Cancés | Periodic Schrödinger operators with local defects and spectral pollution[END_REF].

In this case, as the localized modes are exponentially decaying, this convergence is exponentially fast with respect to the size of the truncated domain. In practice, this approach replaces the eigenvalue problem set in an unbounded domain to an approximated one set in a bounded domain. See [START_REF] Schmidt | Efficient computation of photonic crystal waveguide modes with dispersive material[END_REF], for example, for numerical results. The main drawback of this strategy relies on the increase of the computational cost, especially when a mode is not well confined.

By adapting to eigenvalue problems the construction of Dirichlet-to-Neumann operators originally developed for scattering problems and presented in Section 1.3.4, we want to offer a rigorously justified alternative to existing methods. Compared to the Supercell method, the DtN method allows us to reduce the numerical computation to a small neighborhood of the defect independently from the confinement of the computed guided modes. Moreover, as the method is exact, we improve the accuracy for non well-confined guided modes. Obviously, there is a price to be paid : the reduction of the problem leads to a non linear eigenvalue problem, of 3.2 A Dirichlet-to-Neumann approach for the exact computation of guided modes in photonic crystal waveguides 77 a fixed point nature. However, this difficulty has been already overcome for homogeneous open waveguides for which the DtN approach is well known [START_REF] Joly | A numerical method for the computation of electromagnetic modes in optical fibres[END_REF],Pedreira and Joly, 2001, Pedreira and Joly, 2002]. This is presented in Section 3.2.

3.2 A Dirichlet-to-Neumann approach for the exact computation of guided modes in photonic crystal waveguides

The propagation model we consider is a simple 2D space-(x = (x, y)) time domain scalar problem and we look for w satisfying the wave equation

ρ(x) ∂ 2 w ∂t 2 -w = 0, x ∈ Ω, t ≥ 0 (P)
The domain of propagation Ω is infinite in the two directions, its geometry is periodic in the y-direction with the period L y and periodic in the x-direction with the period L x outside a straight band

Ω 0 = ] -a, a[ × R
The function ρ is periodic as well in

Ω + ∪ Ω -= Ω \ Ω 0 (Ω ± = Ω ∩ ((±a, ±∞) × R))
, with the same periodicity than the geometry (see figure 3.1)

ρ(x, y) = ρ p (x, y) in Ω + ∪ Ω - ρ 0 (x, y) in Ω 0 with ρ p (x ± L x , y ± L y ) = ρ p (x, y) ∀ (x, y) ∈ Ω ρ 0 (x, y ± L y ) = ρ 0 (x, y) ∀ (x, y) ∈ Ω 0
Let us denote B one period of the domain in the y-direction (see figure 3.2):

B = R×] -L y 2 , L y 2 [.
A guided mode of this problem is by definition a solution w to (P) which can be written in the form w(x, y) = v(x, y) e ı(βy-ωt) where the frequency ω ∈ R + , v is periodic in the y-direction with period L y and v B ∈ L 2 (B), the quasi period β can then be considered in ]-π/L y , π/L y [. This corresponds to finding couples (ω 2 , β) such that there exists u ∈ H 1 (B), u = 0 solution of

- 1 ρ u = ω 2 u, in B u Σ + = e ıβLy u Σ -, ∂ y u Σ + = e ıβLy ∂ y u Σ -
where (see Figure 3.2)

Σ ± = R × ± Ly 2 .
The classical approach for solving this eigenvalue problem set in an unbounded domain is the Supercell Method [Soussi, 2005] which is based on the exponentially decaying property of the mode in the x-direction and consists in truncating the band B far enough. Let us mention that the exponential decay rate depends on the distance between the corresponding eigenvalue to the essential spectrum of the operator. We could then enlighten some drawbacks of the method:

1. the essential spectrum of the operator has to be computed initially.

2. for a fixed β, it seems important to have an estimation of the distance between the not yet computed eigenvalue and the essential spectrum of the operator to choose a relevant truncated domain.

3. The size of the truncated domain depends on β and on the distance between the eigenvalue and the essential spectrum. If the eigenvalue approaches the essential spectrum of the operator when β varies, the corresponding eigenvector becomes less and less confined (less and less exponentially de-creasing) and then the truncated domain has to be bigger and bigger.

Here we propose a novel method based on a DtN approach which offers a rigorously justified alternative to the super-cell method.

Spectral theory results

We have reduced the problem of finding the guided modes to the following problem :

Find β ∈ ] - π L y , π L y [, ω 2 ∈ R + , s.t. ∃ u ∈ H 1 (B), u = 0, A(β)u = ω 2 u (E)
where

A(β) = - 1 ρ D(A(β)) = u ∈ H 1 ( , B), u| Σ + = e ıβLy u| Σ - ∂ y u| Σ + = e ıβLy ∂ y u| Σ - .
3.2 A Dirichlet-to-Neumann approach for the exact computation of guided modes in photonic crystal waveguides 79

with H 1 ( , B) = {u ∈ H 1 (B), u ∈ L 2 (B)}.
To solve (E), there are two different approaches: the ω-formulation (which consists in fixing β and looking for ω) and the β-formulation (which consists in fixing ω and looking for β). To simplify the presentation, we choose here the first one but the method extends to the other formulation -which could be more adapted for dispersive media for example: ρ(ω).

Using [Kuchment, 1993] and proving that A(β) is the compact perturbation of an operator with the perfectly periodic coefficient ρ p , we can show that the operator A(β) is selfadjoint in L 2 (B, ρdxdy), positive and its essential spectrum satisfies

σ ess (β) = R \ n∈ 1,N (β) ]a n (β), b n (β)[
where 0 ≤ a n (β) < b n (β) and N (β) ≤ +∞. The intervals ]a n (β), b n (β)[ are the gaps of the essential spectrum.

Remark 3.2.1

The classical characterization of the essential spectrum involves the eigenvalues of a cell problem with quasi periodic conditions. In Proposition 3.2.4, we will give another characterization of the essential spectrum with a by-product of the method.

Let us suppose now that at least one gap exists (see the introduction and Section 3.3 for sufficient conditions on the periodic medium which ensure the existence of gaps). We are interested in characterizing and then computing the eigenvalues (λ m (β)) m which are in the gaps of the essential spectrum (see the introduction and Section 3.3 for sufficient conditions which ensure the existence of eigenvalues). We can choose, using classical arguments that the dispersion curves β → λ m (β) are continuous, 2π/L y -periodic and even. We deduce in particular that it is sufficient to study the dispersive curves for β ∈ [0, π/L y ].

The non linear eigenvalue problem

From this point, we will focus on the eigenvalues ω 2 / ∈ σ ess (β).

We recall the reader the definition of the DtN operators associated to each half-band problems (see 1.3.4): • compute the DtN operators.

for any β ∈ [0, π/L y ], any α ∈ R + , a given ϕ in H 1/2 β (Γ ± ), let u ± ∈ H 1 (B ± ) be the solution of -u ± -ρ p α 2 u ± = 0 in B ± u Γ ± = ϕ
We have, indeed, explained in Section 1.3.4 that when α 2 does not lie in the essential spectrum of the periodic operator, the stationary Riccati equation has only one solution of spectral radius less than one and its spectral radius is strictly less than one. More precisely, we have the following result.

Proposition 3.2.4

The stationary Riccati equation associated to (P ± ) has a unique solution whose spectral radius is strictly less than one if and only if α 2 / ∈ σ ess (β).

The solution of the stationary Riccati equation enables then a characterization of the essential spectrum. If the spectral method described in Section 1.3.3 is used to solve this Riccati equation, it suffices to look at the solutions. More precisely, if the eigenvalues are all either inside or outside the unit circle then a solution of the Ricatti equation having a spectral radius strictly less than one can be constructed and then α 2 / ∈ σ ess (β). If one eigenvalue is on the unit circle, it is not possible to construct such solution and then α 2 ∈ σ ess (β). This characterization of the essential spectrum (whose computation is necessary) is then a by-product of the method. One does not have to compute it, as it is done classically, through the solution of a family of eigenvalue problems set on the periodicity cell.

Numerical results

We consider the periodic medium represented Figure 3.3. We plot on Figure 3.4 the isovalue lines For a fixed β, the white regions corresponds to the essential spectrum σ ess (β). Note that the function µ 1 (β, α)-α 2 vanishes several times in this interval but once per gap. The dispersion curves are given by the blue lines. We can also solve the Newton algorithm for fixed β's in order to compute the ω's such that ∃n, µ n (β, ω) = ω 2 . For such solution, from the associated eigenvector u 0 to µ n (β, ω), we can construct the guided mode by using (3.1). For two values of β, We have represent in Figure 3.5 the guided modes in eight periods from each side of the line defect.

Other related works

In [9], we have performed the numerical analysis of the DtN method and make precise comparisons with the supercell method. In [START_REF] Fliss | Robin-to-robin transparent boundary conditions for the computation of guided modes in photonic crystal wave-guides[END_REF], we have performed the mathematical and numerical analysis by using Robin-to-Robin operators, the motivation was to compute the edge resonances (or the Dirichlet surface modes). Finally in [25], we apply the method to 3D applications.

Trapped and guided modes in topographic periodic domains

In this section, we propose a family of periodic media for which the existence of gaps is ensured and we give also lineic perturbation for which guided modes exist. Compared to other similar results, the conditions on the periodic media and the lineic perturbation are quite simple.

More precisely, we consider a ladder-like periodic structure, namely a thin periodic structure (the thickness being characterized by a small parameter ε > 0) whose limit (as ε tends to 0) is a periodic graph. In other words, the periodic structure is R 2 minus a periodic (in the two directions) set of rectangles, the distance between two consecutive rectangles being ε. We introduce a perturbation consisting in changing the geometry of the reference medium by modifying the distance between two consecutive lines of rectangles (see Figure 3. 6). The question we investigate is whether such a geometrical perturbation is able to produce guided modes. We have investigated this question when the propagation model is the scalar Helmholtz equation with Neumann boundary conditions. This amounts to solving an eigenvalue problem for the Laplace operator in a band with quasi-periodic boundary conditions.

To address this question, we used a standard approach of analysis (used for instance in [Figotin andKuchment, 1996a, Nazarov, 2012a]) that consists of three main steps. We first find the formal limit of the eigenvalue problem as ε tends to 0. In the present case, it corresponds to an eigenvalue problem for a self-adjoint operator defined along the limit periodic graph ( [Rubinstein and Schatzman, 2001, Kuchment and Zeng, 2001, Saito, 2000[START_REF] Post | Spectral convergence of quasi-one-dimensional spaces[END_REF], Panasenko and Perez, 2007]). This limit operator consists of the second order derivative operator on each edge of the graph together with transmission conditions (called Kirchhoff conditions) at its vertices [Exner, 1996, Carlson, 1998, Kuchment and Zeng, 2001]. Then, we proceed to an explicit calculation of the spectrum of the limit operator using a finite difference scheme [START_REF] Avishai | [END_REF]Luck, 1992, Exner, 1995]. Finally, we prove that the spectrum of the initial operator is close to the spectrum of the limit operator. In particular, we prove the existence of localized modes provided that the geometrical perturbation consists in diminishing the width of one rung of the periodic thin structure. Moreover, in that case, it is possible to create as many eigenvalues as one wants, provided that ε is small enough. We obtain an asymptotic expansion at any order of the eigenvalues, which can be used for instance to compute a numerical approximation of these eigenvalues and associated eigenvectors. Numerical experiments illustrate the theoretical results.

This was done in the framework of the Phd of Elizaveta Vasilevskaya with who I have worked, with also her Phd advisors Patrick Joly (POEMS) and Bérangère Delourme (LAGA). We have written two papers and a research report [15] [19,23] on a similar situation -existence of trapped modes in a ladder-like domain infinite in one direction and bounded in the other one with Neumann boundary conditions and not quasi-periodic ones -with the same approach.

Model problem

We consider the propagation of acoustic waves in a particular periodic medium that consists of the plane R 2 minus an infinite set of equispaced perfect conductor rectangular obstacles with Neumann boundary conditions. The parameter ε represents the distance between the obstacles. We introduce a lineic defect in this perfectly periodic domain by changing the distance between two consecutive columns of obstacles from ε to µε, where µ > 0 (cf. Figure 3.6 for µ ∈ (0, 1)).

Our aim is to find guided modes, that is to say solutions of the homogeneous wave equation propagating along the defect. We know from the previous section that this problem can be reformulated as an eigenvalue problem for the Laplacian in the periodicity band Ω µ ε with βquasi-periodic boundary conditions in the y-direction (cf. Figure 3.7). More precisely, we seek the couples

(u ε , λ ε ) ∈ H 1 (Ω µ ε ) × R + satisfying -∆u ε = λ 2 ε u ε in Ω µ ε , (3.2)
together with β-quasi-periodicity boundary conditions on

Σ ± = ∂Ω µ ε ∩ {y = ±L/2} with β ∈ [0, 2π/L], u ε | Σ + = e ıβL u ε | Σ -, ∂ y u ε | Σ + = e ıβL ∂ y u ε | Σ -, (3.3)
and homogeneous Neumann boundary conditions on the remaining part of the boundary: In this section, we restrict ourselves to the case µ ∈ (0, 1). It is known that the spectrum of the operator A µ ε (β) approaches in some sense the spectrum of the operator A µ (β) as ε is small enough (see [Kuchment andZeng, 2001,Post, 2006] for more details). By constructing a so called quasi-mode of the A µ ε (β) from a mode of the limit operator A µ (β), we can obtain a more precise result:

∂ n u ε = 0 on ∂Ω µ ε \ (Σ + ∪ Σ -). ( 3
A µ ε (β) : L 2 (Ω µ ε ) → L 2 (Ω µ ε ), A µ ε (β)v = -∆v, D(A µ ε (β)) = v ∈ H 1 (∆, Ω µ ε ), ∂ n v| ∂Ω µ ε \(Σ + ∪Σ -) = 0, v| Σ + = e ıβL v| Σ -, ∂ y v| Σ + = e ıβL ∂ y v| Σ -.
Theorem 3.3.2 Let µ ∈ (0, 1). For any β ∈ [0, 2π/L], let {(a n 0 (β), b n 0 (β))
, n ∈ N * } be the gaps of the limit operator A µ (β). Then, for each n 0 ∈ N * , there exists ε 0 > 0 such that if ε < ε 0 , the operator It is worth noting that imposing Dirichlet boundary condition on the obstacles (instead of the Neumann one) would lead to a totally different asymptotic as ε goes to 0 ( [START_REF] Nazarov | On the spectrum of the laplace operator on the infinite dirichlet ladder[END_REF]).

A µ ε (β) has at least n 0 gaps {(a n ε (β), b n ε (β)), 1 ≤ n ≤ n 0 } such that a n ε = a n 0 + O(ε), b n ε = b n 0 + O(ε), ε → 0, 1 ≤ n ≤ n 0 .

Asymptotic analysis

The error estimates (3.6) obtained for the eigenvalues is suboptimal. In fact, using matched asymptotic expansions and writing a high order asymptotic expansion of λ ε restore the optimal convergence rate λ ε -λ = O(ε).

Moreover, we can go further. More precisely, we can show that the eigenvalue has the following asymptotic expansion:

∀ n ∈ N, λ ε = λ 0 + n k=1 λ k ε k + O ε n+1 . (3.7)
The coefficients {λ k } k∈N can be computed by an explicit recurrence procedure that involves:

• the solutions, called the far fields, of family of problems set on the graph G with jump conditions and non homogeneous Kirchhoff conditions at each vertices;

• the solutions, called the near fields, of family of elliptic problems set on infinite normalized junctions (two reference geometry (see Figure 3.9)) with polynomial growth at infinity ;

• the two families of problems are coupled through matching conditions linking the behaviour of the far fields on each vertices to the behaviour at each infinity of the near fields.

Numerical results

To illustrate and validate the existence of gaps for ε small enough and the existence of eigenvalues for µ ∈ (0, 1) and for ε small enough (see Theorem 3.3.2), in Figure 3.10, we have fixed µ = 0.25 and ε = 0.1 and have computed a part of the essential spectrum and eigenvalues of A µ ε (β) for different values of β, using the method presented in Section 3.2.

In Figure 3.11, for a fixed ε = 0.1 and β = 1, we have represented the position of the eigenvalues in the second gap of A µ ε (β) when µ varies from 0 to 1. We remark that when µ tends to 1, the eigenvalues tends to the boundary of the gaps. We know that for µ = 1, the operator is perfectly periodic and then cannot have any eigenvalues.

The different terms of the asymptotic expansion were computed for the case described in [15] [19,23], i.e. for the same problem replacing the quasi-periodic boundary conditions by Neumann ones. For a fixed µ, we compare the computations, for different values of ε, of one eigenvalue using the DtN approach (see Section 3.2) with the asymptotic expansion (3.7) of the eigenvalue at different order n = 1 to 5. We have validated with this numerical result represented in Figure 3.12, the computation of each term of the expansion. Surprisingly, the asymptotic expansion approaches the eigenvalue for quite large values of ε. The computation of the asymptotic expansion is much cheaper since it does not require a mesh adapted to the geometry of the ladder. It requires the solution of problem set on graph and the solution of problems set on normalized junctions, which are independent of ε.

Some perspective works

In all of the works that I have mentioned in this chapter, for instance the existence of gaps, the materials are frequency independent. However, the frequency dispersion can in general not be In presence of a boundary, an interface or more generally a lineic perturbation in a periodic medium, we have seen that energy localization can be created. The mathematical analysis and the numerical methods have been developed for frequency independent materials and for the case of periodic perturbations with the same periodicity than the bulk. Here again, a challenging question is to take the dispersion into account. But, more importantly, it seems that an essential property for applications is robustness: is that possible to propose sufficient conditions to ensure the existence of localized guided waves which are stable with respect to local imperfections of the lineic periodic perturbation (see Figure 3.13)? In that case, we speak about the so called topologically protected states, as in condensed matter physics. All the classical mathematical methods based on the Floquet-Bloch theory, used in the fore-mentioned papers for instance, fail to deal with these local imperfections. In the last decades, the most complete mathematical works on topologically protected states are [START_REF] Fefferman | Honeycomb lattice potentials and Dirac points[END_REF], Fefferman et al., 2014,Fefferman and Weinstein, 2014[START_REF] Fefferman | Edge States in Honeycomb Structures[END_REF],Fefferman et al., 2016a]. Until now, they have considered the Schrodinger operator with periodic potentials having a honeycomb structure. In [START_REF] Lee-Thorp | Elliptic operators with honeycomb symmetry: Dirac points, edge states and applications to photonic graphene[END_REF], the authors have begun the analysis for 2D Maxwell's equations. We want to use our expertise on periodic media and perturbations techniques to contribute to that subject, beginning with graph-type problems. This will be done in collaboration with Bérangère Delourme (LAGA, Paris 13).

For the computation of the localized guided waves, in presence of dispersion, since both the Supercell method and the DtN method would require the solution of non linear eigenvalue problems and since it seems more difficult to evaluate the confinement of the wave, the DtN method seems to be the most relevant method. This has to be investigated. Finally, to our knowledge, there exists no numerical method to compute the topological protected guided waves. This is a configuration which is closed to the ones mentioned in Sections 1.5 and 2.6 and for which we want to propose a numerical method.

Introduction

Recent discoveries have shown the possibility of producing weakly dissipative electromagnetic materials whose effective dielectric and magnetic constants have negative real parts. These "metamaterials", of complex multiscale structure, lead to extraordinary phenomena as regards the propagation of electromagnetic waves (negative refraction, resonance of "wavelength" cavities, etc.) and thus arouse great interest in view of many potential applications (super lenses, stealth coating, miniaturization of antennas, ...).

The optimization of devices exploiting or controlling these metamaterials requires the development of appropriate numerical methods. Their structure presenting several scales of very different size, it is very expensive or even impossible to simulate the wave propagation in these media taking into account all their complexity. An attractive alternative is to model the metamaterial by a homogeneous material, with physical constants of negative real part. This approach is now widely used by physicists and is the subject of active mathematical research in the homogenization community. For a good overview on the homogenization theory of elliptic system with periodic coefficients, we refer the reader to [START_REF] Bensoussan | Asymptotic analysis for periodic structures[END_REF], Cioranescu and Donato, 1999, Tartar, 2009, Jikov et al., 2012]. Thus, one can find in the literature that for certain periodic media whose structure has resonance mechanisms (being related to the geometry via Helmholtz resonators for example or to the characteristics of the materials), the dielectric permittivity [Felbacq andBouchitté, 1997, Silveirinha and[START_REF] Silveirinha | [END_REF] or magnetic permeability [START_REF] Pendry | Magnetism from conductors, and enhanced non-linear phenomena[END_REF], Lindell et al., 2001, Bouchitté and Felbacq, 2004, Bouchitté et al., 2017] or even both [Felbacq and Bouchitté, 2005b,Felbacq and Bouchitté, 2005a,Bouchitté et al., 2009,Bourel, 2010,Bouchitté and Bourel, 2012] can become negative for certain frequency ranges. Specific techniques of homogenization (reiterated homogenization techniques for instance) have to be introduced in order to take into account the resonances phenomena and the multi-scale effects. Convergence (of two-scale type) results have also been proved.

94

Chapter 4. Homogenization of transmission problems However, it is well known that classical homogenization process poorly takes into account boundaries or interfaces. This is particularly unfortunate when considering negative materials, because important phenomena arise precisely at their surface (plasmonic waves for instance) and it seems that the effective model may be imprecise or even completely false. Indeed, when we consider an interface between a dielectric and a metamaterial and that the permittivity and/or permeability contrast is equal to -1, it appears at the interface an accumulation of energy that is not compatible with the usual mathematical/physical framework [START_REF] Costabel | A direct boundary integral equation method for transmission problems[END_REF], Ola, 1995, Nguyen, 2015, Nguyen, 2016, Li, 2016]. It seems that these difficulties are due to an insufficiently fine asymptotic description of the propagation phenomena in the vicinity of the interfaces.

This is why we have proposed to revisit the asymptotic process in order to propose a new homogenized model that is simple to implement and more accurate near the boundaries or the interfaces.

Of course this subject is linked to the presence of boundary layers which appear when considering asymptotic model near boundaries or interfaces. It has already been pointed out for instance in [START_REF] Bensoussan | Asymptotic analysis for periodic structures[END_REF] and studied and analysed in [START_REF] Moskow | First order corrections to the homogenized eigenvalues of a periodic composite medium. the case of neumann boundary conditions[END_REF], Moskow and Vogelius, 1997a, Allaire and Amar, 1999, Birman and Suslina, 2006a, Birman and Suslina, 2006b, Gérard-Varet and Masmoudi, 2012, Gérard-Varet and Masmoudi, 2011, Shen and Zhuge, 2016, Armstrong et al., 2017] for elliptic systems with Dirichlet and Neumann conditions and in a more general setting in [START_REF] Blanc | On correctors for linear elliptic homogenization in the presence of local defects[END_REF]. Concerning transmission problems, very few results are available [START_REF] Bakhvalov | Homogenisation averaging processes in periodic media: mathematical problems in the mechanics of composite materials[END_REF] with the notable exception of the recent article [START_REF] Cakoni | On the homogenization of a scalar scattering problem for highly oscillating anisotropic media[END_REF].

Before dealing with metamaterials and Maxwell's equations, in the context of the Valentin Vinoles's PhD, we have began with an intermediate problem: the propagation of acoustic waves in harmonic regime between a homogeneous half-space and a standard periodic half-space (see Section 4.2). Because the periodic medium is standard, classical homogenization tools can be used. But even for this simplified problem, the effective model is less accurate in presence of interfaces and an enriched homogenized model can be proposed. Coupling a classical multiscale expansion in the periodic half-space with matched asymptotic techniques near the interface, we have proposed an appropriate and accurate (at any order) asymptotic model in Section 4.3. We can also derive high order transmission conditions. They are more precise than the classical transmission conditions (corresponding to the continuity of the solution and its normal derivative across the interface) but less standard since they involve differential operators along the interface (see Section 4.4). An error analysis confirms this accuracy and numerical results, shown in Section 4.5, illustrate the efficiency and the accuracy of these new conditions. This research topic is at its beginning and these first results seem very encouraging to continue in this direction. I will give in Section 4.6 my perspectives on this topic.

Let us mention that the results of the chapter are part of Valentin Vinoles's PhD and presented in two forthcoming papers [START_REF] Claeys | Homogenization of transmission problems between homogeneous and periodic halfspaces. part 1 : asymptotic expansions[END_REF][START_REF] Claeys | Homogenization of transmission problems between homogeneous and periodic halfspaces. part 2 : high order approximate model and numerical aspects[END_REF]. • the two functions (characterizing the material properties of the medium) a , ρ : R 2 → (0, +∞) are defined as follows

Model problem and classical results

Ω

- Ω + Γ ε supp f

Let us denote in the following Ω

± := {x = (x 1 , x 2 ) ∈ R 2 | ± x 1 > 0}, and Γ := {0} × R = ∂Ω ± . Let us consider the problem Find u ε ∈ H 1 (R 2 ) such that -∇ • (a ε ∇u ε ) -ω 2 ρ ε u ε = f in R 2
a ε (x) = a 0 1 Ω -(x) + a p ( x ε )1 Ω + (x) and ρ ε (x) = ρ 0 1 Ω -(x) + ρ p ( x ε )1 Ω + (x)
where 1 Ω ± is the characteristic function associated to Ω ± , a 0 , ρ 0 ∈ (0, +∞) are two positive constants, and a p , ρ p ∈ L ∞ (R 2 ) satisfies λ ≤ a p (x) ≤ Λ and λ ≤ ρ p (x) ≤ Λ for all

x ∈ R 2 (for fixed λ, Λ ∈ (0, +∞)) and they are 1-periodic functions in the two directions.

We denote Y = (0, 1) 2 the unit periodicity cell.

• a source functional f ∈ H -1 (R 2 ) such that Supp(f ) ⊂ Ω -(in particular Γ ∩ Supp(f ) = ∅)
This is pictured in Figure 4.1 above. This problem admits a unique solution according to Lax-Milgram's lemma. Besides satisfying the PDE both in Ω + and Ω -, any solution to the above equation must satisfy transmission conditions at the interface Γ, namely

u ε | + Γ = u ε | - Γ and e 1 • ∇u ε | + Γ = e 1 • ∇u ε | - Γ .
We are interested in studying the behaviour of the solution as → 0. This analysis and the study of this limit is given by the homogenization theory (see for instance [START_REF] Bensoussan | Asymptotic analysis for periodic structures[END_REF], Cioranescu and Donato, 1999, Tartar, 2009, Jikov et al., 2012]).

Classical results of the homogenization theory

To review the classical results of the homogenization theory, let us consider the unique solution

u ε in H 1 (R 2 ) of -∇ • (a ε ∇u ε ) -ω 2 ρ ε u ε = f in R 2 (4.2) with a ε (x) = a p ( x ε ) and ρ ε (x) = ρ p ( x ε ).
To study the behaviour of u ε when ε goes to 0, there exist several methods (two-scale methods, compound method, Floquet-Bloch approach). Here, we use the two-scale asymptotic expansion which consists in postulating an ansatz for the solution

u (x) ∞ n=0 n u n (x, x/ ) x ∈ R 2 where ∀n ≥ 0, u n (x, •) ∈ H 1 (Y). (4.3)
where H 1 (Y) is the closed subspace of H 1 (Y) containing the 1-periodic functions. We denote H 1 (Y) its dual space.

Plugging Ansatz (4.3) into (4.2) yields the following cascade of equations

L 0 u n + L 1 (∂ x )u n-1 + L 2 (∂ x )u n-2 = 0 n ≥ 0, x ∈ R 2 , y ∈ Y, where L 0 v := -∇ y • (a(y)∇ y v) where L 1 (∂ x )v := -∇ y • (a(y)∇ x v) -∇ x • (a(y)∇ y v)
where L 2 (∂ x )v := -a(y)∆ x v + ρ(y)v (4.4)

In the equations above, we took the usual convention that u -1 = u -2 = 0. Besides L j , j = 0, 1, 2 should be understood as polynomials in ∂ x , and as variational operators with respect to y. More precisely, for any η = (η 1 , η 2 ) ∈ C 2 , define L 0 , L 1 (η), L 2 (η) regarded as continuous operators mapping H 1 (Y) into H 1 (Y) and defined variationally by

L 0 v, w := Y a(y)∇ y v • ∇ y w dy L 1 (η)v, w := Y a(y)η • (v∇ y w + w∇ y v) dy L 2 (η)v, w := Y (ρ(y) -a(y) |η| 2 )v(y)w(y) dy ∀v, w ∈ H 1 (Y). (4.5)
Let us show how the calculation of the u k (x, y) that are functions depending a priori on both x and y can be reduced to the calculation on functions ûk (x) that depend on x only and on solution of cell problems set in Y .

First of all, observe that the operator L 0 is not invertible on H 1 (Y), its kernel consists in constant functions over the unit cell Y. According to Fredholm alternative, for g ∈ H 1 (Y) there exists u ∈ H 1 (Y) satisfying L 0 u = g if and only if g, 1 = 0. This is the compatibility condition associated to the problem L 0 u = g in Y . If this condition is satisfied, the solution is defined up to an element of the kernel of L 0 , i.e. up to a constant.

The process is then as follows : for each n, there exists a n-th term in the expansion (4.3) if and only if the compatibility condition associated to Problem (4.4) is satisfied. In that case, the solution is unique up to a function ûn (x) that depend on x only. More precisely,

• (4.4) for n = 0 gives that u 0 (x, y) = u 0 (x);

• (4.4) for n = 1 gives that u 1 (x, y) = w(y) • ∇u 0 (x) + û1 (x) where w(y) = (w 1 , w 2 ) and w i is the unique solution in

Ḣ1 (Y) = {w ∈ H 1 (Y), Y w = 0} of the cell problem L 0 w i = L 1 (e i )1
where e 1 = (1, 0) and e 2 = (0, 1) ;

• the compatibility condition of (4.4) for n = 2 gives that

-∇ x • (A * ∇ x u 0 (x)) -ω 2 ρ * u 0 (x) = f in R 2
where

A * = (a * j,k ) ∈ R 2×2 , a * j,k := Y a(y)(∂ y j w k (y) + δ j,k )dy, ρ * = Y ρ(y)dy. Then u 2 (x, y) = ∇ x • (w(y)∇ x u 0 )(x) + w(y) • ∇ x u 1 (x) + û2 (x)
where w = (w ij ) i,j and the w ij 's are solutions of cell problems whose right hand side depends on w ;

• by induction it can be shown that

u n (x, y) = n k=0 Θ k (y, ∂ x )û n-k (x). (4.6)
where

-Θ 0 = 1 and for all k ≥ 1, ∀η = (η 1 , η 2 ) ∈ C 2 , Θ k (y, η) is solution in Ḣ1 (Y )
of the cell problem parametrized by η and defined by induction

L 0 Θ k (y, η) := -( L1 (η)Θ k-1 (y, η) + L2 (η)Θ k-2 (y, η) ) with Li (η) = L i (η) -Y L i (η). One can prove by induction that η → Θ k (y, η) is a polynomial of degree k: Θ k (y, η) = |α|≤k Θ α k (y)η α , where η α = η α 1 1 η α 2 2 (note that we have Θ 1 (y, η) = w • η;
the ûn 's are solution in H 1 (R 2 ) to an homogenized second order elliptic problem with a right-hand side depending on the previous terms,

-∇ x • (A * ∇ x ûn (x)) -ω 2 ρ * ûn (x) = δ n0 f + n k=1 |α|≤k+2 Q α k ∂ α x ûn-k (x) in R 2
where the coefficients Q α k ∈ R admit explicit expressions in terms of the Θ α k 's. Note that, in the case n = 0, the sum vanishes. Although not obvious, the right hand side above also vanishes for n = 1. To be more specific, we have Q α 1 = 0, ∀α ∈ N 2 . This was established in [Moskow andVogelius, 1997a,Allaire et al., 2016]. We deduce that û1 = 0. Note that the right hand side is non-trivial for n ≥ 2.

Observe that (4.6) takes the form of a discrete convolution. Besides, ûn (x) is the mean value in Y of u n (x, y), the remaining part is called the oscillating part of u n or the n-th volume Chapter 4. Homogenization of transmission problems corrector of the expansion and is defined thanks to the derivatives of the previous terms.

Finally, error estimates validate the asymptotic expansion : for all N if f, a p , ρ p are sufficiently smooth

∀N ∈ N * , u ε - N n=0 n n k=0 Θ k (x/ , ∂ x )û n-k (x) L 2 (R 2 ) = O(ε N +1 ), and
∀N ∈ N * , u ε - N n=0 n n k=0 Θ k (x/ , ∂ x )û n-k (x) H 1 (R 2 ) = O(ε N ).

Application to the transmission problem

A natural idea would be to apply this approach to our transmission problem, the ansatz being introduced only in the periodic part of the domain Ω + . We obtain that u ε tends in L 2 to a function

u 0 ∈ H 1 (R 2 ) which is solution of -∇ • (A * 0 ∇u 0 ) -ω 2 ρ * 0 u 0 = f in R 2 (4.7)
where

A * 0 (x) := a 0 1 Ω -(x)Id + 1 Ω + (x)A * and ρ * 0 (x) := ρ 0 1 Ω -(x)Id + 1 Ω + (x)ρ * .
All the other terms of the expansion would live only in Ω + . For instance, the second term would be given by u 1 (x, y) = 1 Ω + (x) [w(y) • ∇u 0 (x) + û1 (x)] with û1 (x) satisfying

-∇ • (A * ∇û 1 ) -ω 2 ρ * û1 = 0 in Ω + .
Unfortunately, plugging (4.3) to the natural transmission conditions of (4.1) leads to boundary condition for û1 that depends on the fast variable x/ε, which is incompatible with the macroscopic nature of û1 . This incompatibility appears at any order. This shows that the classical asymptotic two-scale expansion is not adapted to the presence of an interface (or similarly to a boundary). Moreover, neglecting û1 by assuming that û1 = 0 as in the free space leads to the error estimates [START_REF] Cakoni | On the homogenization of a scalar scattering problem for highly oscillating anisotropic media[END_REF]]

u ε -u 0 L 2 (R 2 ) = O(ε) u ε -u 0 H 1 (Ω -) + u ε -u 0 -εu 1 (•, •/ ) H 1 (Ω + ) = O( √ ε).
Adding the volume correctors of higher order would not improve these convergence rates.

This deteriorated error estimate is well documented in the existing literature on homogenization theory. It is related to the presence of boundary layers at the interface Γ. There are already a lot of contributions on this topic, as already mentioned in the introduction. We propose in our work to use matched asymptotic expansion in order to propose a more accurate asymptotic description of the solution.

Remark 4.2.1

The method that we propose here can also be used for periodic halfspace problem with Dirichlet or Neumann boundary conditions.

The matched asymptotic expansion applied to the transmission problem

One cannot expect a simple asymptotic expansion that would be valid uniformly in the whole space. We propose to use different asymptotic expansions in different areas using matched asymptotic expansions [Ilin, 1992, Maz'Ya et al., 2012] in order to get a more precise homogenized model. For our problem, we are going to use three asymptotic expansions: two on each side of the interfaces and a third near the interface. These three expansions must coincide on intermediate areas (matching areas).

Our work differs from [START_REF] Cakoni | On the homogenization of a scalar scattering problem for highly oscillating anisotropic media[END_REF]: we do not seek to construct bulk correctors but several asymptotic expansions in different areas. Notably there are two novelty in our approach: [START_REF] Joly | Exact boundary conditions for periodic waveguides containing a local perturbation[END_REF] an asymptotic expansion near the interface and ( 2) an asymptotic expansion in the homogeneous medium. For more details on matched asymptotic expansions, we refer to [Ilin, 1992, Maz'Ya et al., 2012, Van Dyke, 1964]. Our work is inspired by previous works on diffraction of thin layers (periodic or not), see for instance [Schmidt, 2008, Delourme et al., 2012, Delourme et al., 2013, Claeys and Delourme, 2013].

Remark 4.3.1

This work relies on the matched asymptotic method, but one also can use multi-scales techniques. These two methods are actually "equivalent", see for instance [START_REF] Tordeux | Matching and multiscale expansions for a model singular perturbation problem[END_REF] for a precise comparison.

In this work, for simplicity we restrict ourselves to the 2D case, to the Helmholtz equation and to functions a ε that are scalar-valued, but all our results easily generalize to higher dimensions, to any other elliptic equation and to the case where a ε is a symmetric tensor.

In this section, we present first the formal steps of the matched asymptotic expansion, then the construction by induction of the different terms of each expansion and finally give the error estimates which validate the asymptotic expansion. From this asymptotic expansion, we can derive high order approximate problem which computes directly an approximation of the first terms of the asymptotic. This will be done in the next section where numerical results will be also shown. As indicated in Figure 4.2, we will distinguish different regions in which we postulate differ-are missing. They will be derived thanks to the matching condition steps.

Near field terms: Consider the different ansatz

u ε (x) ∞ n=0 ε n U n (x 2 , x/ε) x ∈ Ω 0 ε (4.11)
This time we assume for all n that y 2 → U n (x 2 , y 1 , y 2 ) is 1-periodic for all x 2 and y 1 . However, in contrast with Ansatz (4.9), periodicity is not assumed anymore with respect to the y 1 variable. Plugging (4.11) yields a new cascade of equations that we shall refer to as "near field equations"

L 0 U n + L 1 (∂ x 2 )U n-1 + L 2 (∂ x 2 )U n-2 = 0 n ≥ 0, x 2 ∈ R, y ∈ B := R × (0, 1), where L 0 v := -∇ y • (a(y)∇ y v) where L 1 (∂ x 2 )v := -∂ y 2 (a(y)∂ x 2 v) -∂ x 2 (a(y)∂ y 2 v) where L 2 (∂ x 2 )v := -a(y)∂ 2 x 2 v + ρ(y)v
This time, the natural domain for the equations is the domain B that is periodized according to the y 2 -variable only. In these equations, x 2 plays the role of a parameter. We do not have specified yet the functional space in which the near field terms have to be looked for. As we will see in the matching condition step, the near field terms may have a polynomial growth at infinity.

Matching conditions :

The different asymptotic expansions have to coincide in the regions where they coexist, i.e. the overlapping zones {(x 1 , x 2 ), -2η(ε) < x 1 < -η(ε)} and {(x 1 , x 2 ), 2η(ε) < x 1 < 2η(ε)}. We suppose that lim ε→0 η(ε) = 0 and lim ε→0 η(ε) ε = +∞ so that seen from the far field terms the overlapping zones tend respectively to Γ -and Γ + (where Γ ± = {(x 1 , x 2 ), x 1 = 0 ± }) when ε goes to 0 and seen from the near field terms, they tend respectively to -∞ and +∞. In consequence the matching conditions consists in linking Predominant behaviour of

N n=0 ε n U n (x 2 , x/ε) for ε → 0, x 1 /ε → ±∞ = Predominant behaviour of N n=0 ε n u n (x, x/ε) for ε → 0, x 1 → 0 ±
This relies on Taylor expansion of the far field terms ûn for x 1 = 0 ± and the behaviour at ±∞ of the near field terms U n .

Remark 4.3.2

The choice of η(ε) does not really matter for the formal derivation of the equations. It is important only to obtain optimal error estimates. To fix the ideas, η(ε) is typically ε γ for γ ∈ (0, 1).

Of course, as already discussed in the existing abundant literature on matched asymptotic method (see in particular [Ilin, 1992]) the matching principle is, at this stage of our analysis, a vague statement that only offers a guide to further calculus. Genuine justification of the whole analysis can only be obtained through error estimates that will be provided in Section 4.3.4.

To apply this principle, let us derive the predominant behaviour of the far field expansion at the interface. Examining the behaviour of an arbitrary function v(x, y) (supposed to be C ∞ with respect to x) Taylor expansion formally yields v(x 1 , x 2 , y) = +∞ p=0

x p 1 p! (∂ p x 1 v) Γ ± (x 2 , y) for x 1 → 0 ± with the trace operator on Γ ± defined as (∂ p x 1 v) Γ ± := (∂ p x 1 v)(x 1 = 0 ± , x 2 , y). Replacing formally x 1 = εy 1 , we can rewrite the previous relation as

Ξ ε (v) = +∞ p=0 ε p Ξ ± p (v) for εy 1 → 0 ± , ,
where Ξ ε v(x 2 , y) := v(εy 1 , x 2 , y) and Ξ ± p v (x 2 , y) :=

y p 1 p! (∂ p x 1 v) Γ ± .
(4.12)

We first examine the case ε → 0, x 1 < 0. To obtain formally the behaviour of u ε (x) = +∞ n=0 ε n ûn (x) as x 1 → 0, x 1 < 0, we first write u ε (x) = u ε (x 1 , x 2 ) = u ε (εy 1 , x 2 ), which boils down to applying the operator Ξ ε to the far field ansatz n ε n ûn (x). Letting ε → 0, we then simply use (4.12) for functions independent of y. This yields Note that the definition we take above for m - k (x 2 , y) holds for y 1 < 0.

Now we focus on the boundary behaviour of u ε for ε → 0 and x 1 → 0, x 1 > 0. In this part of the problem, we have to consider functions that depend both on x and y. Writing Θ p (resp. ûn ) instead of Θ p (y, ∂ x ) (resp. ûn (x)) for the sake of conciseness, leads to

Ξ ε (u ε ) = Ξ ε +∞ n=0 ε n u n = +∞ p=0 ε p Ξ + p +∞ n=0 ε n n q=0 Θ q ûn-q Ξ ε (u ε ) = Ξ ε +∞ n=0 ε n u n = +∞ k=0 ε k k p=0 Ξ + p • k-p q=0
Θ q ûk-p-q = 4.14) This time, the above definition holds for y 1 > 0. The expression (4.14) for m + k are less explicit compared with (4.13). To get a more concrete idea of these functions, let us give an explicit expression for the first terms. We have m + 0 (x 2 , y) = û+ 0 (0 + , x 2 ) m + 1 (x 2 , y) = û1 (0 + , x 2 ) + y 1 ∂ x 1 û0 (0 + , x 2 ) + w(y) • ∇û 0 (0 + , x 2 ) m + 2 (x 2 , y) = û2 (0 + , x 2 ) + y 1 ∂ x 1 û1 (0 + , x 2 ) + w(y) • ∇û 1 (0 + , x 2 ) m 2 (x 2 , y) = + y 2 1 /2 ∂ 2 x 1 û0 (0 + , x 2 ) + w(y) • ∇∂ x 1 û0 (0 + , x 2 ) + Θ 2 (y, ∂ x ) û0 (0 + , x 2 ).

(4.15)

We have then shown that Predominant behaviour of

N n=0 ε n u n (x, x/ε) for ε → 0, x 1 → 0 ± = N n=0 ε n m ± n (x 2 , y)
The matching conditions reduces to

U n (x 2 , y) ∼ y 1 →±∞ m ± n (x 2 , y).
The functions m ± k shall play a pivotal role in the subsequent analysis. We shall refer to them as "the matching functions".

Equations caracterising the near field terms

The previous computation leads to imposing the following equations for the caracterization of the near field expansion As is easily seen from Formula (4.15), the matching functions m ± n (x 2 , y) are definitely not bounded with respect to y 1 and admits a polynomial growth (combined with more complicated dependency in y inherited from the Θ j 's).

L 0 U n + L 1 (∂ x 2 )U n-1 + L 2 (∂ x 2 )U n-2 = 0 n ≥ 0,
As a consequence, the problems above defining the U n 's do not fit the standard variational framework for equations associated with the operator L 0 in B, since such a standard framework would enforce boundedness on the U n 's. So we have to introduce another framework allowing polynomially growing behaviour of solution to such equations as (4.16).

In these equations, the variable x 2 plays the role of a parameter; that is why the problems in this section do not involve it. Thus we focus on problems of the form: "find u that is 1-periodic with respect to y 2 and such that -∇ • (a∇U ) = g" posed in the infinite strip B = R×(0,1) (see Figure 4.3). Before imposing any behaviour at infinity, we are first interested in the behaviour when y 1 tend to ±∞ of such solutions.

When the function a is constant at infinity as in previous works on thin layers (see for instance [START_REF] Claeys | High order asymptotics for wave propagation across thin periodic interfaces[END_REF],Delourme et al., 2012,Delourme et al., 2013,Schmidt, 2008]), analytic computations using Fourier series are available. In our problem the function a is periodic in the right part of the strip B. One cannot expect explicit representations of the solutions.

where we have used that U 0 (x 2 ) = û0 (x 2 ). We are not in the appropriate framework to apply Theorem 4.3.3 since the r.h.s is not in V 1 β (B) * (it is not vanishing at +∞). But by using a truncature function χ(y 1 ) = 0 for y 1 < 1 and χ(y 1 ) = 1 for y 1 > 2, we can write a problem for Ũ1 = U 1 -χw 2 ∂ x 2 u 0 as follows

L 0 Ũ1 = ∂ y 2 a -L 0 (χw 2 ) ∂ x 2 u 0 x 2 ∈ R, y ∈ B Ũ1 (x 2 , y) ∼ y 1 →-∞ y 1 (∂ x 1 û0 ) Γ -+ û1 Γ - Ũ1 (x 2 , y) ∼ y 1 →+∞ y 1 (∂ x 1 û0 ) Γ + + w 1 (y)(∂ x 1 û0 ) Γ + + û1 Γ + .
(4.17) Applying Theorem 4.3.3 with g = ∂ y 2 a -L 0 (χw 2 ), α + = û1 Γ + , α -= û1 Γ -, β -= a 0 ∂ x 1 û0 Γ -and β + = A * 11 ∂ x 1 û0 Γ + we have that there exists a unique Ũ1 ∈ V 1 β (B) solution of the previous near field problem if and only if (the proof is not straighforward)

[û 1 ] Γ = C 1 ∂ x 1 û0 + C 2 ∂ x 2 û0 Γ and [∂ n A u 0 ] Γ = 0
where C 1 depends on the profile function N and C 2 depends on N and another profile function (defined thanks to Ũ1 ).

• This can be done at any order. As you have noticed for n = 1, a difficulty is that the r.h.s of the near field problem (4.16) is not in V 1 β (B) * for n ≥ 2 and the behaviour at infinity of the near field terms are not only linear but polynomial. But by definition of the matching functions m ± n from the far field terms ûn (solution of the far field equations (4.8) and (4.5)), it is easy to show that the matching functions m ± n , n ≥ 0 satisfy the recurrent system of PDEs satisfying by the near field terms respectively in B ± = B ∩ {±y 1 > 0}. By using truncature functions, it is then possible to introduce a problem which fits in the framework of Theorem 4.3.3 (as for n = 1). The compatibility conditions yield finally to relations linking [∂ nA ûn-1 ] Γ and [û n ] Γ to (normal and tangential) derivatives of the previous far field terms.

Equations on the far field terms

The uniqueness and existence of the near field terms satisfying the cascad of equations and the matching conditions yield to conditions on the jump of each far field term and its conormal derivative across Γ.

• The first far field term û0 is solution of

-∇ • (A * 0 ∇u 0 ) -ω 2 ρ * 0 u 0 = f, in R 2 , [û 0 ] Γ = 0, [∂ n A û0 ] Γ = 0.
This is exactly the classical homogenized transmission problem. Obviously this problem is well-posed in H 1 (R 2 ).

• The second far field term û1 is solution of

-∇ • (A * 0 ∇u 1 ) -ω 2 ρ * 0 u 1 = 0, in Ω -∪ Ω + , [û 1 ] Γ = C 1 ∂ x 1 û0 + C 2 ∂ x 2 û0 , [∂ n A û1 ] Γ = C 3 ∂ 2 x 1 x 2 û0 + C 4 ∂ 2 x 2 û0 + C 5 ω 2 û0 .
This far field term is not continuous across Γ and its co-normal derivative as well. This problem is well-posed in V = {v ∈ L 2 (Ω), v Ω ± ∈ H 1 (Ω ± )}. All the constants C j are defined thanks to cell solutions and solutions of problem set in the band B.

• By induction, the n-th far field terms can be determined thanks to the previous terms through the equation set on R 2 , for n ≥ 2 where in T 1 k (resp. T 2 k ) are differential operator on Γ ± of order at most k (resp. at most k + 1).

Algorithm of construction of the expansion and error estimates

The construction of each far field term and each near field term is done by induction as follows. For all n ≥ 0 1. the far field term ûn is solution of the transmission problem (4.18) (which depends on the previous terms);

2. the transmission conditions of ûn (jump condition) and of ûn-1 (the jump condition on the conormal derivative) are necessary and sufficient condition for existence and uniqueness of the near field term U n .

We can propose a global approximation of the solution u ε . Introducing three smooth cut-off functions ψ ± , χ : R → [0, 1] defined by ψ + (t) := 1 for t ≥ 2, 0 for t ≤ 1, ψ -(t) := 1 for t ≤ -2, 0 for t ≥ -1, and χ := 1 -(ψ -+ ψ + ),

the global expansion at order N of u ε is given by

u ε,N (x) = ψ -( x 1 η(ε) ) N n=0 ε n u n (x) + χ( x 1 η(ε) ) N n=0 ε n U n x 2 , x ε + ψ + ( x 1 η(ε) ) N n=0 ε n u n x, x ε .
Error estimates validate the asymptotic expansion : for all N if f, a p , ρ p are sufficiently smooth

∀N ∈ N * , u ε -u ε,N L 2 (R 2 ) = O(ε N +1 ),
and

∀N ∈ N * , u ε -u ε,N H 1 (R 2 ) = O(ε N ).

A high order approximate homogenized model

We have explained in the previous section how to construct each term of the asymptotic expansion by induction. One can be interested in deriving a problem whose solution is an approximation of the first N far field terms for instance, the near field approximation being computed a posteriori if needed. This problem corresponds to an approximate model of order N. Let us explain how to derive an approximate model of order 1 (the approximate model of order 0 being the classical homogenized transmission problem).

A natural better approximation than û0 would be ûε,1 = û0 (x) + ε û1 (x) in Ω - û0 (x) + ε w( x ε ) • ∇û 0 (x) + û1 (x) in Ω + .

By using the problems satisfied by û0 and û1 , we find

-∇ • (A * 0 ∇û ε,1 ) -ω 2 ρ * 0 ûε,1 = O(ε 2 ), in Ω -∪ Ω + , [û ε,1 ] Γ = C 1 ε < ∂ n A ∂ x 1 ûε,1 > Γ +C 2 ε < ∂ x 2 ûε,1 > Γ +O(ε 2 ), [∂ n A ûε,1 ] Γ = C 3 ε < ∂ 2 x 1 x 2 ûε,1 > Γ +C 4 ε < ∂ 2 x 2 ûε,1 > Γ +C 5 ω 2 ε < ûε,1 > Γ +O(ε 2 )
where < v > Γ = (v + + v -)/2. By neglecting the O(ε 2 ) terms, we find the homogenized transmission problem of order 1

-∇ • (A * 0 ∇û app ε,1 ) -ω 2 ρ * 0 ûapp ε,1 = 0, in Ω -∪ Ω + , [û app ε,1 ] Γ = C 1 ε < ∂ n A ∂ x 1 ûapp ε,1 > Γ +C 2 ε < ∂ x 2 ûapp ε,1 > Γ , [∂ n A ûapp ε,1 ] Γ = C 3 ε < ∂ 2 x 1 x 2 ûapp ε,1 > Γ +C 4 ε < ∂ 2 x 2 ûapp ε,1 > Γ +C 5 ω 2 ε < ûapp ε,1 > Γ + (4.19)
The volume equation is as simple as the homogenized transmission problem of order 0. The problem depends simply on ε without introducing a microscopic scale. In the transmission conditions, differential operator of the interface up to order 2 are involved. The constants C j can be computed by solving cell problem (as in classical homogenization) and strip problems (taking into account the position of the interface).

The appropriate functional framework is

Ṽ = {v ∈ L 2 (R 2 ), v Ω ± ∈ H 1 (Ω ± ), v Γ ± ∈ H 1 (Γ)}.
Unfortunately, this problem is not necessarily well posed in Ṽ (the sign of the constant C 4 is not a priori known). Adapting [START_REF] Delourme | Approximate models for wave propagation across thin periodic interfaces[END_REF], we write the transmission conditions in two separated boundaries Γ ± α = Γ ± ± α ± ε, α ± > 0 on both sides of the interface, for instance

[u] α = u Γ + α -u Γ - α
and use Taylor expansions, for instance

u Γ + α = u Γ + + α + ε ∂ x 1 u Γ + + O(ε 2 ).
The transmission conditions are similar to the ones of (4.19) but the constants C j are replaced by constants C α j depending on α ± . For well-chosen α ± , the sign of the constant C α controlled and the associated problem is well posed.

Up to this change, error estimates validate that we have at hand a better approximation (than û0 ) of u ε . Indeed, let v app ε,1 (x) = ûapp ε,1 (x) + ε w(

x ε ) • ∇û app ε,1 (x) 1 Ω + (x), x ∈ R 2
we can show that for any open set O ⊂ Ω -∪ Ω + , we have

u ε -v app ε,1 L 2 (O)
= O(ε 2 ),

and u ε -v app ε,1 H 1 (O)
= O(ε).

The latter estimate becoming of order 2 if we add the second volume corrector in Ω + : let us define

v app ε,2 = v app ε,1 + ε 2 Θ 2 ( x ε , ∂ x )û app ε,1 (x) 1 Ω + (x)
then we prove that

u ε -v app ε,2 H 1 (O)
= O(ε 2 ).

Near field terms can be computed a posteriori in order to have an approximation of u ε near Γ and by using truncature functions, one can recover a good approximation in the whole domain.

Numerical results and validation

To compute the approximate problem, it suffices then to 1. solve the first and second family of cell problems appearing in the classical homogenization theory;

2. solve problems set on the strip B (for instance the problem given in Proposition 4. For a fixed ε, the exact solution u ε can also be computed using the method described in Section 1.3.5.

We consider two different media characterized by a function a ε represented Figure 4.4 (and ρ ε = 1 in these examples). We have compared for different values of ε

• the classical homogenized transmission problem : we represent the difference x → u ε (x) -û0 (x) (where û0 is solution of (4.7)) and in order to recover a part of the oscillation of u ε (and have a better H 1 -approximation), we compute a posteriori the first volume corrector and plot x → u ε (x) -û0 (x) -εw(x/ε) • ∇û 0 (x) 1 Ω + (x) ;

• the high order homogenized transmission problem : we represent the difference x → u ε (x)ûε,1 (x) (where ûε,1 is solution of (4.19)) and in order to recover a part of the oscillation of u ε (and have a better H 1 -approximation), we compute the first volume corrector and plot x → u ε (x) -ûε,1 (x) -εw(x/ε) • ∇û ε,1 (x) 1 Ω + (x) .

The differences are represented with the same color scales.

The qualitative comparison are clear : we see that the classical homogenized transmission problem is not accurate near the interface and also in the homogeneous medium. The volume corrector helps for a better accuracy in the periodic medium. By contrast, the high order homogenized transmission problem is more accurate near the interface and in the homogeneous medium. Let us mention that even if the error estimates have been obtained only for sufficiently smooth coefficients, the high order homogenized transmission problem gives good results.

Some perspective works

This research topic is at its beginning and the first results seem very encouraging. Clement Beneteau's thesis have begun in October 2017 to deal with the following questions: (1) the extension to Maxwell's equations, [START_REF] Fliss | Exact boundary conditions for time-harmonic wave propagation in locally perturbed periodic media[END_REF] the transmission problem when a half-space is a metamaterial (a resonant periodic medium). We want also to extend the method to more general geometry than halfspaces. This extension is not simple at all. First, it would be necessary to treat the case of two half-spaces but with an interface which does not cut the periodic medium in a direction of periodicity. This situation is a real challenge for different reasons. First of all, the calculation of the exact solution is already interesting, this subject was already mentioned in Section 1.5. Moreover, to obtain the high order homogenization model, it is necessary to analyze accurately the behavior of the near field terms. This analysis was performed recently for the case of Dirichlet boundary conditions [START_REF] Gérard-Varet | Homogenization in polygonal domains[END_REF],Gérard-Varet and Masmoudi, 2012,Armstrong et al., 2017]. From a numerical point of view, until now, no method has been proposed to compute these near field terms. I think that in the recent papers [Gérard-Varet andMasmoudi, 2011, Gérard-Varet andMasmoudi, 2012], we can find a promising lead.

Concerning the application to metamaterials, we will derive new effective models with high order transmission conditions and the mathematical analysis of these models has to be done. Indeed, we have mentioned in the introduction that when we consider an interface between a dielectric and a metamaterial, the problem is not well-posed in a classical mathematical framework when the permittivity and/or permeability contrast is equal to -1. We want to show that these high order transmission conditions regularize the problem.

In addition, our approach can be used to derive boundary or transmission conditions for all the terms involved in the two-scale asymptotic expansion of the solution. For instance, the third term, u 2 (x, y), is really important if one needs an accurate model of the time-domain wave equation at long time scale (see for instance [START_REF] Santosa | A dispersive effective medium for wave propagation in periodic composites[END_REF], Dohnal et al., 2014, Allaire et al., 2016]). This model is obtained by introducing the equation satisfied approximately by û0 + ε 2 û2 . The model is not hyperbolic anymore since an operator of order 4 is involved. For the moment these models were derived and analyzed in infinite domain. With our approach, we can deal with the presence of boundaries or interfaces. The mathematical analysis of the problem has to be done.

As the construction of the transmission conditions requires the a priori solution of cell and band problems, it could be interesting to propose a numerical method which computes directly the more accurate effective solution, by adapting the FE-HMM which computes the classical effective solution (see [14] for the extension of the classical FE-HMM to Maxwell's equations).

Introduction Chapitre 1 Chapitre 2

 12 Scattering problem in waveguides 1.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.2 Classical results for homogeneous acoustic waveguides . . . . . . . . . . . . . . . 1.3 Periodic waveguides . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.3.1 Model problem . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.3.2 Limiting absorption principle and radiation condition for perfectly periodic waveguides . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.3.3 Construction of transparent boundary conditions in presence of dissipation 1.3.4 Transparent boundary conditions in absence of dissipation . . . . . . . . . 1.3.5 Application to the transmission problem between periodic half-spaces . . . 1.4 A modal transparent boundary conditions for general waveguide problems . . . . 1.4.1 The Poynting-to-Neumann map for isotropic homogeneous acoustic waveguides . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.4.2 Construction of the PtN operator for periodic waveguides . . . . . . . . . 1.4.3 The PtN operator for anisotropic elastic waveguide . . . . . . . . . . . . . 1.5 Ongoing works and perspectives . . . . . . . . . . . . . . . . . . . . . . . . . . . . Scattering problem in infinite complex media 2.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

Figure 1 . 1 :

 11 Figure 1.1: (Left) Straight waveguide, (right) general waveguide

Figure 1 . 2 :

 12 Figure 1.2: Geometry and notation

Figure 1 .

 1 Figure 1.3: (Left) Perfectly periodic waveguide, (right) Junction of two periodic waveguides

Figure 1 . 4 :

 14 Figure 1.4: Notation for the locally perturbed waveguide problem.

Figure 1 . 5 :

 15 Figure 1.5: Solutions of the cell problems and associate local DtN operators

Figure 1 . 6 :

 16 Figure 1.6: The solution of the quadratic eigenvalue problem for different values of ε = 0.5, 0.1, 0.01, in red the eigenvalues of the propagation operator.

Figure 1 . 7 :

 17 Figure 1.7: The solution (bottom figure) of the Helmholtz equation with dissipation in a locally perturbed periodic waveguide (whose coefficient is represented in the top figure) for a compactly supported source (represented in the middle figure).

Figure 1 . 8 :

 18 Figure 1.8: The solution of the quadratic eigenvalue problem for different values of ε = 0, in red the eigenvalues of the propagation operator which satisfy the condition (1.30).

  [START_REF] Chandler-Wilde | Existence, uniqueness, and variational methods for scattering by unbounded rough surfaces[END_REF], Chandler-Wilde and Elschner, 2010]) or boundary integral formulations (see for instance [Chandler-Wilde

Figure 1 .

 1 Figure 1.9: The total field (bottom figure) of the diffraction problem between two different periodic half-guides (whose coefficient is represented in the top figure) for the incident field (represented in the middle figure).

Figure 1 . 10 :

 110 Figure 1.10: Transmission problem between two periodic media (same period in the direction of the interface).

Figure 1 . 11 :

 111 Figure 1.11: Negative refraction phenomena. The solution (right figure) of the transmission problem between an homogeneous media and a periodic one (whose coefficient is represented in the left figure) for the incident field (a gaussian beam represented in the middle figure).

Figure 1 . 12 :

 112 Figure 1.12: Superlens phenomena. The solution (right figure) of the transmission problem between an homogeneous media and a periodic one (whose coefficient is represented in the left figure) for the incident field (a gaussian beam represented in the middle figure).

Figure 1 .

 1 Figure 1.13: A reference solution computed using DtN operator (left figure) with the solution computed using PtN transparent boundary conditions (middle figure) and their difference (right figure).

Figure 1 . 14 :

 114 Figure 1.14: Modulus of the computed solution in the 2D anisotropic waveguide for different Ω a . Solution reconstructed in the half-guides Ω ± a

Figure 1 . 15 :

 115 Figure 1.15: Modulus of the computed solution in the 3D anisotropic waveguide for 3 different size of the computational domain.

Figure 1 . 16 :

 116 Figure 1.16: The periodic media is not periodic in the direction of the interface.

Figure 1 . 17 : 2

 1172 Figure 1.17: Examples of junction of open waveguides for which the embedding medium is homogeneous. PML can be used in the transverse directions.

Figure 2 . 1 :

 21 Figure 2.1: The notation defined in (2.3)-(2.4)-(2.5).

Figure 2 . 2 :

 22 Figure 2.2: Second set of notation (2.17) (where the dashed lines are the Σ j a 's defined in (2.4)).

1 ,

 1 is the set of 2 obstacles, Neumann boundary conditions are imposed on it) with the following parameters ω = 10 + 0.001ı, A = 1 and ρ = 1.

Figure 2 .

 2 Figure 2.3: On the left: Real part of the computed trace ϕ 0 on Σ 0 a (the black line) and real part of the exact solution (the dashed red line). On the right: the computed solution in Ω b .

Figure 2 . 4 :

 24 Figure 2.4: Top left: Reconstruction of the solution in the right halfspace using ϕ 0 and Formulae (2.6). Top right: Reconstruction in two halfspaces : the representations coincide. Bottom left: Reconstruction of the solution outside Ω a . Bottom right: the discretization is not fine enough so that the halfspace representations do not coincide in the quarter planes.

Figure 2 . 5 :

 25 Figure 2.5: Real part of the solution of (2.1)-(2.26) in Ω b and reconstruction of the solution in the halfspace Ω j a .

Figure 2 . 6 :

 26 Figure 2.6: The solution in an anisotropic non dissipative case with ω = 10

Figure 2 . 7 :

 27 Figure 2.7: Isotropic case : we represent Reu x (left), Reu y (center) and |u| (right) for a source term with a null rotational (first line, P wave) and with a null divergence (second line, S wave).

Figure 2 . 8 :

 28 Figure 2.8: Anisotropic case : we represent Reu x (left), Reu y (center) and |u| (right) for a source term with a null rotational (first line, Quasi-P wave) and with a null divergence (second line, Quasi-S wave).

Figure 2 .

 2 Figure 2.9: Isotropic case : we represent |u| computed with the HsMM (left) and the PML (right).

Figure 2 . 10 :

 210 Figure 2.10: Anisotropic case : we represent |u| computed with the HsMM (left) and the PML (right).

Figure 2 . 11 :

 211 Figure 2.11: Sketch of the locally perturbed periodic media.

Figure 2 . 12 :

 212 Figure 2.12: Right figure : Solution of (2.30) for A = I and ρ represented in the left figure.

Figure 3 . 1 :

 31 Figure 3.1: Domain of propagation : typically ρ = 1 in the white region, ρ = 2 in the dark grey regions and ρ = 3 in the light grey regions.

Figure 3 . 2 :

 32 Figure 3.2: The band B.

u

  y=Ly/2 = e ıβLy u y=-Ly/2 , ∂ y u y=Ly/2 = e ıβLy ∂ y u y=-Ly/2 . (P ± ) where B ± = B ∩ Ω ± , Γ ± = {±a} × (-L y /2, L y /2) and the space H 1/2 β (Γ ± ) is the set of traces of functions in H 1 (B ± ) satisfying the β-quasi-periodic boundary conditions u| y=Ly/2 = e ıβLy u| y=-Ly/2 .

Figure 3 . 3 :

 33 Figure 3.3: Isovalues of the photonic crystal waveguide

  Figure 3.5(a) corresponds to the eigenvalue in the first gap of A(0.5) with a well confined mode whereas in Figure 3.5(b) 3.2 A Dirichlet-to-Neumann approach for the exact computation of guided modes in photonic crystal waveguides 83 the eigenvalue belongs to the fourth gap of A(1.9)) and the associated guided mode is not well confined. For the position of the eigenvalues, see the dark point in Figure 3.4.

Figure 3 . 4 :

 34 Figure 3.4: Contours of log(µ 1 (β, α) -α 2 ); β ∈ [0, π/L], α 2 ∈ [0, 20]. The two dark points represent the values ω with µ 1 (β, ω) = ω 2 for which the guided modes are represented in Figure 3.5

Figure 3 . 5 :

 35 Figure 3.5: Modes: well confined (left); not well confined (right).

Figure 3 . 6 :

 36 Figure 3.6: Periodic and perturbed domains

. 4 )

 4 Figure 3.7: The periodicity band Ω µ ε

Figure 3 . 8 :

 38 Figure 3.8: Limit graph

  and the operator A µ ε (β) has 1 or 2 eigenvalues inside each(a n ε (β), b n ε (β)). Moreover, if λ ∈ (a n 0 (β), b n 0 (β)), n ≤ n 0 , is an eigenvalue of A µ (β), there exists an eigenvalue λ ε of A µ ε (β) such that λ ε -λ = O

Figure 3

 3 Figure 3.9: The two reference junctions

Figure 3 . 10 :Figure 3 . 11 :

 310311 Figure 3.10: For a ε = 0.1, µ = 0.25, different values of β, (blue) essential spectrum of A µ ε (β) and (red asterisks) eigenvalues of A µ ε (β)

Figure 3 . 12 :

 312 Figure 3.12: Computation of the eigenvalue λ ε for different values of ε = δ and comparison with the asymptotic expansion at different order n = 1 to 5

Figure 3 . 13 :

 313 Figure 3.13: Robustness of the localized guided modes : under which conditions are they stable with respect to local imperfections?

Figure 4 . 1 :

 41 Figure 4.1: Geometry of the transmission problem (4.1).

•

  Im ω 2 > 0;

Figure 4 . 2 :

 42 Figure 4.2: Matched asymptotic expansion for the transmission problem (4.1).

  p • Θ m-p ûk-m for y 1 > 0.

(

  x 2 ∈ R, y ∈ B := R × (0, 1) U n (x 2 , y) ∼

-a 0

 0 ∆û n -ω 2 ρ 0 ûn = 0 in Ω - -∇ • (A * ∇û n ) -ω 2 ρ * ûn =

  3.1 satisfied by the profile function N or Problem (4.17)), the combined difficulties being (1) the polynomial growth of the solution at infinity (2) the half-infinite periodic part of the strip. We use DtN operators (see Chapter 1) to restrict the computation on the interface between the homogeneous half-strip and the periodic half-strip (see Figure4.5 for 2 examples of profile functions for the medium represented in Figure 4.4(b) );3. compute all the constants appearing in the transmission conditions; 4. solve the approximate problem with the non classical transmission conditions involving differential operators on the interface.

Figure 4 . 4 :

 44 Figure 4.4: Level set of two functions a ε .

Figure 4 . 5 :

 45 Figure 4.5: Two examples of profile functions for the medium represented in Figure 4.4(b).

Figure 4 . 6 :

 46 Figure 4.6: For the case (1) (medium represented in Figure 4.4(a)) and ε = 0.5, comparison of classical homogenized transmission problem (bottom figures) and the high order homogenized transmisssion problem (top figures). The same color scales are used.

Figure 4 . 7 :

 47 Figure 4.7: For the case (2) (medium represented in Figure 4.4(a)) and ε = 0.5, comparison of classical homogenized transmission problem (bottom figures) and the high order homogenized transmisssion problem (top figures). The same color scales are used.

Figure 4 . 8 :

 48 Figure 4.8: For the case (2) (medium represented in Figure 4.4(a)) and ε = 0.25, comparison of classical homogenized transmission problem (bottom figures) and the high order homogenized transmisssion problem (top figures). The same color scales are used.

  est en cours de rédaction avec Khac-Long Nguyen(LAMCOS, INSA Lyon). Mes perspectives aujourd'hui sur ce sujet se situent autour de la prise en compte de la dispersion et l'étude des ondes guidées dans des structures périodiques hexagonales (de type graphene) qui seraient, d'après la bibliographie abondante sur le sujet, topologiquement protégées (c'est à dire stables pour des petites perturbations localisées du milieu). C'est un travail réalisé en collaboration avec Bérangère Delourme. Il semble que ces difficultés soient dues à une description asymptotique insuffisamment fine des phénomènes de propagation au voisinage des interfaces. asymptotiques raccordés, nous avons proposé des conditions de transmission plus précises mais moins standards puisqu'elles font intervenir des opérateurs différentiels le long de l'interface. Une analyse d'erreur confirme cette précision et des résultats numériques illustrent l'efficacité de ces nouvelles conditions. Deux articles[START_REF] Claeys | Homogenization of transmission problems between homogeneous and periodic halfspaces. part 1 : asymptotic expansions[END_REF][START_REF] Claeys | Homogenization of transmission problems between homogeneous and periodic halfspaces. part 2 : high order approximate model and numerical aspects[END_REF] sont en cours de rédaction. Dans le cadre de la thèse de Clément Benneteau (débutée en septembre 2017), nous voulons traiter le cas des équations de Maxwell et en présence de métamatériaux. Sur ce sujet, il y a d'autres perspectives qui permettront de traiter des configurations plus générales (le cas d'interface courbe par exemple.

Chapitre 4 -Homogénéisation en présence de bords ou d'interfaces : ce sujet a été motivé par le projet ANR METAMATH (2012-2016) dont j'ai été le porteur qui traitait de la modélisation mathématique et numérique pour la propagation des ondes en présence de métamatériaux. Des découvertes récentes ont montré la possibilité de réaliser des matériaux électromagnétiques faiblement dissipatifs, dont les constantes diélectriques et magnétiques effectives ont des parties réelles négatives. Ces « métamatériaux », de structure multi-échelle complexe, conduisent à des phénomènes extraordinaires en ce qui concerne la propagation des ondes électromagnétiques (réfraction négative, résonance de cavités « sous longueur d'onde»,...) et suscitent donc un grand intérêt en vue de nombreuses applications potentielles (super lentilles, revêtement furtif, miniaturisation des antennes,. . . ). Leur structure présentant plusieurs échelles de taille très différente, il est très coûteux voire impossible de simuler la propagation des ondes dans ces milieux en prenant en compte toute leur complexité. L'alternative séduisante consiste à modéliser le métamatériau par un matériau homogène, à constantes physiques de partie réelle négative. Ainsi, on trouve dans la littérature que pour certains milieux périodiques dont la structure présente des mécanismes de résonance (étant liés à la géométrie ou aux caractéristiques des matériaux), les permittivité diélectrique et/ou perméabilité magnétique effectives pouvai(en)t devenir négatives pour certaines gammes de fréquences. Cependant, le modèle homogénéisé est souvent obtenu en négligeant les effets de bords et c'est d'ailleurs pour cette raison qu'il est beaucoup moins précis aux bords ou aux interfaces des matériaux. Comme un certain nombre de phénomènes intéressants apparaissent à la surface des métamatériaux (comme la propagation des ondes plasmoniques par exemple), il semble que le modèle effectif classique peut s'avérer peu précis voire complètement faux. En effet, lorsqu'on considère une interface entre un diélectrique et un métamatériau et que le contraste de permittivité et/ou de perméabilité est égal à -1, il apparaît à l'interface une accumulation d'énergie qui n'est pas compatible avec le cadre mathématique/physique usuel. ) Ce document n'est pas une présentation chronologique et exhaustive de mes travaux de recherche mais un résumé des sujets que je viens de mentionner pour lesquels je positionne mon travail en faisant référence à des résultats classiques de la littérature, je présente mes contributions et j'ai également inclus des résultats qui viennent d'être soumis ou le seront très bientôt. J'aimerais enfin signaler que la bibliographie est séparée en 2 catégories. Mes papiers sont cités avec des nombres [N] et le reste de la littérature est citée avec le premier auteur et la convention et al.

[First Author et al., YEAR]

.

It is intricate but not impossible. Indeed in that case, the test functions have to be in H -1/2 so (1) the equations are written piece by piece so a particular attention must be paid at the junctions and (2) inf-sup conditions have to be satisfied (and the discrete version of them).

Chapter 3. Trapped and guided modes in periodic structures Theorem 3.2.1 (Well-posedness of the problems (P ± ) [START_REF] Fliss | A dirichlet-to-neumann approach for the exact computation of guided modes in photonic crystal waveguides[END_REF])

∈ σ ess (β), the problem (P ± ) is well-posed in H 1 except for a countable set of frequencies which depends on β (called the edge resonances). If the periodicity cell is symmetric with respect to the axis x = 0 and if ω 2 / ∈ σ ess (β), the problem (P ± ) is always well-posed in H 1 .

Remark 3.2.2

In Remark 1.3.4, we explain how to avoid the edge resonances.

Suppose that the problems (P ± ) are well posed. Then, the DtN operators

where u ± (β, α; ϕ) is the unique H 1 solution to (P ± ). The next theorem is therefore straightforward.

Theorem 3.2.2 (Problem with DtN conditions [START_REF] Fliss | A dirichlet-to-neumann approach for the exact computation of guided modes in photonic crystal waveguides[END_REF])

The problem (E) is equivalent to the problem set on B 0 = B ∩ Ω 0 Find ω 2 / ∈ σ ess (β), s.t. ∃u 0 ∈ H 1 (B 0 ), u 0 = 0, -1 ρ 0 u 0 = ω 2 u 0 , in B 0 (E 0 ) u 0 satisfying the boundary conditions

(BC 0 )

These problems are equivalent in the sense that if (ω, u) is solution of (E) then (ω, u| B 0 ) is solution of (E 0 ). Conversely, if (u 0 , ω) is solution of (E 0 ) then u defined by

associated to the same value ω is solution of (E). Moreover, the multiplicity of ω is the same for the two problems.

Whereas the problem (E) was linear with respect to the eigenvalue ω 2 but defined on an unbounded domain, the problem (E 0 ) is set on a bounded domain but non linear. Note that the problem (E 0 ) is also non linear with respect to β (whereas the problem (E) can be rewritten as a quadratic eigenvalue problem). In other words, this difficulty would be present if we have decided to fix ω and look for β.

We now introduce the solution algorithm of the non linear eigenvalue problem and explain how to compute the DtN operators in the case where ω 2 / ∈ σ ess (β).

3.2 A Dirichlet-to-Neumann approach for the exact computation of guided modes in photonic crystal waveguides 81

Solution algorithm

For α 2 / ∈ σ ess (β), we denote by A 0 (β, α) the operator

where we have replaced in (BC 0 ) ω by α. By showing a Gärding inequality for Λ(β, α), we can show that Theorem 3. 2.3 ( [7])

The operator A 0 (β, ω) is selfadjoint and with compact resolvant.

Its spectrum is then a pure point one and consists of a sequence of eigenvalues (µ n (β, ω)) n of finite multiplicity tending to +∞. The explicit expression of these eigenvalues using the Min-Max principle and the norm-continuity of the DtN operator Λ ± (β, α) with respect to α yield some regularity properties of each eigenvalue with respect of ω.

Consequently, the solutions of the non-linear problem (E 0 ) are the roots of the equations :

We then infer the iterative algorithm for the computation of the guided modes and associated eigenvalues with two nested loops:

• the outer loop consists in a fixed point algorithm to solve the non linear equation:

• each iteration of this fixed point algorithm requires the computation of the m-th eigenvalue µ m (β, α) of the operator A 0 (β, α) (and possibly the derivative of µ m (β, α) with respect to α if a Newton method is used to solve the fixed point problem).

This algorithm is quite classical for the computation of guided modes in open waveguides (see [Bonnet-Ben Dhia and Starling, 1994]).

Here the novelty comes from the fact that the eigenvalues ω 2 could belong to any gap of the spectrum and moreover that the operators Λ ± (β, α) have no analytical expression, but they only can be computed numerically.

For the construction of the DtN operators Λ ± (β, α) for any β and α, we refer to Section 1.3.4 which has to be adapted to take into account the β-quasi periodic boundary conditions. Let me remind the algorithm, written for α 2 / ∈ σ ess (β). Actually, for α 2 / ∈ σ ess (β), Problem (P ± ) being well-posed in H 1 , the algorithm is the same as the one for the case with dissipation:

• solve the two cell problems associated to (P ± );

• compute the local DtN operators;

• compute the unique solution of spectral radius strictly less than 1 of the associated Riccati equation;

Chapter 3. Trapped and guided modes in periodic structures

) turns out to be an eigenvalue problem for A µ ε (β) in L 2 (Ω µ ε ). For µ = 1, there is no eigenvalue. We investigate the possibility of creating eigenvalues by playing with the parameter µ.

Limit problem

The investigation of the spectral problem (3.2)-(3.3)-(3.4) is based on its asymptotic analysis as ε tends to zero. First, we identify the limit spectral problem. This problem is posed on the graph G obtained by taking the geometrical limit of the domain Ω µ ε when its thickness tends to zero (cf. Figure 3.8). More precisely, we look for the eigenpairs (u, λ), where, for any edge e of G, the restriction u e of u to e is solution of -u e = λ u e , u is β-quasi periodic in the y-direction, and, at each "interior" vertex M of the graph, u is continuous and satisfies the so-called Kirchhoff transmission conditions

u e (M ) being defined outward. In (3.5), E(M ) denotes the set of the edges sharing M as a common vertex, w µ (e) = 1 for any unperturbed edge and w µ (e) = µ for the two perturbed edges.

In order to describe the limit operator defined on the graph G we need to introduce the following function spaces.

. Its spectrum can be characterized explicitly through the Floquet-Bloch theory and using a finite difference scheme and we have in particular the following result: (β) < a n 0 (β). Moreover, if µ ≥ 1, the discrete spectrum of A µ (β) is empty, while if µ < 1 it contains one or two eigenvalue(s) in each gap.

Chapter 4. Homogenization of transmission problems

ent ansatz: two far field zones (regions far from the interface Ω ± ε = {(x 1 , x 2 ), ±x 1 ≥ η(ε)) and one near field zone (region in a neighborhood of the interface

. The regions overlap and the asymptotic expansions have to coincide in the overlapping zone.

Asymptotic expansion : ansatz and equations

Far field terms: we start by describing a suitable asymptotic expansion of u in the left halfplane which is the homogeneous part of the propagation medium. This is the easiest part of our analysis. We start with an ansatz for the expansion taking the form

As usual in asymptotic analysis, the series above is not a priori convergent and should be understood in the sense of asymptotic series, see e.g the introduction of [Ilin, 1992]. Plugging this ansatz in the equations of (4.1) yields a PDE to be satisfied by each term of the expansion,

Next we focus on the derivation of an expansion in the other half plane, namely Ω + ε . We can use the two-scale asymptotic expansion to describe u in the periodic part of the medium:

where for all n ≥ 0, u n (x, •) ∈ H 1 (Y) and as explained in Section 4.2.1

with the Θ k 's defined in the previous section and the u n 's satisfy

all the coefficients A * , ρ * and Q α k being defined also in the previous section.

Let us sum up : for all n, if the preceding terms ûn-1 , ûn-2 , . . . , û0 are assumed already known and sufficiently regular, the macroscopic far field term ûn satisfies (4.8) in Ω -and (4.10) in Ω + . The oscillating part of the u n 's in Ω + can be computed then a posteriori. To make the problem well-posed, transmission conditions involving the jump and the jump in the conormal derivative associated to the homogenised operator ∇ • (A * 0 ∇•) defined by

Chapter 4. Homogenization of transmission problems Thus we shall use more theoretical tools.

As already mentioned, the variational framework is too restrictive because it does not allow unbounded solution at y 1 = ±∞. We shall use Kondratiev theory [START_REF] Kozlov | Elliptic boundary value problems in domains with point singularities[END_REF]] that requires weighted Sobolev spaces. So let us introduce for β > 0 the space V k β (B) defined as

This space contains functions that are allowed to grow when y 1 tends to ±∞ but not faster than e 2πβ|y 1 | . We also need to consider the dual space

equipped with the dual norm. Intuitively, the space V 1 β (B) * contains objects that have the same behaviour at infinity that the functions of V 1 -β (B) but less regular. An important example of elements of V 1 β (B) * are L 2 (B) functions with compact support. We will look for the near field terms in the space V 1 β (B). The choice of β has to be chosen carefully and we show (without giving any detail here) that there exists a β for which all the following results hold.

Recall that we want to solve problems of form: "find U ∈ V 1 β (B) such that -∇ • (a∇U ) = g". First observe that there is no uniqueness because constants functions (which belongs to V 1 β (B)) are non-trivial solutions of the homogeneous problem. We can build another function which is linearly independent of the constants, noting that for y 1 < 0, the equation -∇ • (a∇u) is simply -a 0 ∆u with a non-trivial solution y 1 (which belongs to in V 1 β (B)); and for y 1 > 0, that a is periodic and a non-trivial solution (in V 1 β (B)) of -∇ • (a∇u) is given by y 1 + w 1 (see Section 4.2.1). By coupling Floquet-Bloch Transform and the Kondratiev theory [START_REF] Kozlov | Elliptic boundary value problems in domains with point singularities[END_REF], we have shown the two following results.

Moreover, there exists N ∞ ∈ R such that N satisfies the conditions at infinity

The function N is called a profile function and depends not only on the periodic medium but also on the position of the interface.

4.3

The matched asymptotic expansion applied to the transmission problem

β (B) * . Then there exist α ± , β ± ∈ R such that U has the following behaviour at infinity

Since the kernel of L 0 is of dimension 2, the constants α ± , β ± ∈ R and the r.h.s g have to satisfy necessarily two compatibility conditions (obtained using generalized Green's formulas with the two functions 1 and N ). These necessary conditions are also sufficient to ensure existence and uniqueness of strip problems, as stated by the following theorem.

Theorem 4.3.3 (Existence and uniqueness of the strip problems)

if and only if

We have now all the ingredients at hand to formulate a precise definition for the terms U n , n = 0 of the ansatz (4.16). We suppose that the far fields terms ûn exist and necessary conditions on this far field terms are derived that ensure the existence and uniqueness of the near field terms. We use the definition (4.13) and (4.14) of the matching functions m ± n . Let us remind that x 2 plays a role of a parameter in the near field problems.

• We look for a near field term U 0 such that

Applying Theorem 4.3.3 with g = 0, α -= û0 Γ -, α + = û0 Γ + , β ± = 0 we have that there exists a unique U 0 ∈ V 1 β (B) solution of this near field problem if and only if

In that case, U 0 (x 2 , y) = û0 Γ (x 2 ) for all x 2 .

• We look for a near field term U 1 such that