Projet Ondes

previous up next contents
Précédent : Logiciels Remonter : Projet ONDES, Modélisation et Simulation Suivant : Actions industrielles

Résultats nouveaux

Éléments finis triangulaires non-conformes d'ordre élevé pour l'équation des ondes



Participants : G. Cochard , G. Cohen , P. Monk


Mots-clés : Éléments finis non-conformes, équation des ondes, condensation de masse


Résumé : Recherche d'éléments finis triangulaires non-conformes d'ordre élevé avec condensation de masse pour résoudre l'équation des ondes en régime transitoire en dimension 2.


Les éléments finis conformes triangulaires avec condensation de masse pour l'équation des ondes obtenus dans le cadre de la thèse de Nathalie Tordjman [18] présentent l'inconvénient de mal se prêter à l'utilisation de conditions aux limites absorbantes. En dimension 2, cela est dû au fait que les noeuds apparaissant sur les arêtes ne sont pas ceux de formules de quadrature 1D. Pour résoudre ce problème, nous avons essayé de construire de nouveaux éléments finis triangulaires dont les points des arêtes sont ceux de formules de quadrature de Gauss. Deux éléments ont été exhibés, l'un contenant 10 degrés de liberté qui donne une dispersion d'ordre 4 et l'autre contenant 12 degrés de liberté qui donne une dispersion d'ordre 6. Une étude numérique des méthodes obtenues à partir de ces éléments, menée par Gary Cochard, dans le cadre de son stage de deuxième année du CUST, a permis de montrer que, sur des maillages réguliers, les premiers éléments étaient d'ordre 3 et les seconds d'ordre 4.

Éléments finis mixtes quadrilatéraux ou héxaédriques pour les équations de Maxwell en milieu anisotrope



Participants : G. Cohen , P. Monk


Mots-clés : Éléments finis mixtes, équations de Maxwell, condensation de masse, milieux anisotropes


Résumé : Construction d'une méthode d'éléments finis précise et rapide pour résoudre les équations de Maxwell en régime transitoire pour des géométries complexes et en milieu anisotrope.


Résoudre à peu de frais les équations de Maxwell pour des géométries complexes en milieu anisotrope était un défi difficile à relever. Après plusieurs tentatives infructueuses, nous avons abouti cette année à une méthode qui réunit d'excellentes qualités de précision, de rapidité et de stockage. Les ingrédients de base de cette méthode ne sont pas nouveaux. Ce sont des éléments finis mixtes, ou d'arête dans le langage industriel, qui constituent la ``nouvelle famille d'éléments finis d'arête de Nédélec'', introduite aussi par G. Mur. Ces éléments ont leurs degrés de liberté aux points des degrés de liberté des éléments finis de Lagrange, ce qui permet, pour obtenir une condensation de la matrice de masse du système, d'utiliser des formules de Gauss-Lobatto, comme l'avait fait Nathalie Tordjman dans sa thèse dans les années précédentes. Le caractère vectoriel des équations conduit ici, pour des maillages non orthogonaux et des milieux anisotropes, à une matrice de masse diagonale par blocs de petites tailles. Après condensation de la masse, le coût principal de l'algorithme tient dans le calcul des termes de rigidité qui présentent des molécules volumineuses dans le cas de géométries complexes. La solution à ce problème repose notamment sur l'introduction d'un nouvel espace d'approximation pour $L^2$ qui coïncide localement avec les espaces habituellement utilisés pour $H(rot)$.L'étude de dispersion numérique des éléments finis de types $Q_{1}$ à $Q_{3}$sur des maillages déformés périodiques a montré la robustesse de la méthode pour des géométries complexes. La programmation d'une méthode $Q_{3}$ et d'ordre 2 en temps en dimension 2 a donné d'excellents résultats tant du point de vue de la précision que de celui du stockage et du temps calcul. À précision égale, la méthode est aussi rapide qu'un schéma d'ordre 2 en temps et en espace (schéma de Yee) pour une place mémoire six fois moins importante. Un résultat numérique est donné à la figure [*]. Cette étude fera l'objet d'une publication dans [20].


   Figure: Diffraction d'une onde magnétique par une ogive recouverte d'un matériau anisotrope.

\begin{figure} \begin{center} \leavevmode {\includegraphics[width=.46\linewidth]{grandban104.ps}} \end{center}\end{figure}


Eléments finis d'arête et condensation de masse pour les équations de Maxwell en maillage triangulaire et tétraédriques



Participants : A. Elmkies , P. Joly , D. Martin , P. Monk


Mots-clés : Eléments finis d'arête, équations de Maxwell, condensation de masse, milieux anisotropes


Résumé : Construction et implémentation de nouvelles familles d'éléments finis d'arête d'ordre élevé avec condensation de masse pour les équations de Maxwell.


Les éléments finis d'arête introduits par Nédélec, s'ils s'imposent du point de vue mathématique et physique comme un outil naturel pour la résolution des équations de Maxwell [19], introduisent après discrétisation en espace une matrice dite ``de masse'' qui doit être inversée à chaque pas de temps. Pour garantir l'efficacité de la méthode, notre objectif est de construire de nouvelles familles d'éléments finis d'arêtes permettant de ``condenser'' la matrice de masse, c'est-à-dire l'approcher par une matrice diagonale sans perturber la précision de la méthode. Cette année a été consacrée à la généralisation des éléments finis triangulaires conçus l'année dernière à la dimension 3 (maillages tétraédriques). L'idée de base repose sur le théorème de trace dans $H(\mbox{\bf rot})$ et consiste à construire de nouveaux espaces de discrétisation dont la trace tangentielle, sur chaque face, décrit les espaces 2D construits l'année dernière. En réincorporant certaines composantes tangentielles sur chaque faces, on obtient, par intégration au moyen d'une formule de quadrature adaptée, une matrice de masse diagonale par blocs dont la taille est donnée par le nombre de faces du maillage concourant en chaque arête, ceci même en milieu anisotrope. Ces nouveaux éléments ont été analysés à l'aide d'une étude de dispersion numérique ([29], [30]). D'autre part, en s'inspirant de la méthode développée par G. Cohen et P. Monk pour les quadrilatères, le nouvel élément fini d'arête triangulaire 2D du second ordre a été implémenté et est actuellement en cours d'étude. Mentionnons le fait que dans une formulation des équations de Maxwell comme un système différentiel du premier ordre (en $(E,H)$), un choix adéquat de l'espace de discrétisation du champ magnétique permet de diminuer le coût calcul de 20${\%}$. Ce code doit être introduit dans le code MELINA développé à l'ENSTA. Enfin, en collaboration avec P.Monk, une analyse d'erreur plus théorique est en cours.

