# Equations intégrales pour l'équation des ondes

Eliane Bécache

► To cite this version:

Eliane Bécache. Equations intégrales pour l'équation des ondes. École d'ingénieur. France. 2001. cel-01529244

HAL Id: cel-01529244

https://cel.hal.science/cel-01529244v1

Submitted on 30 May 2017

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
# Equations intégrales pour l'équation des ondes

E. Bécache*

<footnote>*, INRIA, Domaine de Voluceau-Rocquencourt, BP 105, F-78153 Le Chesnay Cédex</footnote>
Table des matières

<table><tr><td>1 Introduction</td><td>2</td></tr><tr><td>2 Le Problème de diffraction. Notations et équations.</td><td>4</td></tr><tr><td>3 Quelques outils.</td><td>5</td></tr><tr><td>3.1 Solution fondamentale.</td><td>6</td></tr><tr><td>3.2 Formule de Kirchhoff pour l'onde incidente</td><td>8</td></tr><tr><td>3.3 Premier théorème de représentation</td><td>9</td></tr><tr><td>3.4 Définition des potentiels retardés et deuxième théorème de représentation</td><td>11</td></tr><tr><td>3.5 Traces des potentiels retardés</td><td>14</td></tr><tr><td>3.5.1 Traces d'un potentiel de simple couche</td><td>14</td></tr><tr><td>3.5.2 Traces d'un potentiel de double couche</td><td>15</td></tr><tr><td>4 Quelques exemples d'équations intégrales</td><td>17</td></tr><tr><td>4.1 Problème de Dirichlet.</td><td>18</td></tr><tr><td>4.2 Problème de Neumann</td><td>20</td></tr><tr><td>5 Analyse mathématique d'un problème modèle: le problème de Neumann double couche</td><td>23</td></tr><tr><td>5.1 Transformation de Fourier-Laplace du problème</td><td>24</td></tr><tr><td>5.2 Analyse du problème transformé.</td><td>29</td></tr><tr><td>5.3 Retour en temps</td><td>38</td></tr><tr><td>6 Mise en oeuvre et résultats numériques</td><td>42</td></tr><tr><td>6.1 Discrétisation</td><td>42</td></tr><tr><td>6.2 Méthodes numériques. Mise en oeuvre</td><td>45</td></tr><tr><td>6.3 Autres problèmes traités et quelques unes de leurs difficultés respectives...</td><td>50</td></tr><tr><td>7 Bibliographie</td><td>53</td></tr></table>
# 1 Introduction

On s'intéresse à des problèmes de diffraction des ondes par un obstacle dans un milieu homogène. Ces problèmes se situent souvent dans un domaine extérieur (infini). Il peut s'agir d'ondes acoustiques, élastiques ou électromagnétiques. Typiquement, on va considérer un domaine $\Omega_e = \tilde{\Omega}_i^c$, où $\Omega_i$ est un ouvert borné de $R^N$ en dimension $N = 2$ ou 3, dans lequel une onde se propage. La question est de savoir comment se propage cette onde après s'être diffractée sur l'obstacle. La solution vérifie une Equation aux Dérivées Partielles dans $\Omega_e$

$$Lu = f \quad \text{dans } \Omega_e \times R$$

avec des conditions initiales et des conditions aux limites sur le bord du domaine.

Il existe plusieurs méthodes permettant de résoudre ce problème numériquement qui sont essentiellement :

## • Les méthodes de Différences Finies (DF)/Eléments Finis (EF)

Elles conduisent à mailler tout le domaine de calcul. Le domaine $\Omega_e$ étant infini, elles doivent introduire des frontières artificielles sur lesquelles on pose des Conditions aux Limites Absorbantes (CLA) (voir cours de F. Collino). Ces conditions ne sont pas exactes et produisent donc des réflexions. Pour éviter de trop grandes perturbations, on doit donc placer les frontières artificielles assez "loin" de l'obstacle (loin par rapport à la taille de l'obstacle et de la longueur d'onde) et utiliser des CLA d'ordre suffisamment élevées, ce qui conduit à des systèmes à résoudre assez coûteux (de la taille du nombre de noeuds dans le domaine). De plus, le maillage introduit une anisotropie numérique (dispersion), qui est atténuée en se fixant un nombre minimum de noeuds par longueur d'onde ce qui limite également la fréquence de l'onde incidente si on ne veut pas avoir des systèmes trop grands à résoudre. Cependant, l'avantage des méthodes d'éléments finis est que les matrices intervenant, bien que grandes, sont en général très creuses. S'agissant de problèmes d'évolution, on peut utiliser, pour éviter d'inverser ces systèmes, des schémas explicites en temps. La stabilité de ces schémas impose alors une condition de type CFL sur les paramètres de discrétisation du type $c\Delta t/\Delta x < Cte$, ce qui impose un paramètre de discrétisation en temps $\Delta t$ très petit lorsque la vitesse des ondes est élevée. Pour toutes ces raisons, il est souhaitable (et même crucial pour certains problèmes 3D) de pouvoir utiliser des méthodes numériques moins coûteuses. Cependant, ce sont souvent les seules dont nous disposons pour des problèmes posés en milieu non homogène (sauf cas particuliers, par exemple les approximations paraxiales lorsque la propagation se fait suivant une
direction privilégiée ...).

• Les méthodes d'Equations Intégrales (EI)

Lorsque le milieu est homogène, il existe des méthodes plus performantes qui sont basées sur la résolution d'équations intégrales. Ces équations sont posées sur la frontière de l'obstacle qui est dans la plupart des cas un domaine borné. De plus, on gagne une dimension d'espace. Elles conduisent à la résolution de systèmes beaucoup plus petits, mais pleins.

Ces méthodes sont maintenant très bien étudiées (voir la partie bibliographie 7) et couramment utilisées pour les problèmes stationnaires et les problèmes harmoniques (voir partie 7), beaucoup moins pour les problèmes transitoires pour lesquels apparaissent des difficultés supplémentaires, comme nous allons le voir.

Elles reposent sur la connaissance de la solution fondamentale (ou fonction de Green) ce qui explique leur limitation aux milieux homogènes. Cependant, si une partie seulement du milieu est non homogène, on peut alors coupler les méthodes EI(en milieu homogène)-EF(en milieu non homogène), ce qui permet de ne mailler que la partie non homogène (voir les travaux de Johnson et Nédélec, [35] et l'article de revue sur le couplage de Hsiao, [31]). Elles ont également été utilisées pour résoudre des problèmes posés dans des milieux non homogènes particuliers (milieux stratifiés, voir les travaux de Aubry et Clouteau, [2], [13]).

En ce qui concerne les ondes transitoires, avant les travaux de Ha Duong, seules des applications des équations intégrales espace-temps (potentiels retardés) existaient dans la littérature mais étant donnée leur nature (opérateurs intégro-différentiels) aucune analyse mathématique n'en avait été faite. On trouvera les références de quelques unes de ces applications dans la partie bibliographie 7. Dans [22] (et par la suite dans [4], [5], [23]) Ha Duong a montré comment analyser les potentiels retardés par l'intermédiaire de leur transformée de Fourier-Laplace en temps. Grâce à cette étude, des résultats d'existence, d'unicité et de stabilité ont pu être obtenus dans un cadre fonctionnel de type espaces de Sobolev.

Par la suite, ces résultats ont été appliqués et étendus et sont encore actuellement l'objet de recherches. Certains problèmes concernant l'analyse numérique de ces méthodes restent ouverts. Citons enfin les travaux récents de Lubich, [40], qui propose une autre méthode numérique (discrétisation par un schéma de Galerkin en espace et par une méthode multipas en temps) et en fait l'analyse numérique (étude d'erreurs, ...).

L'objet de ce cours est de présenter l'analyse mathématique des potentiels retardés et leur approximation par des schémas de type Galerkin (en espace et en temps) pour les ondes transitoires.

Il y aura une présentation de l'analyse théorique des potentiels retardés. Pour
cela, on aura préalablement rappelé quelques résultats de base sur l'équation des ondes, et introduit un certain nombre d'outils. Puis, nous présenterons les aspects pratiques liés à cette méthode: mise en oeuvre numérique, les difficultés liées à la méthode (singularités des noyaux), au caractère temporel (géométrie, causalité des noyaux), difficultés plus spécifiques (hypersingularité, différences entre acoustique, élastique, Maxwell...)

Nous présenterons enfin quelques résultats numériques.

## 2 Le Problème de diffraction. Notations et équations.

On considère un ouvert borné $\Omega_i$ de $R^N$, $N = 2$ ou $3$, de frontière $\Gamma = \partial\Omega_i$ et le domaine extérieur $\Omega_e = \bar{\Omega}_e^c$. On suppose dans toute la suite (sauf mention du contraire) que $\Gamma$ est une surface fermée au moins $C^1$ par morceaux ou bien un écran infini et que $\Omega_e$ est situé localement d'un seul côté de $\Gamma$. On notera $n$ la normale extérieure à $\Omega_i$ (orientée de $\Omega_i$ vers $\Omega_e$).

**Remarque 2.1** : on peut traiter de manière analogue le cas de la diffraction d'une onde par une ou plusieurs fissures. Dans ce cas $\Gamma$ n'est pas une courbe fermée et $\Omega_e$ n'est pas situé localement d'un seul côté de $\Gamma$. Nous n'exposons pas ce cas ici, par soucis de lisibilité, car les espaces fonctionnels qui interviennent sont un peu différents et font intervenir le comportement des fonctions au voisinage des extrêmités de ces fissures. Nous renvoyons aux travaux [23], [7], [8] sur ce sujet.

On veut calculer pour $t > 0$ la solution du problème mixte d'ondes acoustiques suivant :

