Précédent : Logiciels Remonter : Rapport activite 1997 Suivant : Actions industrielles
Précédent : Logiciels Remonter : Rapport activite 1997 Suivant :
Actions industrielles
Résumé : L'étude des méthodes DIMISIM, SIRK et Radau5 s'est pousuivie, et a conduit en particulier à l'élaboration du code Radau5M. Par ailleurs, le projet a introduit la notion de pseudo-symplecticité, qui semble posséder d'intéressantes applications pratiques. Enfin, l'application d'une variante de la méthode de A. Bellen et M. Zennaro au calcul de trajectoires a donné des résultats encourageants.
Participants : Anne Aubry , Philippe Chartier
Nos travaux ont porté sur une technique de
<<smoothing>> permettant de rétablir l'ordre de
convergence habituel [15].
L'idée consiste à utiliser les approximations de la composante
z disponibles aux pas précédents et de les combiner
de façon à éliminer les termes d'erreur indésirables. Cette
technique occasionne un surcoût quasi nul et a été incorporée au
code RADAU5 devenu RADAU5M et mis à disposition publique sur les
sites Web
http://www.irisa.fr/aladin/perso/aubry.html ,
http://www.irisa.fr/aladin/perso/chartier.html
Participant : Philippe Chartier
L'implémentation des méthodes DIMSIM a nécessité la définition d'une stratégie de contrôle du pas fiable. Pour ce faire, nous avons eu recours à une réécriture adéquate des formules d'intégration, permettant une estimation rigoureuse (asymptotiquement correcte) des erreurs locales et autorisant des changements de pas d'un coût négligeable. L'étude de la stabilité à pas variable s'en trouve en outre simplifiée. Ces travaux font l'objet d'une collaboration avec J.C. Butcher, de l'université d'Auckland, Nouvelle-Zélande et Z. Jackiewicz de l'université d'Arizona, Etats-Unis et ont d'ores et déjà donné lieu à publication [18].
Participant : Philippe Chartier
Les méthodes SIRK classiques
souffrent d'un handicap important qui les a jusqu'alors confinées
à des applications particulières. Il est lié à la localisation
des approximations intermédiaires. L'extension proposée permet de
lever cette contrainte tout en conservant la propriété
fondamentale des méthodes classiques, et ce en ayant recours à la
notion d'``ordre effectif''. L'étude de la stabilité de ces
méthodes et de leur implantation est en cours.
La méthode de Newton modifiée (voir section ![[*]](../icons/cross_ref_motif.gif){kind=icon}
) reste au coeur de l'intégration
numérique aussi bien pour les méthodes SIRK que pour les autres méthodes de Runge-Kutta. La
possibilité de découpler le système non-linéaire à résoudre en
sous-systèmes est conditionnée à l'application de certaines
transformations linéaires, de coût négligeable pour les systèmes
de grande taille, mais important pour ceux de faible dimension.
Même dans les cas favorables, la résolution de ces différents
sous-systèmes reste séquentielle. L'introduction d'un niveau
supplémentaire d'itérations (relaxation non-linéaire), permet au
contraire un découplage complet et l'utilisation du parallélisme.
Le choix de la matrice d'itération sous-jacente est crucial et a
fait l'objet d'une étude de convergence dans le cas des méthodes
de Runge-Kutta ``classiques'' [20].
Participants : Anne Aubry , Philippe Chartier
Les méthodes symplectiques ont leur talon d'Achille : elles sont implicites et donc coûteuses. C'est pour cette raison que nous avons introduit le concept de pseudo-symplecticité [16]. Une méthode numérique d'ordre p est dite pseudo-symplectique d'ordre q si son flot vérifie
![\begin{eqnarray*}{\varphi '}^T \left[\begin{array}{cc}0 & I_n \\ -I_n & 0\e... ...ray}{cc}0 & I_n \\ -I_n & 0\end{array}\right] + {\cal O}(h^q). \end{eqnarray*}](img220.gif){kind=small}
On peut alors montrer qu'une méthode d'ordre p d'ordre de pseudo-symplecticité 2p a une erreur croissant linéairement elle-aussi. L'avantage majeur est que de telles méthodes sont explicites [10].
Participants : Philippe Chartier , Jocelyne Erhel ,
Bernard Philippe , Stéphanie Rault
La trajectoire d'un satellite satisfait une équation
différentielle que l'on cherche à résoudre sur calculateur
parallèle. Nous utilisons une partition de l'intervalle
d'intégration en sous-intervalles sur lesquels nous pouvons mener
des intégrations indépendantes. Il faut alors résoudre un
problème non linéaire de point fixe, à l'aide d'une méthode de
Newton modifiée (voir section ). Pour garantir un faible nombre
d'itérations, il faut d'une part partir d'une approximation
initiale proche de la solution et d'autre part utiliser une
approximation du Jacobien assez précise.
Nous avons étudié trois approximations par des modèles analytiques, de plus en plus complexes. Le dernier modèle améliore nettement la précision et par conséquent la convergence de Newton. Pour calculer le Jacobien de ce modèle, nous avons utilisé le logiciel de différentiation automatique Odyssée développé par le projet Safir ; en effet, il est difficile de dériver formellement la fonction parce qu'elle est définie par un schéma itératif.
Nous avons d'abord programmé une parallélisation semi-automatique, qui s'est avérée peu performante. Nous avons donc choisi une autre technique de parallélisation, avec un meilleur équilibrage des charges et une meilleure convergence de Newton.
Enfin, nous avons constaté que le découpage influait fortement sur la convergence. Nous avons donc défini un découpage automatique, au prix parfois d'un léger déséquilibre des tâches. Le nombre d'itérations peut ainsi être réduit, souvent divisé par deux. L'accélération actuelle est environ 3 pour 8 processeurs, ce qui est déjà bien pour une équation différentielle [36,28].
Résumé : Le projet poursuit l'étude des méthodes itératives pour résoudre des systèmes linéaires. Cette année, le projet a démontré un résultat théorique important permettant d'accélérer la résolution successive de plusieurs systèmes linéaires symétriques. Les travaux sur GMRES(m) avec déflation ont abouti à la réalisation d'un logiciel parallèle prototype. Le projet a également obtenu un résultat marquant pour définir une stratégie du choix de m dans GMRES(m).
Quant aux problèmes non linéaires, ils sont abordés cette année par le biais d'une application en mécanique des matériaux composites. Le projet a adapté une méthode de résolution permettant de découpler les équations linéaires et les équations non linéaires du problème global. Les performances obtenues montrent clairement une convergence quadratique qui réduit fortement les coûts de calcul.
Participants : Jocelyne Erhel , Frédéric Guyomarc'h
Il s'agit ici de résoudre une succession de systèmes linéaires
avec différents seconds membres, la matrice étant symétrique
définie positive. Nous avons étudié une approche de type
sous-espace augmenté (voir section ). L'intérêt dans ce contexte est
de préserver les relations de récurrences courtes lors des
itérations du gradient conjugué. Nous avons montré que la vitesse
de convergence est effectivement plus grande dans le gradient
conjugué augmenté que dans le gradient conjugué classique. Les
essais numériques ont confirmé les résultats théoriques [34]. Nous appliquons maintenant
cette méthode à un problème d'électromagnétisme, en coopération
avec le projet M3N.
Participants : Jocelyne Erhel , Bernard Philippe
Nous avons poursuivi nos travaux sur l'algorithme GMRES(m) où
est le paramètre de redémarrage.
Nous avons étudié et comparé deux principales méthodes de
déflation, par sous-espace augmenté ou par préconditionnement
(voir section ). Si le sous-espace qui sert à la
déflation est exactement invariant, les deux méthodes sont
semblables. Dans les deux cas, il apparaît important de calculer
à chaque redémarrage une approximation de plus en plus précise du
sous-espace invariant. La méthode DEFLATION que nous avons
définie est basée sur un préconditionnement. Elle est aussi
efficace et plus robuste qu'une méthode de type sous-espace
augmenté. Ces travaux font l'objet d'une collaboration avec K.
Burrage, de l'université du Queensland, Australie [17].
Par contre, l'approche par sous-espace augmenté est plus aisément parallélisable. Nous l'avons incluse dans notre version parallèle de GMRES(m). Le logiciel est écrit en MPI et les tests sur la machine Paragon ont montré son efficacité. Ces travaux sont décrits dans un rapport de stage.
La surveillance de la convergence des méthodes itératives est en général limitée à la détection du passage sous un seuil donné d'un résidu. Nous avons voulu examiner des modes plus élaborés de surveillance. A cette fin, nous avons appliqué des méthodes de traitement du signal sur la suite des résidus afin de segmenter l'évolution en des phases homogènes et d'identifier les paramètres du signal sur un modèle particulier. Cette technique a été appliquée au cas de la résolution de systèmes par la méthode GMRES(m) ce qui a permis de définir une stratégie de démarrage à taille variable : l'utilisateur peut indiquer en entrée la valeur du paramètre qui correspond à la limite autorisée en espace mémoire, et la résolution adapte automatiquement à une taille plus efficace si m a été choisi trop grand. Ce travail a été réalisé avec l'appui du projet Sigma2 (Bernard Delyon). Il est décrit dans un rapport de stage.
Participants : Mathias Brieu , Jocelyne Erhel
Un des soucis principaux des industriels lors de la conception puis de l'utilisation de plus en plus fréquente de matériaux composites est la détermination de leur comportement à moindre coût numérique. Pour atteindre ce comportement, des méthodes d'homogénéisation sont classiquement mises en oeuvre.
Dans le cas de matériaux composites hyper-élastiques, l'homogénéisation induit des problèmes fortement non linéaires qui sont difficiles à résoudre à cause de la nature des composants.
Les méthodes jusqu'à maintenant proposées sont des techniques de résolution incrémentale. Nous avons mis en oeuvre une méthode de résolution non incrémentale adaptée au cas particulier des problèmes d'homogénéisation considérés. Cette méthode repose sur une décomposition du problème initial en une succession de problèmes vérifiant alternativement les équations linéaires et non linéaires du problème d'origine.
Cette nouvelle méthode de résolution nous a permis dès à présent de réduire les coûts de calcul d'un facteur 20. Nous étudions désormais le couplage de cette méthode avec une méthode de décomposition en sous domaines afin de pouvoir étudier des structures et des problèmes de plus grande taille.
Ce travail s'effectue en collaboration avec F. Devries de l'université de Paris VI.
Résumé : Les activités suivent deux axes : amélioration d'algorithmes de calcul de valeurs propres et fiabilité des calculs en tenant compte de l'imprécision des données. Dans le premier cas, on traite le difficile cas du calcul des plus petites valeurs propres d'une très grande matrice issue d'une discrétisation. Dans le deuxième, on qualifie la confiance à accorder à la localisation des valeurs propres calculées.
Participants : Mickaël Robbe , Miloud Sadkane
Nous nous sommes intéressés au calcul du sous-espace invariant correspondant aux valeurs propres de plus petit module d'une matrice A non singulière, issue de la discrétisation d'un opérateur différentiel elliptique obtenu par un schéma aux différences finies. Nous avons travaillé à l'élaboration d'un algorithme permettant de résoudre ce problème, en s'inspirant de la méthode de descente et des méthodes de projection de type Davidson. Le coeur de l'algorithme est un opérateur K de lissage basé sur les méthodes multigrilles de résolution de systèmes linéaires. Ces méthodes sont particulièrement performantes pour les problèmes elliptiques. De plus, leur efficacité est indépendante de la dimension de A, et sont donc bien adaptées aux problèmes de très grande taille [30].
Participants : Olivier Bertrand , Bernard Philippe ,
Sophie Robert
L'activité dans ce domaine se poursuit dans deux directions en particulier sous l'impulsion du projet européen STABLE. Un nouvel algorithme parallèle a été mis au point ; il est fondé sur un algorithme de Lanczos utilisant les inverses de matrices de Hessenberg dérivées de la matrice dont on recherche le portrait spectral. Le portrait est construit sur une grille du plan complexe avec un procédé ingénieux de choix de vecteur initial qui supprime des dépendances nuisibles au parallélisme. Cet algorithme se présente comme le plus rapide de tous pour les matrices pleines [13].
La deuxième direction consiste à suivre les contours des pseudo-spectres.
L'algorithme déjà connu pour cela présente deux
inconvénients : d'une part il ne dispose pas de contrôle du
pas et d'autre part il ne permet pas de calculer dans le même
temps une intégrale qui permette de donner le nombre de valeurs
propres entourées par la courbe construite. Nous avons d'abord
déterminé un contrôle de pas qui assure un calcul fiable du
nombre de valeurs propres entourées par la courbe fermée
:
Nous améliorons actuellement le procédé par un changement d'intégrateur qui utilise de plus grands pas et donc moins de calcul de valeurs singulières.
Participants : Pierre-François Lavallée , Miloud
Sadkane
Une alternative au calcul du 
-spectre
standard a été étudiée. Basée sur les méthodes de dichotomies
spectrales, cette approche fiable et complète permet dans un
premier temps de localiser les différentes taches spectrales dans
le plan complexe, puis dans un second temps de calculer le
sous-espace invariant associé à chacune de ces taches. Pour
diminuer les coûts de calcul, une version parallèle (MPI) de la
détermination récursive des taches spectrales a été développée
[12,32].
Résumé : Dans le cas où la matrice définissant le système linéaire ou le problème aux valeurs propres est une matrice d'intervalles, il est nécessaire de définir des méthodes qui construisent un encadrement de tous les résultats possibles. Le projet propose de nouveaux algorithmes pour le faire.
Participants : Olivier Beaumont , Bernard Philippe
Nous nous sommes intéressés à deux problèmes liés à la résolution de systèmes linéaires. Considérons le système
![\begin{displaymath}[A]x = [b].\end{displaymath}](img233.gif){kind=small}
Une première approche consiste à déterminer l'ensemble des
qui peuvent être solution d'un système 
pour ![$A \in [A]$](img236.gif){kind=small}
et ![$b \in [b]$](img237.gif){kind=small}
.
(<<outer problem>>)
Une seconde approche consiste à déterminer dans quel intervalle
peuvent varier les entrées 
du système pour assurer que
la sortie sera dans ![$[b]$](img239.gif){kind=small}
quel que soit ![$A \in [A]$](img240.gif){kind=small}
.
(<<tolerance problem>>)
Dans les deux cas, de nouveaux algorithmes ont été proposés, qui
donnent des résultats de meilleure qualité que les méthodes
itératives avec un coût du même ordre (
opérations),
dans un domaine d'application plus large que celui des méthodes
directes.
Nous nous sommes également intéressés à la résolution de
problèmes aux valeurs propres. C'est un domaine assez peu étudié,
bien que d'une grande importance dans le domaine du contrôle.
Seule une méthode directe aux conditions d'applications très
restrictives est connue.
Nous avons développé un algorithme dans le cas symétrique, qui
donne des résultats très satisfaisants en un temps raisonnable.
Malheureusement, ses conditions d'application sont encore trop
restrictives pour pouvoir traiter le cas de matrices aux valeurs
propres trop voisines.
Résumé : Les logiciels de simulation numérique doivent gérer de très grandes tables de données qui dans un contexte parallèle sont distribuées dans les différentes mémoires des processeurs. Le projet a modélisé et simulé deux stratégies d'accès à ces tables. Ce travail s'est effectué en collaboration avec le CEA.
Le projet utilise fréquemment le logiciel Matlab pour ses développements algorithmiques. Toutefois, le passage à une version parallèle performante nécessite un nouveau développement dans un langage de programmation tel que Fortran 90. C'est pourquoi le projet, en coopération avec le projet Caps, a développé un compilateur du langage Matlab qui génère un code parallèle.
Participants : Nicolas Malléjac , Bernard Philippe
Dans de nombreuses simulations numériques, l'algorithme de
résolution utilise des fonctions de contantes physiques définies
sur le domaine de l'étude. Elles dépendent souvent de plusieurs
variables comme par exemple les variables d'espace, la
température, la pression... Ces fonctions étant discrétisées sur
l'espace des paramètres, le calcul d'une constante suppose une
localisation de maille et une interpolation. Lorsque la table, de
très grande dimension, est répartie sur la mémoire distribuée
d'un calculateur parallèle il est nécessaire de définir une
stratégie d'accès. Le travail a consisté à modéliser deux
statégies pour satisfaire les requêtes à des parties de la table
situées dans un autre processeur : envoi de la requête au
processeur qui possède les données ou demande d'envoi des
données. Deux modélisations reposent sur des modèles markoviens.
Deux autres simulations utilisent le logiciel QNAP. L'étude
permet de délimiter les conditions qui entraînent la supériorité
d'une approche sur l'autre.
Participants : François Bodin (Caps), Stéphane Chauveau ,
Bernard Philippe
L'exécution des programmes Matlab par un interpréteur de code permet une grande souplesse d'utilisation mais réduit les performances par rapport à l'exécution de programmes compilés. Notre approche consiste à exploiter le haut niveau de description du langage Matlab pour le compiler et le paralléliser efficacement. Ce projet est issu d'une collaboration entre l'équipe Caps pour la partie compilation et l'équipe Aladin pour l'utilisation des bibliothèques numériques parallèles. Une deuxième version du compilateur permettant de compiler et d'exécuter la majorité des codes Matlab est maintenant disponible.
Précédent : Logiciels Remonter : Rapport activite 1997 Suivant :
Actions industrielles