Equations intégrales de précision élevée



Participants : Marc Lenoir , Eric Masson


Mots-clés : diffraction, équations intégrales


La mise en oeuvre des méthodes numériques pour la résolution des problèmes de diffraction se heurte toujours aux difficultés liées à la singularité du noyau. A cet égard le choix de méthodes variationnelles constitue un progrès notable ; au prix d'une intégration par parties le long du bord, elles permettent en effet d'exprimer la forme bilinéaire associée à la résolution du problème de Neumann à l'aide du noyau du potentiel de simple couche. Dans le cas où on souhaite s'affranchir des fréquences irrégulières, l'utilisation d'un potentiel mixte amène des intégrales faisant intervenir un potentiel de double couche.

L'utilisation d'un code d'éléments finis relativement général tel que MELINA permet sans difficulté de réaliser des discrétisations d'ordre élevé, tant pour les inconnues que pour la géométrie. Il est alors possible de respecter l'estimation théorique, qui invite à conserver un ordre de précision supplémentaire pour la discrétisation de la géométrie par rapport à celle de l'inconnue. Il reste qu'alors, le calcul des intégrales singulières mérite d'être réalisé avec soin pour ne pas perdre l'avantage apporté par une méthode de discrétisation d'ordre élevé.

Le choix que nous avons effectué consiste à donner de la surface une représentation par des éléments finis hyperparamétriques, et à développer l'ensemble de l'intégrand (fonction de Green, fonction de base, carte locale et élément de surface) au voisinage du plan tangent au point courant dans l'intégrale externe ; compte tenu de notre ambition modeste consistant à considérer comme régulière une fonction continue, ce calcul reste d'une complexité modérée, même dans le cas d'un potentiel de double couche.

Remarquons bien que ce calcul ne dépend pas du problème résolu, il est en fait le même pour tous les opérateurs du second ordre elliptiques, le choix de la fonction de Green précise à utiliser n'intervient que dans la phase d'intégration numérique de la partie régulière, ce qui confère une grande souplesse à la méthode.

L'opération précédente permet de se ramener au calcul de l'intégrale d'une ou plusieurs fonctions homogènes sur un polygone dans le plan tangent (triangle ou quadrilatère à bords droits selon le type d'éléments finis choisi). On utilise alors une formule spécifique pour l'intégration des fonctions homogènes (dite formule d'Euler) qui permet de se ramener à deux quadratures successives, l'une le long de la normale à l'élément et l'autre le long de son bord ; il s'agit en fait d'une généralisation de l'équivalence bien connue entre répartition de sources surfaciques et de tourbillons linéiques.

On aboutit ainsi à des formules explicites ne faisant intervenir que des fonctions connues calculées aux sommets du triangle projeté sur le plan tangent, le calcul de la partie régulière s'opérant à l'aide d'une quadrature numérique classique.

La méthode décrite ci-dessus a été mise en oeuvre en bidimensionnel, où le bord du domaine est constitué d'une ligne courbe, les résultats obtenus corroborent parfaitement les estimations d'erreur théoriques obtenues par Nédélec, et permettent d'obtenir une précision jusqu'ici inégalée avec un faible nombre d'éléments.

Dans le cas tridimensionnel, le travail analytique a été réalisé et l'ensemble des formules obtenu, dans une situation d'ailleurs plus générale dont les équations intégrales ne forment qu'un cas particulier ; les modifications générales nécessaires ont été apportées au code MELINA, restent à réaliser l'implémentation et les tests de convergence.

Singularités et équations de Maxwell



Participants : Patrick Ciarlet , Christophe Hazard , Stephanie Lohrengel


Mots-clés : équations de Maxwell, singularités


Nous développons depuis quelques années une méthode de régularisation des équations de Maxwell destinée à simuler des phénomènes de propagation d'ondes électromagnétiques à l'aide d'une discrétisation classique (par éléments finis nodaux, et non par éléments finis d'arête). L'idée de base de cette méthode est assez ancienne et revient formellement à remplacer les équations de Maxwell par les équations de l'élasticité: aux ondes électromagnétiques qui s'apparentent aux ondes dites de cisaillement (à divergence nulle), on ajoute ainsi les ondes de pression (qui dérivent d'un potentiel). Par contre, les conditions aux limites restent celles de l'électromagnétisme, auxquelles il faut adjoindre une condition dont le rôle est précisément de supprimer les ondes de pression. En résumé, on autorise dans le milieu les ondes de pression mais les conditions aux limites les interdisent, de sorte que la seule solution possible est bien l'onde électromagnétique.

Ces équations régularisées, bien que mathématiquement équivalentes aux équations de Maxwell, ne peuvent être directement discrétisées si le champ électromagnétique présente un comportement singulier, comme par exemple au voisinage des coins et arêtes d'un conducteur parfait polyédrique. La difficulté provient du fait que les éléments finis habituels ne permettent pas d'approcher un tel comportement. L'idée que nous avons cherché à développer consiste simplement à décomposer le champ électromagnétique en deux parties: une partie régulière, dont l'approximation numérique peut être obtenue par une discrétisation standard des équations de Maxwell régularisées, et une partie singulière qui sera traitée de façon explicite ; il faut pour cela connaître a priori la nature du comportement singulier du champ électromagnétique au voisinage des arêtes et des coins. Nous avons montré comment ces singularités peuvent se déduire des singularités scalaires liées à l'opérateur Laplacien dans le même domaine géométrique.

Ces travaux nous ont conduit à développer une méthode de champs singuliers dont la mise en oeuvre numérique (actuellement implémentée dans le cas bidimensionnel sur le code MELINA) revient à résoudre un problème couplé entre les parties régulière et singulière du champ électromagnétique.