$$ (1) \quad \left\{ \begin{array}{ll} \square u \equiv \frac{\partial^2 u}{\partial t^2} - \Delta u = f & \text{dans } \Omega_e \times R_+ \quad (a) \\ u(x, 0) & = u_0(x) \quad \text{dans } \Omega_e \quad (b) \\ \frac{\partial u}{\partial t}(x, 0) & = u_1(x) \quad \text{dans } \Omega_e \quad (c) \\ Bu & = 0 \quad \text{sur } \Gamma \times R_+ \quad (d) \end{array} \right. $$

où $f, u_0$ et $u_1$ ont leurs supports respectivement dans $\Omega_e \times R_+$ et $\Omega_e$, et où la vitesse des ondes acoustiques a été choisie égale à $c = 1$. La condition aux limites $B$ sera par la suite soit une condition de Dirichlet soit une condition de Neumann.

De manière classique, on décompose l'onde comme une superposition d'un champ incident $u^I$ représentant l'onde qui se propage dans tout l'espace libre $R^N$ et d'un
champ diffracté $u^D$. Par définition, le champ incident vérifie le problème :

$$ (2) \qquad \left\{ \begin{array}{lcl} \Box u^I & = & f \quad \text{dans } R^N \times R_+ \quad (a) \\ u^I(x, 0) & = & u_0(x) \quad \text{dans } R^N \quad (b) \\ \frac{\partial u^I}{\partial t}(x, 0) & = & u_1(x) \quad \text{dans } R^N \quad (c) \end{array} \right. $$

et le champ diffracté s'obtient comme différence entre l'onde totale et l'onde incidente :

$$ (3) \qquad u^D = u - u^I $$

et vérifie donc le problème diffracté suivant :

$$ (4) \qquad \left\{ \begin{array}{lcl} \Box u^D & = & 0 \quad \text{dans } \Omega_e \times R_+ \quad (a) \\ u^D(x, 0) & = & 0 \quad \text{dans } \Omega_e \quad (b) \\ \frac{\partial u^D}{\partial t}(x, 0) & = & 0 \quad \text{dans } \Omega_e \quad (c) \\ Bu^D & = & -Bu^I \quad \text{sur } \Gamma \times R_+ \quad (d) \end{array} \right. $$

Il s'agit du problème de diffraction (ou "scattering") que nous cherchons à ré-soudre.

L'onde incidente est connue explicitement grâce à la formule de Kirchhoff et l'onde diffractée admet une représentation intégrale en fonction de ses traces sur $\Gamma$. C'est en écrivant la condition aux limites satisfaite par l'onde diffractée en fonction de cette représentation intégrale qu'on obtient une équation intégrale sur $\Gamma$.

Avant de déterminer les équations intégrales associées à quelques problèmes de diffraction, nous devons nous munir d'un certain nombre d'outils...

### 3 Quelques outils.

Le point essentiel des méthodes d'équations intégrales est le théorème de représentation intégrale de la solution (formule de Kichhoff, pour les ondes acoustiques). Il permet d'obtenir une expression analytique de l'onde incidente et une expression de l'onde diffractée en fonction d'inconnues définies sur la surface de l'obstacle. Ce théorème de représentation nécessite l'utilisation de la formule de Green (ou formule analogue d'intégration par parties, pour un autre opérateur que celui des ondes acoustiques) ainsi que la connaissance de la solution fondamentale.
## 3.1 Solution fondamentale.

La solution fondamentale de l'équation des ondes $G(x, t)$ est par définition le champ de déplacement au point $x$ et à l'instant $t$ dû à une impulsion ponctuelle appliquée à l'origine et à l'instant $0$ et à support $t > 0$:

$$ (5) \qquad \frac{\partial^2 G}{\partial t^2} - \Delta G = \delta(x)\delta(t) $$

### • Solution fondamentale en dimension 2

$$ (6) \qquad G(x, t) = \frac{1}{2\pi} \frac{H(t - |x|)}{(t^2 - |x|^2)^{1/2}} $$

où $H$ représente la fonction de Heaviside.

### • Solution fondamentale en dimension 3

$$ (7) \qquad G(x, t) = \frac{1}{4\pi} \frac{\delta(t - |x|)}{|x|} $$

### Remarque 3.1 : différence fondamentale entre le 2D et le 3D

On voit apparaître sur les expressions des solutions fondamentales 2D et 3D une différence essentielle entre ces deux dimensions, qui est le support de ces fonctions.

En dimension 2, la solution fondamentale remplit le cône de lumière, c'est à dire le cône d'équation $t = |x|$, alors qu'en dimension 3 elle est localisée au bord du cône. Ce phénomène se comprend bien physiquement. En 2D, on peut penser à une onde se propageant à la suite d'une perturbation dans l'eau (par exemple, si on lance une pierre dans l'eau); un front d'onde se propage, mais des vaguelettes continuent à exister à l'intérieur du cercle (elles s'atténuent avec le temps, en $1/t$ mais ne disparaissent pas en théorie). En 3D, on peut penser à un son qui se propage: une fois qu'il est passé, heureusement, il ne subsiste pas éternellement!...

Cette différence va avoir des conséquences importantes sur la résolution numérique: la représentation intégrale de la solution fait intervenir, comme nous allons le voir, un produit de convolution de la solution fondamentale avec d'autres fonctions. La présence de la masse de Dirac, en 3D, va permettre d'éliminer l'intégration en temps et de ne prendre en compte qu'un seul instant, correspondant au passage de l'onde. En 2D, par contre, la fonction d'Heaviside impose d'intégrer en temps et de prendre en compte non seulement l'instant correspondant au passage du front d'onde, mais aussi tout son passé.
On peut déjà noter qu’en élasticité, la situation se complique du fait de la pré-
sence des ondes P et S (pression et cisaillement). En 3D, il y a deux fronts d’onde et
il faut tenir compte de tous les instants compris entre les instants de passage des deux
ondes (donc une intégration en temps supplémentaire par rapport à l’acoustique). En
2D, il faut considérer le passé des deux fronts d’ondes.

Calcul de la solution fondamentale

Nous indiquons brièvement la méthode permettant de calculer la solution fondamentale (formellement). De manière générale, si on considère une EDP

$Lu = f$

(où $L$ est un polynôme en les dérivées partielles), la solution fondamentale vérifie

$$
LG = \delta \quad \text{dans } \mathbb{R}^N \quad (i)
$$

On calcule $G$ en prenant la transformée de Fourier de (i) dans $\mathbb{R}^N$, $L$ étant un polynôme par rapport aux dérivées partielles on a

$$
\widehat{LG} = \hat{L}\hat{G} = 1 \quad (i)
$$

Il suffit alors de savoir inverser $\hat{L}$ pour en déduire $\hat{G}$. L'étape la plus délicate est alors
de revenir aux variables initiales.

En ce qui concerne l'équation des ondes (scalaire), il est très facile d'obtenir (formellement)

$$
(8) \qquad \hat{G}(\xi, \omega) = \frac{1}{|\xi|^2 - \omega^2}
$$

dont on connait la transformée de Fourier inverse donnée en (6) et (7).

**Remarque:** en toute rigueur, la transformée de Fourier en temps diverge et pour justifier ces calculs, on doit introduire une partie imaginaire $\varepsilon$ dans $\omega$ et faire tendre $\varepsilon$ vers 0 (i.e., $\hat{G}(\xi, \omega) = \lim_{\varepsilon\to 0} \hat{G}(\xi, \omega - i\varepsilon) = \lim_{\varepsilon\to 0} \frac{1}{|\xi|^2 - (\omega - i\varepsilon)^2}$. Le calcul rigoureux de la solution fondamentale (en particulier de la transformée de Fourier inverse de (8) est fait dans [55].

**Remarque :** à partir de l’expression (8), on obtient également la solution fondamentale de l’équation de Helmholtz en prenant la transformée de Fourier inverse seulement par rapport aux variables d’espace :

$$
(9) \qquad \mathcal{F}_{\xi}^{-1}(\hat{G}(\xi, \omega)) = \begin{cases} \frac{i}{4} H_0^{(1)}(\omega |x|) & \text{en } 2D \\ \frac{1}{4\pi |x|} e^{i\omega|x|} & \text{en } 3D \end{cases}
$$

**Remarque 3.2** : la solution fondamentale de l'équation des ondes avec une vitesse $c \neq 1$ s'obtient en remplaçant $\omega$ par $\frac{\omega}{c}$.
## 3.2 Formule de Kirchhoff pour l'onde incidente

Si les données $f, u_0$ et $u_1$ sont assez régulières et à support dans $\Omega_e \times R_+, \Omega_e$ et $\Omega_e$, l'onde incidente peut s'écrire sous forme intégrale :

$$ (10) \qquad \left\{ \begin{array}{l} u^I(x,t) = (G^{(x,t)} * f)(x,t) + (G(\cdot, t) * u_1)(x) + (\dot{G}(\cdot, t) * u_0)(x) \\ x \in \mathbb{R}^N, \ t > 0 \end{array} \right. $$

où $\dot{G} = \frac{\partial G}{\partial t}$, ce qui peut encore s'écrire, en explicitant les produits de convolution :

$$ (11) \qquad \left\{ \begin{array}{l} u^I(x,t) = \displaystyle\int_{R_+} \int_{\mathbb{R}^N} G(x-y, t-\tau) f(y, \tau) dy d\tau \\ \phantom{u^I(x,t) = \displaystyle\int_{R_+} \int_{\mathbb{R}^N} } + \displaystyle\int_{\mathbb{R}^N} G(x-y, t) u_1(y) dy + \displaystyle\int_{\mathbb{R}^N} \dot{G}(x-y, t) u_0(y) dy \\ x \in \mathbb{R}^N, \ t > 0 \end{array} \right. $$

### Application

En remplaçant $G$ par son expression donnée en (6) et (7) on obtient :

* en dimension $N = 3$

$$ (12) \qquad \left\{ \begin{aligned} u^I(x,t) = & \frac{1}{4\pi} \int_{|x-y| \le t} \frac{f(y, t - |x-y|)}{|x-y|} dy \\ & + \frac{1}{4\pi t} \int_{|x-y|=t} u_1(y) d\gamma_y \\ & + \frac{1}{4\pi} \frac{\partial}{\partial t} \left( \frac{1}{t} \int_{|x-y|=t} u_0(y) d\gamma_y \right) \end{aligned} \right. $$

* en dimension $N = 2$

$$ (13) \qquad \left\{ \begin{aligned} u^I(x,t) = & \frac{1}{2\pi} \int_0^t \int_{|x-y| \le t-\tau} \frac{f(y, \tau)}{((t-\tau)^2 - |x-y|^2)^{1/2}} dy d\tau \\ & + \frac{1}{2\pi} \int_{|x-y| \le t} \frac{u_1(y)}{(t^2 - |x-y|^2)^{1/2}} d\gamma_y \\ & + \frac{1}{2\pi} \frac{\partial}{\partial t} \left( \int_{|x-y| \le t} \frac{u_0(y)}{(t^2 - |x-y|^2)^{1/2}} d\gamma_y \right) \end{aligned} \right. $$
**Remarque 3.3** : de nouveau on voit la différence sur ces expressions entre le cas 2D et 3D. En 3D, le temps caractéristique qui apparaît est $\tau = t- |x-y|$ indiquant le retard (ce qui se passe au point $x$ et à l'instant $t$ dépend de ce qui se passait au point $y$ à l'instant $\tau = t- |x-y|$). En 2D, il apparaît des intégrales supplémentaires $(\int_0^t d\tau$ et on intègre sur tout le disque $|x-y|\le t$ au lieu de la sphère $|x-y|=t$ en 3D). La solution en $(x,t)$ dépend de tout ce qui se passait en $y$ aux instants précédents $t- |x-y|$. Dans les deux cas, ce terme de retard explique le nom des *potentiels retardés*. Une difficulté supplémentaire apparaît en 2D liée à la forme de la fonction de Green : les variables espace et temps ne sont pas “découplées” comme en 3D, ce qui donne des noyaux plus compliqués à intégrer. D'autre part, le noyau possède une singularité (intégrable) sur le bord du cône. Le cas 2D peut en fait être obtenu simplement à partir du 3D en intégrant par rapport à la troisième variable (on se ramène à un problème 2D lorsqu'il y a une invariance par rapport à $z$) et cette intégration fait perdre le découplage en $x$ et $t$.

## 3.3 Premier théorème de représentation

La solution du problème mixte admet elle aussi une représentation intégrale qui fait intervenir ses traces sur la frontière de l'obstacle $\Gamma$.

**Notations:** soit $v \in C^1(\bar{\Omega}_e)$. On note

$$ (14) \qquad \left\{ \begin{array}{ll} v_e(x) \equiv v_+(x) = \displaystyle\lim_{x' \in \Omega_e \setminus x} v(x'), & x \in \Gamma \\ \\ \frac{\partial v_e}{\partial n}(x) \equiv \frac{\partial v_+}{\partial n}(x) = \displaystyle\lim_{x' \in \Omega_e \setminus x} n_x \cdot \nabla v(x'), & x \in \Gamma \end{array} \right. $$

les traces extérieures de $v$ sur $\Gamma$, où $n$ est la normale extérieure à $\Omega_i$. De manière analogue, on note avec un indice “i” ou “-” les traces intérieures et $[f] = f_i - f_e$ le saut de $f$ à travers $\Gamma$.

**Théorème 3.1**

Si $v$ est une solution classique de $\square v = 0$ dans $\Omega_e \times R$ à support dans $t \ge t_0$ et possédant des traces régulières sur $\Gamma \times R$, on a :

$$ (15) \qquad \left\{ \begin{array}{l} v(x,t) = \displaystyle\int_{\Gamma \times R} \frac{\partial G}{\partial n_y}(x-y,t-\tau)v_e(y,\tau)d\gamma_y d\tau \\ \\ \qquad \qquad \qquad - \displaystyle\int_{\Gamma \times R} G(x-y,t-\tau)\frac{\partial v_e}{\partial n}(y,\tau)d\gamma_y d\tau, \quad x \in \Omega_e, t \in R \end{array} \right. $$

Le second membre garde un sens pour $x \in \Omega_i$ et vaut 0.
**Remarque 3.4** : la condition de support dans $t \ge t_0$ est essentielle. Si on ramène l'origine des temps à $t = t_0$, cela revient à dire qu'on considère des conditions initiales nulles ainsi que toutes les dérivées en temps initiales.

**Remarque 3.5** : cette formule de représentation est la base de la méthode des équations intégrales. La remarque importante est que le champ $v$ peut être reconstitué dans tout l'espace à partir uniquement de ses traces sur la frontière de l'obstacle. Le but est alors d'obtenir une équation (bien posée !) sur ces traces afin de les déterminer et de pouvoir calculer la solution partout.

**Démonstration** : on définit le prolongement de $v$ à tout l'espace

$$ (i) \quad \tilde{v}(x,t) = \begin{cases} v(x,t) & \text{dans } \Omega_e \times R \\ 0 & \text{dans } \Omega_e^c \times R \end{cases} $$

1ère étape: on montre que la distribution $\tilde{v}$ vérifie :

$$ (ii) \quad -\Box \tilde{v} = \frac{\partial v_e}{\partial n}(x,t)\delta_{\Gamma} + v_e(x,t)\delta'_{\Gamma} $$

où $\delta'_{\Gamma}$ est la dérivée normale de $\delta$. Soit $\Phi \in D(R^N \times R)$, on a :

$$ \langle \Box \tilde{v}, \Phi \rangle = \langle \tilde{v}, \Box \Phi \rangle = \int_{\Omega_e \times R} v(x,t) \left( \frac{\partial^2 \Phi}{\partial t^2} - \Delta \Phi \right) (x,t) dx dt $$

On applique la formule de Green (avec $n$= normale entrante dans $\Omega_e$):

$$ \langle \Box \tilde{v}, \Phi \rangle = \int_{\Omega_e \times R} \left( \frac{\partial^2 v}{\partial t^2} - \Delta v \right) (x,t) \Phi(x,t) dx dt \\ + \int_{\Omega_e} \left[ v \frac{\partial \Phi}{\partial t} - \frac{\partial v}{\partial t} \Phi \right]^{+}_{-\infty} dx \\ + \int_{\Gamma \times R} \left( v_e \frac{\partial \Phi}{\partial n} - \frac{\partial v_e}{\partial n} \Phi \right) d\gamma_x dt $$

On en déduit que

$$ -\Box \tilde{v} = -\Box v1_{\Omega_e} + \frac{\partial v_e}{\partial n}(x,t)\delta_{\Gamma} + v_e(x,t)\delta'_{\Gamma} $$

où la deuxième distribution est définie par

$$ \langle v_e(x,t)\delta'_{\Gamma}, \Phi \rangle = - \int_{\Gamma \times R} v_e(x,t) \frac{\partial \Phi}{\partial n} d\gamma_x dt $$

Par conséquent si $v$ est solution de $\Box v = 0$ dans $\Omega_e$, on obtient l'identité (ii).

2ème étape:

- La distribution $\Box \tilde{v}$ est à support compact en espace. On peut donc la convoler en espace avec la solution fondamentale. L'hypothèse sur le support de $v$ en temps permet de plus de vérifier que $G$ et $v$ ont aussi des supports convolutifs en temps (car $G$ est elle-même à support $t > 0$).
- On calcule $\tilde{v} * G$ en remplaçant $\tilde{v}$ par son expression (ii) et on identifie : $\tilde{v} * G$ à $\tilde{v} * \Box G$ ce qui donne le résultat annoncé.

<figure><img src="image_2.png" /></figure>

## Application

On utilise de nouveau les expressions spécifiques de la solution fondamentale dans chacun des cas 3D et 2D. Le point délicat vient du calcul du terme contenant $\frac{\partial G}{\partial n_y}$. Cette dérivée ainsi que le produit de convolution sont à prendre au sens des distributions. On réarrange les termes de façon à obtenir des noyaux intégrables.

### en dimension 3

$$ (16) \quad \left\{ \begin{aligned} v(x,t) &= \frac{1}{4\pi} \int_{\Gamma} \frac{(x-y) \cdot n_y}{|x-y|} \left( \frac{v_e(y, t-|x-y|)}{|x-y|^2} + \frac{\dot{v}_e(y, t-|x-y|)}{|x-y|} \right) d\gamma_y \\ &\quad - \frac{1}{4\pi} \int_{\Gamma} \frac{1}{|x-y|} \frac{\partial v_e}{\partial n}(y, t-|x-y|) d\gamma_y \\ &= \frac{1}{4\pi} \int_{\Gamma} n_y \cdot \nabla_x \left( \frac{v_e(y, t-|x-y|)}{|x-y|} \right) d\gamma_y \\ &\quad - \frac{1}{4\pi} \int_{\Gamma} \frac{1}{|x-y|} \frac{\partial v_e}{\partial n}(y, t-|x-y|) d\gamma_y \end{aligned} \right. $$

### en dimension 2

$$ (17) \quad \left\{ \begin{aligned} v(x,t) &= \frac{1}{2\pi} \int_{\Gamma} \int_{0}^{t-|x-y|} \frac{1}{((t-\tau)^2 - |x-y|^2)^{1/2}} \times \\ & \quad \left( \frac{(x-y) \cdot n_y}{|x-y|} \left[ \frac{v_e(y, \tau)}{t-\tau+|x-y|} + \dot{v}_e(y, \tau) \right] - \frac{\partial v_e}{\partial n}(y, \tau) \right) d\tau d\gamma_y \end{aligned} \right. $$

Le théorème 3.1 fait apparaître deux termes intégraux (deux potentiels) dans la représentation de $v$. Nous allons maintenant étudier plus précisément chacun de ces potentiels.

## 3.4 Définition des potentiels retardés et deuxième théorème de représentation

Dans toute cette partie on suppose que les fonctions ont suffisamment de régularité, et on définit les potentiels classiquement. Nous verrons par la suite les
espaces fonctionnels de type Sobolev dans lesquels on peut se placer pour l'analyse mathématique.

Pour $p$ et $\phi$ assez régulières, on appelle potentiels retardés de simple et double couche sur $\Gamma$ les fonctions définies sur $(R^N \setminus \Gamma) \times R$ par :

• Potentiel de simple couche

$$ (18) \qquad \left\{ \begin{aligned} Lp(x,t) &= \int_{\Gamma} G(x-y, \cdot) \stackrel{(t)}{*} p(y, \cdot) d\gamma_y \\ &= \int_{\Gamma} \int_{\Gamma} G(x-y, t-\tau) p(y, \tau) d\gamma_y d\tau \end{aligned} \right. $$

• Potentiel de double couche

$$ (19) \qquad \left\{ \begin{aligned} M\phi(x,t) &= \int_{\Gamma} \frac{\partial G}{\partial n_y}(x-y, \cdot) \stackrel{(t)}{*} \phi(y, \cdot) d\gamma_y \\ &= \int_{\Gamma} \int_{\Gamma} \frac{\partial G}{\partial n_y}(x-y, t-\tau) \phi(y, \tau) d\gamma_y d\tau \end{aligned} \right. $$

Le théorème de représentation peut alors se réécrire et être complété :

Théorème 3.2

a) Une solution assez régulière de (1) peut s'écrire

$$ (20) \qquad \left\{ \begin{array}{l} u(x,t) = u^I(x,t) + M u_e^D(x,t) - L \frac{\partial u_e^D}{\partial n}(x,t) \\ t \ge 0, x \in \Omega_e \end{array} \right. $$

où $u^I$ est donnée par (10).

b) Si de plus $u^I$ est assez régulière dans un voisinage de $\Gamma$ pour $t > 0$, on a

$$ (21) \qquad \left\{ \begin{array}{l} u(x,t) = u^I(x,t) + M u_e(x,t) - L \frac{\partial u_e}{\partial n}(x,t) \\ t \ge 0, x \in \Omega_e \end{array} \right. $$

Démonstration : la partie a) du théorème découle directement du théorème 3.1 appliqué à $v = u^D H(t)$ (i.e., on prolonge $u^D$ par 0 pour $t \le 0$).
Etant données les hypothèses sur les supports des données, l'onde incidente vérifie le problème intérieur suivant :

$$
\left\{
\begin{array}{ll}
\square u^I & = 0 \text{ dans } \Omega_i \times R \\
u^I(x, 0) & = 0 \text{ dans } \Omega_i \\
\frac{\partial u^I}{\partial t}(x, 0) & = 0 \text{ dans } \Omega_i
\end{array}
\right.
$$

On applique le théorème 3.1 à $w = H(t)u^I/\Omega_i \times R$ d'où :

$$
w = M u_i^I - L \frac{\partial u_i^I}{\partial n}
$$

et pour $x \in \Omega_e$ on a donc

$$
0 = M u_i^I - L \frac{\partial u_i^I}{\partial n}
$$

On somme cette égalité avec l'identité du a) pour obtenir le b), sous réserve que $u^I$
vérifie $u_e^I = u_i^I$ et $\frac{\partial u_i^I}{\partial n} = \frac{\partial u_e^I}{\partial n}$ ce qui est en particulier vrai lorsque les données $f, u_0$ et $u_1$
sont de classe $C^2$. ■

On peut en fait établir d'autres représentations intégrales de la solution en lui
associant un prolongement particulier dans le domaine intérieur :

**Corollaire 3.1** Si $v$ est une solution classique de $\square v = 0$ dans $(\Omega_e \cup \Omega_i) \times R$ à support dans $t \ge t_0$ et possédant des traces régulières sur $\Gamma \times R$, on a :

$$
(22) \qquad v = L \left[ \frac{\partial v}{\partial n} \right] - M[v] \quad \text{dans } (\Omega_e \cup \Omega_i) \times R
$$

**Démonstration :** ce résultat est une conséquence immédiate du théorème de représentation 3.1 appliqué successivement à $v_{/\Omega_e}$ et à $v_{/\Omega_i}$. ■

**Remarque 3.6** : d'après le corollaire 3.1, on peut choisir a priori de représenter la solution par exemple :

- par un potentiel de simple couche de densité *p*, ce qui implique que le prolongement choisi est tel que *v* est continue à travers Γ et *p* coïncide alors avec le saut de la dérivée normale de *v*.
- par un potentiel de double couche de densité $\varphi$, ce qui implique que le prolongement choisi est tel que $\frac{\partial v}{\partial n}$ est continue à travers $\Gamma$ et $\varphi$ coïncide alors avec le saut de $v$.

Il faut cependant être prudent au moment du choix du prolongement (ou de manière équivalente du choix de la représentation intégrale) car chacun mène à une équation intégrale différente et un mauvais choix pourrait conduire à un problème mal posé.

Afin d'établir les équations intégrales, on doit passer à la limite sur $\Gamma$ dans les expressions intégrales (20), (21)... de la solution, c'est à dire en prendre les traces.

## 3.5 Traces des potentiels retardés

### 3.5.1 Traces d'un potentiel de simple couche

Pour $p$ assez régulière, on pose

$$ (23) \qquad \begin{cases} v(x,t) \equiv Lp(x,t) = \int_{\Gamma \times R} G(x-y, t-\tau) p(y, \tau) d\gamma_y d\tau \\ x \in \mathbb{R}^N \setminus \Gamma, t \in R \end{cases} $$

$v$ est alors solution classique de l'équation des ondes homogène dans $\Omega_e \times R$ et $\Omega_i \times R$ et possède les traces suivantes pour $(x, t) \in \Gamma \times R$:

$$ (24) \qquad \left\{ \begin{aligned} v_+(x,t) &= v_-(x,t) = \int_{\Gamma \times R} G(x-y, t-\tau) p(y, \tau) d\gamma_y d\tau \\ \frac{\partial v_+}{\partial n}(x,t) &= -\frac{1}{2} p(x,t) + \int_{\Gamma \times R} \frac{\partial G}{\partial n_x}(x-y, t-\tau) p(y, \tau) d\gamma_y d\tau \\ \frac{\partial v_-}{\partial n}(x,t) &= \frac{1}{2} p(x,t) + \int_{\Gamma \times R} \frac{\partial G}{\partial n_x}(x-y, t-\tau) p(y, \tau) d\gamma_y d\tau \end{aligned} \right. $$

**Remarque 3.7** : le noyau $G$ est faiblement singulier donc l’intégrale dans (24) est bien définie pour $x \in \Gamma$. Dans le cas des ondes acoustiques, la limite quand $x' \to x \in \Gamma$ de $\nabla v(x', t) \cdot n_x$ ne pose pas de problème non plus car, grâce au terme $n_x \cdot (x - y) = O(|x - y|^2)$ si $x, y \in \Gamma$, le noyau est encore faiblement singulier. C'est une différence importante avec l'élasticité, car dans ce cas, on n'a plus $n_x \cdot (x - y)$ en facteur et le noyau $\frac{\partial G}{\partial n_y}$ devient singulier. L'intégrale dans (24) est cependant encore définie au sens d'une valeur principale de Cauchy. Rappelons brièvement que par définition,
l'intégrale sur $\Gamma$ en valeur principale de Cauchy d'une fonction $f(x, y)$ singulière en $y = x$ est (cf [39]):

$$ (25) \qquad \int_{\Gamma} f(x, y) d\gamma_y = \lim_{\epsilon \to 0} \int_{\Gamma \setminus \{|x-y| \le \epsilon\}} f(x, y) d\gamma_y $$

Ces intégrales sont une difficulté en ce qui concerne la mise en oeuvre numérique et certains auteurs préfèrent utiliser des techniques de régularisation qui transforment ces intégrales en intégrales à noyaux faiblement singuliers (par des intégrations par parties, pour l'élastodynamique voir par exemple [11], [51] ...).

On note $S$ et $K$ les opérateurs définis sur $\Gamma \times R$ par

$$ (26) \qquad \left\{ \begin{aligned} Sp(x,t) &= \int_{\Gamma \times R} G(x-y, t-\tau) p(y, \tau) d\gamma_y d\tau \\ Kp(x,t) &= \int_{\Gamma \times R} \frac{\partial G}{\partial n_x}(x-y, t-\tau) p(y, \tau) d\gamma_y d\tau \end{aligned} \right. $$

On a donc pour $v = Lp$

$$ (27) \qquad \left\{ \begin{aligned} v_+ = v_- &= Sp \\ \frac{\partial v_+}{\partial n} &= \left(-\frac{1}{2}I + K\right)p \\ \frac{\partial v_-}{\partial n} &= \left(\frac{1}{2}I + K\right)p \end{aligned} \right. $$

et $v$ est solution de

$$ (28) \qquad \left\{ \begin{aligned} \Box v &= 0 \quad \text{dans } (\Omega_e \cup \Omega_i) \times R_+ \\ [v] &= 0 \\ \left[\frac{\partial v}{\partial n}\right] &= p \end{aligned} \right. $$

### 3.5.2 Traces d'un potentiel de double couche

Pour $\varphi$ assez régulière, on pose

$$ (29) \qquad \left\{ \begin{array}{l} v(x,t) = -M\varphi(x,t) = -\int_{\Gamma \times R} \frac{\partial G}{\partial n_y}(x-y, t-\tau)\varphi(y, \tau) d\gamma_y d\tau \\ x \in \mathbb{R}^N \setminus \Gamma, t \in \mathbb{R} \end{array} \right. $$
$v$ est alors solution classique de l'équation des ondes homogène dans $\Omega_e \times R$ et
$\Omega_i \times R$ et posséde les traces suivantes pour $(x, t) \in \Gamma \times R$:

$$
(30) \qquad \left\{
\begin{aligned}
v_+(x,t) &= -\frac{1}{2}\varphi(x,t) - \int_{\Gamma \times R} \frac{\partial G}{\partial n_y}(x-y,t-\tau)\varphi(y,\tau)d\gamma_y d\tau \\
v_-(x,t) &= \frac{1}{2}\varphi(x,t) - \int_{\Gamma \times R} \frac{\partial G}{\partial n_y}(x-y,t-\tau)\varphi(y,\tau)d\gamma_y d\tau
\end{aligned}
\right.
$$

et

$$
(31) \qquad \left\{
\begin{aligned}
\frac{\partial v_+}{\partial n}(x,t) &= \frac{\partial v_-}{\partial n}(x,t) = \lim_{\substack{x'=x+\epsilon n_x \\ \epsilon>0, \epsilon \to 0}} \nabla_{x'} v(x',t) \cdot n_x \\
&\equiv D\varphi(x,t)
\end{aligned}
\right.
$$

On note $K'$ l' opérateur défini sur $\Gamma \times R$ par

$$
(32) \qquad K'\varphi(x,t) = \int_{\Gamma \times R} \frac{\partial G}{\partial n_y}(x-y,t-\tau)\varphi(y,\tau)d\gamma_y d\tau
$$

L'opérateur $K'$ est en quelque sorte le dual de l'opérateur $K$ défini précédemment (pas tout à fait vrai en hyperbolique car $t$ et $\tau$ ne jouent pas le même rôle). Encore une fois, l'intégrale définissant $K'$ a un noyau faiblement singulier dans le cas de l'équation des ondes acoustiques, mais est définie au sens d'une valeur principale de Cauchy pour l'élastodynamique.

L'opérateur $D$ défini en (31) comme la dérivée normale du potentiel de double
couche n'est pas connu explicitement en fonction de $\varphi$ car le noyau $\nabla_{x'}(\frac{\partial G}{\partial n_y})$ est
hypersingulier. On peut cependant lui donner un sens, au sens des distributions (voir
[29] pour le 3D harmonique, [22] pour le 3D temporel, [6] pour le 2D temporel). Soit
$\phi \in \mathcal{D}(\Gamma)$, on a :

- en dimension 3

$$
(33) \qquad \begin{aligned}
\langle D\varphi(\cdot, t), \phi \rangle &= \frac{1}{4\pi} \int_{\Gamma \times \Gamma} \frac{n_x \cdot n_y}{|x-y|} \ddot{\varphi}(y, t-|x-y|) \phi(x) d\gamma_x d\gamma_y \\
&\quad + \frac{1}{4\pi} \int_{\Gamma \times \Gamma} \frac{\vec{rot}_{\Gamma} \varphi(y, t-|x-y|) \cdot \vec{rot}_{\Gamma} \phi(x)}{|x-y|} d\gamma_x d\gamma_y
\end{aligned}
$$

où $\vec{rot}_{\Gamma}$ désigne le rotationnel surfacique.

- en dimension 2

$$
(34) \qquad \begin{aligned}
\langle D\varphi(\cdot, t), \phi \rangle &= \frac{1}{2\pi} \int_{\Gamma \times \Gamma} \int_0^{t-|x-y|} \frac{1}{((t-\tau)^2 - |x-y|^2)^{1/2}} \times \\
&\quad \left( n_x \cdot n_y \ddot{\varphi}(y, \tau) \phi(x) + \frac{d\varphi}{ds}(y, \tau) \frac{d\phi}{ds}(x) \right) d\tau d\gamma_x d\gamma_y
\end{aligned}
$$
On a donc pour $v = -M\varphi$

$$ (35) \qquad \left\{ \begin{aligned} v_+ &= \left(-\frac{1}{2}I - K'\right)\varphi \\ v_- &= \left(\frac{1}{2}I - K'\right)\varphi \\ \frac{\partial v_+}{\partial n} &= \frac{\partial v_-}{\partial n} = D\varphi \end{aligned} \right. $$

et $v$ est solution de

$$ (36) \qquad \left\{ \begin{aligned} \Box v &= 0 && \text{dans } (\Omega_e \cup \Omega_i) \times R_+ \\ [v] &= \varphi \\ \left[\frac{\partial v}{\partial n}\right] &= 0 \end{aligned} \right. $$

**Remarque 3.8** : des formules analogues à (33), (34) peuvent être établies dans le cas des ondes élastiques, mais cela demande un peu plus de sueur!... Pour les ondes acoustiques, elles s'obtiennent facilement à partir des formules établies pour Helmholtz, en prenant une transformée de Fourier inverse (en temps). Pour l'élasticité, on procède de la même manière, mais la formulation en fréquences n'est pas unique et il s'agit de choisir la bonne formulation pour le passage en temps, i.e., celle qui conserve à la formulation son caractère causal. Cette formulation repose sur une décomposition d'une fonction de Green en deux termes et cette décomposition doit être faite de telle sorte que les deux noyaux vérifient la causalité (voir [9],[8]).

**Remarque 3.9** : on peut aussi donner un sens à $D\varphi$ “directement” (sans passer par les distributions), en utilisant des méthodes de régularisation. Plusieurs techniques existent et il existe une littérature très fournie dans ce domaine (notamment en élasticité fréquentielle, voir par exemple [11], [38],[49]...).

# 4 Quelques exemples d'équations intégrales

Nous disposons maintenant des outils nécessaires pour pouvoir écrire les équations intégrales à résoudre. Nous allons voir qu'il existe en fait plusieurs équations intégrales associées à un même problème. D'après le théorème 3.2 et le corollaire 3.1, on peut utiliser plusieurs représentations intégrales, que nous notons ici :

$$ (i) \quad u^D = M u_+^D - L \frac{\partial u_+^D}{\partial n} $$
$$
(ii) \quad u = u^I + M u_+ - L \frac{\partial u_+}{\partial n}
$$

$$
(iii) \quad u^D = -M[u^D] + L \left[ \frac{\partial u^D}{\partial n} \right]
$$

4.1 Problème de Dirichlet

La condition aux limites s'écrit alors

$$
(37) \qquad Bu^D = u_+^D = g \text{ sur } \Gamma \times R_+
$$

avec $g = -u^I|_{\Gamma \times R_+}$, ou encore

$$
(38) \qquad Bu = u_+ = 0 \text{ sur } \Gamma \times R_+
$$

• Equations intégrales de première espèce

★ Représentation (i)

La formule de représentation de l'onde diffractée (i) nous donne

$$
(39) \qquad u^D = Mg - L \frac{\partial u_+^D}{\partial n} \text{ dans } \Omega_e \times R_+
$$

On exprime la condition aux limites (37) en prenant la trace de (39), ce qui,
compte tenu des résultats précédents sur les traces des potentiels, donne l'équation
intégrale de première espèce suivante :

$$
(40) \qquad Sp^D = \left(-\frac{1}{2}I + K'\right)g
$$

$$
\text{avec } p^D = \frac{\partial u_+^D}{\partial n}.
$$

★ Représentation (ii)

Si on écrit cette fois la formule de représentation de l'onde totale (ii), on
obtient, compte tenu de la condition aux limites (38), une représentation de la solution
comme somme de l'onde incidente et d'un potentiel de simple couche:

$$
(41) \qquad u = u^I - L \frac{\partial u_+}{\partial n} \quad \text{dans } \Omega_e \times R_+
$$
On exprime la condition aux limites (38) en prenant la trace de (41), et on obtient l'équation intégrale de première espèce d'inconnue $p = \frac{\partial u_+}{\partial n}$ suivante:

$$ (42) \qquad Sp = -g $$

Les équations intégrales (42) et (40) ont la même structure, elles reviennent toutes les deux à inverser l'opérateur $S$.

## ★ Représentation (iii)

Les équations précédentes proviennent d'une représentation de la solution à l'aide de ses traces sur $\Gamma$, c'est à dire en considérant un prolongement par 0 à l'intérieur. Il est également possible de choisir le type de représentation que l'on veut, en utilisant (iii). Si, par exemple, on veut l'onde diffractée a priori sous forme d'un potentiel de simple couche de densité $q$:

$$ u^D = Lq \quad \text{dans } (\Omega_i \cup \Omega_e) \times R_+ $$

ceci revient à prolonger l'onde diffractée à l'intérieur de telle sorte que $u^D$ soit continue à la traversée de $\Gamma$. La densité représente alors le saut de la dérivée normale, $q = \left[ \frac{\partial u^D}{\partial n} \right]$, et est solution d'une équation de première espèce de même type que les précédentes :

$$ (43) \qquad Sq = g $$

## • Equations intégrales de deuxième espèce

Cherchons la solution $u^D$ sous forme d'un potentiel de double couche de densité $\varphi$:

$$ (44) \qquad u^D = -M\varphi \quad \text{dans } (\Omega_i \cup \Omega_e) \times R_+ $$

D'après (iii), ceci revient à prolonger l'onde diffractée à l'intérieur de telle sorte que sa dérivée normale soit continue à la traversée de $\Gamma$, i.e., $\frac{\partial u_+^D}{\partial n} = \frac{\partial u_-^D}{\partial n}$, et la densité représente alors le saut de $u^D$, $\varphi = [u^D]$.

En prenant la trace de (44), on obtient l'équation intégrale de deuxième espèce d'inconnue $\varphi$:

$$ (45) \qquad (\frac{1}{2}I + K')\varphi = -g $$
## 4.2 Problème de Neumann

On procède exactement de la même manière que pour le problème de Dirichlet. La condition aux limites s'écrit ici

$$ (46) \qquad Bu^D = \frac{\partial u_+^D}{\partial n} = g \text{ sur } \Gamma \times R_+ $$

avec $g = -\frac{\partial u_+^I}{\partial n}$, ou encore

$$ (47) \qquad Bu = \frac{\partial u_+}{\partial n} = 0 \text{ sur } \Gamma \times R_+ $$

### • Equations intégrales de première espèce

#### ★ Représentation (ii)

La formule de représentation de l'onde totale (ii) nous donne

$$ (48) \qquad u = u^I + M u_+ \quad \text{dans } \Omega_e \times R_+ $$

On exprime la condition aux limites (47) en prenant la dérivée normale de (48), ce qui donne l'équation intégrale de première espèce, d'inconnue $u_+$, suivante :

$$ (49) \qquad Du_+ = -g $$

**Remarque 4.1** : on peut vérifier que la représentation (48) sur l'onde totale revient à chercher l'onde diffractée sous la forme d'un potentiel de double couche de densité $\varphi$, $u^D = -M\varphi$. On a alors $\varphi = [u^D] = [u] = -u_+$, ce qui signifie que le prolongement intérieur $u_i^D$ tel que $[\frac{\partial u^D}{\partial n}] = 0$ correspond à un prolongement par 0 sur l'onde totale.

#### ★ Représentation (i)

Si on procède de la même manière avec la représentation (i) de l'onde diffractée :

$$ (50) \qquad u^D = M u_+^D - Lg \quad \text{dans } \Omega_e \times R_+ $$

on obtient l' équation de première espèce d'inconnue $u_+^D$ :

$$ (51) \qquad Du_+^D = (\frac{1}{2}I + K)g $$
Les équations intégrales (49) et (51) ont la même structure, elles reviennent toutes les deux à inverser l'opérateur $D$. Le second membre de (51) étant plus compliqué à évaluer, elle n'est pas utilisée en pratique.

• Equations intégrales de deuxième espèce

★ Représentation (ii)

On peut prendre la trace de (48) au lieu de sa dérivée normale et on obtient alors l'équation de deuxième espèce d'inconnue $u_+$

$$ (52) \qquad \left(\frac{1}{2}I - K'\right)u_+ = u_+^I $$

★ Représentation (iii)

Cherchons la solution $u^D$ sous forme d'un potentiel de simple couche de densité $p$:

$$ (53) \qquad u^D = Lp \quad \text{dans } (\Omega_i \cup \Omega_e) \times R_+ $$

D'après (iii), ceci revient à prolonger l'onde diffractée à l'intérieur de telle sorte que $u^D$ soit continue à la traversée de $\Gamma$, i.e., $u_+^D = u_-^D$, et la densité représente alors le saut de la dérivée normale, $p = \left[\frac{\partial u^D}{\partial n}\right]$. La dérivée normale de (53) donne l'équation

$$ (54) \qquad \left(\frac{1}{2}I - K\right)p = -g $$

en quelque sorte duale de (52).

**Remarque 4.2** : on utilise la dénomination abusive “première espèce” et “deuxième espèce” par analogie avec les équations intégrales obtenues pour les problèmes elliptiques. Dans ces cas-là, ce sont des équations de première et deuxième espèce au sens de Fredholm, mais ce n'est plus vrai pour les potentiels retardés. Ce n'est donc pas à l'aide de la théorie de Fredholm qu'on va les étudier. Avant d'aller plus loin, faisons quelques remarques sur les équations intégrales de première espèce et deuxième espèce "habituelles", c'est à dire celles qui entrent dans le cadre d'étude de la théorie de Fredholm.

Traditionnellement, les équations de deuxième espèce sont les plus utilisées numériquement, celles de première espèce ayant mauvaise réputation. Ceci provient du fait qu'on dispose de la théorie de Fredholm pour analyser une équation du type $(I + N)f = g$ lorsque $N$ est un opérateur compact. Examinons la situation en ce qui concerne l'opérateur de Laplace. Les opérateurs $K$ et $K'$ sont alors compacts dans $L^2$ et on peut montrer que les équations de deuxième espèce associées sont bien
posées dans $L^2$. Ces équations sont en pratique approchées par collocation et une analyse numérique permet d'obtenir des résultats de stabilité et de montrer que les matrices à inverser ont un conditionnement en $O(h)$. Cette théorie ne s'applique plus aux équations de première espèce. De plus, ces équations reviennent à inverser une convolution, ce qui est en général très instable car l'intégration est une opération régularisante. Cependant ces raisonnements ne sont plus vrais en ce qui concerne les équations intégrales car les noyaux possèdent une singularité en $x = y$ ce qui gomme le lissage de l'intégration et on comprend intuitivement que cela va permettre de les inverser. D'un point de vue mathématique, les opérateurs $S$ et $D$ ne sont pas bijectifs dans $L^2$. Cependant, on peut montrer par d'autres techniques que les équations intégrales de première espèce sont bien posées lorsqu'on cherche à les résoudre dans le bon cadre fonctionnel (voir par exemple [28], [32], [45]). Pour l'opérateur de Laplace, par exemple, l'opérateur $S$ est bijectif de $H^s$ dans $H^{s+1}$ et l'opérateur $D$ de $H^s$ dans $H^{s-1}$. L'analyse est faite par des techniques variationnelles et conduit donc à des approximations par des schémas de Galerkin. La perte ou le gain d'un cran de régularité de l'opérateur se traduit par un conditionnement des matrices en $O(1/h)$, ce qui est a priori moins bon que pour les équations de seconde espèce. En fait, les schémas numériques ont de très bonnes propriétés: en effet, la formulation variationnelle découle directement de celle du problème volumique associé et préserve donc la coercivité et la symétrie. En pratique, ces méthodes ont été mises en oeuvre dans de nombreuses situations (ondes acoustiques, électromagnétiques, élastiques, ...) et se sont révélées très stables et performantes. Il est à noter que les méthodes de collocation approchant les équations de seconde espèce nécessitent en moyenne dix points par longueur d'onde pour obtenir de bons résultats, alors que les méthodes variationnelles approchant les équations de première espèce n'en nécessitent que cinq.

Il semble donc que la préférence pour l'utilisation systématique des équations de deuxième espèce, d'un point de vue numérique, plutôt que celles de première espèce ne se justifie pas vraiment.

En ce qui concerne les potentiels retardés, les équations intégrales appelées ici de deuxième espèce n'entrent plus dans le cadre de la théorie de Fredholm et en particulier, Ha Duong a montré (cf [22]) que dans des espaces de Sobolev naturels (liés à l'énergie) les opérateurs $K$ et $K'$ ne sont pas compacts. L'analyse qui va être menée dans la suite, montre comment obtenir des résultats sur ces équations dans des espaces fonctionnels appropriés.

**Conclusion** : on a vu 4 types d'équations intégrales, avec des seconds membres plus ou moins compliqués à évaluer :

- l'équation intégrale de première espèce 1 (Dirichlet simple couche) $Sp = g$

- l'équation intégrale de première espèce 2 (Neumann double couche) $D\varphi = g$

- les équations de deuxième espèce $(\frac{1}{2}I \pm K)p = g$ et $(\frac{1}{2}I \pm K')\varphi = g$
Pour chacune de ces équations, il est possible de déterminer un cadre fonctionnel (de type espaces de Sobolev) dans lequel une analyse mathématique peut être menée (existence, unicité, stabilité ...). La plus utilisée numériquement est l'équation de deuxième espèce. Cependant, ce choix ne se justifie pas ici, car comme nous allons le voir, les équations de première espèce conduisent à des problèmes bien posés.

Nous allons étudier dans le paragraphe suivant l'une de ces équations en détail.

## 5 Analyse mathématique d'un problème modèle : le problème de Neumann double couche

Nous avons choisi de traiter ici l'équation intégrale de première espèce as-sociée au problème de Neumann double couche. La démarche pour le problème de Dirichlet simple couche est analogue. Pour plus de détails sur ce sujet, ou pour voir comment traiter un autre cas, on peut se référer à Ha Duong, [22]. Rappelons d'abord le problème.

On cherche la solution $u^D$ du problème diffracté :

$$ (55) \qquad \left\{ \begin{array}{lll} \Box u^D & = & 0 \quad \text{dans } \Omega_e \times R_+ \\ u^D(x, t \le 0) & = & 0 \quad \text{dans } \Omega_e \\ \frac{\partial u^D}{\partial n} & = & g \quad \text{sur } \Gamma \times R_+ \end{array} \right. $$

sous forme d'un potentiel de double couche de densité $\varphi$

$$ (56) \qquad u^D = -M\varphi $$

La densité représente alors le saut du champ à travers $\Gamma$, $\varphi = [u^D]$, et le champ a été prolongé à l'intérieur par la solution du problème suivant :

$$ (57) \qquad \left\{ \begin{array}{lll} \Box u^D & = & 0 \quad \text{dans } \Omega_i \times R_+ \\ u^D(x, t \le 0) & = & 0 \quad \text{dans } \Omega_i \\ \frac{\partial u^D}{\partial n} & = & g \quad \text{sur } \Gamma \times R_+ \end{array} \right. $$

Le champ $u^D$ est donc défini dans $\Omega = \Omega_i \cup \Omega_e$ et sa dérivée normale est continue à la traversée de $\Gamma$. La condition de Neumann conduit alors à résoudre l'équation intégrale de première espèce suivante :

$$ (58) \qquad D\varphi = g \quad \text{sur } \Gamma \times R_+ $$
Rappelons que l'énergie de l'onde a pour expression :

$$ (59) \qquad E(t) = \frac{1}{2} \| \frac{\partial u^D}{\partial t} \|_{0,\Omega}^2 + \frac{1}{2} \| \nabla u^D \|_{0,\Omega}^2 $$

On se place dans $\mathbb{R}^N$ avec $N = 2$ ou $3$. Nous allons présenter l'étude de l'opérateur $D$ et de l'équation (58) en passant par une transformée de Fourier-Laplace en temps, ce qui est possible car on travaille avec des fonctions causales. Le problème (55), (57) est alors transformé en un problème de Helmholtz avec une fréquence $\omega$ complexe et on peut étudier l'équation intégrale associée. La forme bilinéaire de la formulation variationnelle associée à cette équation intégrale vérifie une inégalité de coercivité, ce qui permet d'obtenir l'existence et l'unicité de la solution ainsi que des estimations par rapport aux données. Pour pouvoir transposer ces résultats en temps, il faut connaître les dépendances par rapport à $\omega$, ce qui n'est pas classique pour l'étude de l'équation de Helmholtz. On peut alors faire une transformée de Fourier-Laplace inverse et obtenir des résultats sur l'équation espace-temps (58). L'analyse du problème en fréquence peut être faite dans deux types d'espaces : les espaces de Sobolev munis des normes usuelles et les espaces de Sobolev munis de normes liées à l'énergie. Ils conduisent à deux types d'espaces en temps, et on indiquera les différences entre les résultats obtenus dans chacun des cas. Les espaces liés à l'énergie permettront d'avoir des estimations plus fines et de perdre moins de régularité en temps sur la solution. Cependant, ils ont l'inconvénient de coupler les variables d'espace et de temps, ce qui rend plus délicate leur approximation. Pour l'analyse numérique, on se placera donc dans les espaces découplés (c'est à dire ceux qui sont liés aux normes usuelles).

## 5.1 Transformation de Fourier-Laplace du problème

### Définition

Faisons d'abord quelques rappels. Si $E$ est un espace de Banach, on peut définir la transformée de Fourier sur l'ensemble $S'(E)$ des distributions tempérées de $R$ à valeurs dans $E$. On étend alors cette transformée à un demi-plan complexe pour des distributions causales. On considère une distribution $f$ de $R$ à valeurs dans $E$, causale - i.e. supp($f$) $\subset R_+$ - où supp($f$) représente le support de $f$ en temps et telle que

$$ \exists \omega_I^0(f) > 0 \text{ tel que } e^{-\omega_I^0 t} f(t) \in S'(E). $$

La distribution $e^{-\omega_I^0 t} f(t)$ est alors encore dans $S'(E)$ pour tout $\omega_I > \omega_I^0 > 0$ et on définit la transformée de Fourier-Laplace de $f$ comme étant la transformée de Fourier de $e^{-\omega_I^0 t} f(t)$, c'est à dire par :

$$ \hat{f}(\omega) = \int_{\mathbb{R}} f(t)e^{i\omega t}dt = \mathcal{F}_{(t)}(e^{-\omega_I^0 t}f(t))(\omega_R) $$
avec $\omega \in \mathbb{C}$, $\omega = \omega_R + i\omega_I$. On note $\mathcal{L}'(E)$ l'ensemble formé de telles distributions, c'est à dire des distributions Fourier-Laplace transformables de $R$ à valeurs dans $E$:

$$ \mathcal{L}'(E) = \{ f \in \mathcal{D}'(E), \text{supp}(f) \subset R_+, \exists \omega_I^0(f) \text{ tel que } e^{-\omega_I^0 t} f(t) \in \mathcal{S}'(E) \} $$

## Théorème de Paley-Wiener

Rappelons le théorème de Paley-Wiener donnant des conditions nécessaires et suffisantes pour qu'une distribution soit la transformée de Fourier-Laplace d'une distribution de $\mathcal{L}'(E)$. D'après le théorème de Paley-Wiener, on a l'équivalence entre les assertions (i) et (ii) suivantes :

(i) $h(\omega) = \hat{f}(\omega)$ avec $f \in \mathcal{L}'(E)$

(ii) (a) $h$ est holomorphe dans un demi-plan complexe $\{\Im(m(\omega)) > \omega_I^0\}$ à valeurs dans $E$ et

(b)

$$ \exists \sigma_1 > \omega_I^0, C > 0 \text{ et } k \ge 0 \text{ tels que} $$

$$ \| h(\omega) \|_{E} \le C(1 + |\omega|)^k \quad \forall \omega \text{ tel que } \Im(m(\omega)) \ge \sigma_1 $$

## Plan de travail

A quoi va nous servir le théorème de Paley-Wiener ? au risque de nous répéter (!...), donnons les étapes de l'analyse (étapes qui reviendront d'ailleurs dans les démonstrations qui vont suivre) : rappelons que le but de cette partie est l'étude de l'équation intégrale (58) (ou encore de l'opérateur $D$) associée à un problème volumique de Neumann transitoire. Des méthodes classiques ne permettent pas de l'étudier directement. Moyennant des hypothèses de régularité sur la donnée $g$ et en cherchant la solution dans un espace approprié, on peut grâce à l'implication (i) $\implies$ (ii) définir la transformée de Fourier-Laplace de cette équation pour toute fréquence appartenant à un demi-plan complexe.

L'étude à fréquence $\omega$ fixée permet d'obtenir l'existence et l'unicité d'une solution $\Phi(\omega)$ de l'équation intégrale transformée.

**L'implication inverse (ii) $\implies$ (i)** va nous servir pour revenir au problème transitoire: en effet, si on peut vérifier que $\Phi(\omega)$ satisfait les conditions (ii) (a) et (b) on pourra en déduire l'existence d'une distribution $\varphi = \mathcal{F}^{-1}(\Phi)$ dans un espace approprié. On vérifiera que cette distribution est bien solution de l'équation intégrale en temps (58).

En fait, l'analyse du problème en fréquences va fournir des estimations sur $\Phi$ plus précises que (ii), (b) demandée par Paley-Wiener. On va exploiter ces estimations et être amenés à déterminer des espaces fonctionnels plus fins que $\mathcal{L}(H^{1/2}(\Gamma))$, ce qui donnera le bon cadre fonctionnel pour chercher une solution de l'équation en temps
(58). L'analyse de l'équation intégrale en fréquences va découler de celle du problème volumique associé. On peut récapituler cela par le diagramme formel suivant :

$$
\begin{array}{rcl@{\quad}rcl}
\varphi & \text{solution de} & D\varphi = g, & \xrightarrow{\mathcal{F}} & \Phi(\omega) & \text{solution de} & \hat{D}\Phi = \hat{g}, \\
       & \varphi = [u] & \xleftarrow{\mathcal{F}^{-1}} & & \Phi = [U(\omega)] & \xleftarrow{\mathcal{F}^{-1}} & & \\
       & \downarrow & & & & & \downarrow \\
\multicolumn{1}{c}{u} & \text{solution de} & (55), (57), & \xrightarrow{\mathcal{F}} & U(\omega) & \text{solution de} & (\tilde{55}), (\tilde{57}), \\
       & u = -M\varphi & \xleftarrow{\mathcal{F}^{-1}} & & U = -\hat{M}\Phi & \xleftarrow{\mathcal{F}^{-1}} & &
\end{array}
$$

$$
u \quad \text{solution de} \quad (55), (57), \quad \xrightarrow{\mathcal{F}} \quad U(\omega) \quad \text{solution de} \quad (\tilde{55}), (\tilde{57}), \\
u = -M\varphi \quad \xleftarrow{\mathcal{F}^{-1}} \quad U = -\hat{M}\Phi
$$

Transformation de Fourier-Laplace du problème

On suppose que $g$ est une distribution Fourier-Laplace transformable à valeurs
dans $H^{-1/2}(\Gamma)$ (i.e., $g \in \mathcal{L}'(H^{-1/2}(\Gamma)))$ et on suppose que le problème (55), (57) admet
une solution $u^D$ dans $\mathcal{L}'(H^1(\Omega))$. La transformée de Fourier-Laplace de la solution,
$\hat{u}^D$, vérifie alors le problème transformé suivant :

$$
(60) \qquad \left\{
\begin{array}{ll}
(-\Delta - \omega^2) \hat{u}^D & = 0 \quad \text{dans } \Omega \\
\displaystyle\frac{\partial \hat{u}^D}{\partial n} & = \hat{g} \quad \text{sur } \Gamma \\
\hat{u}^D & \in H^1(\Omega)
\end{array}
\right.
$$

où $\omega \in \mathfrak{S}m(\omega) > \omega_I^0 > 0$. La dernière condition exprime le fait qu'on cherche les solutions d'énergie finie et grâce à la partie imaginaire de $\omega$ cela suffit à définir un problème bien posé, sans avoir besoin de condition de radiation à l'infini.

De même que pour le problème en temps, on peut définir les potentiels de simple
et double couche et montrer qu'ils correspondent aux transformées de Fourier-Laplace
de ceux définis en temps (le produit de convolution est remplacé par un simple produit)

• Potentiel de simple couche

$$
\begin{equation} \tag{61}
\begin{split}
L_{\omega} \hat{p}(x, \omega) &= \int_{\Gamma} \hat{G}(x-y, \omega) \hat{p}(y, \omega) d\gamma_y, \quad x \in \mathbb{R}^N \setminus \bar{\Gamma} \\
&= \hat{L} p(x, \omega)
\end{split}
\end{equation}
$$

• Potentiel de double couche

$$
\begin{equation} \tag{62}
\begin{split}
M_{\omega} \hat{\phi}(x, \omega) &= \int_{\Gamma} \frac{\partial \hat{G}}{\partial n_y}(x - y, \omega) \hat{\phi}(y, \omega) d\gamma_y, \quad x \in \mathbb{R}^N \setminus \bar{\Gamma} \\
&= \hat{M} \hat{\phi}(x, \omega)
\end{split}
\end{equation}
$$
ainsi que leurs traces.

**Remarque 5.1** : l'égalité $L_\omega \hat{p} = \hat{L}p$ ou encore formellement $L_\omega = \hat{L}$ exprime que le potentiel en fréquences est simplement la transformée de Fourier-Laplace du potentiel correspondant en temps. Ceci est une conséquence du fait que la solution fondamentale du problème transformé en fréquences coïncide avec la transformée de Fourier-Laplace de la solution fondamentale du problème en temps (ce qui est clair car $\hat{\delta}(t) = 1$). Cette propriété est donc vérifiée par tous les potentiels et leurs traces.

On peut alors montrer que le problème (60) est bien posé dès que $\hat{g} \in \mathcal{L}'(H^{-1/2}(\Gamma))$. La solution de (60) peut se représenter sous forme d'un potentiel de double couche de densité $\hat{\varphi} = [\hat{u}^D]$ solution de l'équation intégrale

$$ (63) \qquad D_{\omega} \hat{\varphi} = \hat{g} \quad \text{sur } \Gamma $$

où l'opérateur $D_\omega$ est défini par

$$ (64) \qquad \left\{ \begin{aligned} D_{\omega} \hat{\varphi}(x, \omega) &= \lim_{\substack{x' = x + \varepsilon n_x \\ \varepsilon > 0, \varepsilon \to 0}} \nabla_{x'} (-M_{\omega} \hat{\varphi})(x', \omega) \cdot n_x \\ &= D \widehat{\varphi(x, \cdot)}(\omega) \end{aligned} \right. $$

et encore une fois, on remarque que l'opérateur $D_\omega$ coïncide avec $\hat{D}$.

Nous allons montrer un certain nombre de propriétés de $D_\omega$ ce qui permettra d'analyser l'opérateur en temps $D$. Donnons d'abord quelques définitions. L'énergie transformée s'exprime par :

$$ (65) \qquad \hat{E}(\hat{u}^D, \omega) = \frac{1}{2} \| \omega \hat{u}^D \|_{0,\Omega}^2 + \frac{1}{2} \| \nabla \hat{u}^D \|_{0,\Omega}^2 $$

On remarque que pour $\omega \in \{\Im m(\omega) > \omega_I^0\}$, l'énergie définit sur $H^1(\Omega)$ une norme équivalente à la norme usuelle. Nous la noterons :

$$ (66) \qquad \| u \|_{1,\omega,\Omega}^2 = \| \omega u \|_{0,\Omega}^2 + \| \nabla u \|_{0,\Omega}^2 = 2 \hat{E}(u, \omega) $$

On peut de même définir une norme liée à l'énergie sur les espaces $H^s(\Gamma)$. Considérons le cas particulier où $\Gamma$ est soit l'axe $R$ si $N = 2$, soit le plan $R^2$ si $N = 3$ (i.e. $\Gamma = R^{N-1}$). Les normes usuelles sur $H^s(\Gamma)$, pour $s > 0$, peuvent alors être définies à partir de la transformée de Fourier :

$$ \| \varphi \|_{s, \Gamma}^2 = \int_{\mathbb{R}^{N-1}} (1 + |\xi|^2)^s |\tilde{\varphi}|^2 (\xi) d\xi $$
où $\tilde{\varphi}$ désigne la transformée de Fourier par rapport à la variable d'espace $(x_1, \dots, x_{N-1})$.
Les normes liées à l'énergie sont alors définies par :

$$ \| \varphi \|_{s, \omega, \Gamma}^2 = \int_{\mathbb{R}^{N-1}} (|\omega|^2 + |\xi|^2)^s |\tilde{\varphi}|^2(\xi) d\xi $$

Pour $s < 0$, on définit la norme par dualité. Dans le cas général d'un obstacle borné quelconque, on procède de la même manière en définissant des cartes locales.

Les équivalences entre les normes usuelles et les normes liées à l'énergie s'expriment par les inégalités suivantes :

• Pour tout $u \in H^1(\Omega)$

$$ (67) \qquad C_m(\omega_I^0) \| u \|_{1,\Omega} \le \| u \|_{1,\omega,\Omega} \le C_M(\omega) \| u \|_{1,\Omega} \le \tilde{C}_M(\omega_I^0) \| \omega u \|_{1,\Omega} $$

avec

$$ (68) \qquad \left\{ \begin{array}{l} C_m(\omega_I^0) = \min(1, \omega_I^0) = C_m \\ C_M(\omega) = \max(1, |\omega|) = C_M \\ \tilde{C}_M(\omega_I^0) = \left( \frac{1 + \omega_I^{02}}{\omega_I^{02}} \right)^{1/2} = \tilde{C}_M \end{array} \right. $$

• Pour tout $s > 0$ et tout $\varphi \in H^s(\Gamma)$

$$ (69) \qquad C_m^s \| \varphi \|_{s,\Gamma} \le \| \varphi \|_{s,\omega,\Gamma} \le C_M^s \| \varphi \|_{s,\Gamma} \le \tilde{C}_M^s |\omega|^s \| \varphi \|_{s,\Gamma} $$

L'analyse du problème sera présentée dans les espaces de Sobolev munis des normes usuelles. Nous donnerons à titre indicatif les résultats obtenus dans les espaces liés à l'énergie, ces résultats s'obtenant de la même façon.

Avant d'étudier le problème (60) ainsi que l'équation intégrale (63) qui lui est associée, donnons quelques résultats préliminaires.

**Lemme 5.1 (de trace)**

L'application $\gamma = [\dots]$ est une application linéaire continue et surjective de $H^1(\Omega)$ dans $H^{1/2}(\Gamma)$. Il existe donc une constante $C = C(\Gamma)$ telle que

$$ (70) \qquad \| \gamma v \|_{1/2, \Gamma} \le C \| v \|_{1, \Omega} \quad \forall v \in H^1(\Omega) $$

**Démonstration :** ce lemme est une conséquence immédiate du lemme de trace habituel appliqué successivement dans $\Omega_i$ et dans $\Omega_e$. ■

**Lemme 5.2 (de relèvement).**

Pour tout $\psi \in H^{1/2}(\Gamma)$ on peut trouver un “relèvement”, $\Lambda\psi$ dans $H^1(\Omega)$, i.e. une fonction $v = \Lambda\psi \in H^1(\Omega)$ telle que $\gamma v = \psi$, qui vérifie :
$$
(71) \quad \| v \|_{1,\omega,\Omega}^2 \le \frac{C}{C_m} |\omega| \| \psi \|_{1/2,\Gamma}^2 = C \max(1, \frac{1}{\omega_I^0}) |\omega| \| \psi \|_{1/2,\Gamma}^2
$$

où $C = C(\Gamma)$ ne dépend pas de $\omega$.

Démonstration : voir [4].

5.2 Analyse du problème transformé

L'analyse de l'opérateur $D_\omega$ (et de l'équation (63)) découle de celle du
problème volumique :

$$
(72) \qquad \left\{
\begin{array}{ll}
(-\Delta - \omega^2)U = 0 & \text{dans } \Omega \\
\left(\frac{\partial U}{\partial n}\right)_{\pm} = \hat{g}(\omega) & \text{sur } \Gamma \\
U \in H^1(\Omega)
\end{array}
\right.
$$

Chercher une solution du problème (72) , $U \in H^1(\Omega)$, revient à chercher la
solution du problème variationnel suivant :

$$
(73) \qquad \left\{
\begin{array}{ll}
\text{Trouver } U \in H^1(\Omega) \text{ telle que} & \\
a(U,v) = L(v) & \quad \forall v \in H^1(\Omega)
\end{array}
\right.
$$

avec :

$$
(74) \qquad \left\{
\begin{array}{l}
a(u, v) = \int_{\Omega} (\nabla u \cdot \nabla \bar{v} - \omega^2 u \bar{v}) dx \\
L(v) = \int_{\Gamma} \hat{g}[\bar{v}] d\gamma
\end{array}
\right.
$$

Nous regroupons un certain nombre de résultats concernant ce problème dans
le :

Théorème 5.1

a) La forme bilinéaire $a(\cdot, \cdot)$ vérifie l'inégalité de coercivité suivante :

$$
(75) \qquad \left\{
\begin{array}{l}
\Re e(a(v, -i\omega v)) = \omega_I \|v\|_{1,\omega,\Omega}^2 = 2\omega_I \hat{E}(v, \omega) \\
\ge \omega_I^0 C_m^2 \|v\|_{1,\Omega}^2 \quad \forall v \in H^1(\Omega)
\end{array}
\right.
$$
b) Si $\hat{g} \in H^{-1/2}(\Gamma)$, le problème (73) admet une unique solution $U$ dans $H^1(\Omega)$.
Cette solution est continue par rapport à la donnée $\hat{g}$ :

$$
(76) \qquad \omega_I^0 C_m^2 \| U \|_{1,\Omega} \le C \mid \omega \mid \| \hat{g} \|_{-1/2,\Gamma}
$$

c) L'opérateur :

$$
\hat{g} \in H^{-1/2}(\Gamma) \longrightarrow U \in H^1(\Omega) \quad \text{solution de (73)}
$$

est un isomorphisme et la constante de continuité de son inverse est donnée par :

$$
(77) \quad \| \hat{g} \|_{-1/2, \Gamma} \le C \mid \omega \mid^{1/2} \max\left(1, \frac{1}{\sqrt{\omega_I^0}}\right) \| U \|_{1, \omega, \Omega}
$$

**Remarque 5.2** : notons $X_\omega$ l'espace :

$$
X_{\omega} = \{v \in H^{1}(\Omega); (-\Delta - \omega^{2})v = 0 \text{ dans } \Omega; \left[\frac{\partial v}{\partial n}\right]_{\partial \Omega} = 0\}
$$

On peut montrer que $X_\omega$ est un sous-espace fermé de $H^1(\Omega)$. Pour tout $v \in X_\omega$,
on peut définir sa dérivée normale $\partial v / \partial n$ par le crochet de dualité entre $H^{-1/2}(\Gamma)$ et
$H^{1/2}(\Gamma)$ (formule de Green) :

$$
\left\langle \frac{\partial v}{\partial n}, \psi \right\rangle_{-1/2, 1/2, \Gamma} = a(v, \Lambda \psi) \quad \forall \psi \in H^{1/2}(\Gamma)
$$

où $\Lambda$ est l'opérateur de relèvement défini au lemme 5.2.

Il est alors facile de montrer que $U$ est solution du problème (73) si et seulement
si $U \in X_{\omega}$ et $\frac{\partial U}{\partial n} = \hat{g}$.

**Démonstration** : (du théorème 5.1) **a)** On a :

$$
a(v, -i\omega v) = \int_{\Omega} (|\nabla v|^2 i\bar{\omega} - i\omega |\omega v|^2)
$$

et en prenant la partie réelle, on obtient (75).

b) La forme linéaire $L$ vérifie :

$$
(78) \qquad |L(v)| \le \|\hat{g}\|_{-1/2, \Gamma} \|\gamma v\|_{1/2, \Gamma}
$$

Le terme de droite peut être majoré grâce à l'estimation (70) du lemme de trace, ce qui montre la continuité de $L$ sur $H^1$ :

$$
(79) \qquad |L(v)| \le \|\hat{g}\|_{-1/2, \Gamma} \|v\|_{1, \Omega}
$$
En appliquant Cauchy-Schwartz, on obtient aussi la continuité de la forme bilinéaire $a(.,.)$:

$$ (80) \qquad |a(u, v)| \le \|u\|_{1,\omega,\Omega} \|v\|_{1,\omega,\Omega} \le \tilde{C}_M^2 |\omega|^2 \|u\|_{1,\Omega} \|v\|_{1,\Omega} $$

Encore une fois, dans l'estimation (80), on voit que les normes naturelles sont les normes liées à l'énergie.

Etant donnée la coercivité de $a(.,.)$ établie dans le a), on peut appliquer Lax-Milgram et en déduire l'existence et l'unicité de la solution de (73).

L'estimation de stabilité (76) provient de la coercivité de $a(.,.)$ et de la continuité de $L(\cdot)$. En effet, pour $v = -i\omega U$, on a

$$ |a(U, -i\omega U)| = |L(-i\omega U)| $$

On utilise (75) pour minorer le terme de gauche et (78) pour majorer celui de droite

$$ (81) \qquad \omega_I^0 \| U \|_{1,\omega,\Omega}^2 \le |\omega| \|\hat{g}\|_{-1/2,\Gamma} \|\gamma U\|_{1/2,\Gamma} $$

et en utilisant à nouveau le théorème de trace, on obtient

$$ \omega_I^0 \| U \|_{1,\omega,\Omega}^2 \le C |\omega| \|\hat{g}\|_{-1/2,\Gamma} \| U \|_{1,\Omega} $$

Le terme de gauche de l'inégalité (67) permet finalement d'obtenir (76).

c) Pour montrer l'estimation de continuité de l'opérateur inverse, on utilise le lemme de relèvement : pour tout $\psi \in H^{1/2}(\Gamma)$, on définit $v = \Lambda\psi$ tel que

$$ [v] = \psi $$

et

$$ \| v \|_{1,\omega,\Omega} \le \left( \frac{C}{C_m} |\omega| \right)^{1/2} \| \psi \|_{1/2,\Gamma} $$

et pour tout $\psi$ on a

$$ \langle \hat{g}, \psi \rangle_{-1/2,1/2} = \langle \hat{g}, [v] \rangle_{-1/2,1/2} = a(U, v) $$

D'après la continuité (80) de $a(.,.)$, on a donc

$$ | \langle \hat{g}, \psi \rangle_{-1/2,1/2} | \le \| U \|_{1,\omega,\Omega} \| v \|_{1,\omega,\Omega} $$

On utilise la majoration de $\| v \|$ en fonction de la norme de $\psi$ :

$$ | \langle \hat{g}, \psi \rangle_{-1/2,1/2} | \le \| U \|_{1,\omega,\Omega} \left( \frac{C}{C_m} |\omega| \right)^{1/2} \| \psi \|_{1/2,\Gamma} $$
ce qui, par définition de la norme dans $H^{-1/2}$, donne (77). ■

De même que pour le problème en temps, il est facile de montrer, à l'aide des théorèmes de représentation, que la solution $U$ de (73) peut se représenter sous forme d'un potentiel de double couche de densité $\Phi$,

$$ (82) \qquad U = -M_{\omega}\Phi $$

où, par définition

$$ (83) \qquad \Phi = [U] $$

Les estimations (81) et (77) donnent alors une majoration de la solution en fonction de son saut, i.e. du potentiel de double couche en fonction de sa densité:

$$ (84) \qquad \| U \|_{1,\omega,\Omega} \le \frac{C}{\omega_I^0} \max(1, \frac{1}{\omega_I^{0.5}}) |\omega|^{3/2} \| \Phi \|_{1/2,\Gamma} $$

ou encore, en utilisant (75),

$$ (85) \qquad \| U \|_{1,\Omega} \le \frac{C}{C_m^{5/2}} |\omega|^{3/2} \| \Phi \|_{1/2,\Gamma} $$

En définissant l'opérateur $D_\omega$ par (64), on en déduit que $\Phi$ est solution de (63). En particulier, l'existence d'une solution du problème variationnel (73) entraine donc l'existence d'une solution de l'équation intégrale (63). Réciproquement, pour $\Phi$ donnée dans $H^{1/2}(\Gamma)$, si on définit $U(\Phi)$ dans $H^1(\Omega)$ comme le potentiel de double couche de densité $\Phi$, d'après les propriétés du potentiel de double couche, $U(\Phi)$ vérifie :

$$ (86) \qquad \left\{ \begin{array}{l} \Delta U + \omega^2 U = 0 \text{ dans } \Omega \\ \left[ \frac{\partial U}{\partial n} \right] = 0 \\ [U] = \Phi \end{array} \right. $$

c'est à dire que $U(\Phi) \in X_\omega$ et $[U] = \Phi$. Si $\Phi$ vérifie l'équation (63), on en déduit, par définition de $D_\omega$, que $U(\Phi)$ vérifie la condition aux limites

$$ \frac{\partial U}{\partial n}(\Phi) = \hat{g} \text{ sur } \Gamma $$

ce qui signifie que $U(\Phi)$ est solution de (73).

En particulier, si $\Phi$ vérifie l'équation (63) homogène, on en déduit que $U(\Phi)$ est solution de (73) homogène, ce qui implique d'après l'unicité de la solution de (73) que $U(\Phi) \equiv 0$ d'où $\Phi = [U(\Phi)] \equiv 0$. On a donc unicité de la solution de (63).
L'équivalence entre (73) et (63) peut être résumée dans le :

**Lemme 5.3**

a) Les deux assertions suivantes sont équivalentes :

(i) $U \in H^1(\Omega)$ est solution de (73) et $\Phi = [U]$.

(ii) $\Phi \in H^{1/2}(\Gamma)$ est solution de (63) et $U = -M_{\omega}\Phi$.

b) Pour tout $\hat{g} \in H^{-1/2}(\Gamma)$, l'équation (63) admet une unique solution $\Phi \in H^{1/2}(\Gamma)$.

Afin d'étudier l'opérateur $D_{\omega}$, on introduit la formulation variationnelle intégrale de (63), qui grâce à la formule de Green est reliée à la formulation variationnelle volumique (73):

$$ (87) \qquad \left\{ \begin{array}{l} \text{Trouver } \Phi \in H^{1/2}(\Gamma) \text{ telle que} \\ \langle D_{\omega} \Phi, \psi \rangle_{-1/2,1/2,\Gamma} = \langle \hat{g}, \psi \rangle_{-1/2,1/2,\Gamma} \\ \forall \psi \in H^{1/2}(\Gamma) \end{array} \right. $$

Par la suite, on notera

$$ b_{\omega}(\Phi, \psi) = \langle D_{\omega}\Phi, \psi \rangle_{-1/2,1/2,\Gamma} $$

Donnons une "formule de passage" entre la formulation variationnelle volumique et la formulation variationnelle intégrale. On se donne $\Phi \in H^{1/2}(\Gamma)$ et on pose $U(\Phi) = -M_{\omega}\Phi \in H^1(\Omega)$. Pour tout $\psi \in H^{1/2}(\Gamma)$, d'après la surjectivité du saut, on peut définir $v(\psi) \in H^1(\Omega)$ tel que $[v] = \psi$. Et pour tout $\psi \in H^{1/2}(\Gamma)$, on a par définition de $D_{\omega}$:

$$ b_{\omega}(\Phi, \psi) = \left\langle \frac{\partial U(\Phi)}{\partial n}, [v] \right\rangle_{-1/2,1/2,\Gamma} $$

En appliquant la formule de Green, on en déduit que, $\forall \Phi \in H^{1/2}(\Gamma)$ et $\forall \psi \in H^{1/2}(\Gamma)$, on a :

$$ (88) \qquad b_{\omega}(\Phi, \psi) = a(U(\Phi), v(\psi)) $$

où $U(\Phi) = -M_{\omega}\Phi$ et $v(\psi)$ est n'importe quel relèvement de $\psi$ dans $H^1(\Omega)$.

On a l'équivalence, entre les formulations variationnelles volumiques et intégrales, analogue à celle établie au lemme 5.3 :

**Lemme 5.4**

Les deux assertions suivantes sont équivalentes :

(i) $U \in H^1(\Omega)$ est solution de (73) et $\Phi = [U]$.
(ii) $\Phi \in H^{1/2}(\Gamma)$ est solution de (87) et $U = -M_{\omega}\Phi$.

### Démonstration :

(i) ⇒ (ii) Si $U \in H^1(\Omega)$ est solution de (73) et $\Phi = [U]$, d'après le lemme 5.3, $\Phi \in H^{1/2}(\Gamma)$ est solution de (63) et $U = -M_{\omega}\Phi$. Il est clair que alors $\Phi$ vérifie aussi (87).

(ii) ⇒ (i) Si $\Phi \in H^{1/2}(\Gamma)$ est solution de (87) et $U = -M_{\omega}\Phi$ on a pour tout $\psi \in H^{1/2}(\Gamma)$:

$$b_{\omega}(\Phi, \psi) = \langle \hat{g}, \psi \rangle_{-1/2, 1/2, \Gamma}$$

Soit $v \in H^1(\Omega)$, posons $\psi = [v]$. La formule de passage (88) permet de conclure que $U$ vérifie

$$a(U, v) = b_{\omega}(\Phi, \psi) = \langle \hat{g}, [v] \rangle$$

donc (73).

■

### Théorème 5.2

a) Pour $\omega \in \{\Im m(\omega) > \omega_I^0\}$, l'opérateur $D_\omega$ est un isomorphisme entre $H^{1/2}(\Gamma)$ et $H^{-1/2}(\Gamma)$ et sa norme vérifie la majoration suivante

$$(89) \qquad \| D_{\omega} \Phi \|_{-1/2, \Gamma} \le \frac{C}{\omega_I^0 C_m} |\omega|^2 \| \Phi \|_{1/2, \Gamma} \quad \forall \Phi \in H^{1/2}(\Gamma)$$

b) Pour $\omega \in \{\Im m(\omega) > \omega_I^0\}$, la forme bilinéaire de (87) vérifie l'inégalité de coercivité suivante :

$$(90) \qquad \Re \epsilon(b_{\omega}(\Phi, -i\omega\Phi)) \ge C_m^2 \omega_I^0 \|\Phi\|_{1/2, \Gamma}^2 \quad \forall \Phi \in H^{1/2}(\Gamma)$$

c) Pour tout $\hat{g} \in H^{-1/2}(\Gamma)$, le problème (87) admet une unique solution $\Phi \in H^{1/2}(\Gamma)$, également unique solution de (63) et vérifiant :

$$(91) \qquad C_m^2 \omega_I^0 \| \Phi \|_{1/2, \Gamma} \le |\omega| \| \hat{g} \|_{-1/2, \Gamma}$$

### Démonstration : • Pour montrer l'estimation de continuité de $D_\omega$, on utilise la formule de passage. Soient $\Phi$ et $\psi$ dans $H^{1/2}(\Gamma)$. On pose $u(\Phi) = -M_\omega\Phi$ et $v(\psi)$ le relèvement de $\psi$ défini par le lemme 5.2. La formule de passage 88 et la continuité de la forme bilinéaire $a(\cdot, \cdot)$ donnent alors

$$| \langle D_{\omega} \Phi, \psi \rangle | = | a(u(\Phi), v(\psi)) | \le \| u \|_{1, \omega} \| v \|_{1, \omega}$$

On utilise l'estimation (84) majorant la norme du potentiel de double couche en fonction de la norme de sa densité et l'estimation (71) du lemme de relèvement, ce qui donne :

$$| \langle D_{\omega} \Phi, \psi \rangle | \le \frac{C}{\omega_I^0 C_m} |\omega|^2 \| \Phi \|_{1/2, \Gamma} \| \psi \|_{1/2, \Gamma}$$
L'opérateur $D_\omega$ est donc continu et sa norme vérifie (89).

• Si on utilise de nouveau 88 appliquée à $v = -i\omega u$, on obtient l'inégalité de coercivité sur $b_\omega(.,.)$ à partir de celle de $a(.,.)$ et du lemme de trace.

• La forme linéaire

$$ \psi \in H^{1/2}(\Gamma) \longrightarrow \langle \hat{g}, \psi \rangle_{-1/2, 1/2, \Gamma} $$

est bien sûr continue sur $H^{1/2}(\Gamma)$.

• On peut donc appliquer Lax-Milgram et en déduire le théorème. ■

<figure><img src="image_3.png" /></figure>

Afin de faire une transformation de Fourier-Laplace des résultats en fréquence et d'obtenir ainsi des résultats analogues sur le problème en temps, montrons le

**Lemme 5.5 (i) L'opérateur $D_\omega$ est holomorphe dans $\mathbb{C}$**.

(ii) Si $g \in \mathcal{L}'(H^{-1/2}(\Gamma))$, les distributions $U(\omega)$ et $\Phi(\omega)$ définies comme les solutions des problèmes (72) et (63) dans $H^1(\Omega)$ et $H^{1/2}(\Gamma)$ sont holomorphes dans le demi-plan complexe $\{\Im(\omega) > \omega_I^0\}$.

(iii) Si $g \in \mathcal{L}'(H^{-1/2}(\Gamma))$, il existe une distribution $u^D \in \mathcal{L}'(H^1(\Omega))$ telle que $U(\omega) = \hat{u}^D(\omega)$ et une distribution $\varphi \in \mathcal{L}'(H^{1/2}(\Gamma))$ telle que $\Phi(\omega) = \hat{\varphi}(\omega)$.

**Démonstration :**

(i) L'analycité de la solution fondamentale $G$ sur $\mathbb{C}$ implique l'analycité de la forme bilinéaire $b_\omega(\Phi, \psi) = \langle D_\omega \Phi, \psi \rangle$. Or l'analycité “faible” implique l'analycité “forte” (cf Kato, [36]) donc le résultat est encore vrai pour $D_\omega$.

(ii) La distribution $\Phi(\omega) \in H^{1/2}(\Gamma)$ est solution de l'équation intégrale

$$ D_{\omega} \Phi = \hat{g} $$

Or l'opérateur $D_\omega$ est un isomorphisme pour tout $\omega$ tel que $\Im(\omega) > 0$. On en déduit que $D_\omega^{-1}$ est holomorphe dans $\{\omega; \Im(\omega) > 0\}$ (voir [36]). On a supposé que $g \in \mathcal{L}'(H^{-1/2}(\Gamma))$, ce qui signifie en particulier que il existe $\omega_I^0$ tel que $\hat{g}$ est holomorphe dans $\{\Im(\omega) > \omega_I^0\}$ et

$$ \Phi = D_{\omega}^{-1} \hat{g} $$

donc $\Phi$ vérifie la même propriété. La solution $U(\omega)$ est obtenue comme potentiel de double couche de densité $\Phi$, $U(\omega) = -M_\omega \Phi$, et toujours grâce à l'analycité de la solution fondamentale, on peut montrer l'analycité de l'opérateur $M_\omega$.

(iii) Pour montrer ce dernier point, on utilise la caractérisation donnée par le théorème de Paley-Wiener. On a supposé que $g \in \mathcal{L}'(H^{-1/2}(\Gamma))$ ce qui implique l'existence d'une valeur $\omega_I^0(g)$ telle que $\hat{g}$ est holomorphe dans $\{\Im(\omega) > \omega_I^0\}$ et telle que

$$ \begin{gathered} \exists \sigma_1(g) > \omega_I^0(g), \quad C(g) > 0 \quad \text{et} \quad k(g) \ge 0 \quad \text{t. q.} \\ \| \hat{g}(\omega) \|_{-1/2, \Gamma} \le C(1 + |\omega|)^k \quad \forall \omega \text{ t. q. } \Im(\omega) \ge \sigma_1 \end{gathered} $$
L'estimation (91) montre alors que

$$C_m^2 \omega_I^0 \| \Phi \|_{1/2, \Gamma} \le | \omega | \| \hat{g} \|_{-1/2, \Gamma} \le C (1 + | \omega |)^{k+1}$$

ce qui montre le résultat pour $\Phi$. De même, l'estimation (76) permet de conclure pour $U$.

<figure><img src="image_4.png" /></figure>

Avant de revenir au problème en temps, nous donnons un tableau indiquant une comparaison entre les résultats obtenus avec les normes usuelles et les normes liées à l'énergie :
<table><thead><tr><th></th><th>Normes usuelles</th><th>Normes de l'énergie</th></tr></thead><tbody><tr><td>Continuité de l'opérateur :<br>ĝ ∈ H<sup>-1/2</sup>(Γ) → U ∈ H<sup>1</sup>(Ω)</td><td>ω<sub>I</sub><sup>0</sup>C<sub>m</sub><sup>2</sup>(ω<sub>I</sub><sup>0</sup>) || U ||<sub>1,Ω</sub> ≤ C |ω||| ĝ ||<sub>-1/2,Γ</sub></td><td>ω<sub>I</sub><sup>0</sup> || U ||<sub>1,ω,Ω</sub> ≤ C |ω||| ĝ ||<sub>-1/2,ω,Γ</sub></td></tr><tr><td>Continuité de l'opérateur inverse</td><td>|| ĝ ||<sub>-1/2,Γ</sub> ≤ (<sup>C</sup>/<sub>C<sub>m</sub></sub>)<sup>1/2</sup> || U ||<sub>1,ω,Ω</sub></td><td>|| ĝ ||<sub>-1/2,ω,Γ</sub> ≤ C || U ||<sub>1,ω,Ω</sub></td></tr><tr><td>Continuité de l'opérateur :<br>ĝ ∈ H<sup>-1/2</sup>(Γ) → Φ ∈ H<sup>1/2</sup>(Γ)</td><td>ω<sub>I</sub><sup>0</sup>C<sub>m</sub><sup>2</sup>(ω<sub>I</sub><sup>0</sup>) || Φ ||<sub>1/2,Γ</sub> ≤ C |ω||| ĝ ||<sub>-1/2,Γ</sub></td><td>ω<sub>I</sub><sup>0</sup> || Φ ||<sub>1/2,ω,Γ</sub> ≤ C |ω||| ĝ ||<sub>-1/2,ω,Γ</sub></td></tr><tr><td>Continuité du potentiel<br>de double couche<br>U = -M<sub>ω</sub>Φ</td><td>|| U ||<sub>1,ω,Ω</sub> ≤ <sup>C</sup>/<sub>ω<sub>I</sub><sup>0</sup>C<sub>m</sub><sup>1/2</sup></sub> |ω|<sup>3/2</sup> || Φ ||<sub>1/2,Γ</sub></td><td>|| U ||<sub>1,ω,Ω</sub> ≤ <sup>C</sup>/<sub>ω<sub>I</sub><sup>0</sup></sub> |ω||| Φ ||<sub>1/2,ω,Γ</sub></td></tr><tr><td>Continuité de D<sub>ω</sub></td><td>|| D<sub>ω</sub>φ ||<sub>-1/2,Γ</sub> ≤ <sup>C</sup>/<sub>ω<sub>I</sub><sup>0</sup>C<sub>m</sub></sub> |ω|<sup>2</sup> || φ ||<sub>1/2,Γ</sub></td><td>|| D<sub>ω</sub>φ ||<sub>-1/2,ω,Γ</sub> ≤ <sup>C</sup>/<sub>ω<sub>I</sub><sup>0</sup></sub> |ω||| φ ||<sub>1/2,ω,Γ</sub></td></tr><tr><td>Coercivité de b<sub>ω</sub></td><td>Re (D<sub>ω</sub>φ, -iωφ)<sub>-1/2,1/2</sub> ≥ ω<sub>I</sub><sup>0</sup>C<sub>m</sub><sup>2</sup> || φ ||<sub>1/2,Γ</sub></td><td>Re (D<sub>ω</sub>φ, -iωφ)<sub>-1/2,1/2</sub> ≥ Cω<sub>I</sub><sup>0</sup> || φ ||<sub>1/2,ω,Γ</sub></td></tr></tbody></table>
On remarque en particulier que les puissances de $|\omega|$ ne sont pas les mêmes, dans l'estimation de continuité de l'opérateur $D_\omega$, ce qui se traduit par une perte de régularité en temps si on utilise les normes usuelles.

## 5.3 Retour en temps

Les résultats en temps sont obtenus en appliquant une transformée de Fourier-Laplace inverse aux résultats en fréquence. Rappelons l'égalité de Parseval, pour deux distributions de $\mathcal{L}'(E)$:

$$ (92) \qquad \frac{1}{2\pi} \int_{R+i\omega_I} (\hat{f}(\omega), \hat{g}(\omega))_E d\omega = \int_R e^{-2\omega_I t} (f(t), g(t))_E dt $$

Revenons maintenant à l'analyse de l'équation intégrale espace-temps (58). Si la donnée $g$ est dans $\mathcal{L}'(H^{-1/2}(\Gamma))$, d'après l'implication (i) $\implies$ (ii) du théorème de Paley-Wiener, sa transformée $\hat{g}$ vérifie (ii) (a) et (b). L'étude précédente montre que pour tout $\omega$ fixé dans un demi-plan complexe $(\{\Im m(\omega) > \omega_I^0\})$, il existe une unique distribution $\Phi(\omega) \in H^{1/2}(\Gamma)$ solution de l'équation intégrale en fréquences (63). De plus le lemme 5.5 montre que $\Phi(\omega)$ vérifie les conditions (ii) (a) et (b) du théorème de Paley-Wiener. L'implication inverse (ii) $\implies$ (i) permet alors de définir une distribution $\varphi = \mathcal{F}^{-1}(\Phi) \in \mathcal{L}'(H^{1/2}(\Gamma))$ et, d'après la remarque 5.1, $\varphi$ est solution de (58).

Etant donné le type d'estimations obtenues entre $U(\omega)$, $\Phi(\omega)$ et $\hat{g}(\omega)$, on introduit des espaces fonctionnels indiquant un peu plus précisément le comportement par rapport à $\omega$ des solutions. On va supposer non seulement que $g \in \mathcal{L}'(H^{-1/2}(\Gamma))$ mais qu'elle vérifie en plus

$$ (93) \qquad \int_{R+i\omega_I} |\omega|^{2s} \| \hat{g}(\omega) \|_{-1/2, \Gamma}^2 d\omega < +\infty $$

De manière plus générale, pour un espace de Hilbert $E$, on introduit les espaces

$$ \mathcal{H}_{\omega_I}^s(R_+, E) = \{ f \in \mathcal{L}'(E) ; \int_{R+i\omega_I} |\omega|^{2s} \| \hat{f}(\omega) \|_E^2 d\omega < \infty \} $$

et on note

$$ \| f \|_{s,E,\omega_I}^2 = \int_{R+i\omega_I} |\omega|^{2s} \| \hat{f}(\omega) \|_E^2 d\omega $$

la norme dans ces espaces. Dans le cas particulier où $E$ est un espace de Sobolev $E = H^r(\Gamma)$, on notera $\| f \|_{s,r,\omega_I}^2$ cette norme.
On peut remarquer que pour $s = k \in N$, le terme $|\omega|^{2k}$ représente alors une
dérivation en temps d'ordre $k$ et d'après Parseval on a:

$$
\frac{1}{2\pi} \int_{R+i\omega_I} |\omega|^{2k} \| \hat{f}(\omega) \|_{E}^{2} d\omega = \int_{R} e^{-2\omega_I t} \| f^{(k)}(t) \|_{E}^{2} dt
$$

ce qui signifie qu'on peut caractériser cet espace par :

$$
\mathcal{H}_{\omega_I}^k (R_+, E) = \{ f \in \mathcal{L}'(E) ; e^{-\omega_I t} f^{(k)}(t) \in L^2(R, E) \}
$$

On introduit la formulation variationnelle espace-temps, obtenue en intégrant
celle en fréquence :

$$
(94) \qquad \left\{
\begin{array}{ll}
\text{Trouver } \varphi & \text{dans un espace à préciser, telle que} \\
& \\
\quad b(\varphi, \psi) = L(\psi) & \\
\text{pour tout } \psi & \text{appartenant au même espace}
\end{array}
\right.
$$

avec

$$
(95) \qquad \left\{
\begin{aligned}
b(\varphi, \psi) &= \frac{1}{2\pi} \int_{R+i\omega_I} <D_\omega \hat{\varphi}, -i\omega \hat{\psi}> d\omega \\
&= \int_R e^{-2\omega_I t} <D\varphi, \dot{\psi}> dt \\[\medskipamount]
L(\psi) &= \frac{1}{2\pi} \int_{R+i\omega_I} <\hat{g}, -i\omega \hat{\psi}> d\omega \\
&= \int_R e^{-2\omega_I t} <\hat{g}, \dot{\psi}> dt
\end{aligned}
\right.
$$

On précise le cadre fonctionnel dans lequel est posée cette formulation variation-
nelle dans le

Théorème 5.3

a) L'opérateur $D$ est linéaire, continu de $\mathcal{H}_{\omega_I}^s(R_+, H^{1/2}(\Gamma))$ dans $\mathcal{H}_{\omega_I}^{s-2}(R_+, H^{-1/2}(\Gamma))$ et sa norme vérifie la majoration suivante

$$
(96) \quad \| D\varphi \|_{s-2, -1/2, \omega_I} \le C(\omega_I^0) \| \varphi \|_{s, 1/2, \omega_I} \quad \forall \varphi \in \mathcal{H}_{\omega_I}^s (R_+, H^{1/2}(\Gamma))
$$

b) La forme bilinéaire du problème variationnel espace-temps vérifie l'inégalité
de coercivité suivante :

$$
(97) \quad b(\varphi, \varphi) \ge C(\omega_I^0) \| \varphi \|_{0,1/2,\omega_I}^2 \quad \forall \varphi \in \mathcal{H}_{\omega_I}^2 (R_+, H^{1/2}(\Gamma))
$$
## Démonstration :

a) Soit $\varphi \in \mathcal{H}_{\omega_I}^s(R_+, H^{1/2}(\Gamma))$. Par définition $\hat{\varphi}$ est holomorphe dans un demi-plan complexe $\{\Im m(\omega) > \omega_I^0\}$ et vérifie

$$ (i) \quad \int |\omega|^{2s} \| \hat{\varphi} \|_{1/2, \Gamma}^2 d\omega < \infty $$

et d'après le lemme 5.5, $D_\omega$ est holomorphe dans $\mathbb{C}$. On en déduit que $D_\omega \hat{\varphi}$ est holomorphe dans $\{\Im m(\omega) > \omega_I^0\}$. De plus, l'estimation (89) nous donne

$$ \| D_{\omega} \hat{\varphi} \|_{-1/2, \Gamma} \leq \frac{C}{\omega_I^0 C_m} |\omega|^2 \| \hat{\varphi} \|_{1/2, \Gamma} $$

ce qui permet de conclure que d'une part $D_\omega \hat{\varphi}$ est bien la transformée de Fourier-Laplace d'une distribution $(D_\omega \hat{\varphi} = \tilde{D}\varphi)$ de $\mathcal{L}'(H^{-1/2}(\Gamma))$, et d'autre part que

$$ \int |\omega|^{2(s-2)} \| D_{\omega} \hat{\varphi} \|_{-1/2, \Gamma}^2 d\omega \leq \left( \frac{C}{\omega_I^0 C_m} \right)^2 \int |\omega|^{2s} \| \hat{\varphi} \|_{1/2, \Gamma}^2 d\omega < \infty $$

l'inégalité de droite provenant de (i). Ceci montre que $D\varphi \in \mathcal{H}_{\omega_I}^{(s-2)}(R_+, H^{-1/2}(\Gamma))$ et qu'on a l'estimation (96) avec $C(\omega_I^0) = \frac{C}{\omega_I^0 C_m}$.

b) Pour $\varphi \in \mathcal{H}_{\omega_I}^2(R_+, H^{1/2}(\Gamma))$, on a :

$$ \begin{aligned} b(\varphi, \varphi) &= \int_R e^{-2\omega_I t} \langle D\varphi, \dot{\varphi} \rangle dt \\ &= \frac{1}{2\pi} \int_{R+i\omega_I} \langle D_\omega \hat{\varphi}, -i\omega \hat{\varphi} \rangle d\omega \end{aligned} $$

On utilise l'inégalité de coercivité sur $b_\omega$, (90) :

$$ \begin{aligned} b(\varphi, \varphi) &\ge \frac{1}{2\pi} \int_{R+i\omega_I} C_m^2 \omega_I^0 \| \hat{\varphi} \|_{1/2, \Gamma}^2 d\omega \\ &\ge C_m^2 \omega_I^0 \| \varphi \|_{0, 1/2, \omega_I}^2 \end{aligned} $$

<figure><img src="image_5.png" /></figure>

**Remarque 5.3** : on peut remarquer que l'inégalité (97) n'est pas réellement une inégalité de coercivité car le terme de droite est seulement une semi-norme sur $\mathcal{H}_{\omega_I}^2(R_+, H^{1/2}(\Gamma))$.

Ces résultats permettent d'obtenir le

## Théorème 5.4

Si $g \in \mathcal{H}_{\omega_I}^3(R_+, H^{-1/2}(\Gamma))$ alors

a) le problème (94) admet une unique solution $\varphi \in \mathcal{H}_{\omega_I}^2(R_+, H^{1/2}(\Gamma))$ également unique solution de (58)et vérifiant :

$$ (98) \qquad C_m^2 \omega_I^0 \| \varphi \|_{2, 1/2, \omega_I} \le \| g \|_{3, -1/2, \omega_I} $$
b) Le problème (55), (57) admet une unique solution $u^D \in \mathcal{H}_{\omega_I}^2(R_+, H^1(\Omega))$ qui vérifie

$$ (99) \qquad C_m^2 \omega_I^0 \| u^D \|_{2, H^1(\Omega), \omega_I} \le C \| g \|_{3, -1/2, \omega_I} $$

et qui s'exprime par

$$ (100) \qquad u^D = -M\varphi $$

où $\varphi$ est la solution définie au a).

## Démonstration :

a) On part de $g \in \mathcal{H}_{\omega_I}^3(R_+, H^{-1/2}(\Gamma))$. On peut donc prendre sa transformée de Fourier-Laplace et définir la solution $\Phi$ du problème en fréquence associé. D'après le lemme (5.5), $\Phi$ est la transformée de Fourier-Laplace d'une distribution $\varphi \in L'(H^{1/2}(\Gamma))$. On a même un peu plus, grâce à l'estimation (91):

$$ C_m^2 \omega_I^0 \| \Phi \|_{1/2, \Gamma} \le |\omega| \| \hat{g} \|_{-1/2, \Gamma} $$

On en déduit que $\varphi \in \mathcal{H}_{\omega_I}^2(R_+, H^{1/2}(\Gamma))$ et vérifie (98). Il est enfin immédiat de vérifier que $\varphi$ est solution de (58) puisque $\Phi$ est solution de l'équation intégrale en fréquence. On a donc montré l'existence d'une solution de l'équation intégrale ainsi que de sa formulation variationnelle.

L'unicité découle très clairement de l'unicité de la solution du problème en fréquence.

b) On définit de nouveau $u^D$ à partir de $U(\omega)$, solution du problème en fréquence. On définit ainsi une solution de (55), (57) dans $\mathcal{H}_{\omega_I}^2(R_+, H^1(\Omega))$ grâce à l'estimation (76). ■

L'étude en fréquence menée avec les normes liées à l'énergie conduit aux espaces fonctionnels suivants :

$$ H_{+}^{s, s_0, \omega_I}(\Gamma \times R) = \{f \text{ t. q. } \int_{R+i\omega_I} |\omega|^{2s_0} \| \hat{f}(\omega) \|_{s, \omega, \Gamma}^2 d\omega < +\infty\} $$

On peut retrouver tous les résultats énoncés dans ces espaces en utilisant le tableau de comparaison donné en fin de partie 5.2. En particulier, l'opérateur $D$ est alors continu de $H_+^{1/2, s_0, \omega_I}(\Gamma \times R)$ dans $H_+^{-1/2, s_0-1, \omega_I}(\Gamma \times R)$, ce qui constitue une différence essentielle avec les espaces usuels car on remarque qu'on ne perd qu'un cran de régularité en temps au lieu de deux.
# 6 Mise en oeuvre et résultats numériques

## 6.1 Discrétisation

Nous traitons ici le problème modèle étudié précédemment, c'est à dire le problème de Neumann. La formulation variationnelle de ce problème s'écrit :

$$ (101) \quad \left\{ \begin{array}{l} \text{Pour } g \in \mathcal{H}_{\omega_I}^3(R_+, H^{-1/2}(\Gamma)) \text{ trouver } \varphi \in \mathcal{H}_{\omega_I}^2(R_+, H^{1/2}(\Gamma)) \text{ telle que} \\ \qquad \int_R e^{-2\omega_I t} \langle D\varphi, \dot{\psi} \rangle dt = \int_R e^{-2\omega_I t} \langle g, \dot{\psi} \rangle dt \\ \text{et } \forall \psi \in \mathcal{H}_{\omega_I}^2(R_+, H^{1/2}(\Gamma)) \end{array} \right. $$

On approche alors cette formulation en utilisant des éléments finis frontières, suivant la démarche proposée par Nédélec ([45]). Cette démarche étant classique, on se contente de la rappeler brièvement et on se limite au cas où $\Gamma$ est un domaine polyédrique. Dans le cas général, il faudrait définir une application transformant chaque morceau en un domaine polyédrique et tenir compte de l'erreur d'approximation de la surface. On découpe la surface $\Gamma$ en $p$ morceaux fermés $(\Gamma_i)_{i=1,...,p}$ tels que

$$ (a) \quad \cup_{i=1}^p \Gamma_i = \Gamma $$

$$ (b) \quad \Gamma_i \cap \Gamma_j = \emptyset \text{ ou est une courbe portée par } \Gamma \text{ (en 3D).} $$

### Semi-discrétisation en espace

On introduit un espace $V_h$ de dimension finie construit à l'aide d'éléments finis approchant $H^{1/2}(\Gamma)$. Le problème semi-discrétisé consiste à

$$ (102) \quad \left\{ \begin{array}{l} \text{trouver } \varphi_h \in \mathcal{H}_{\omega_I}^2(R_+, V_h) \text{ telle que} \\ b(\varphi_h, \psi_h) = L_h(\psi_h) \quad \forall \psi_h \in \mathcal{H}_{\omega_I}^2(R_+, V_h) \end{array} \right. $$

$$ (103) \quad \text{où} \quad L_h(\psi_h) = \int_R e^{-2\omega_I t} \langle g_h, \dot{\psi}_h \rangle dt $$

et $g_h$ est une approximation de $g$.

La construction de $V_h$ est décrite par Nédélec dans [45]. L'espace $V_h$ est constitué de fonctions continues et polynômiales par morceaux de degré $m_1 \ge 1$ ce qui assure l'appartenance à $H^{1/2}(\Gamma)$. Le choix le plus simple est celui de l'élément fini $P_1$ constitué de fonctions affines par élément et continues sur $\Gamma$. Soit $(\varphi_h^j)_{j=1,...,N_h}$ une base de $V_h$. La solution approchée $\varphi_h(t)$ se décompose sur cette base :

$$ (104) \quad \varphi_h(x, t) = \sum_{j=1}^{N_h} \alpha_j(t) \varphi_h^j(x) \quad \forall x \in \Gamma $$
où les fonctions $\alpha_j(t)$ doivent appartenir à $\mathcal{H}_{\omega_I}^2(R_+, R)$ pour que $\varphi_h \in \mathcal{H}_{\omega_I}^2(R_+, V_h)$.
Le problème approché consiste à trouver $\alpha_j(t) \in \mathcal{H}_{\omega_I}^2(R_+, R)$ telles que

$$ (105) \qquad \sum_{j=1}^{N_h} b(\alpha_j(t)\varphi_h^j(x), \psi_h) = L(\psi_h) \quad \forall \psi_h \in \mathcal{H}_{\omega_I}^2(R_+, V_h) $$

## Discrétisation totale

On doit approcher les fonctions en temps $\alpha_j(t) \in \mathcal{H}_{\omega_I}^2(R_+, R)$. L'espace $\mathcal{H}_{\omega_I}^2(R_+, R)$ est approché par le sous-espace $H^{m_2}(\Delta t, R)$ constitué de fonctions continues sur $R$ et polynômiales de degré $m_2 \ge 2$ dans chaque intervalle de temps $I_n = [t_n, t_{n+1}]$ où $t_n = n\Delta t$.

**Remarque 6.1** : $\alpha \in \mathcal{H}_{\omega_I}^2(R_+, R)$ est équivalent à dire que $\alpha \in \mathcal{L}'(R_+)$ et $e^{-\omega_I t}\alpha'' \in L^2(R_+, R)$. Si on approche $\alpha$ par des polynômes par morceaux de degré inférieur à 2 (c'est à dire par des fonctions affines par intervalle) sa dérivée seconde devient une somme de dirac et n'est donc pas $L^2$. Ce qui explique que le plus bas degré permis est $m_2 \ge 2$. On verra cependant que, en pratique, on utilise des fonctions affines par morceaux. On peut encore donner un sens aux intégrales. La méthode ainsi obtenue donne de bons résultats. Mais on ne peut pas lui appliquer l'analyse numérique qui va suivre.

Dans ce qui suit, nous allons traiter en détail le cas $m_2 = 2$. La fonction $\varphi_h$ est donc approchée par

$$ (106) \qquad \varphi_{h,\Delta t}(x,t) = \sum_{j=1}^{N_h} \alpha_{j,\Delta t}(t)\varphi_h^j(x) \quad \forall x \in \Gamma, \forall t $$

avec

$$ \alpha_{j,\Delta t/I_n} \in P_2, \forall n $$

La dérivée seconde de $\alpha_{j,\Delta t}$ est alors une constante par intervalles et on peut montrer facilement que, en notant $a_j^n$ cette constante, on a

$$ (107) \qquad \alpha_{j,\Delta t} = \frac{1}{2}(t-t_n)^2 a_j^n + \Delta t \sum_{m=0}^{n-1} (t-t_m + \frac{1}{2}\Delta t) a_j^m \quad \forall t \in I_n $$

Rappelons la forme bilinéaire du problème :

### – en dimension 3

la forme bilinéaire est donnée par (33) :

$$ \left\{ \begin{aligned} b(\varphi, \psi) = & \int_R e^{-2\omega_I t} \frac{1}{4\pi} \int_{\Gamma \times \Gamma} \frac{n_x \cdot n_y}{|x-y|} \ddot{\varphi}(y, t-|x-y|) \dot{\psi}(x,t) d\gamma_x d\gamma_y dt \\ & + \int_R e^{-2\omega_I t} \frac{1}{4\pi} \int_{\Gamma \times \Gamma} \frac{\vec{rot}_\Gamma \varphi(y, t-|x-y|) \cdot \vec{rot}_\Gamma \dot{\psi}(x,t)}{|x-y|} d\gamma_x d\gamma_y dt \end{aligned} \right. $$
- en dimension 2

la forme bilinéaire est donnée par (34) :

$$
\left\{
\begin{array}{l}
b(\varphi, \psi) = \int_R e^{-2\omega_I t} \int_{\Gamma \times \Gamma} \frac{H(t - |x - y|)}{2\pi} \int_0^{t - |x - y|} \frac{1}{((t - \tau)^2 - |x - y|^2)^{1/2}} \times \\
\qquad \left( n_x \cdot n_y \ddot{\varphi}(y, \tau) \dot{\psi}(x, t) + \frac{d\varphi}{ds}(y, \tau) \frac{d\dot{\psi}}{ds}(x, t) \right) d\tau d\gamma_x d\gamma_y dt
\end{array}
\right.
$$

La fonction test $\psi_{h,\Delta t}$ n'intervient dans la formulation variationnelle que par l'intermédiaire de sa dérivée en temps. On la choisit comme le produit d'une fonction test en temps $\beta_{\Delta t}^n(t)$ par une fonction test en espace $\varphi_h^i(x)$. Il suffit, pour que les fonctions tests décrivent les fonctions $P_2$ par morceaux, que leurs dérivées décrivent les fonctions continues et $P_1$ par morceaux, c'est à dire qu'on choisit les $\beta_{\Delta t}^n(t)$ telles que

$$
(108) \qquad (\beta_{\Delta t}^n)'(t) = \begin{cases} t - t_{n-1} & \text{si } t_{n-1} < t < t_n \\ t_{n+1} - t & \text{si } t_n < t < t_{n+1} \\ 0 & \text{sinon} \end{cases}
$$

En reportant dans la formulation variationnelle, on obtient un système de la
forme :

$$
(109) \qquad \sum_{k=0}^{n} M^{n-k} A^k = B^n
$$

où les vecteurs $A^k$ sont les inconnues $(a_1^k, \dots, a_{N_h}^k)^t$ et les matrices $M^{n-k}$ ne dépendent que de $n-k$.

**Remarque 6.2** : on peut montrer que, en 3D, les matrices $M^{n-k}$ s'annulent dès que $n-k > \left[ \frac{d(\Gamma)}{\Delta t} \right]$. Pour $n$ grand, on n'a donc à prendre en compte dans la somme que les termes compris entre $n - \left[ \frac{d(\Gamma)}{\Delta t} \right]$ et $n$. Ceci n'est plus vrai en 2D, où on doit garder tous les termes depuis le début, ce qui peut être une limitation en ce qui concerne le stockage, si on s'intéresse à des temps élevés.

On peut montrer un certain nombre de propriétés du schéma (109).

Analyse du schéma

L'analyse de la formulation variationnelle du problème continu permet d'obtenir
une étude du schéma, des estimations d'erreur et de stabilité. On se contente ici
d'énoncer les résultats :

**Lemme 6.1** La matrice $M^0$ est symétrique, définie, positive pour $\Delta t$ assez petit.
**Démonstration :** voir [5] pour le 3D.

**Théorème 6.1**

a) **Stabilité.** Si $g_{h,\Delta t}$ est une approximation consistante de $g$ dans $\mathcal{H}_{\omega_I}^1(R_+, H^{-1/2}(\Gamma))$ le schéma est stable au sens que

$$ (110) \qquad \| \varphi_{h,\Delta t} \|_{0,1/2,\omega_I} \le Cte \\ \text{quand } h, \Delta t \to 0. $$

b) **Estimations d'erreur.** Si $g$ est suffisamment régulière, le schéma converge quand $m_1 \ge 1$ et $m_2 \ge 3$ et on a

$$ (111) \qquad \| \varphi - \varphi_{h,\Delta t} \|_{0,1/2,\omega_I} \le C \left\{ \begin{aligned} & \| g - g_{h,\Delta t} \|_{1,-1/2,\omega_I} \\ & + \frac{h^{m_1+1/2}}{\Delta t} \| \varphi \|_{2,m_1+1,\omega_I} + \Delta t^{m_2-2} \| \varphi \|_{m_2+1,1/2,\omega_I} \end{aligned} \right\} $$

## 6.2 Méthodes numériques. Mise en oeuvre

Les estimations d'erreur obtenues précédemment ne sont pas optimales, en particulier en ce qui concerne la régularité en temps. Il est donc tentant de prendre des fonctions tests et l'inconnue moins régulières en temps. D'autant plus que les calculs se trouveront ainsi simplifiés. En pratique, on va donc considérer une discrétisation à l'aide d'éléments $P^1$ en temps. On aborde dans cette partie les calculs concernant la mise en oeuvre de la méthode de façon un peu plus détaillée et on indiquera les problèmes d'ordre numériques qui peuvent être rencontrés.

On reprend les notations de la partie précédente.

La première simplification consiste à choisir $\omega_I = 0$. Dans ce cas, on ne dispose plus des estimations obtenues précédemment. Cependant, le paramètre $\omega_I$ permet de contrôler le comportement à l'infini de la solution, or on s'intéresse à cette solution sur des intervalles de temps fini, et on s'attend à ce que ce paramètre n'ait alors pas trop d'influence.

Nous notons les fonctions de base en temps :

$$ (112) \qquad \gamma^m(t) = \begin{cases} t - t_{m-1} & \text{si } t_{m-1} < t < t_m \\ t_{m+1} - t & \text{si } t_m < t < t_{m+1} \\ 0 & \text{sinon} \end{cases} $$

c'est à dire $\gamma^m(t) = (\beta_{\Delta t}^m)'(t)$ définie en (108) et $\chi^m$ la fonction indicatrice de l'intervalle $I_m = [t_m, t_{m+1}]$ :

$$ (113) \qquad \chi^m(t) = \begin{cases} 1 & \text{si } t_m < t < t_{m+1} \\ 0 & \text{sinon} \end{cases} $$
Les dérivées de ces fonctions de base s’expriment alors comme des combinaisons linéaires des fonctions indicatrices :

$$ (114) \qquad (\gamma^m)' = \chi^{m-1} - \chi^m $$

Les fonctions tests n’interviennent que par l’intermédiaire de leurs dérivées en temps. On peut donc prendre des fonctions tests dont la dérivée en temps soit le produit d’une fonction de base en espace par une fonction indicatrice :

$$ (115) \qquad \psi_i^n(x, t) = \chi^n(t)\varphi_i(x) $$

où les $\varphi_i$ ont été définies dans la partie précédente comme les fonctions de base de l’élément fini $P^1$ en espace. La fonction $\varphi_h$ est approchée par

$$ (116) \qquad \varphi_{h,\Delta t}(x,t) = \sum_{j=1}^{N_h} \sum_{m \ge 1} a_j^m \gamma^m(t) \varphi_j(x) \quad \forall x \in \Gamma, \forall t $$

les $a^m, m \le 0$ étant nuls d’après la causalité de $\varphi$.

La formulation variationnelle approchée consiste donc à chercher les $a_j^m$ vérifiant :

$$ (117) \qquad \sum_{j=1}^{N_h} \sum_{m \ge 1} a_j^m b(\gamma^m(t)\varphi_j(x), \psi_i^{n-1}) = L(\psi_i^{n-1}) \quad \forall i = 1, \dots, N_h \text{ et } \forall n \ge 1 $$

ce qui peut encore s’écrire sous forme vectorielle :

$$ (118) \qquad \sum_{m \ge 1} M^{nm} A^m = B^n \quad \forall n \ge 1 $$

où $A^m = (a_1^m \dots a_{N_h}^m)^t$ est l’inconnue à l’étape $m$ et

$$ (119) \qquad \left\{ \begin{array}{l} M_{ij}^{nm} = b(\gamma^m(t)\varphi_j(x), \psi_i^{n-1}) \\ B_i^n = L(\psi_i^{n-1}) \end{array} \right. $$

De même que dans la partie précédente, on peut réécrire ce système sous forme d’une convolution discrète :

$$ (120) \qquad \sum_{m=1}^{n} M^{n-m} A^m = B^n \quad \forall n \ge 1 $$

et dans ce cas simple, nous allons le démontrer. Nous allons en fait indiquer en détail comment on traite les différents termes de la forme bilinéaire qui contribuent à la construction de la matrice $M$. Remarquons d’abord que le choix des fonctions de
base $P^1$ en temps n'est pas complètement naturel, car le terme contenant $\ddot{\varphi}$ fait alors
apparaître des masses de dirac en temps:

$$
(121) \quad \ddot{\varphi}_{h,\Delta t}(x,t) = \sum_{j=1}^{N_h} \sum_{m \ge 1} a_j^m (\delta^{m-1} - 2\delta^m + \delta^{m+1})(t) \varphi_j(x) \quad \forall x \in \Gamma, \forall t
$$

Il est donc nécessaire de donner un sens à ces intégrales, ce qui est possible
grâce à la deuxième intégration en temps (par rapport à *t*) qui a un effet régularisant.
Distinguons les cas, dans chaque dimension:

• en dimension 3

la forme bilinéaire est donnée par (33), ce qui donne comme matrice :

$$
\left\{
\begin{aligned}
M_{ij}^{nm} &= \int_R \int_{\Gamma \times \Gamma} \left( \frac{n_x \cdot n_y}{4\pi |x-y|} \varphi_i(y) \varphi_j(x) (\delta^{m-1} - 2\delta^m + \delta^{m+1}) (t-|x-y|) \chi^{n-1}(t) \right. \\
&\qquad \left. + \frac{\vec{rot}_\Gamma \varphi_i(y) \cdot \vec{rot}_\Gamma \varphi_j(x)}{4\pi |x-y|} \gamma^m (t-|x-y|) \chi^{n-1}(t) \right) d\gamma_x d\gamma_y dt
\end{aligned}
\right.
$$

Pour le premier terme, on doit évaluer des intégrales contenant un dirac, du
type :

$$
\int_R \delta^k (t - |x - y|) \chi^{n-1}(t) dt
$$

qu'on interprète comme le produit au sens des distributions :

$$
< \delta(t - |x - y| - t_k), \chi^{n-1}(t) > = \chi^{n-1}(t_k + |x - y|)
$$

Les termes de la matrices $M$ peuvent alors se réécrire sous la forme :

$$
\left\{
\begin{aligned}
M^{nm} &= \int_R \int_{\Gamma \times \Gamma} R(x, y) \gamma^m(t - |x - y|) \chi^{n-1}(t) d\gamma_x d\gamma_y dt \\
&\quad + \int_{\Gamma \times \Gamma} F(x, y) \left( \begin{aligned}[t]
&\chi^{n-1}(t_{m-1} + |x - y|) \\
&-2\chi^{n-1}(t_m + |x - y|) + \chi^{n-1}(t_{m+1} + |x - y|)
\end{aligned} \right) d\gamma_x d\gamma_y
\end{aligned}
\right.
$$

Il est alors très facile de vérifier que :

* Pour $m \ge n + 1$ les intégrales sont nulles pour des raisons de supports.

• En utilisant que

$$
\chi^n(t) = \chi^0(t - t_n)
$$

$$
\gamma^m(\tau) = \gamma^1(\tau - t_{m-1})
$$
et en faisant un changement de variable dans la seconde intégrale, ces termes ne dépendent que de la différence $n - m$. On notera donc dorénavant $M^{nm} = M^{n-m}$, et on a:

$$ \left\{ \begin{array}{l} M^k = \int_0^{\Delta t} \int_{\Gamma \times \Gamma} R(x, y)\gamma^1(t + t_k - |x - y|)d\gamma_x d\gamma_y dt \\ \\ + \int_{\Gamma \times \Gamma} F(x, y)(\chi^0(|x - y| - t_k) - 2\chi^0(|x - y| - t_{k-1}) + \chi^0(|x - y| - t_{k-2}))d\gamma_x d\gamma_y \end{array} \right. $$

Regardons plus attentivement les domaines sur lesquels on va devoir intégrer (en espace). La surface $\Gamma$ est divisée en éléments $\Gamma_l$. L'intégrale double sur $\Gamma \times \Gamma$ va donc donner des intégrales sur $\Gamma_l \times \Gamma_p$. Si on suppose qu'on fait une première intégration de manière numérique, cela revient à se fixer un point $x$ dans l'élément $\Gamma_l$. Il reste donc à évaluer sur $\Gamma_p$:

• une première intégrale :

$$ \int_{\Gamma_p} F_{ij}(x, y)(\chi^0(|x - y| - t_k) - 2\chi^0(|x - y| - t_{k-1}) + \chi^0(|x - y| - t_{k-2}))d\gamma_y $$

et on voit sur cette expression qu'on va intégrer sur l'intersection de $\Gamma_p$ avec des fronts d'onde, c'est à dire sur des domaines du type :

$$ \{y \in \Gamma_p \text{ t. q. } t_m \le |x - y| \le t_{m+1}\} $$

pour $m = k, k-1, k-2$.

• une deuxième intégrale :

$$ \int_{\Gamma_p} R_{ij}(x, y) \int_0^{\Delta t} \gamma^1(t + t_k - |x - y|) dtd\gamma_y $$

cette fois-ci l'intégrale en temps ne disparait pas ce qui complique un peu la situation. Etant donné le support de $\gamma^1$, cette intégrale s'exprime encore :

$$ \int_{\Gamma_p} R_{ij}(x, y) \int_{\mathcal{O}^k(x,y)} \gamma^1(t + t_k - |x - y|) dtd\gamma_y $$

où $\mathcal{O}^k(x, y) = [0, \Delta t] \cap [|x - y| - t_k, |x - y| - t_{k-2}]$. Le domaine $\mathcal{O}^k(x, y)$ est donc vide pour $x, y$ tels que $|x - y| \le t_{k-2}$ ou $|x - y| \ge t_{k+1}$. Il reste donc à intégrer

$$ \int_{y \in \Gamma_p} \text{ et } t_{k-2} \le |x-y| \le t_{k+1} R_{ij}(x, y) \int_{\mathcal{O}^k(x,y)} \gamma^1(t + t_k - |x - y|) dtd\gamma_y $$

et les domaines d'intégration sont du même type que précédemment.

On peut remarquer que, en 3D, le calcul de l'intégrale en temps est très simple car il reste seulement la fonction de base $\gamma^1$ à intégrer c'est à dire une fonction $P^1$ par morceaux. Ceci provient du fait que le noyau contient les variables temps et espace de façon découplée. Ceci ne sera plus vrai en 2D. La première intégration sur $\Gamma$ est
calculée numériquement (par une formule de Gauss par exemple). En ce qui concerne
la deuxième intégration, il faut être prudent, car les noyaux ont une singularité faible
en $x = y$ du type $\frac{1}{|x - y|}$. Il faut donc distinguer les calculs selon que les éléments
sont disjoints ou qu'ils sont confondus ou adjacents. Dans le premier cas, il n'y a pas de
problème de singularité et on peut donc de nouveau utiliser une formule d'intégration
numérique sur $\Gamma_p$. Dans le second cas, on doit calculer l'intégrale de manière exacte.

• en dimension 2

la forme bilinéaire est donnée par (34), ce qui donne comme matrice :

$$
\left\{
\begin{aligned}
M_{ij}^{nm} &= \int_R \int_{\Gamma \times \Gamma} \frac{H(t - |x - y|)}{2\pi} \int_0^{t - |x - y|} \frac{1}{((t - \tau)^2 - |x - y|^2)^{1/2}} \times \\
&\quad (n_x \cdot n_y(\delta^{m-1} - 2\delta^m + \delta^{m+1})(\tau)\chi^{n-1}(t)\varphi_i(y)\varphi_j(x) \\
&\quad + \gamma^m(\tau)\chi^{n-1}(t)\frac{d\varphi_i}{ds}(y)\frac{d\varphi_j}{ds}(x)) \, d\tau \, d\gamma_x \, d\gamma_y \, dt
\end{aligned}
\right.
$$

Pour le premier terme, on doit évaluer des intégrales contenant un dirac et qu’on définit par :

$$
\int_0^{t-|x-y|} \frac{\delta^k(\tau)}{((t-\tau)^2 - |x-y|^2)^{1/2}} d\tau = \frac{H(t-t_k - |x-y|)}{((t-t_k)^2 - |x-y|^2)^{1/2}}
$$

pour $k = m - 1, m, m + 1$. On obtient donc :

$$
\left\{
\begin{aligned}
M_{ij}^{nm} &= \int_R \int_{\Gamma \times \Gamma} \frac{H(t - |x - y|)}{2\pi} \int_0^{t - |x - y|} \frac{\gamma^m(\tau)}{((t - \tau)^2 - |x - y|^2)^{1/2}} d\tau \times \\
&\quad \chi^{n-1}(t) \frac{d\varphi_i}{ds}(y) \frac{d\varphi_j}{ds}(x) d\gamma_x d\gamma_y dt \\
&\quad + \int_R \int_{\Gamma \times \Gamma} \frac{H(t - |x - y|)}{2\pi} n_x \cdot n_y \chi^{n-1}(t) \varphi_i(y) \varphi_j(x) \times \\
&\quad \left( \frac{H(t - t_{m-1} - |x - y|)}{((t - t_{m-1})^2 - |x - y|^2)^{1/2}} - 2 \frac{H(t - t_m - |x - y|)}{((t - t_m)^2 - |x - y|^2)^{1/2}} + \frac{H(t - t_{m+1} - |x - y|)}{((t - t_{m+1})^2 - |x - y|^2)^{1/2}} \right) d\gamma_x d\gamma_y dt
\end{aligned}
\right.
$$

De nouveau, on peut vérifier que :

$$
\bullet M^{nm} = 0 \text{ pour } m \ge n + 1.
$$
• En faisant des changements de variables dans les intégrales en temps, ces termes ne dépendent que de la différence $n - m$. On notera encore $M^{nm} = M^{n-m}$, et on a:

$$ M_{ij}^k = \int_{\Gamma \times \Gamma} \int_0^{\Delta t} \frac{H(t + t_k - |x - y|)}{2\pi} \chi^0(t) \frac{d\varphi_i}{ds}(y) \frac{d\varphi_j}{ds}(x) \times \\ \qquad \int_0^{t+t_k-|x-y|} \frac{\gamma^1(\tau)}{((t - \tau + t_k)^2 - |x-y|^2)^{1/2}} d\tau dt d\gamma_x d\gamma_y \\ \qquad + \int_{\Gamma \times \Gamma} \int_0^{\Delta t} \frac{n_x \cdot n_y}{2\pi} \chi^0(t) dt \varphi_i(y) \varphi_j(x) \left( \frac{H(t + t_k - |x - y|)}{((t + t_k)^2 - |x-y|^2)^{1/2}} \right. \\ \qquad \left. - 2 \frac{H(t + t_{k-1} - |x - y|)}{((t + t_{k-1})^2 - |x-y|^2)^{1/2}} + \frac{H(t + t_{k-2} - |x - y|)}{((t + t_{k-2})^2 - |x-y|^2)^{1/2}} \right) dt d\gamma_x d\gamma_y $$

Les deux termes sont donc de la forme :

$$ \int_{\Gamma \times \Gamma} \int_0^{\Delta t} H(t + t_m - |x - y|) N(x, y, t) dt d\gamma_x d\gamma_y $$

avec $m = k, k-1, k-2$ et avec un noyau $N(x, y, t)$ différent pour chaque terme. Encore une fois, regardons plus attentivement les domaines sur lesquels on va devoir intégrer (en espace). La surface $\Gamma$ est divisée en éléments $\Gamma_l$. L'intégrale double sur $\Gamma \times \Gamma$ va donc donner des intégrales sur $\Gamma_l \times \Gamma_p$. Si on suppose qu'on fait une première intégration de manière numérique, cela revient à se fixer un point $x$ dans l'élément $\Gamma_l$. Il reste donc une intégrale à évaluer sur $\Gamma_p$ et il est facile de voir que cette intégrale est nulle si $|x-y| \ge t_{m+1}$. Il reste donc

$$ \int_{y \in \Gamma_p; |x-y| \le t_{m+1}} \int_{sup(0, |x-y|-t_m)}^{\Delta t} N(x, y, t) dt d\gamma_y $$

Cette fois-ci les domaines d'intégration sont donc les intersections de $\Gamma_p$ avec l'intérieur du disque de centre $x$ et de rayon $t_{m+1}$ pour $m = k, k-1, k-2$ (au lieu des couronnes dans le cas 3D). Les noyaux sont plus compliqués qu'en 3D à cause du couplage entre la variable temps et la variable espace du terme $((t - \tau)^2 - |x - y|^2)^{-1/2}$. En 2D, il y a également une singularité au bord du cône (c'est à dire pour $\tau = t - |x - y|$). En pratique, on calcule analytiquement les intégrales en temps. Il reste alors à intégrer sur $\Gamma_p \cap \{y; |x - y| \le t_{m+1}\}$ un noyau $K(x, y)$ qui a une singularité faible en $Log |x - y|$. Il faut de nouveau distinguer les cas où les éléments sont disjoints et où ils sont confondus et adjacents.

## 6.3 Autres problèmes traités et quelques unes de leurs difficultés respectives...

La méthode présentée a été développée initialement par T. Ha Duong (voir par exemple [22], [5], [4], [23]). Cette démarche a été reprise et appliquée depuis
dans plusieurs autres problèmes, et nous en signalons un certain nombre (par ordre chronologique!). Rappelons que pour l'étude mathématique, deux cadres fonctionnels peuvent être envisagés : les espaces de Sobolev munis des normes usuelles et les espaces de Sobolev munis de normes liées à l'énergie et indexés par $\omega_I$.

Ding, [20], a le premier mis en oeuvre cette méthode pour un problème de diffraction d'ondes acoustiques par une surface rigide. Il a traité l'équation intégrale de deuxième espèce liée au problème de Neumann. Ce sont les premiers résultats qui ont permis de comparer la méthode variationnelle des potentiels retardés avec les méthodes de collocation jusqu'alors utilisées. Elle s'est révélée beaucoup plus robuste, stable et précise que les méthodes de collocation mais bien sûr... plus coûteuse aussi.

E. Bécache, [6], [8], [7], a traité un problème de diffraction d'ondes élastiques avec condition de Neumann. L'étude théorique concerne la diffraction soit par un obstacle soit par une fissure dans un domaine 2D ou 3D. Il s'agit alors du problème de Neumann et la solution est représentée par un potentiel de double couche. Une des difficultés essentielles de ce problème provient de l'hypersingularité du noyau. En ce qui concerne l'équation des ondes acoustiques, on a présenté dans la partie 5 une formulation qui traitait cette difficulté. Pour le problème élastique, on peut obtenir une formulation analogue en utilisant une méthode de régularisation décrite dans [9]. Une autre difficulté, liée à l'élasticité, est la présence de plusieurs types d'ondes. Contrairement aux ondes acoustiques scalaires, il y a dans un milieu élastique des ondes de pression et des ondes de cisaillement. En 3D, la formulation fait intervenir tous les instants compris entre le passage des deux types d'ondes au lieu d'un seul instant et en 2D elle doit prendre en compte le passé des deux ondes. La méthode a été implémentée pour la diffraction par une fissure en 2D.

A. Pujols, [50], a traité la diffraction d'ondes électromagnétiques par un objet. L'étude théorique et l'analyse numérique ont été faites, dans les espaces munis des normes usuelles, pour :

- un obstacle cylindrique infini : dans ce cas, les équations sont ramenées à des équations d'ondes 2D. Elle a traité deux types de conditions aux limites : un obstacle réfléchissant (condition de Dirichlet) et un obstacle absorbant (condition dissipative).

- un obstacle borné parfaitement conducteur. Les équations sont alors les équations de Maxwell 3D.

et l'implémentation a été faite pour le problème de Dirichlet 2D. Des comparaisons avec des méthodes en fréquences ont été présentées afin de valider la méthode.

I. Terrasse, [54], a également traité un problème de diffraction d'ondes électromagnétiques par un obstacle parfaitement conducteur en domaine 3D. Les équations sont alors les équations de Maxwell 3D vectorielles. Pour l'étude théorique, elle s'est placée dans des espaces munis de normes indexées par la fréquence qui permettent
d'obtenir des résultats optimaux. Ces espaces couplent le comportement en espace et celui en temps. Pour comprendre la différence entre les travaux de A. Pujols et ceux de I. Terrasse examinons la relation de conservation de la charge (reliant la divergence surfacique du courant et la dérivée en temps de la charge):

$$div_{\Gamma} \vec{j} + \varepsilon \frac{\partial \rho}{\partial t} = 0$$

Notons (formellement) $E^s$ et $E_\omega^s$ respectivement les espaces munis des normes usuelles et les espaces couplant le temps et l'espace. Dans les espaces $E_\omega^s$, une dérivation en espace s'accompagne d'une perte de régularité en temps. Cela implique que si on cherche $\vec{j}$ dans $E_\omega^s$, la relation de conservation étant vérifiée pour le problème continu, $\rho$ doit appartenir au même type d'espace (car $div_{\Gamma} \vec{j}$ est alors dans $E_\omega^{s-1}$, donc $\frac{\partial \rho}{\partial t}$ est dans $E_\omega^{s-1}$). Supposons maintenant qu'on approche $E_\omega^s$ par des fonctions $P^1$ en temps: on choisit $\vec{j}_{h,\Delta t} \in P^1$ en temps et $\rho_{h,\Delta t} \in P^1$ en temps. Or les espaces discrets découplent l'espace et le temps ce qui signifie que $div_{\Gamma} \vec{j}_{h,\Delta t} \in P^1$ et $\frac{\partial \rho_{h,\Delta t}}{\partial t} \in P^0$. La relation de conservation ne peut donc pas être satisfaite numériquement de manière exacte. Pour quand même prendre cette relation en compte de manière numérique, I. Terrasse impose à la charge de vérifier la relation de conservation seulement en moyenne (et non pas de manière exacte) ceci à l'aide d'un multiplicateur de Lagrange.

Cette formulation a été implémentée, en 3D. Une des difficultés numériques essentielles vient alors de la géométrie. Les intégrations portent en effet sur des domaines qui représentent l'intersection des fronts d'onde avec les éléments du maillage de la surface. Il faut donc considérer tous les cas d'intersection d'une couronne avec un élément et calculer les intégrales sur ces domaines. Cependant les calculs en 3D ont deux avantages: ils permettent de ne pas prendre en compte tout le passé de l'onde, comme on l'a vu pour les ondes acoustiques et les noyaux sont plus simples qu'en 2D et conduisent donc à calculer des intégrales très simples.

M. Filipe-Rocha, [24], [25], s'est intéressée à deux problèmes:

- Un problème de couplage fluide-structure: il s'agit d'un problème de diffraction d'une onde acoustique par un obstacle élastique. Elle a proposé une formulation couplant la méthode des éléments finis dans la partie structure avec la méthode des potentiels retardés dans le domaine extérieur fluide. L'analyse mathématique de cette formulation couplée a été faite (coercivité, existence et unicité ...).

- un problème de diffraction d'une onde acoustique par un obstacle absorbant, c'est à dire avec une condition aux limites absorbante (de type $\frac{\partial u}{\partial n} + \gamma \frac{\partial u}{\partial t} = 0$). Elle a montré qu'on ne dispose plus de coercivité, mais que l'opérateur est de type Fredholm ce qui permet en utilisant la théorie de Fredholm de montrer quand-même l'existence et l'unicité de la solution.

La programmation de ce second problème est actuellement en cours.

V. Lubet, [3], étudie le couplage entre la méthode des potentiels retardés et celle
des éléments finis pour la diffraction d'une onde acoustique par un obstacle ou des hétérogénéités. La condition aux limites sur la frontière de l'obstacle est une condition de transmission parfaite. Elle a fait l'analyse mathématique de cette formulation et a obtenu des résultats d'existence et d'unicité de la solution (encore une fois sans coercivité). Elle travaille actuellement sur l'implémentation de cette formulation.

V. Lange a implémenté la formulation proposée par A. Pujols pour résoudre la diffraction d'une onde électromagnétique par un obstacle parfaitement conducteur dans un domaine 3D. Elle a également obtenu une analyse de ce schéma en utilisant les espaces fonctionnels introduits par I. Terrasse. Elle travaille actuellement sur la diffraction par un obstacle sur lequel est posée une condition d'impédance.

D. Barbier travaille actuellement sur la diffraction d'une onde élastique par un obstacle ou une fissure. La programmation doit être effectuée pour le problème 3D. Il doit également étudier le problème du couplage entre la méthode des potentiels retardés et celle des éléments finis dans le cas où une partie du milieu est anisotrope ou hétérogène.

Dans tous les cas, en pratique, et pour des raisons de simplicité, les formulations ont été implémentées avec $\omega_I = 0$, valeur pour laquelle on ne dispose plus d'aucun résultat... Cependant, le paramètre $\omega_I$ permet en fait de contrôler le comportement de la solution pour des temps très grands et on s'attend à ce qu'il ait une faible influence sur des temps finis. Dans pratiquement tous les cas, des problèmes de stabilité ont été observés. Les formulations utilisant les éléments de plus bas degré en temps semblent donner des résultats plus stables. Pour ces éléments, on n'est plus dans le cadre d'application de l'analyse numérique théorique et l'étude de ces phénomènes reste à faire.

# 7 Bibliographie

Dans cette partie, nous citons une liste (non exhaustive) de travaux relatifs aux domaines abordés dans ce cours.

## • Analyse mathématique des équations intégrales pour les problèmes stationnaires et harmoniques

Colton et Kress [15], Dautray-Lions [19] (chap XI, Partie B), Nédélec [45], Giroire (Laplacien, bilaplacien, Helmholtz), [27], [26], [28], Costabel, Stephan, Wendland, [17], [18], Hsiao, Mac Camy, [34], [33], [31], [32], Kleinman, [37]

## • Application des équations intégrales pour les problèmes stationnaires
## et harmoniques

* méthodes variationnelles:

Bonnemay (élasticité stationnaire), Cortey-Dumont [16] (élastodynamique harmonique), Djaoua [21], Bendali(Maxwell harmonique) [10], Mazari [43](Maxwell harmonique: antennes), Starling (Maxwell harmonique: réseaux) [46], [53], Morelot (Maxwell harmonique: réseaux) [?], ...

* méthodes de collocation:

Bonnet [11] (élasticité statique et dynamique harmonique), Bui et al [12], Nishimura [49], [47] (élasticité statique et dynamique harmonique), Martin, Rizzo, Whickham [41], [42] (ondes élastiques harmoniques),...

* Analyse mathématique des équations intégrales pour les problèmes transitoires

Bamberger-Ha Duong [22], [5], [4], [23], Lubich [40]

* Applications des potentiels retardés ; résultats numériques

* méthodes variationnelles:

Y. Ding [20], E. Bécache [6](élastodynamique), [7], [8], I. Terrasse [54] (Maxwell), A. Pujols [50] (Maxwell), M. Rocha-Filipe [24], [25] (ondes acoustiques et élastiques), V. Lubet [3] (acoustique), V. Lange (ondes électromagnétiques), D. Barbier (ondes élastiques)

* méthodes de collocation:

Antes [1], Cole et al [14], Hirose et Achenbach [30], Mitzner [44], Nishimura [48], Shaw [52]...
Références

[1] H. Antes. A Boundary Element Procedure for Transient Wave Propagation in Two0Dimensional Isotropic Elastic Media. *Fin. El. in Anal. and Des.*, 1:313–322, 1985.

[2] D. Aubry and D. Clouteau. A regularised Boundary Element Method for Stratified Media. In G. Cohen, L. Halpern, and P. Joly, editors, *Mathematical and Numerical Aspects of Wave Propagation Phenomena*, pages 660–668. SIAM, 1991.

[3] A. Bachelot and V. Lubet. Couplage Equation Intégrale-Eléments Finis dans le Domaine Temporel. In *Proc. 25ème Congés National d'Analyse Numérique*, page 165, 1993.

[4] A. Bamberger and T. Ha Duong. Formulation variationnelle espace-temps pour le calcul par potentiel retardé d'une onde acoustique. *Math. Methods Appl. Sci.*, 8:405-435, 1986.

[5] A. Bamberger and T. Ha Duong. Formulation variationnelle espace-temps pour le calcul par potentiel retardé d'une onde acoustique; Problème de Neumann. *Math. Methods Appl. Sci.*, 8:598-608, 1986.

[6] E. Bécache. *Résolution par une méthode d'équations intégrales d'un problème de diffraction d'ondes élastiques transitoires par une fissure*. PhD thesis, Université de Paris 6, 1991. Thèse.

[7] E. Bécache. A Variational Boundary Integral Equation Method for an Elastodynamic Antiplane Crack. *Int. J. for Numerical Meth. in Eng.*, 36:969–984, 1993.

[8] E. Bécache and T. Ha Duong. A Space-Time Variational Formulation for the Boundary Integral Equation in a 2D Elastic Crack Problem. *RAIRO, M2AN*, 28, 94, to appear.

[9] E. Bécache, J. C. Nédélec, and N. Nishimura. Regularization in 3D for Anisotropic Elastodynamic Crack and Obstacle Problems. *J. of Elasticity*, 31:25–46, 1993.

[10] A. Bendali. *Approximation par éléments finis de surface de problèmes de diffrraction des ondes électromagnétiques*. PhD thesis, Paris VI, 1984. Thèse de doctorat d'état.

[11] M. Bonnet. *méthode des équations intégrales régularisées en élastodynamique*. PhD thesis, ENPC, 1986. Thèse.

[12] H. D. Bui, B. Loret, and M. Bonnet. Régularisation des équations intégrales de lélastodynamique et de l'élastostatique. *C.R. Acad. Sc.*, 300, 1985.

[13] D. Clouteau. *Propagation d'ondes en milieux hétérogènes, application à la tenue des ouvrages sous séismes*. PhD thesis, Ecole Centrale Paris, 1986. Thèse.

[14] D. M. Cole, D. D. Kosloff, and J. B. Minster. A numerical boundary integral equation method for elastodynamics I. *Bull. Seism. Soc. Amer.*, 68:1331–1357, 1978.
[15] D. Colton and R. Kress. *Integral Equation Methods in Scattering Theory*. John Wiley and Sons, 1983.

[16] P. Cortey-Dumont. *Simulation Numérique de Problèmes de Diffraction d'Ondes par une Fissure*. PhD thesis, Université Paris VI, 1984. Thèse d'Etat.

[17] M. Costabel and E. Stephan. A direct boundary integral equation method for transmission problems. *J. Math. Anal. Appl.*, 106:367-413, 1985.

[18] M. Costabel and W. L. Wendland. Strong Ellipticity of boundary integral operators. *J. fur die reine und angewandte Mathematik*, 372:34-63, 1986.

[19] R. Dautray and J. L. Lions. *Analyse Mathématique et Calcul Numérique pour les Sciences et les Techniques*, volume 2. Masson, 1985.

[20] Yu Ding. *Méthodes numériques pour l'équation intégrale de bord dans le problème de diffraction d'ondes acoustiques par une surface rigide*. PhD thesis, Université d'Orsay, 1989. Thèse.

[21] M. Djaoua. *Equations intégrales pour un problème singulier dans le plan*. PhD thesis, Université Paris VI, 1977. Thèse de 3ème cycle.

[22] T. Ha Duong. *Equations intégrales pour la résolution numérique de problèmes de diffraction d'ondes acoustiques dans $\mathbb{R}^3$*. PhD thesis, Université de Paris 6, 1987. Thèse de doctorat d'état.

[23] T. Ha Duong. On the transient acoustic scattering by a flat object. *Japan J. Appl. Math.*, 7:489-513, 1990.

[24] M. Filipe-Rocha. PhD thesis, à paraître. Thèse.

[25] A. Forestier, T. Ha Duong, and M. Rocha. Résolution d'un problème instationnaire d'interaction fluide-structure par la méthode des potentiels retardés. In *Proc. 25ème Congés National d'Analyse Numérique*, page 209, 1993.

[26] J. Giroire. *Etude de quelques problèmes aux limites extérieurs et résolution par équations intégrales*. PhD thesis, Université de Paris 6, 1987. Thèse de doctorat d'état.

[27] J. Giroire. Mise en oeuvre de méthodes d'éléments finis de frontière. Technical report, Université Paris 6, 1987. Cours de DEA d'Analyse Numérique.

[28] J. Giroire and J. C. Nédélec. Numerical Solution of an Exterior Neumann Problem using a Double Layer Potential. *Math. Comp.*, 32:973-990, 1978.

[29] M. A. Hamdi. *Formulation Variationnelle par Equations Intégrales pour le Calcul de Champs Acoustiques Linéaires Proches et Lointains*. PhD thesis, Université de Compiègne, 1982. Thèse d'Etat.

[30] S. Hirose and J. D. Achenbach. BEM Method to analyze the interaction of an acoustic pulse with a rigid circular disk. *Wave motion*, 10:267-275, 1988.

[31] G. C. Hsiao. The coupling of BEM and FEM - A brief review. In *Brebbia, editor, Proc. Boundary Elements X*, volume 1, pages 431-445. Springer-Verlag, 1988.
[32] G. C. Hsiao. On Boundary Integral Equations of the First Kind. *J. of Comp. Math.*, 7, n° 2:121-131, 1989.

[33] G. C. Hsiao and R. C. Mac Camy. Solution of boundary value problems by integral equations of the first kind. *SIAM Review*, 15:687-705, 1973.

[34] G. C. Hsiao and W. L. Wendland. On a boundary integral method for some exterior problems in elasticity. In *Proc. Tbilisi University 257 Ser. Mat. Mech. Astron.*, volume 18, pages 31-60, 1985.

[35] C. Johnson and J. C. Nédélec. On the coupling of Boundary Integral and Finite Element Methods. *Math. of Comp.*, 35:1063-1079, 1980.

[36] T. KATO. *Perturbation theory of linear operators*. Springer Verlag, 1980.

[37] R. E. Kleinman and G. F. Roach. Boundary Integral Equations for the 3-dimensional Helmholtz Equation. *SIAM Review*, 16:214-236, 1974.

[38] G. Krishnasamy, F. J. Rizzo, and T. J. Rudolphi. Hypersingular boundary integral equations: Their occurence interpretation, regularization and computation. In P. K. Banerjee and S. Kobayashi, editors, *Developments in Boundary Element Methods, Vol. 7:Advanced Dynamic Analysis*. Elsevier Applied Science Publishers, 1991.

[39] V. B. Kupradze. *Dynamical problems in elasticity*. Progress in Solid Mechanics, vol III.

[40] Ch. Lubich. On multistep time discretization of linear initial-boundary value problems and their boundary integral equations. Preprint n° 192, submitted to numerische mathematik, Mathematische Institute der Julius-Maximilians-Universitat, Würzburg, 1993.

[41] P. A. Martin and F. J. Rizzo. On boundary integral equations for crack problems. *Proc. Roy. Soc. London A*, 421:341-355, 1989.

[42] P. A. Martin and G. R. Whickham. Diffraction of elastic waves by a penny-shaped crack: analytical and numerical results. *Proc. Roy. Soc. London A*, 390:91-129, 1983.

[43] A. Mazari. *Détermination par une méthode d'équations intégrales du champ électromagnétique rayonné par une structure filiforme*. PhD thesis, Université de Paris 6, 1991. Thèse.

[44] K. M. Mitzner. Numerical solution for transient scattering from a hard surface of arbitrary shape- retarded potential technique. *J. Acoust. Soc. Amer.*, 42:391-397, 1967.

[45] J. C. Nédélec. Approximation des Equations Intégrales en Mécanique et en Physique. Technical report, Ecole Polytechnique, CMAP, 1977. Cours de l'École d'Eté CEA-EDF-INRIA.

[46] J. C. Nédélec and F. Starling. Integral equation methods in a quasi-periodic diffraction problem for the time harmonic Maxwell's equations. Technical report, Ecole Polytechnique, CMAP, 1988. Rapport Interne n° 179.
[47] N. Nishimura. *Applications of the Boundary Integral Equation Methods to Solid Mechanics*. PhD thesis, 1988. Thesis.

[48] N. Nishimura, Q. C. Guo, and S. Kobayashi. Boundary Integral Equation Methods in Elastodynamic Crack Problems. In Brebbia, Wendland, and Kuhn, editors, *Proc. 9th Int. Conf. BEM*, volume 2: Stress Analysis Applications, pages 279-291. Springer-Verlag, 1987.

[49] N. Nishimura and S. Kobayashi. A regularized boundary integral equation method for elastodynamic crack problems. *Computat. Mech.*, 4:319-328, 1989.

[50] A. Pujols. *Equations intégrales Espace-Temps pour le système de Maxwell - Application au calcul de la Surface Equivalente Radar*. PhD thesis, Université de Bordeaux I, 1991. Thèse.

[51] F. J. Rizzo, D. J. Shippy, and M. Rezayat. A boundary integral equation method for radiation and scattering of elastic waves in three dimensions. *Int. J. Num. Meth. Eng.*, 21:115-129, 1985.

[52] R. P. Shaw. Integral equation formulation of dynamic acoustic fluid-elastic solid interaction problems. *J. Acoust. Soc. Am.*, 53:514-520, 1973.

[53] F. Starling. *Etude mathématique de quelques problèmes de diffraction en électromagnétisme*. PhD thesis, Ecole Polytechnique, 1991. Thèse.

[54] I. Terrasse. *Résolution mathématique et numérique des équations de Maxwell instationnaires par une méthode de potentiels retardés*. PhD thesis, Ecole Polytechnique, 1993. Thèse.

[55] Treves. *Basic Linear Partial Differential Equations*. Academic Press, 1975.