Login| Sign Up| Help| Contact|

Patent Searching and Data


Title:
ELECTROMAGNETIC SIMULATION ALGORITHM, IN PARTICULAR FOR AN ELECTROMAGNETIC ANTENNA PERFORMANCES
Document Type and Number:
WIPO Patent Application WO/2001/095159
Kind Code:
A1
Abstract:
The invention concerns an electromagnetic simulation algorithm, which enables to calculate the electromagnetic wave diffracted by a conductor operating in single-frequency. The invention concerns an electromagnetic simulation algorithm based on iterative resolution of a system of integral equations comprising a preconditioner Said preconditioner is in particular derived from the adaptation of Calderon formulae to integral equations of the electromagnetic boundary, also known as Electric Field Integral Equation (EFIE). It further consists in using an original representation of the residue of the calculations at each iteration. Said representation, as well as a projection and a composition are involved in the expression of said precondition factor. The invention is in particular applicable to simulation instruments used for designing receiving and transmitting antennae such as antennae of cellular telephones, of anti-collision radar, of electronic countermeasure system, of surveillance or tracking radar, or of satellite. The invention is also useful for calculating radar reflective area of objects whereof the geometrical properties are known.

Inventors:
CHRISTIANSEN SNORRE HARALD (FR)
BEREUX FRANCOIS (FR)
NEDELEC JEAN-CLAUDE (FR)
MARTINAUD JEAN-PAUL (FR)
Application Number:
PCT/FR2001/001715
Publication Date:
December 13, 2001
Filing Date:
June 01, 2001
Export Citation:
Click for automatic bibliography generation   Help
Assignee:
THALES SA (FR)
CHRISTIANSEN SNORRE HARALD (FR)
BEREUX FRANCOIS (FR)
NEDELEC JEAN CLAUDE (FR)
MARTINAUD JEAN PAUL (FR)
International Classes:
G06F17/50; (IPC1-7): G06F17/50; G06F17/16
Other References:
YAGHJIAN A D: "Banded-matrix preconditioning for electric-field integral equations", IEEE ANTENNAS AND PROPAGATION SOCIETY INTERNATIONAL SYMPOSIUM 1997. DIGEST (CAT. NO.97CH36122), IEEE ANTENNAS AND PROPAGATION SOCIETY INTERNATIONAL SYMPOSIUM 1997. DIGEST, MONTREAL, QUE., CANADA, 13-18 JULY 1997, 1997, New York, NY, USA, IEEE, USA, pages 1806 - 1809 vol.3, XP002162413, ISBN: 0-7803-4178-3
DENG H ET AL: "Preconditioning of electromagnetic integral equations using pre-defined wavelet packet basis", ELECTRONICS LETTERS, 8 JULY 1999, IEE, UK, vol. 35, no. 14, pages 1144 - 1146, XP002162414, ISSN: 0013-5194
BOWEN J M ET AL: "Use of preconditioners to ease disparate grid size problem", IEEE ANTENNAS AND PROPAGATION SOCIETY INTERNATIONAL SYMPOSIUM. 1998 DIGEST. ANTENNAS: GATEWAYS TO THE GLOBAL NETWORK. HELD IN CONJUNCTION WITH: USNC/URSI NATIONAL RADIO SCIENCE MEETING (CAT. NO.98CH36194), IEEE ANTENNAS AND PROPAGATION SOCIETY INTERN, 1998, New York, NY, USA, IEEE, USA, pages 1770 - 1773 vol.3, XP002162415, ISBN: 0-7803-4478-2
BANNELIER P ET AL: "Mixed finite elements method for antenna design", ANNALES DES TELECOMMUNICATIONS, SEPT.-OCT. 1989, FRANCE, vol. 44, no. 9-10, pages 464 - 474, XP000989716, ISSN: 0003-4347
Attorney, Agent or Firm:
Lucas, Laurent (av. du Prés. Salvador Allende Arcueil Cédex, FR)
Download PDF:
Claims:
REVENDICATIONS
1. Algorithme de simulation électromagnétique, pour déterminer l'onde électomagnétique diffractée par un corps en régime monofréquentiel, à partir d'un maillage dudit corps et de l'excitation électromagnétique, caractérisé en ce qu'il comporte au moins : (a) une détermination d'une matrice M, dite matrice d'interaction, dont les coefficients sont déterminés à partir du maillage dudit corps ; (b) une détermination d'un préconditionneur Z de la matrice M, ce préconditionneur étant la traduction matricielle de l'opérateur 1* y, où est l'opérateur associé à M, 7 l'opérateur produit vectoriel avec la normale à la surface dudit corps, et 1* l'opérateur adjoint de « 1 ; (b) une détermination des courants qui circulent à la surface dudit corps, par un algorithme itératif de type gradient conjugué utilisant ledit préconditionneur Z, I'algorithme itératif permettant de résoudre une équation, dite d'équation intégrale de frontière de l'électromagnétisme, s'écrivant sous forme matricielle de la manière suivante : MU=L où U est un vecteur que l'on cherche à déterminer, dont les coefficients représentent les courants de surface, et L est un vecteur connu, dont les coefficients représentent l'excitation électromagnétique ; (c) une détermination de l'onde diffractée par ledit corps, à partir desdits courants de surface.
2. Algorithme de simulation électromagnétique selon la revendication 1, caractérisé en ce que les coefficients du vecteur U sont exprimés dans la base usuelle de ('espace de RaviartThomas.
3. Algorithme de simulation électromagnétique selon l'une des revendications précédentes, caractérisé en ce que le préconditionneur Z est déterminé de manière implicite.
4. Algorithme de simulation électromagnétique selon l'une des revendications précédentes, caractérisé en ce que le préconditionneur Z est défini par la relation suivante : Z=tJ M J où J est une traduction matricielle de l'opérateur J, et tJ la matrice transposée de J.
5. Algorithme de simulation électromagnétique selon la revendication 4, caractérisé en ce que la matrice J est définie par la relation suivante : (pj est la base usuelle de l'espace de RaviartThomas ; Oi est la base des fonctions affines par triangles ; ri est la base des fonctions constantes par triangles ; A un point appartenant à la surface r dudit corps ; z (A) est un vecteur de norme unitaire, normal à la surface dudit coprs au point A, et orienté vers l'extérieur.
6. Algorithme de simulation électromagnétique selon l'une des revendications précédentes, caractérisé en ce que l'algorithme itératif utilisé est un algorithme rapide, du type méthode multiples multiniveaux.
7. Algorithme de simulation électromagnétique selon l'une des revendications 1 à 5, caractérisé en ce que l'algorithme itératif utilisé est un algorithme rapide, du type méthode méthode intégrale adaptative.
8. Algorithme de simulation électromagnétique selon l'une des revendications précédentes, caractérisé en ce que le corps est une antenne dont on cherche à déterminer une forme optimale, en utilisant l'algorithme de simulation dans un outil de conception d'antennes.
9. Algorithme de simulation électromagnétique selon l'une des revendications 1 à 7, caractérisé en ce que le corps est un objet de forme connue dont on cherche à déterminer la surface équivalente radar (SER).
Description:
ALGORITHME DE SIMULATION ELECTROMAGNETIQUE, NOTAMMENT DES PERFORMANCES D'UNE ANTENNE.

La présente invention se rapporte à un algorithme de simulation électromagnétique, notamment des performances d'une antenne, qui permet de calculer l'onde électomagnétique diffractée par un conducteur en régime monofréquentiel. Elle s'applique notamment aux outils de simulation utilisés lors de la conception d'antennes de réception ou d'émission telle que les antennes de téléphone cellulaire, de radar anti-collision, de système de contre mesure électronique (CME), de radar de veille ou de poursuite, ou de satellite. L'invention s'applique aussi au calcul des surfaces équivalentes radar (SER) d'objets dont on connaît les propriétés géométriques.

On utilise les simulations d'antennes pour limiter le nombre de maquettes et de prototypes lors de la conception desdites antennes. Ces simulations permettent notamment de calculer le diagramme de rayonnement en champ lointain des antennes et d'adapter les antennes en émission ou en réception, en présence ou non d'une structure environnante.

Elles utilisent comme donnée d'entrée un maillage de l'antenne dont on souhaite évaluer les performances, ainsi que les caractéristiques de l'excitation électromagnétique à laquelle elle est soumise. L'invention ne se limite pas aux simulations d'antennes. Elle s'applique aussi par exemple aux calculs de SER de cibles. On décrit maintenant, à titre d'illustration, l'application de l'invention dans les simulations d'antennes, en réception.

On distingue deux principales méthodes dans les simulations couramment employées. Une première méthode est basée sur le calcul par différences finies, encore connue sous le nom de méthode d'éléments finis de volume. Selon cette méthode, on utilise un maillage d'un volume entourant l'antenne. Un inconvénient de cette méthode est que le maillage est nécessairement borné, alors qu'on s'intéresse au diagramme de rayonnement à l'infini. On doit alors réaliser un compromis entre la dimension du volume maillé, c'est à dire la précision de calcul, et le temps de calcul.

Pour pallier cet inconvénient, on utilise une seconde méthode basée sur des équations intégrales dans le domaine fréquentiel. Selon cette méthode, on utilise un maillage en surface de l'antenne uniquement. On calcule directement le diagramme de rayonnement à l'infini à partir de courants électriques et magnétiques sur la surface de l'antenne.

Certaines techniques connues utilisant les équations intégrâtes . calculent les courants électriques et magnétiques (desquels on déduit le

diagramme de champ rayonné à l'infini) grâce à une factorisation d'une matrice. Cette matrice est connue sous le nom de matrice d'interaction ou encore matrice d'impédance. Cette factorisation permet un calcul direct, c'est à dire non itératif des courants de surface. Un inconvénient de ces techniques est que le temps de calcul est long. Si on note N le nombre de points intervenant dans le maillage de l'antenne (encore appelé triangulation de surface), le temps de calcul selon ces techniques varie comme N3. Or le nombre de points N est lui-même lié à la longueur d'onde (et par conséquent la fréquence) de l'onde rayonnée par l'antenne. Supposons que l'on réalise une simulation à 10GHz, en utilisant N points dans le maillage de l'antenne, et que le temps de calcul est. Pour transposer cette simulation à 20GHz, il sera nécessaire d'utiliser un maillage comportant 4xN points, ce qui représentera un temps de calcul de l'ordre de 43x. Un problème de temps de calcul se pose aussi lorsqu'on cherche à simuler des géométries complexes d'antennes, telles que les réseaux de petites tailles. Ceci rend ces techniques inutilisables notamment dans les outils de conception qui requièrent un temps de calcul réduit pour permettre aux concepteurs de réaliser plusieurs tests.

D'autres techniques connues utilisant les équations intégrales permettent de diminuer le temps de calcul grâce à une méthode de résolution itérative. Si on note IT le nombre d'itérations, le temps de calcul selon ces techniques varie comme ITxN2. Un problème de ces techniques est que rien ne garantit la convergence des calculs. En d'autres termes, il existe des formes d'antennes pour lesquelles on ne peut pas calculer le diagramme de champ rayonné avec ces techniques.

Un but de l'invention est de pallier les inconvénients précités, et notamment de réduire les temps de calcul.

A cet effet, l'invention concerne un algorithme de simulation des performances d'une antenne basé sur une résolution itérative d'un système d'équations intégrales comprenant un préconditionneur. Ce préconditionneur est notamment issu de l'adaptation des formules de Calderon aux équations intégrales de frontière de l'électromagnétisme. En particulier dans le cas d'une antenne entièrement métallique, on propose un préconditionneur pour l'équation appelée « Electric Field Integral Equation » (EFIE) dans la littérature anglo-saxonne. On utilise aussi une représentation originale du

résidu des calculs lors de chaque itération. Cette représentation, ainsi qu'une projection et une composition interviennent dans l'expression dudit préconditionneur.

L'invention a les principaux avantages suivants : elle converge rapidement ; elle permet de simuler des géométries arbitraires et des excitations quelconques ; son conditionnement est indépendant de la finesse du maillage ; elle s'adapte aux algorithmes basés sur le calcul d'une matrice d'impédance, en réutilisant ladite matrice d'impédance ; elle permet de traiter des antennes comportant, en plus du métal, des matériaux diélectriques.

D'autres caractéristiques et avantages de l'invention apparaîtront plus clairement dans la description qui va suivre et dans les figures annexées qui représentent : -la figure 1, un maillage d'une antenne ; -la figure 2, une vue en coupe du maillage de la figure 1 ; -la figure 3, un détail du maillage de la figure 1 sur lequel est représenté un champ de vecteurs ; -la figure 4, un diagramme fonctionnel d'un algorithme itératif ; -la figure 5, une fonction de base constante par triangles sur un maillage ; -la figure 6, une fonction de base affine par triangles et continue sur un maillage ; -les figures 7 et 8, deux illustrations des performances de l'algorithme selon l'invention comparé aux techniques connues.

On se réfère maintenant aux figures 1 et 2 qui représentent un exemple de forme d'antenne dont on cherche à déterminer le champ diffracté lorsque l'antenne est éclairée par une onde incidente. En d'autres termes, on cherche à simuler une antenne en réception. L'antenne prise dans cet exemple est une cavité sphérique de demi-angle d'ouverture a/4. Le rayon intérieur est 7/8, le rayon extérieur est 9/8 (unité de longueur arbitraire). La surface de l'antenne est notée r. Cette surface r est maillée par des

triangles. La surface r est supposée être celle d'un conducteur parfait (I'antenne) Q-plongé dans le vide Q+.

Dans cet exemple, l'onde électromagnétique incidente qui illumine l'antenne est une onde plane monofréquentielle. Cette onde électromagnétique incidente, de nombre d'onde kin, connu, est représentée par deux champs de vecteurs notés inc et ci ine correspondant respectivement aux champ électrique et au champ magnétique. Bien entendu, l'invention ne s'applique pas aux seulement aux ondes planes. II est possible de substituer cette onde incidente par le champ émis par un dipôle rayonnant par exemple (mode émetteur de l'antenne).

On cherche à déterminer l'onde électromagnétique diffractée par l'antenne à l'infini, c'est à dire le diagramme de rayonnement en champ lointain. Cette onde électromagnétique diffractée est représentée par deux champs de vecteurs notés Edif et Hdif correspondant respectivement aux champ électrique et au champ magnétique.

A partir du champ de courant circulant à la surface de ladite antenne, encore appelé courants électriques et magnétiques de surface, il est possible de calculer le champ électromagnétique rayonné en tout point de l'espace. L'expression asymptotique à l'infini du champ électromagnétique rayonné est le diagramme de rayonnement que l'on cherche à déterminer.

Ce calcul bien connu en électromagnétisme est rappelé dans le document Integral Equation Methods In Scattering Theory)) de D. Colton et R. Kress -John Wiley & Sons, New-York, 1983.