Modélisation numérique de la propagation d'ondes en milieu élastique fissuré



Participants : E. Bécache , F. Collino , P. Joly , C. Tsogka


Mots-clés : éléments finis mixtes, élastodynamique, fissure


 Nous avons poursuivi nos travaux sur la diffraction d'une onde élastique par une fissure. Rappelons que pour prendre en compte la géométrie de la fissure, tout en utilisant un maillage structuré, nous avons choisi d'appliquer la méthode des domaines fictifs. Ceci nous a conduit à considérer la formulation vitesse-contraintes des équations de l'élastodynamique, cette formulation étant par ailleurs bien adaptée à l'utilisation de couches absorbantes Perfectly Matched Layers (voir des exemples dans http://www-rocq.inria.fr/~collino). Nous avons construit un nouvel élément fini mixte permettant de faire la condensation de masse et donc d'obtenir des schémas en temps explicites [6]. Une des caractéristiques de cet élément est qu'il permet de prendre en compte la symétrie du tenseur des contraintes de façon forte. Une fois effectuée la semi-discrétisation en espace, on s'aperçoit qu'on peut éliminer les inconnues relatives aux contraintes et se ramener à une équation différentielle d'ordre 2 en temps ne portant que sur les vitesses, de même type que les schémas classiques obtenus directement sur la formulation en vitesses (différences finies ou éléments finis Q1...). Deux types d'analyse ont été menées.

La première analyse est une étude de la dispersion et de la stabilité du schéma dans un milieu homogène isotrope. Cette étude a été menée pour l'élément d'ordre le plus bas $Q_1-Q_0$ ainsi que, dans le cadre d'un stage effectué par F. Imedni, pour l'élément $Q_2-Q_1$. Dans les deux cas, nous avons observé un phénomène de superconvergence, dû à l'utilisation de maillages réguliers, c'est à dire une dispersion d'ordre deux pour le schéma $Q_1-Q_0$ et d'ordre quatre pour le schéma $Q_2-Q_1$. On a pu déterminer la condition de stabilité par une étude théorique pour le schéma $Q_1-Q_0$ et par une étude numérique pour le schéma $Q_2-Q_1$.

La deuxième analyse est une analyse type éléments finis sur l'erreur d'approximation. Celle ci revêt essentiellement deux difficultés. La première est due au fait que les hypothèses classiques (voir Brezzi-Fortin) pour obtenir des estimations d'erreur ne sont pas satisfaites et ceci parce que, par rapport aux approximations classiques (Raviart-Thomas), on a enrichi l'espace d'approximation des contraintes sans enrichir celui des vitesses. La deuxième difficulté est liée à la prise en compte de la symétrie du tenseur des contraintes.

Pour traiter la première difficulté, nous avons envisagé un modèle simplifié d'ondes anisotropes scalaire. Rappelons que l'obtention d'estimations d'erreur classiques pour les problèmes mixtes, qui utilisent deux espaces fonctionnels, repose essentiellement sur deux points : un résultat de coercivité et une condition inf-sup. Nos deux espaces d'approximation ne vérifient pas ces deux hypothèses. Néanmoins, nous avons montré que moyennant l'introduction d'un troisième espace, une hypothèse de coercivité ``affaiblie'' et une hypothèse de condition inf-sup ``renforcée'', il est encore possible d'obtenir des estimations d'erreur (non classiques). [7].

Le deuxième point concernant la prise en compte de la symétrie et l'analyse d'erreur sur le problème initial d'élastodynamique est à l'étude actuellement.

Nous avons implémenté le schéma $Q_1-Q_0$ dans des milieux homogènes anisotropes comportant des fissures. La figure [*] représente un exemple de propagation dans un demi-espace élastique isotrope. Notons que l'utilisation des couches PML permet une bonne absorption non seulement des ondes S et P, mais aussi des ondes de Rayleigh. Sur la figure [*], on présente un exemple de diffraction par une fissure dans un milieu élastique anisotrope. Nous avons également réalisé un film présentant la diffraction d'une onde élastique par une fissure dans un milieu isotrope .
   Figure: Propagation dans un demi-espace élastique

\begin{figure} \begin{center} \leavevmode {\includegraphics[width=.46\linewidth]{rayleigh.ps}} \end{center}\end{figure}



   Figure: Diffraction par une fissure en milieu anisotrope

\begin{figure} \begin{center} \leavevmode {\includegraphics[width=.46\linewidth]{enr120.ps}} \end{center}\end{figure}


Raffinement de maillage pour une méthode en différences finies de résolution de l'équation des ondes



Participants : F. Collino , T. Fouquet , P. Joly


Mots-clés : raffinement de maillage, FDTD, stabilité, précision


Résumé : Étude théorique et numérique du raffinement de maillage en différences finies pour la résolution de l'équation des ondes.


Afin de tenir compte de géométries complexes dans des méthodes de différences finies, il est intéressant de pouvoir raffiner une partie de la grille de calcul en espace. Pour des raisons de condition de stabilité et de coût de calcul, il est souhaitable de raffiner également le pas de temps. Le point délicat est le traitement de la frontière entre les grilles. Celui-ci est en effet susceptible de détruire la stabilité et d'altérer la précision notamment par des phénomènes de réflexion à l'interface.

Nous avons commencé l'étude de raffinement de maillages pour l'équation des ondes en une dimension par l'analyse d'un schéma d'interpolation temporelle pour le raccordement des grilles. Nous avons choisi de discrétiser l'équation de propagation des ondes à l'aide d'un schéma saute-mouton qui est d'ordre 2 en temps et en espace. Notre objectif est d'obtenir une analyse théorique complète de notre schéma en terme de stabilité et de précision puis d'implémenter un code de calcul.

En adaptant la technique de Kreiss à notre schéma, nous avons démontré la stabilité dans le cas d'une seule interface. L'étude de la précision nous montre que malgré un schéma dans chaque grille d'ordre 2 et une interpolation temporelle à l'interface d'ordre 2, notre schéma n'est que d'ordre 1.

Les calculs théoriques de coefficients de réflexion et de transmission à l'interface des deux grilles sont corroborés par les expériences numériques et montrent que la réflexion est imperceptible sur la plages de fréquences couramment utilisée.

Dans le cas de plusieurs interfaces, par exemple un maillage composé d'une grille grossière d'une grille fine puis à nouveau d'une grille grossière, des expériences numériques font apparaître des phénomènes de type instabilité qui ne sont pour le moment pas très bien compris. Une analyse de ces phénomènes est actuellement à l'étude.

Par la suite, nous tenterons d'étendre notre étude aux équations de Maxwell à une dimension puis à deux dimensions.

Décomposition de Domaine et Contrôle



Participants : J.-D. Benamou


La technique classique pour la résolution de problèmes de contrôle optimal de systèmes gouvernés par des équations linéaires d'évolution repose sur la reformulation en un système d'optimalité comprenant une équation pour l'état direct dans laquelle intervient le contrôle, une équation pour l'état adjoint couplée à la précédente par l'intermédiaire du critère à minimiser et les conditions d'optimalités faisant intervenir le gradient du critère.

Une des difficultés principales pour la résolution numérique de tels problèmes vient du caractère rétrograde en temps de l'équation adjointe : nous sommes en présence d'un problème ``au deux bouts''. Dans le cas ou le critère est quadratique, on peut utiliser soit un algorithme de descente soit la technique du découplage, aussi appelée synthèse quadratique. Dans le premier cas, on est amené à calculer une suite de points dans l'espace des paramètres de contrôle avec pour chacun de ces points la résolution de l'état direct puis de l'état adjoint pour calculer le gradient du critère qui donne une direction de descente. Dans le deuxième cas, on remarque que l'état adjoint peut pour tout temps s'exprimer comme une loi affine de l'état direct. On peut alors dériver une équation de type Riccati non linéaire pour les coefficients de cette loi, dite de ``feed-back''. Il ne reste plus qu'à résoudre l'équation de Riccati et à ré-injecter la loi de ``feed back'' dans l'équation de l'état direct pour résoudre le problème.

La technique de synthèse quadratique peut poser des problèmes théoriques et numériques, en particulier pour les problèmes hyperboliques. Du point de vue numérique, on est ramené à la résolution d'un système d'équations non-linéaires de grande taille. Cette résolution peut être difficile ou coûteuse.

La méthode de décomposition de domaine mise au point consiste après décomposition du domaine à ramener la résolution du problème ``global'' à une suite de résolutions de sous-problèmes locaux posés chacun sur les sous-domaines de cette décomposition.

L'utilisation d'une technique de ce type se justifie, soit au niveau informatique pour l'utilisation d'ordinateurs parallèles afin de résoudre des problèmes de grande taille autrement inaccessibles, soit au niveau de la modélisation lorsque la décomposition du domaine permet de réduire le caractère hétérogène de la géométrie ou de la nature du milieu. Elle permet également d'utiliser des lois de feed back locales indépendantes du processus itératif. Ceci ramène le coût de ce type de méthodes dans les ordres de grandeurs donnés pour les méthodes de descente. Un problème de contrôle pour l'équation des ondes a ainsi été résolu numériquement.

Notons par ailleurs que la méthode de décomposition de domaine a par ailleurs été étendue au contrôle de réseaux de poutres élastiques par G. Leugeuring (Université de Bayreuth). Finalement l'extension de la même méthode à la contrôllabilité exacte de l'équation des ondes, via la méthode HUM, est l'objet de recherche communes avec F. Bourquin (LCPC).

Conditions aux limites absorbantes



Participants : E. Bécache , G. Cohen , F. Collino , P. Joly , P. Monk , C. Tsogka


Ce thème est transversal aux divers domaines d'application que nous traitons. Il s'agit d'abord de faire progresser les conditions absorbantes en même temps que les méthodes numériques et de les adapter aux nouveaux modèles que nous abordons. Nous devons également intégrer et adopter les progrès récents dans ce domaine. C'est l'objet de nos travaux sur les couches absorbantes. Le travail sur les équations paraxiales est maintenant achevé, [26]. Cette année a été consacrée à l'étude et la conception de couches absorbantes pour l'électromagnétisme en coordonnées cylindriques ainsi que pour l'élastodynamique en milieux isotrope et anisotrope. Dans les deux cas il s'agit de prolonger la technique des couches parfaitement absorbantes de Bérenger (PML) initialement conçue pour le système de Maxwell en coordonnées cartésiennes. Ces extensions sont basées sur une ré-interprétation du modèle de couches à l'aide d'un changement de variable complexe appliqué à la variable latérale, [24]. Une étude par ondes planes permet de souligner l'importance de la discrétisation dans le choix optimal des paramètres du modèle, [25]. Leur optimisation permet d'obtenir des résultats spectaculaires avec un sur-coût marginal et une grande facilité de mise en oeuvre. G. Cohen et P. Monk ont implémenté, pour les équations de Maxwell, la version de ces couches absorbantes proposée par L. Zhao et A. C. Cangellaris. Cette nouvelle formulation est parfaitement adaptée aux formulations variationnelles. Reste un point négatif : malgré nos efforts , le modèle de couches reste mal compris du point de vue mathématique. Si des résultats partiels ont été obtenus dans le cas stationnaire (existence de solutions ... sauf et éventuellement pour un nombre fini de fréquences exceptionnelles dans le cas stationnaire), ceux-ci ne sont pas encore suffisants pour bien comprendre les limites de la méthode.

Couches minces et ferromagnétisme



Participants : P. Joly , H. Haddar , C. Poirier


Mots-clés : Ferromagnétisme, Couches minces, Conditions équivalentes


Résumé : Étude asymptotique de la diffraction d'une onde électromagnétique par une couche mince ferromagnétique. Mise en oeuvre numérique du cas monodimensionnel.


Les matériaux absorbants, tels les matériaux ferromagnétiques, sont souvent utilisés comme revêtements pour rendre des objets furtifs. Dans bon nombre de situations, l'épaisseur de la couche absorbante est très petite devant la longueur d'onde. Dans ce cas, modéliser numériquement la propagation des ondes à l'intérieur du matériau ferromagnétique devient très coûteux en temps de calcul. Il est alors intéressant de pouvoir approcher l'action de cette couche très fine par des conditions aux limites équivalentes. Cela a déjà fait l'objet des travaux de Bouchitté, d'Engquist-Nédelec ou de Bendali-Lemrabet, dans le cas des matériaux linéaires. En ce qui nous concerne, nous avons abordé l'étude du cas non linéaire introduit par les matériaux ferromagnétiques. Nous avons considéré que les caractéristiques de la couche absorbante sont invariants, ou du moins lentement variables, dans les directions tangentiellles à la couche. Le seul petit paramètre est alors l'épaisseur de la couche. Notre travail s'est composé de deux parties :

1.
Construction des conditions approchées, dans le cas de frontières planes ou courbes, à partir de développements asymptotiques formels par rapport à l'épaisseur de la couche absorbante. Étude de la stabilité des nouveaux problèmes [42]. Cette étude a notamment révélé que la condition d'ordre 4 est instable dans le cas monodimensionnel. Cependant son approximation par un développement de Padé admet les bonnes propriétés de stabilité, ce qui a été confirmé par l'expérimentation numérique.
2.
Conception d'un schéma d'approximation numérique pour la prise en compte des conditions équivalentes d'ordre 2 (qui sont aussi d'ordre 3) et d'ordre 4 dans le cas monodimensionnel. Analyse de la stabilité des différents schémas obtenus. Mise en oeuvre numérique de ces différentes conditions. Validation de la démarche à l'aide de comparaisons avec des simulations numériques dans lesquelles la propagation de l'onde à l'intérieur de la couche absorbante est intégralement prise en compte.

