Méthodes itératives

Dans tnfem les systèmes linéaires sont r/'esolus par factorisation de Gauss bande. Il est trivial de rempacer le stockage bande par un stockage profil (skyline) mais les méthodes directes ne sont plus compétitives dès que le nombre d'inconnues est supérieure a 5000 environ.

Pour résoudre Ax=b, la méthode du gradient conjuguée préconditionnée par une matrice C repose sur l'algorithme suivant qui démarre de x quelconque:

D\'emarrage:    r\'esoudre Cg = (b-Ax).  
                Poser oldgnorm = 1e30, h = g.
It\'erations:   r\'esoudre Cr = Ah,
                poser rho = g.Ch / h.r,
                x <- x + rho h,
                g <- g - rho r.
                Poser gnorm = g.Cg, 
                      gamma = gnorm / oldgnorm,  
                      oldgnorm = gnorm
                h <- g + gamma h

Notations:

la flèche vers la gauche signifie "écriture du résultat dans la mémoire occupée par la variable de gauche".

Les matrices A et C doivent être définies positives.

Notez que le produit Cg s'obtient par Cg = (Cg)old - rho Ah. De même Ch =Cg + gamma (Ch)old. On aura sans doute intéret a stocker 2 vecteurs de plus pour (Cg)old et (Ch)old. Ainsi chaque itération coute un produit par A et une résolution par C.

Cet algorithme se généralise très bien aux produits par blocs et pourra donc s'appliquer au système issu de l'élasticité.

dans un premier temps on conservera le stockage bande. Ensuite on adpotera le stockage compact (Morse) pour utiliser pleinement les performances de cet algorithme.