Le champ électromagnétique satisfait les équations de Maxwell dans le vide #+ qui s'écrivent : <BR> <BR> <BR> <BR> <BR> <BR> #(#) = i#µ# (1)<BR> <BR> <BR> <BR> <BR> <BR> <BR> <BR> <BR> <BR> <BR> #(#) = -i### (2) Les termes utilisés dans ces équations représentent : rot l'opérateur rotationnel ; # # = #inc + Edif le champ électrique total en notation complexe ; # # = #inc + Hdif le champ magnétique total en notation complexe ;

w la pulsation de l'onde électromagnétique ; # µ la perméabilité magnétique ; c la permittivité électrique.

On rappelle que le nombre d'onde est relié à la pulsation simplement par la relation suivante : Le champ électrique de l'onde électromagnétique diffractée Edif s'exprime en tout point B de #+ à partir des courants de surface notés û par les relations suivantes : Edif(B) =k jG, (A, B) û (A) dSgradJ jG, (A, B) div (û) dsJ (4) AE=r k er Les termes des relations (4) et (5) représentent : diVA (u) la divergence prise en un point A du champ de vecteur grade le gradient pris en un point B ; # dSA un élément différentiel de surface ; Gk la fonction de Green standard ; AB la distance entre les points A et B ; # k la norme du nombre d'onde défini par la relation (3).