Structures filaires



Participants : E. Bécache , P. Joly


Mots-clés : fils, domaines fictifs


Nous avons poursuivi l'étude d'une nouvelle approche pour la diffraction par des structures filaires en partant de la méthode des domaines fictifs analysée par C. Atamian et P. Joly. Rappelons que cette méthode a, par rapport à la théorie classique des potentiels, notamment les avantages de conduire à une équation de caractère auto-adjoint et de supprimer le problème des fréquences irrégulières (pour lesquelles le problème n'est pas bien posé). Dans le cas particulier d'un obstacle filaire, on s'intéresse à l'équation obtenue en effectuant un développement asymptotique sur l'équation intégrale issue de la méthode des domaines fictifs, le petit paramètre étant naturellement le rayon du fil. Pour le problème scalaire (équation de Helmholtz), cette équation est à comparer à l'équation intégrale dérivant de la théorie classique du potentiel étudiée par F. Rogier. Notre objectif est d'obtenir une approximation aussi précise et conservant les deux avantages précédents. L'étude menée a conduit à une équation qui ne satisfait pas encore cet objectif, l'analyse d'erreur indiquant qu'on a perdu de la précision. Dans le cas statique cependant, il a été possible de déterminer le terme suivant du développement asymptotique qui permettrait de regagner la précision perdue. Il reste bien sûr à voir si c'est toujours possible dans le cas qui nous intéresse, c'est à dire dans le cas dynamique.

Propagation Haute Fréquence en milieux inhomogène



Participants : J.-D. Benamou


A la suite des travaux sur la méthode des ``gros rayons'', un lien a pu être établi, à l'aide de la théorie du contrôle optimal, entre solution géométrique multivaluée de l'équation Eikonale et solutions de viscosité univaluées. La solution multivaluée, que l'on veut calculer, peut être décomposée, au moins dans des cas simples, en une succession de solutions de viscosités exploitant les solutions de viscosités précédentes et les données géométriques qui en découlent naturellement (caustiques et shocks). On s'affranchit ici complètement de la technique de lancer de rayons.

Les premiers tests numériques exploitant ces idées sont encourageants. Notons également que cette décomposition de la solution géométrique est naturellement reliée à la théorie géométrique de la diffraction (Keller, Ludwig). Ceci devrait permettre d'utiliser les corrections proposées par cette théorie aux voisinages des caustiques pour l'équation de transport.

Ce travail donne lieu a de multiples collaborations. Il fait l'objet d'une action intégrée Franco-Tchèque avec L.Klimes et P. Bulant (Charles University, Prague) sur les aspects lancer de rayons. L'aspect géométrie différentielle et application à l'acoustique sous-marine est abordé par G. Kossioris (Forth, Heraklion). L'utilisation de solveur Eikonal sur maillage non structuré est abordé en collaboration avec R. Abgrall (Université de Bordeaux). Finalement, des contacts informels se poursuivent sur ce sujet avec W. Symes (Rice University).

Vers une méthode particulaire pour l'équation de Helmholtz



Participants : J.-D. Benamou


L'utilisation de la transformée de Wigner permet, comme dans le cas de l'équation de Shrodinger (P.L. Lions et Paul), d'établir la limite semi-classique de l'équation de Helmholtz vers le système Hamiltonien des rayons lorsque la fréquence $\omega$ tends vers l'infini. Cette technique produit également une suite, indicée par $\omega$, d'équations de type Liouville pour une densité qui est précisément la transformée de Wigner de la solution de l'équation de Helmholtz.

Avec T. Katsaounis (ENS-Ulm) et B. Perthame (ENS-Ulm/M3N), nous commençons l'étude de la discrétisation particulaire de cette densité. On a ainsi pu écrire divers modèles particulaires de l'équation de Helmholtz consistant à la limite avec le modèle asymptotique des rayons. Il font clairement apparaître, à travers la présence d'un noyau de convolution, le caractère diffusif de la propagation des ondes (principe de Huygens) qui existe avant la limite en fréquence. C'est la source principale des difficultés liées à l'implémentation de cette méthode car les rayons n'existent pas à fréquences finie.

Nous testons un algorithme pour le problème unidimensionel.

Couplage non-linéaire de l'équation de transport et de l'équation Eikonale avant la limite en fréquence



Participants : J.-D. Benamou


L'idée principale de l'optique géométrique est, pour les calculs de propagation d'ondes à haute fréquence $\omega$ , de faire l'hypothèse d'un comportement oscillatoire de la solution du type $A(x)e^{i \omega \phi(x)}$ ou $A$ l'amplitude et $\phi$ la phase sont réels. Après introduction de cet ansatz dans l'équation de Helmholtz et séparation des partie réelle et imaginaires on trouve d'une part l'équation classique de transport pour $A$ (la même que pour l'optique géométrique) et d'autre part l'équation Eikonale pour $\phi$ additionnée d'un terme non-linéaire de couplage avec l'équation de transport. C'est précisément ce terme, en $1/\omega^2$, que l'on néglige pour obtenir l'optique géométrique classique.

