Le solveur GMRES

Comme indiqué dans la section 2.3.1, la résolution des équations de Navier- Stokes linéarisées implique la résolution d’un système linéaire. Ce système présente plusieurs particularités. Premièrement, il est de très grande taille. Un maillage Navier-Stokes 3D standard possède une dizaine de millions de nœuds. Comme cinq inconnues (cf. 2.1) sont présentes à chacun de ces nœuds, le système linéaire possède plus de cinquante millions d’inconnues. Ensuite, ce système est issu d’une discrétisation des équations par une méthode éléments finis. Pour des éléments finis de Lagrange P1, le graphe de remplissage de la matrice [135] est identique au graphe du maillage. Le nombre de sommets voisins d’un autre définit le nombre de blocs non-nuls de chaque ligne et colonne de la matrice, et est de l’ordre de plusieurs dizaines. Cela est à comparer à la dizaine de millions de nœuds du maillage. La matrice est donc extrêmement creuse.Il existe deux grandes catégories de méthodes pour résoudre des sys- tèmes linéaires de grandes tailles : les méthodes itératives et les méthodes directes. Les méthodes itératives, comme leur nom l’indique, créent une suite d’approximations de la solution. L’utilisateur définit une tolérance de résolution, qui est la norme maximale du résidu r = b − Ax à partir de laquelle l’utilisateur est satisfait de la solution. La méthode directe calcule la solution exacte (aux erreurs numériques près) x = Ab.

Les avantages et les inconvénients de ces deux classes de méthodes sont faciles à comprendre. Les méthodes itératives se basent sur la création de suite d’approximations soit par itérés de la matrice (comme pour les méthodes basées sur des espaces de Krylov), soit par multiplication par des parties de la matrice (méthodes de relaxation type Jacobi, Gauss-Seidel et SOR). Ces méthodes utilisent donc des produits matrices-vecteurs, pour lesquels les formats de stockage peuvent être optimisés [135, chapitre 3]. Le coût en mémoire de ces algorithmes est uniquement celui de la matrice creuse et d’un certain nombre de vecteurs. Pour les problèmes abordés dans cette thèse, le nombre de vecteurs composant l’espace de Krylov est de l’ordre de la centaine. à un vecteur. En général, cela passe par une décomposition de la matrice A, par exemple une décomposition LU [10, 39]. Les facteurs de cette décomposition sont facilement inversibles. Une fois cette décomposition calculée, on peut résoudre quasi-instantanément de nombreux seconds membres différents pour une même matrice. Le temps de résolution ne dépend pas du conditionnement de la matrice κ (A) = kAkkAk. À l’inverse, certaines méthodes itératives, comme le gradient conjugué pour des matrices symétriques définies positives, ont un taux de convergence dépendant directement de κ(A). Les méthodes directes ont un important désavantage de mémoire requise. Pour un système plein, la décomposition LU coûte autant en mémoire que la matrice. Une méthode directe est alors très adaptée. Par contre, la décomposition LU d’une matrice creuse donne des facteurs L et U beaucoup plus remplis que la matrice A initiale, comme cela est présenté dans la section 4.5.1. Il existe des solveurs directs comme MUMPS [7, 9], adaptés aux matrices creuses. Le coût de stockage des facteurs pleins est alors identique à celui de la matrice. Pour conclure, la parallélisation des méthodes directes est délicate, comme le montre la parallélisation d’une méthode de préconditionnement ILU [135]. Les méthodes itératives sont très simples à paralléliser, comme cela est montré dans la section 3.2.2.

Cette section présente l’algorithme GMRES avec déflation des petites valeurs propres utilisé par Dassault Aviation pour la résolution des systèmes linéaires Navier-Stokes linéarisées.. On peut montrer [135] [63, sec. 5.2.8] qu’en arithmétique inexacte, cette méthode est sensible aux erreurs numériques qui peuvent conduire à une perte d’orthogonalité entre les vecteurs. La forme modifiée de la méthode d’Arnoldi, qui effectue la projection sur l’orthogonal de chacun des vecteurs de la base V Pour encore améliorer la précision de la méthode, on peut effectuer une deuxième passe d’orthonormalisation après la première [135]. Enfin, une dernière méthode d’orthonormalisation d’une famille de vecteurs est basée sur les réflexions d’Householder, et est encore plus précise en arithmétique inexacte que la méthode d’Arnoldi, au prix d’un coût de calcul plus élevé [63].L’algorithme 3 présente cette version du GMRES avec m itérations. Si, à la fin des ces itérations, la tolérance sur la norme du résidu n’est pas atteinte, l’algorithme est redémarré en prenant comme nouveau point de départ la solution approchée obtenue à la fin des itérations. Cela consiste à recommencer l’algorithme 3 en remplaçant les deux premières lignes par x .