L'homme de l'art pourra aller chercher d'autres éléments techniques relatifs à ce calcul dans le document susmentionné dans la mesure où celui-ci fait partie intégrante de la description.

Le calcul des courants de surface, c'est à dire du champ de vecteurs û, est déterminé à partir de l'équation variationnelle suivante : Vv m(#,#) = l(#) (6)

dans laquelle Cette équation variationnelle (6) est connue sous le nom de d'équation intégrale de frontière de l'électromagnétisme, ou encore « Electric Field Integral Equation » (EFIE) dans la littérature anglo-saxonne.

Afin de résoudre cette équation variationnelle (6) dans une simulation numérique d'antenne, on doit approcher la solution û dans un espace de dimension finie, dit espace de discrétisation. Cet espace contient des champs de vecteurs qui représentent des courants de surface. La dimension cet espace est le nombre de composantes nécessaires pour décrire entièrement ledit champ de vecteurs, en tout point de la surface r.

Cette surface étant maillée, le nombre de composantes servant à décrire Û sera fonction notamment du nombre de points du maillage N, ainsi que de la nature des fonctions de base servant à décrire le champ de vecteurs (par exemple des fonctions linéaires ou de degré 2). Dans la suite de la description, on prendra à titre d'illustration l'espace de Raviart-Thomas de plus bas degré sur le maillage de r. Cet espace de discrétisation, qui représente des courants de surface, est noté Dh. Le symbole h représente la longueur caractéristique, ou encore précision, du maillage. En effet, le nombre de dimension de cet espace Dh dépend du nombre de point N du maillage, qui lui-même dépend de la précision h du maillage.

On se réfère maintenant à la figure 3 pour décrire une base de l'espace Dh. Dans cette base, le champ de courants de surface Û est représenté par les coefficients d'un vecteur noté U.

La base de l'espace des courants de surface contient des éléments notés #i où i est un indice entier associé à une arête du maillage de la surface r. Ces éléments sont des champs de vecteurs définis sur le maillage de la surface r. Le champ de vecteurs (pj, représenté par des flèches sur la figure 3, a un support borné à deux triangles T1 et T2 du maillage. Ces triangles T1 et T2 ont en commun l'arête d'indice i et de longueur ,. Cette arrête est orientée par un vecteur de norme unitaire noté yj. On note Pi) le sommet du triangle T1 non contenu sur l'arête i ; on note P2 le sommet du triangle T2 non contenu sur l'arête i. On note Si et S2 les surfaces des triangles T1 et T2. Soient Z1 et z2les vecteurs de norme unitaire, ayant une direction normale à la surface des triangles T1 et T2, et orientés de l'intérieur Q-vers l'extérieur Q+. On défini les vecteurs #1i et x21es vecteurs de normes unitaires tels que le triplet (X ; yj ; z1) et le triplet (xi2 ; yj ; Z2) soient des trièdres directs. On définit le champ de vecteurs pour tout point A appartenant à la surface r avec les relations suivantes : sinon, (pj (A) = 0. (11) Une telle base de l'espace des courants est connue sous le nom de base usuelle de l'espace Raviart-Thomas ou encore courants élémentaires de Rao-Wilton-Glisson. L'homme de l'art pourra aller chercher d'autres éléments techniques dans « Electromagnetic scattering by surfaces of arbitrary shape » de S. S. M. Rao, D. R. Wilton et A. W. Glisson-IEEE Trans. Ant. Prop. AP-30, pp. 409-418,1982-dans la mesure où ce document fait partie intégrante de la description.

On note Uj les coefficients du vecteur U. Le champ de vecteurs û (A) utilisé dans la relation (4) se décompose sur la base usuelle de l'espace de Raviart-Thomas de la manière suivante : On transpose maintenant en notations matricielles les relations (6), (7) et (8) en utilisant la base ((p ;) décrite ci-dessus. L'équation variationnelle (6) s'écrit sous la forme du système d'équations linéaires suivant : MU=L (13) Les termes de la relation (13) sont les suivants : M est une matrice connue, dite matrice d'interaction ou encore matrice d'impédance ; L est un vecteur connu, dont les coefficients représentent l'onde incidente, c'est à dire l'excitation électromagnétique ; U est le vecteur que l'on cherche à déterminer, dont les coefficients représentent les courants de surface.

On définit les coefficients Mjj de la matrice d'interaction M et les coefficients Li du vecteur L représentant l'onde incidente par les relations suivantes :

L'homme de l'art pourra aller chercher d'autres éléments techniques relatifs au calcul de ces coefficients dans « Approximation par éléments finis de surface de problèmes de diffraction des ondes électromagnétiques » de A. Bendali-Thèse de l'Université Paris VI, 1984- dans la mesure où ce document fait partie intégrante de la description.

On se réfère maintenant à la figure 4 où est illustré un algorithme itératif de résolution du système d'équations linéaires (13) basé sur la technique de gradient conjugué. II est à noter que la technique de préconditionnement que nous illustrons dans le cas de l'algorithme du gradient conjugué s'applique aussi bien d'autres algorithmes itératifs. On peut citer notamment les algorithmes du Résidu Minimum Généralisé, Gradient Bi-Conjugué, Résidu Quasi-Minimum et du Gradient Bi-Conjugué Stabilisé encore connus dans la littérature anglo-saxonne sous les noms de « Generalized Minimum Residual » (GMRES), « BiConjuguate Gradient » (BiCG), « Quasi-Minimal Residual » (QMR), et « BiConjugate Gradient Stabilized » (Bi-CGSTAB). L'homme de l'art pourra aller chercher d'autres éléments techniques sur les méthodes itératives dans « Templates for the Solution of Linear Systems : Building Blocks for Iterative Methods » de R.

Barrett, M. Berry, T. F. Chan, J. Demmel, J. Donate, J. Dongarra, V. Eijkhout, R. Pozi, C. Romine et H. Van der Vorst-SIAM (1994), Philadelphia, PA- dans la mesure où ce document fait partie intégrante de la description.

La figure 4 illustre un algorithme 40 qui prend en entrée la matrice M et le vecteur L de l'équation (13) et donne en sortie le vecteur U. L'homme de l'art pourra aller chercher d'autres éléments techniques relatifs à la résolution de systèmes d'équations linéaires dans « Iterative Methods for Linear and Nonlinear Equations » de C. T. Kelley-SIAM Frontiers in Applied Mathematics, Philadelphia, 1995-dans la mesure où ce document fait partie intégrante de la description.

Une première étape d'initialisation 41 permet d'initialiser quatre suites de vecteurs notés U [n], R [n], S [n] et P [n] où n est un entier. La première suite U [n] est une solution approchée qui converge vers la solution recherchée U. La seconde suite R [n], connue sous le nom de résidu, converge vers le vecteur nul. La troisième suite S [n], que l'on nomme résidu

préconditionné, converge aussi vers le vecteur nul. La dernière suite P [n] est connue sous le nom de direction de recherche. Les premières valeurs de ces suites sont définies par les relations suivantes : UIOI = 0 (16) R[0] = L - M U[0] (17) S[0] = Z R[0] (18) P[0] = S[0] (19) La relation (16) est utilisée par défaut lorsqu'on ne connaît pas de solution approchée. Une variante de cette relation consiste à prendre pour U [0] le résultat d'un calcul de courant de surface réalisé pour une même antenne mais à une autre fréquence.

Le terme Z de la relation (18) est une matrice. Cette matrice est un préconditionneur pour la matrice M selon l'invention. On rappelle qu'un préconditionneur pour M est une matrice approchant l'inverse de M. Selon une variante, on peut remplacer le préconditionneur Z par l'identité. En d'autres termes, on peut initialiser S [0] avec R [0].

Une seconde étape d'itération 42 permet de calculer les valeurs des suites précitées à un rang n+1 à partir des termes de rang n. Cette étape d'itération utilise les relations suivantes : U [n + 1] = U [n] + (x P [n] (20) R [n + 1] = R [n]-a M P [n] (21) S [n+1] =ZR [n+1] (22) (24)

où (,) représente le produit scalaire complexe en notation matricielle.

Une dernière étape 43 réalise un test de convergence. Ce test peut s'exprimer par exemple par l'inégalité suivante : où P est un seuil prédéterminé. En d'autres termes, on compare le <BR> <BR> résidu préconditionné normalisé par rapport au un seuil prédéterminé<BR> Lorsque l'inégalité (25) est vérifiée, le calcul s'achève et on prend U = Une comme solution. Dans le cas contraire, on incrément la valeur de n et on retourne à l'étape 42.

Bien entendu, il est possible d'utiliser un autre test de convergence pour arrêter les itérations. On peut par exemple remplacer l'inégalité (25) par l'inégalité suivante : L'algorithme du gradient conjugué dans les techniques actuelles est utilisé sans préconditionneur, c'est à dire avec la matrice Z égale à l'identité.

L'invention consiste à utiliser un préconditionneur basé sur une généralisation d'une formule de Calderon. L'homme de l'art pourra aller chercher des éléments techniques sur cette formule de Calderon dans « (Mathematical Methods in Electromagnetism, Linear Theory and Applications » de M. Cessenat-World Scientific Publishing Co., page 89, 1996-dans la mesure où ce document fait partie intégrante de la description. Cette formule s'écrit :

(27) Dans cette relation (27), SX, I et t sont trois opérateurs sur les champs de vecteurs tangents, c'est à dire des applications qui à un champ de vecteurs tangents associent un autre champ de vecteurs tangents, et I est l'application identité. 7est l'opérateur de produit scalaire avec la normale à la surface de l'antenne, et M est l'opérateur associé à la matrice d'interaction dont on cherche un préconditionneur. Les opérateurs M, J et R sont définis formellement par les relations suivantes : (u) (B) = k f Gk (A, B) û (A) dSA AEr + kgradB f Gk (A, B divA (û dSA (28) axer (J #)(B) = # (B) Az (B) (29) (tu) (B) = jgrad, (G, (A, B)) Aû (A) dS (30) AeF Dans les relations (28), (29) et (30), X û, y û, û et û sont des champs de vecteurs tangents. L'indice t de la relation (28) représente la composante tangentielle du vecteur entre crochets. B est un point de la surface r. z (B) est un vecteur unitaire, normal à la surface r en B orienté vers l'extérieur.

La Demanderesse a constaté que l'opérateurs, étant compact, c'est à dire négligeable, l'opérateur 4 est approximativement un inverse de l'opérateur M. Ensuite, les préconditionneur étant définis à une constante multiplicative près, il est possible d'éliminer la constante 4. Sachant que zu =-Y,) e préconditionneur selon l'invention est défini à partir de J* X J. On rappelle que 1* est l'opérateur adjoint de J, c'est à dire l'opérateur vérifiant la relation suivante :

# J*### = ##### (31) Le préconditionneur selon invention est une traduction matricielle de l'opérateur J# M J. C'est dans cette traduction matricielle qu'apparaît l'avantage d'utiliser l'opérateur 1*X1 plutôt que l'opérateur j M-7. En effet, la traduction matricielle Z de l'opérateur 1*91d1 est une matrice symétrique, ce qui est essentiel pour les algorithmes itératifs, tel que le gradient conjugué.

On décrit maintenant un exemple de préconditionneur selon l'invention qui permet d'accélérer la convergence de l'algorithme et aussi de rendre cet algorithme plus stable (on converge toujours vers la solution quelles que soient les conditions initiales). Ce préconditionneur est adapté aux problèmes électromagnétiques et tire parti de la structure du problème à résoudre.

On définit d'abord des espaces et leurs bases associées qui serviront dans la suite de la description. On a déjà défini la base ((pu).

L'espace engendré par cette base, Dh, est un espace de champs de vecteurs tangents. En d'autres termes, Dh est un ensemble de fonctions qui à tout point de r associent un vecteur tangent à la surface maillée r. On définit aussi deux espaces de fonctions Ch et Sh. Les fonctions de ces espaces associent une valeur scalaire à chaque point de r.

Ch est l'ensemble des fonctions constantes par triangles, c'est à dire dont la valeur est constante pour tout point appartenant à un triangle donné. La base de Ch que l'on utilise par la suite est l'ensemble des fonctions notées vj, qui valent 1 sur le triangle d'indice i et 0 en dehors. On se réfère à la figure 5 sur laquelle est représentée une fonction de base yr ;.

On a représenté par soucis de clarté un maillage 50 plan et régulier. Les valeurs de la fonction vj sont représentées sur un axe 51, perpendiculaire au maillage 50. La fonction vj représentée sur cette figure vaut 1 sur le triangle i et 0 en dehors.

Sh est l'ensemble des fonctions affines par triangles et continues.

Ces fonctions ont un gradient constant sur tout triangle donné. La base de Sh que l'on utilise par la suite est l'ensemble des fonctions notées Oj, qui valent 1 au noeud d'indice i et 0 sur les autres noeuds. On se réfère à la

figure 6 sur laquelle est représentée une fonction de base Oj. On a représenté par soucis de clarté un maillage 60 plan et régulier. Les valeurs de la fonction Oj sont représentées sur un axe 61, perpendiculaire au maillage 60. La fonction 8 ; représentée sur cette figure est affine par triangles, et vaut 1 au noeud i. Cette fonction a un support (des valeurs non nulles) borné aux triangles 62,63,64,65,66,67 qui ont un sommet confondu avec le noeud i.

L'homme de l'art pourra aller chercher des éléments techniques complémentaires dans « Handbook of Numerical Analysis Vol. II, Finite Elements Methods (Part 1) » de P. G. Ciarlet-Ed. J. L. Lions, North-Holland, 1991-dans la mesure où ce document fait partie intégrante de la description.

On définit maintenant des matrices qui serviront dans la suite de la description du préconditionneur. Ces matrices sont basées sur les bases définies ci-avant. Elles s'expriment par les relations suivantes : où z (A) est un vecteur de norme unitaire, normal à la surface r au point A, et orienté vers l'extérieur.

On décrit maintenant les étapes qui permettent de calculer le résidu préconditionné à partir du résidu dans les relations (18) et (22). Ces relations utilisent le préconditionneur Z que l'on décrit maintenant. Dans la description qui va suivre, on donne une décomposition de cette matrice Z en produit de matrices. Les matrices intervenant dans ce produit sont des matrices creuses ou des inverses de matrices creuses. Ainsi, selon une variante avantageuse de l'invention, on utilisera préférentiellement les matrices creuses plutôt que directement la matrice Z, ce qui permet notamment de réduire la mémoire utilisée et le temps de calcul.

On note maintenant R le résidu, et S le résidu préconditionné. R correspond respectivement à R [O] et à R [n+1] dans les relations (18) et (22).

S correspond respectivement à S [0] et à S [n+1] dans les relations (18) et (22).

Une première étape consiste à utiliser une représentation mixte.

Cette étape sera mieux comprise à l'aide des notations vectorielles qui seront ensuite transposées en notations matricielles. Le résidu R correspond à une forme bilinéaire sur l'espace Dh défini ci-avant. On note p cette forme bilinéaire. On représente p par un champ de vecteurs r1 dans Dh de divergence nulle et une fonction q1 dans Ch d'intégrale nulle. Le champ r1 et la fonction q1 sont définis par la relation suivante pour tout champ de vecteurs v dans Dh : Le champ de vecteurs r1 ayant une divergence nulle, on peut écrire pour toute fonction f dans Ch et d'intégrale nulle : Dans la suite, pour ne pas compliquer inutilement la présentation, on ne rappelle pas que les champs scalaires considérés sont d'intégrale nulle.

Les relations (38) et (39) se transposent de façon matriciel comme suit : Dans la relation (40), on a utilisé une notation matricielle par blocs, dans laquelle : -*M2 est) a matrice transposée de M2 ; 0 est un bloc nul (vecteur ou matrice) ; R1 est la représentation matricielle de r1 ; Q1 est la représentation matricielle de q1.

On rappelle qu'il n'est pas nécessaire d'inverser la matrice définie par blocs pour déterminer R1 et Q1 à partir de R. En effet, cette matrice est une matrice creuse, c'est à dire qui contient beaucoup de termes nuls. La résolution d'un tel système est bien connue de l'homme du métier. Elle est rappelée dans « Handbook of Numerical Analysis Vol.//, Mixed and hybrid methods » pp 523-640 de J. E. Roberts et J.-M. Thomas-Ed. J. L. Lions, North-Holland, 1991. Ce document fait partie intégrante de la description.

Une seconde étape consiste à projeter r1 et q1. Ces projections se traduisent en notations vectorielles par les relations suivantes : r2= (rlAz) (41) q2 = #Sh(q1) où : pDh est l'opérateur de projection dans Dh ; pSh est l'opérateur de projection dans Sh ; est un vecteur unitaire, normal à r et orienté vers l'extérieur.

Ces projections définies par les relations (41) et (42) se traduisent de façon matricielle par :

La résolution de la relation (43) est similaire à la résolution de la relation (40) dans la mesure où les matrices M1 et M5 sont des matrices creuses. On peut noter que cette résolution est encore plus facile que celle de la relation (40) dans la mesure où les matrices M1 et M5 sont de plus symétriques, définies positives et bien conditionnées. On détermine ainsi R2 et Q2 à partir de R1 et Q1.

Une troisième étape consiste à combiner r2 et q2 en utilisant la relation suivante : r3 = r2-rot (q2) (44) Cette troisième étape se traduit de façon matricielle par : R3 = R2-M6 Q2 (45) On appelle j l'application vectorielle qui résulte de la composition des trois étapes décrites ci-avant. Cette application est définie par : j (r) = r3 (46) Cette relation (46) se transpose en notations matricielles par : J R = R3 (47) où on note J la matrice correspondant à l'application vectorielle j, la matrice J étant définie par le produit suivant :

II n'est bien entendu pas nécessaire de calculer directement J, car la décomposition (48) en produit de matrices creuses et d'inverses de matrices creuses est plus simple à utiliser. En d'autres termes, la matrice J est calculée de manière implicite lors de son utilisation.

On note z l'application vectorielle qui correspond à la matrice du préconditionneur Z. Le préconditionneur z est défini à partir de l'application vectorielle j par : .j* offno j (49) Z = J ° m ° J (49) j est l'adjoint de j ; o est l'opérateur de composition d'applications ; lfi est l'application linéaire Dh o Dh associée à la forme bilinéaire m.

La relation (42) se traduit de façon matricielle par : Z=tJ M J (50) On se réfère maintenant aux figures 7 et 8 sur lesquelles sont illustrées les performances d'un algorithme avec un préconditionneur selon l'invention par rapport aux techniques connues n'utilisant pas de préconditionneur.

Sur la figure 7, la courbe 70 représente la fonction en l'absence de préconditionneur. La courbe 71 représente la fonction en présence du préconditionneur décritci-dessus.

Sur la figure 8, la courbe 80 représente la fonction en l'absence d'un préconditionneur. La courbe 81 représente la fonction en présence du préconditionneur

décrit ci-dessus.

On constate que si on prend un critère de convergence classique ii = 104 dans la relation (25), 50 itérations suffisent avec le préconditionneur. Avec un algorithme non préconditionné, on n'atteint pas une précision suffisante en 200 itérations.

Ces simulations numériques on permis de montrer que l'utilisation d'un préconditionneur selon l'invention permet d'accélérer et de stabiliser les algorithmes itératifs.

Bien entendu, l'invention ne se limite pas à l'exemple utilisé pour la décrire. II est possible notamment d'utiliser d'autres fonctions de base ou un maillage différent de ceux pris en exemple.

L'invention se généralise à tout autre discrétisation de l'espace Dh. Si on prend pour Dh un espace différent de l'espace de Raviart-Thomas, on remplace les espaces Ch et Sh respectivement par : Ch = {div(#)## # Dh}+1 (51) Sh = rot (p) e Dh (52) En d'autres termes, Ch est un espace d'éléments finis minimal tel que la divergence des éléments de Dh soit comprise dans Ch ; Sh est un espace d'éléments finis maximal tel que le rotationnel des éléments de Sh soit compris dans Dh.

Une application principale de l'invention se trouve dans les outils de conception d'antenne, mais l'invention ne se limite pas à cette seule application. L'invention s'applique bien entendu aussi à tout outil de simulation basé sur le calcul du champ rayonné par un conducteur. On peut citer notamment le calcul de surfaces équivalentes radar (SER) d'objets dont on connaît les propriétés géométriques.

II est à noter aussi que la technique de preconditionnement décrite dans les cas d'une méthode itérative s'applique aussi à d'autres méthodes numériques rapides. Ces méthodes rapides sont basées sur une résolution itérative, mais seuls les termes utiles aux produits matrice-vecteur sont calculés. Ainsi, on calcule de l'ordre de N x log (N) éléments de la matrice d'impédance, au lieu de N2 éléments selon les techniques itératives classiques. En d'autres termes, on calcule la matrice d'impédance de manière implicite. L'utilisation du préconditionneur selon l'invention dans ces méthodes rapides se fait sans difficulté. Ces méthodes sont intéressantes pour les objets de grandes taille, c'est à dire pour N grand. C'est le cas notamment pour les simulations d'antennes dites sur structure. On peut citer notamment les Méthodes Multiples Multiniveaux et les Méthodes Intégrales Adaptatives, connues aussi dans la littérature anglo-saxonne sous les noms de « Fast Multilevel Multipol Methods » (FMM) et de « Adaptative Integral Method » (AIM). L'homme de l'art pourra aller chercher des éléments techniques sur : les méthodes rapides en générât dans « Fast Solution Methods In Electromagnetics » de W. C. Chew, J.-M. Jin, C.-C. Lu, E. Michielssen, J. M. Song-IEEE Trans. on Antennas and Propagation, 45 (3) : 533-543, Mars 1997 ; les Méthodes Multipôles Multiniveaux dans « Multilevel Fast Multipole Algorithm For Electromagnetic Scattering By Large Complex Objects » de J. M. Song, C.-C. Lu, W. C. Chew, S. W. Lee-IEEE Trans on Antennas and Propagation, 45 (10) : 1488-1493, Octobre 1997 ; les Méthodes Intégrales Adaptatives dans « AIM : Adaptative Integral Method For Solving Large Scale Electromagnetic Scattering And Radiation Problems » de E. Bleszynski, M. Bleszynski, T. Jaroszewicz- Radio Science, 31 (5) : 1225-1251,1996 ; dans la mesure où ces documents font partie intégrante de la description.

L'invention s'étend sans difficulté aux objets (antennes ou cibles) comportant des matériaux diélectriques. Dans ce cas, on cherche des courants électriques et magnétiques équivalents sur chaque interface. La matrice d'interaction entre ces courants comporte des blocs diagonaux du même type que la matrice M décrite ci-dessus. Un préconditionneur pour la matrice d'interaction est donc obtenu en considérant la matrice diagonale par blocs, dont les blocs sont du type de la matrice Z décrite ci-dessus.