L'objet de ce travail en collaboration avec B.Engquist et S.Osher est l'étude théorique et numérique du modèle conservant le terme de couplage. On aimerait plus particulièrement évaluer si la correction induite par le terme de couplage permet de traiter des situations pour lesquelles l'optique géométrique classiques échoue (caustiques par exemple).

Modélisation numérique de la timbale



Participants : F. Collino , P. Joly , L. Rhaouti , C. Hazard


Mots-clés : domaines fictifs, non-linéarité, interaction fluide-structure


 Ce sujet est développé en collaboration avec le Groupe Acoustique de l'Ecole Nationale Supérieure des Télécommunications (A. Chaigne, D. Matignon).

Nous nous sommes intéressés à l'analyse mathématique d'un modèle mécanique simplifié pour l'interaction entre un maillet et une membrane élastique. Mathématiquement, ce modèle est un couplage entre une loi non-linéaire et une équation des ondes. Dans un premier temps, à l'aide de techniques mathématiques élémentaires, nous avons mené une étude assez fine du cas où le maillet est en contact avec un corps rigide. Cette analyse permet notamment de comprendre l'influence des paramètres physiques sur le mouvement du maillet. Dans un deuxième temps, nous avons montré que les résultats de cette étude peuvent s'étendre au cas où le maillet est en contact avec une membrane [48].

Concernant la simulation numérique de l'instrument, à cause de la forme particulière de la cavité, il s'est avéré pertinent d'introduire des méthodes de domaines fictifs. Une nouvelle formulation a été écrite dans laquelle la pression est modélisée en vitesse-pression. Les autres inconnues sont le déplacement de la membrane, le déplacement du maillet, et enfin l'inconnue supplémentaire introduite par la méthode de domaines fictifs : le saut de pression sur le bord de l'instrument (membrane + cavité). Un code en 3D a été implémenté et testé numériquement. Il reste à confronter les résultats obtenus avec des résultats empiriques, sur une vraie timbale. Cette partie du travail sera effectuée à l'ENST. Dans la figure [*], on peut visualiser le saut de pression sur le bord de l'instrument ainsi que le déplacement de la membrane à peine quelques millisecondes après que le maillet soit entré en contact avec la membrane.

Nous réfléchissons à une méthode numérique pour le calcul des fréquences de résonance du dispositif afin, dans un premier temps, de confronter les résultats obtenus à une analyse fréquentielle des signaux temporels calculés par le code transitoire.


   Figure: saut de pression (à gauche) et déplacement de la membrane (à droite) à $t=5$ ms

\begin{figure}\begin{center}\includegraphics[angle=0,width=4cm]{f15tp34l.coupj... ...\includegraphics[angle=0,width=4cm]{f15tp34l.coupm5.ps}\end{center}\end{figure}


Rayonnement d'un guide fermé semi-infini



Participants : E. Lunéville , P.M. Cutzach


Mots-clés : Guide d'onde semi-infini, Méthode de Wiener-Hopf, Développements asymptotiques


A la suite d'une question relative à l'optimisation des fenêtres de sortie latérale des gyrotrons, posée il a plusieurs années par Thomson/CSF/DTE, le calcul du champ rayonné par l'embout d'un guide d'onde fermé semi-infini a été réalisé dans le cadre de la thèse de Benoît Aublin. Dans le cas de l'acoustique bidimensionnelle, ce dernier a mis en oeuvre une méthode numérique reposant sur le calcul (par la méthode de Wiener-Hopf) de la fonction de Green d'un guide fermé semi-infini. Il a ainsi résolu par une méthode de type équation intégrale le problème du rayonnement d'un guide semi-infini localement perturbé.

Suite à ce travail, nous avons réalisé une étude théorique des problèmes de diffraction par une telle structure. A partir de la connaissance de la fonction de Green obtenue par la méthode de Wiener-Hopf, nous avons déterminé son comportement asymptotique à l'aide d'une technique de ``Steepest Descent'' à l'extérieur du guide et à l'aide d'un calcul de résidus à l'intérieur du guide. Ces comportements asymptotiques montrent que le champ diffracté vérifie la condition de rayonnement de Sommerfeld à l'extérieur du guide et une représentation en modes sortants dans le guide. Dans un second temps, nous avons été capables de prouver un résultat d'unicité pour l'équation de Helmholtz munie des conditions de rayonnement précédentes. Cette démonstration s'appuie sur des estimations fines du comportement asymptotique d'équations différentielles de type Bessel perturbées, provenant d'une démarche développée par Kato dans un cadre différent. Ce résultat d'unicité nous permet, via l'utilisation de l'alternative de Fredholm et de la représentation intégrale des solutions, d'établir un résultat d'existence.

La démarche suivie dans ce cas permet également de traiter d'autres situations géométriques, telles le cas d'une fissure semi-infinie dans un milieu stratifié (une seule interface). La généralisation au cas d'un guide cylindrique semi-infini est en cours et fait l'objet d'une collaboration avec Nabil Gmati, chercheur à l'Ecole Nationale des Ingénieurs de Tunis, dans le cadre d'un accord franco-tunisien entre le CNRS et la DGRST.

Problèmes de diffraction au second ordre



Participants : Isabelle Champy , Christophe Hazard , Marc Lenoir , Jean-Luc Bellier


Mots-clés : acoustique sous-marine, condition de rayonnement, amplitude limite


Le traitement numérique de problèmes de diffraction au second ordre (en acoustique et en hydrodynamique navale) a fait l'objet depuis quelques années de recherches effectuées sous la direction de Daniel Euvrard, mort accidentellement en Juillet 1994.

Les travaux poursuivis actuellement concernent la question de l'interaction entre une onde acoustique émise par un sous-marin et une onde de gravité (houle), lorsque les longueurs d'ondes des signaux acoustiques et de la houle sont du même ordre de grandeur. La surface libre a la forme d'une tôle ondulée et joue alors le rôle d'un réseau : l'approche que nous développons consiste à effectuer un développement asymptotique du champ de pression suivant les puissances successives d'un petit paramètre représentant l'amplitude de la houle. Le premier terme de ce développement correspond à la réflexion d'une onde acoustique sur une surface plane (problème bien connu); le problème au second ordre fait apparaître l'interaction entre l'onde acoustique et la tôle ondulée.

La particularité essentielle de ce problème est qu'il comprend une condition aux limites dont la donnée possède un support s'étendant à l'ensemble de la surface libre. Cette situation soulève une difficulté relative à la modélisation même d'un problème au second ordre : on ne sait pas en effet écrire explicitement une condition à l'infini qui permette de sélectionner les ondes sortantes (jouant le rôle de la condition de rayonnement écrite au premier ordre). Cette sélection est en fait réalisée en choisissant une solution particulière du problème au second ordre qui s'exprime à l'aide de la fonction de Green sortante du problème du premier ordre. Nous avons justifié par un principe d'amplitude limite que la solution ainsi sélectionnée est bien la solution physique du problème.

Mais la difficulté principale de ce problème au second ordre est sa résolution numérique : la donnée aux limites qui apparaît sur la surface libre a non seulement un support non borné, mais elle est surtout peu décroissante à l'infini, ce qui interdit l'utilisation des méthodes numériques usuelles (telles que la transformation de Fourier rapide). L'algorithme retenu consiste à prendre en compte analytiquement le comportement asymptotique de cette donnée non homogène, la partie résiduelle relevant quant à elle des méthodes numériques usuelles.

Un problème analogue se pose lors de l'étude des mouvements d'un navire soumis à la houle lorsque l'on prend en compte les termes du second ordre dans le développement de la non-linéarité de la condition de surface libre. La méthode retenue par Jean-Luc Bellier dans sa thèse consiste à décomposer la solution du problème au premier ordre en série de Fourier vis-à-vis de la variable angulaire dans le plan horizontal, à en déduire le second membre de la condition de surface libre du problème au second ordre, puis à en calculer la transformée de Fourier au sens des distributions, grâce à une connaissance détaillée de son comportement asymptotique.

Le produit par la transformée de Fourier de la solution élémentaire et un retour au plan physique par une combinaison de calculs analytiques et numériques parmet de relever le second membre de la condition de surface libre et de se ramener à la solution d'un problème de diffraction classique, imposant par là même implicitement une condition de rayonnement.

Les calculs numériques ont été réalisés dans le cas d'une structure oscillant à une fréquence unique mais les développements théoriques permettent de prendre en compte le cas d'une excitation multichromatique.

Modélisation de la propagation d'ondes en milieu ferromagnétique

 

Participants : P. Joly , P. Monk , O. Vacus


La thèse d'Olivier Vacus intitulée ``Propagation d'ondes en milieu ferromagnétique'', a été soutenue le 25 juin 1997 à l'Ecole Centrale Paris. Rappelons que cette thèse portait sur l'analyse du modèle constitué par le couplage des équations de Maxwell avec l'équation de Landau-Lifschitz-Gilbert (LLG) et sur la conception et l'implémentation d'une méthode numérique pour l'approximation de ce modèle [32]. Les résultats de la dernière année sont:

Obtention de nouveaux résultats numériques dans le cas bidimensionnel. Obtention d'estimation d'erreurs pour une classe d'éléments finis (collaboration avec Peter Monk).

Par ailleurs, d'un point de vue plus théorique, nous avons étudié le comportement asymptotique aux temps longs des solutions du modèle.

Enfin, au cours de son séjour à l'Université du Delaware, Olivier Vacus s'est intéressé, toujours en collaboration avec Peter Monk, à l'introduction du terme d'échange, jusqu'ici négligé, dans le modèle (LLG). Ce terme modifie assez sensiblement la nature du problème à résoudre. Une formulation variationnelle bien adaptée a été proposée, ainsi que des calculs d'estimations d'erreurs justifiant de son intérêt. Des expérimentations numériques sont actuellement en cours, mettant en oeuvre en dimensions 1 et 2 des éléments finis quadratiques (avec condensation de masse).

Guides ondes electromagnétiques en Optique



Participants : P. Joly , D. Pedreira , C. Poirier


Fibres optiques

Il s'agit de la suite des travaux de P. Joly et C. Poirier sur le calcul des électromagnétiques guidés dans une fibre optique (voir référence [31]). Nous nous sommes orientés vers le problème plus spécifique du calcul des fréquences de coupure. Ceci a constitué le sujet de séminaire de Recherche de MM. G. Gacon et E. Halère (Ecole Centrale Paris). Le but poursuivi a été de trouver une formulation en domaine borné pour l'équation des seuils mise en évidence par C. Poirier dans sa thèse. Cette entreprise s'est soldée par un échec partiel : une indétermination non encore comprise apparaît dans la détermination de l'opérateur frontière de type Dirichlet-Neumann que l'on est amené à introduire.

Optique intégrée

P. Joly et D. Pedreira ont achevé l'étude, entamée l'an dernier, consistant à analyser l'erreur commise par la méthode qu'ils ont mise au point pour le calcul de guides d'onde en optique intégrée dans le cadre de l'approximation du faible guidage. Par ailleurs, plusieurs tests numériques ont permis de valider la méthode. Enfin, une variante a été proposée pour contourner le problème des "fréquences critiques" inhérent à la méthode originale.

Etude des lignes de transmission supraconductrices



Participants : Anne-Sophie Bonnet-Ben Dhia , Karim Ramdani , Patrick Joly , Eliane Bécache


Mots-clés : ligne microruban, supraconducteur, modèle de London, mode guidé


Ce thème est développé dans le cadre d'un contrat avec la DGA et fait l'objet d'une collaboration avec le Laboratoire Central de Recherches de Thomson à Orsay, l'IRCOM à Limoges et le laboratoire d'électronique de l'ENSEEIHT à Toulouse. Le but poursuivi est de développer puis de valider par comparaison à des mesures des méthodes numériques pour simuler la propagation dans les lignes microruban ou coplanaires réalisées à l'aide de matériaux supraconducteurs.

Nous nous sommes intéressés à l'étude mathématique des modes guidés dans une ligne microruban supraconductrice. En utilisant le modèle de London, on vérifie que tout se passe comme si le matériau supraconducteur avait une permittivité diélectrique négative à basse fréquence. Le changement de signe de la permittivité (négatif dans le supraconducteur et positif ailleurs) a des répercussions importantes sur l'approche mathématique. L'opérateur étant non borné inférieurement, la stabilité du problème ne peut se déduire du théorème de Lax-Milgram mais repose sur l'écriture d'une équation intégrale sur la frontière du supraconducteur. Des difficultés non résolues subsistent lorsque cette frontière présente des points anguleux.

Une fois la stabilité établie, le Principe du Min-Max, appliqué à l'inverse de l'opérateur, permet d'établir l'existence de modes guidés, le nombre de ces modes étant fonction de la largeur du ruban.

Nous avons également étudié formellement un modèle "de couche mince" qui permet de tenir compte de la présence d'un ruban supraconducteur d'épaisseur négligeable à l'aide d'une condition d'impédance appropriée. L'analyse mathématique de ce modèle est menée en collaboration avec Stéphanie Lohrengel : nous avons établi diverses propriétés de l'espace fonctionnel sur lequel est écrite la formulation variationnelle, et les résultats obtenus sont étroitement reliés à la connaissance des singularités du champ électromagnétique aux extrémités du ruban.

Nous envisageons maintenant de résoudre numériquement ce dernier problème par la méthode de régularisation des équations de Maxwell évoquée dans le paragraphe relatif à la diffraction des ondes électromagnétiques.

Ondes de surface en élasticité



Participants : A.S. Bonnet-Ben Dhia , L. Dahi , P. Joly , C. Poirier


Les travaux sur le cas du demi-espace perturbé avec condition de surface libre (sujet de la thèse de J. Duterte) sont maintenant achevés. Cette année, nos efforts ont été consacrés à la rédaction de deux articles : un sur la partie théorique [14], l'autre sur la partie numérique [28].

Avec la thèse de L. Dahi, nous avons commencé à étendre ces travaux au cas où le demi-espace est immergé dans un fluide : il s'agit alors d'étudier les ondes guidées par l'interface fluide-solide.

Comportement asymptotique aux temps long des solutions du système Maxwell-Landau-Lifschitz

 

Participants : P. Joly , A. Komech , O. Vacus


Le problème que nous avons considéré couple les équations de Maxwell avec une loi de comportement magnétique dissipative et non linéaire, la loi de Landau-Lifchitz-Gilbert. Nous nous sommes intéressés dans le cas monodimensionnel au comportement asymptotique pour les temps longs des solutions du problème de Cauchy ainsi défini. Nous avons montré la décroissance locale vers du champ électromagnétique, résultat directement lié à la décroissance locale de l'énergie. Nous avons ensuite introduit l'ensemble des états d'équilibre pour l'équation de Landau-Lifchitz-Gilbert, et nous avons montré qu'il s'agissait en fait d'un ensemble attracteur pour la distribution d'aimantation. Cette analyse a donné lieu à un rapport de recherche ainsi qu'à un article.

Comportement asymptotique aux temps long des solutions de l´équation de Klein-Gordon



Participants : F. Collino , A. Komech


On s'intéresse au comportement asymptotique aux temps longs des solutions de l'équation de Klein Gordon $\partial_{tt} u - \partial_{xx} u+ u^3 - u = 0$.Cette équation admet comme ondes stationnaires, les solitons de la forme $\tanh( \frac{x- vt}{\sqrt{2-2v^2}})$.Après perturbation, on s'attend à ce qu'un second état stationnaire, éventuellement identique soit progressivement atteint tandis qu'une partie de l'énergie s'évacue vers l'infini. Les questions relatives au taux de décroissance de l'énergie locale, des lois régissant les caractéristiques des états asymptotiques sont encore largement débattues. La situation est encore plus complexe lorsque l'on combine plusieurs états asymptotiques. Afin de disposer d'un outil numérique, nous avons construit et mis en oeuvre deux schémas qui permettent d'expérimenter par ordinateur différentes hypothèses. Pour l'instant, il n'y a pas de résultats si ce n'est quelques beaux dessins, telle la figure [*] où est modélisé le choc de deux solitons se percutant à la même vitesse.


   Figure: Visualisation du choc de deux solitons dans le plan $(x-t)$

\begin{figure} \begin{center} \leavevmode {\includegraphics[width=.46\linewidth]{bo_soliton2.ps}} \end{center}\end{figure}


Méthode de Lagrangien augmenté pour la factorisation polaire des champs de vecteurs et problèmes associés



Participants : J.-D. Benamou


La reconstruction du champ de vitesse d'un écoulement dont on connaît certaines données lagrangiennes est un problème fréquemment rencontré dans de nombreuses applications (reconstitution du vent à partir de cartes de concentration en météorologie, du courant marin à partir de trajectoires de bouées en océanographie, etc...). Le plus souvent l'écoulement est régi par des équations de mécanique des fluides plus ou moins complexes, mais en définitive issues des équations de Navier-Stokes ou d'Euler et, assez souvent, en régime incompressible.

Avec Y. Brenier (P6-LAN), un modèle simplifié consistant a reconstruire le flot incompressible transportant une densité initiale vers une densité finale toutes deux prescrites a été considéré. Ce modèle consiste à minimiser de l'énergie cinétique, où les équations d'Euler apparaissent naturellement comme conditions d'optimalité formelles. Ce problème peut également être reformulé à l'aide de la factorisation polaire des champs de vecteurs ou encore de l'équation de Monge-Ampère.

La recherche d'un algorithme numérique performant pour ce type de problème est largement ouvert (c'était en autre un des aspects abordé dans ma thèse). Nous avons appliqué une méthode de Lagrangien augmenté (voir Glowinski-Fortin) avec succès. Les premiers résultats confirment la robustesse de l'algorithme. Bien que nous ne disposions pas d'estimations sur la vitesse de convergence de la méthode, elle semble plus performante que les méthodes employées jusqu'à présent.


previous up next contents Précédent : Logiciels Remonter : Projet ONDES, Modélisation et Simulation Suivant : Actions industrielles