Discussion :

[membreComplexe12]

Bonjour tous,

j'aimerai une bonne fois pour toute comprendre l'intérêt d'un algo itératif plutôt qu'un algo direct (dans le cas général).

Voici où j'en suis de ce côté là :

1°) algo direct (je connais pas mal):
- Gauss très lent O(n^3)
- LU et cholesky : décomposition lente mais ensuite autres résolution rapides
- pour un algo direct il faut stocker des matrices : tout d'abord ma matrice
à transformer mais aussi la matrice de gauss qui permet d'éliminer les termes que l'on ne souhaite pas.
- on a une solution direct, c.a.d sans critère de tolérance sur un résidu ce qui est bien mais du coup on ne maitrise pas la précision du résultat...? (cette précision est dépendante du conditionnement)

2°) algo itératif (je connais pratiquement pas):
- rapides (plus que Gauss en tout cas)
- on fixe la tolérance sur le résultat et cette précision est indépendante du conditionnement
- il n'y a pas de stockage de matrice à effectuer ??? C'est quelque chose que j'ai déjà entendu dans pas mal de cours mais ça peut parait pas possible : il faut bien au moins stocker la matrice de départ ?

Conclusion :
A part cette histoire de stockage mémoire que je n'ai pas compris il y a à mon avis aucun intérêt à stocker des solveur direct car ils sont lent et en plus on ne maitrise pas forcement la précision comme avec un algo itératif ?

pourriez vous m'éclairez un peu sur tout ceci ? merci

[Aleph69]

Bonjour,

- Gauss très lent O(n^3)

Ca dépend de la matrice!

- LU et cholesky : décomposition lente mais ensuite autres résolution rapides
- pour un algo direct il faut stocker des matrices : tout d'abord ma matrice
à transformer mais aussi la matrice de gauss qui permet d'éliminer les termes que l'on ne souhaite pas.

Il y a 2 cas :
- soit on stocke la factorisation pour résoudre plusieurs systèmes linéaires avec la même matrice,
- soit on ne stocke pas la factorisation et on applique simultanément l'élimination et la remontée-descente pour calculer la solution d'un seul système linéaire.

- on a une solution direct, c.a.d sans critère de tolérance sur un résidu ce qui est bien mais du coup on ne maitrise pas la précision du résultat...? (cette précision est dépendante du conditionnement)

Non, c'est beaucoup plus compliqué. Il faut te renseigner sur les notions de "backward stability" et "backward analysis", sur le raffinement itératif.

2°) algo itératif (je connais pratiquement pas):
- rapides (plus que Gauss en tout cas)

Non, pourquoi penses-tu cela?

- on fixe la tolérance sur le résultat et cette précision est indépendante du conditionnement

on fixe la tolérance sur le résidu calculé mais l'erreur sur la solution dépend de cette tolérance et du conditionnement.

- il n'y a pas de stockage de matrice à effectuer ??? C'est quelque chose que j'ai déjà entendu dans pas mal de cours mais ça peut parait pas possible : il faut bien au moins stocker la matrice de départ ?

Tu n'es pas obligé de stocker la matrice tant que tu sais calculer les produits entre cette matrice et des vecteurs sans la construire. Mais tu peux tout de même avoir à stocker des matrices (arnoldi). Voir le théorème de Faber-Manteuffel pour des précisions théoriques.

A part cette histoire de stockage mémoire que je n'ai pas compris il y a à mon avis aucun intérêt à stocker des solveur direct car ils sont lent et en plus on ne maitrise pas forcement la précision comme avec un algo itératif ?

Non, pas du tout. On ne peut jamais maîtriser la précision sur la solution sans connaître le conditionnement du système. On n'est pas obligé de stocker la factorisation. La complexité d'un solveur dépend du système à résoudre. Il existe des systèmes favorables aux solveurs directs pour les temps de calcul ou la précision des calculs. Il en existe aussi qui sont favorables aux solveurs itératifs. Comparer des solveurs c'est comparer leurs complexité temporelle, complexité spatiale et stabilité numérique. Se restreindre à un seul de ces paramètres n'a pas de sens.

[membreComplexe12]

salut et merci de ton aide!

Ca dépend de la matrice!

1) ouahou, je suis très étonné là....
le nombre d'opérations de la triangulation ne dépend pas de la forme de la matrice donc je ne vois pas pourquoi tu dis ceci
ps: dans ma discussion je considère que le cas où on a une matrice pleine. Je ne considère par le cas de matrice triadiagonal ou diagonale car dans ces cas je comprends bien que ça sera plus rapide.

Non, c'est beaucoup plus compliqué. Il faut te renseigner sur les notions de "backward stability" et "backward analysis", sur le raffinement itératif.

2) je suis un peu étonné car d'après les cours de méthode numérique que j'ai suivi c'est bien le conditionnement qui défini la precision que l'on aura sur la réponse...

Non, pourquoi penses-tu cela?

3) en fait ce qui me fait pensez à ceci c'est que Gauss c'est en O(n^3) et par exemple le gradient conjugué c'est en O(n.log(n)) du coup la compléxité est beaucoup plus faible et donc on va aller plus vite pour trouver la solution puisqu'il y a moins d'opérations à faire.... ?

on fixe la tolérance sur le résidu calculé mais l'erreur sur la solution dépend de cette tolérance et du conditionnement.

4) je ne comprends pas.
Si j'ai comme résidu :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

{R}=[A]{x}-{b}

et que je trouve à la fin de calcul que chaque composantes du résidu est < à 1e-9 ça veut dire que les composantes de x sont précises à 1e-9 près ?
je ne comprends pas pourquoi le conditionnement interviendrais ici...

Tu n'es pas obligé de stocker la matrice tant que tu sais calculer les produits entre cette matrice et des vecteurs sans la construire. Mais tu peux tout de même avoir à stocker des matrices (arnoldi). Voir le théorème de Faber-Manteuffel pour des précisions théoriques.

5) OK, merci pour les infos.

On ne peut jamais maîtriser la précision sur la solution sans connaître le conditionnement du système.

6) our un algo direct je comprends mais pour un algo itératif je ne comprends pas pourquoi le conditionnement à une influence sur le résultat (cf. ma réponse n°4)

La complexité d'un solveur dépend du système à résoudre. Il existe des systèmes favorables aux solveurs directs pour les temps de calcul ou la précision des calculs. Il en existe aussi qui sont favorables aux solveurs itératifs.

7) en fait comme je l'ai précisé dans la réponse n°1 de ce message je cherche à comprendre la différences lorsqu'on a une matrice pleine.

Si la matrice est tridiagonal ou diagonal je comprends qu'un solveur direct sera plus rapide mais je pense que dans tout les autres cas les solveurs itératifs tel que le gradient conjugué seront plus rapides ? (sauf dans le cas où on a un très mauvais conditionnement mais dans ce cas il suffit d'utiliser un préconditionneur et le problème est réglé?)

Comparer des solveurs c'est comparer leurs complexité temporelle, complexité spatiale et stabilité numérique.

aïe
je ne connais pas la différence entre complexité temporelle et spaciale... pour moi il y a qu'une seule complexité : le nombre d'opérations à effectuer

pour la stabilité numérique je n'y avais pas pensé.....

Question générale :
En fait ce que je comprends de cette discussion c'est que si on ne fait pas avant le calcul de la solution :
- le calcul du conditionnement
- la vérification de la stabilité du système

alors on ne s'aura jamais si notre solution sera correcte ou non ?

du coup, on préconditionneur est quelque chose qui devrait être obligatoire pour tout solveur ?

[Aleph69]

1) ouahou, je suis très étonné là....
le nombre d'opérations de la triangulation ne dépend pas de la forme de la matrice donc je ne vois pas pourquoi tu dis ceci
ps: dans ma discussion je considère que le cas où on a une matrice pleine. Je ne considère par le cas de matrice triadiagonal ou diagonale car dans ces cas je comprends bien que ça sera plus rapide.

Ca dépend comment tu vois les choses. Toi, tu as l'air de fixer un algorithme ne comportant aucune instruction conditionnelle (autre qu'un test de division par zéro) que tu appelles triangulation. Dans ce cas, effectivement, le nombre d'opérations est indépendant de la matrice. C'est comme cela qu'on présente les méthodes de factorisation dans les livres.

En revanche, on ne les code jamais ainsi. On ajoute plein d'instructions conditionnelles sur la structure et les valeurs des matrices. Par exemple, dans l'élimination de Gauss, tu peux ajouter le test idiot consistant à regarder si le pivot vaut 1 pour t'épargner plusieurs divisions. Avec un simple ajout comme celui, la complexité de ton algorithme devient alors dépendante de tes données.

Il faut bien comprendre que ce qu'on appelle la méthode de Gauss était à l'origine effectué à la main sur des petits systèmes et tu imagines bien que les gens réfléchissaient un minimum en l'appliquant : ils ne vont pas faire des divisions par 1!!! Dès lors, dire que la version académique naïve que tu considères correspond à la méthode de Gauss est une affirmation contestable et de facto contestée par les implémentations effectuées en pratique.

2) je suis un peu étonné car d'après les cours de méthode numérique que j'ai suivi c'est bien le conditionnement qui défini la precision que l'on aura sur la réponse...

A quel résultat penses-tu?

3) en fait ce qui me fait pensez à ceci c'est que Gauss c'est en O(n^3) et par exemple le gradient conjugué c'est en O(n.log(n)) du coup la compléxité est beaucoup plus faible et donc on va aller plus vite pour trouver la solution puisqu'il y a moins d'opérations à faire.... ?

Je pense que tes deux résultats de complexité sont faux et par conséquent la conclusion que tu en tires aussi. Je suis très curieux de savoir où tu as lu cette complexité pour les gradients conjugués...

4) je ne comprends pas.
Si j'ai comme résidu :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

{R}=[A]{x}-{b}

et que je trouve à la fin de calcul que chaque composantes du résidu est < à 1e-9 ça veut dire que les composantes de x sont précises à 1e-9 près ?
je ne comprends pas pourquoi le conditionnement interviendrais ici...

Non, pas du tout. On a Ax=b. On note cond(A)=||A||*||A^{-1}|| le conditionnement de A, y la solution calculée, e=x-y l'erreur et r=Ae=b-Ay le résidu. Alors ||x||>=||b||/||A|| et ||e||<=||A^{-1}||*||r||. Par conséquent,

||e||/||x||<=cond(A)*||r||/||b||,

c'est-à-dire que l'erreur relative sur la solution est majorée par le produit du conditionnement de A et la norme relative du résidu calculé.

6) our un algo direct je comprends mais pour un algo itératif je ne comprends pas pourquoi le conditionnement à une influence sur le résultat (cf. ma réponse n°4)

A cause de la réponse précédente qui est indépendante de la méthode de résolution choisie.

7) en fait comme je l'ai précisé dans la réponse n°1 de ce message je cherche à comprendre la différences lorsqu'on a une matrice pleine. Si la matrice est tridiagonal ou diagonal je comprends qu'un solveur direct sera plus rapide mais je pense que dans tout les autres cas les solveurs itératifs tel que le gradient conjugué seront plus rapides ? (sauf dans le cas où on a un très mauvais conditionnement mais dans ce cas il suffit d'utiliser un préconditionneur et le problème est réglé?)

Oui mais même dans ce cas les gradients conjugués peuvent converger en plus de n itérations à cause des erreurs d'arrondis et ainsi nécessiter plus de calculs que l'algorithme de Cholesky par exemple. On peut utiliser un préconditionneur mais il faut alors aussi compter le nombre d'opérations nécessaire à sa construction et le nombre d'opérations nécessaires à son application au système. A la fin, tu n'es pas forcément gagnant.

je ne connais pas la différence entre complexité temporelle et spaciale... pour moi il y a qu'une seule complexité : le nombre d'opérations à effectuer

La complexité temporelle correspond au nombre d'opérations et la complexité spatiale au nombre maximal d'informations stockées simultanément durant l'exécution. En l'absence d'instructions conditionnelles (hors tests d'erreur), ces deux complexités sont constantes. Sinon, il faut distinguer les complexités moyenne, dans le pire cas et dans le meilleur cas. C'est comme pour les algorithmes de tri.

Question générale :
En fait ce que je comprends de cette discussion c'est que si on ne fait pas avant le calcul de la solution :
- le calcul du conditionnement
- la vérification de la stabilité du système

alors on ne s'aura jamais si notre solution sera correcte ou non ?

C'est la stabilité du solveur, pas du système. Il faut effectivement connaître le conditionnement et la backward error (mesure de la stabilité) pour pouvoir contrôler une borne supérieure d'une erreur relative. Mais, cela ne te donne pas l'erreur : tu sais juste que cette erreur est majoré par un certain nombre que tu as estimé numériquement. Pour connaître l'erreur, il faut connaître la solution.

du coup, on préconditionneur est quelque chose qui devrait être obligatoire pour tout solveur ?

Ce n'est pas obligatoire mais fortement conseillé et ça dépend bien sûr de la précision de ta machine.

[membreComplexe12]

Par exemple, dans l'élimination de Gauss, tu peux ajouter le test idiot consistant à regarder si le pivot vaut 1 pour t'épargner plusieurs divisions.

Excellent
je n'y avais jamais pensé
as tu d'autres astuces dans ce genre qui permettent de gagner du temps en général ?

A quel résultat penses-tu?

je pensais à |dx/x|<=cond.|db/b|

Je pense que tes deux résultats de complexité sont faux et par conséquent la conclusion que tu en tires aussi. Je suis très curieux de savoir où tu as lu cette complexité pour les gradients conjugués...

j'avais vu ceci sur un lien internet (peut etre pas super correct....)
c'etait pour le gradient conjugué preconditionné

Non, pas du tout. On a Ax=b. On note cond(A)=||A||*||A^{-1}|| le conditionnement de A, y la solution calculée, e=x-y l'erreur et r=Ae=b-Ay le résidu. Alors ||x||>=||b||/||A|| et ||e||<=||A^{-1}||*||r||. Par conséquent,

ouahou, super !!!!!!!!!!!
je n'avais jamais compris ceci, je suis super content d'avoir compris à présent !!

Oui mais même dans ce cas les gradients conjugués peuvent converger en plus de n itérations à cause des erreurs d'arrondis et ainsi nécessiter plus de calculs que l'algorithme de Cholesky par exemple.

Oui, je comprends.
j'avais lu une fois un document qui disais que les directions peuvent ne plus être conjugués à cause des erreurs d'arrondi et du coup il faut de temps en temps (toute les x itérations) réinitialiser la direction de recherche où un truc dans ce genre....

On peut utiliser un préconditionneur mais il faut alors aussi compter le nombre d'opérations nécessaire à sa construction et le nombre d'opérations nécessaires à son application au système. A la fin, tu n'es pas forcément gagnant.

je comprends, tu as tout à fait raison

La complexité temporelle correspond au nombre d'opérations et la complexité spatiale au nombre maximal d'informations stockées simultanément durant l'exécution. En l'absence d'instructions conditionnelles (hors tests d'erreur), ces deux complexités sont constantes. Sinon, il faut distinguer les complexités moyenne, dans le pire cas et dans le meilleur cas. C'est comme pour les algorithmes de tri.

j'ai un peu de mal avec ceci, je fouillerai un peu ceci un de ces 4

C'est la stabilité du solveur, pas du système.

aie, j'ai un peu de mal là car pour moi les deux sont liés...
par exemple si je prends des equadiff (EDO) la stabilité du schéma dépend du solveur mais aussi de l'EDO :
- si j'ai une EDO pas du tout raide il se peut que mon schéma soit stable
- mais si je prends le meme solveur mais une EDO raide alors il se peut que ça ne soit pas stable

Pour connaître l'erreur, il faut connaître la solution.

exact

Ce n'est pas obligatoire mais fortement conseillé et ça dépend bien sûr de la précision de ta machine.

c'est vrai pour la précision je n'y avais pas pensé....
en fait ce n'est pas si critique que ceci le conditionnement :

si j'ai des variables avec une precision de 1e-14 (double) et que j'ai conditionnement énorme 1e6 alors j'aurais quand même une précision relative de 1e-8 ce qui est super.

ps: j'ai pas fais une erreur de raisonnement là ?

[Aleph69]

as tu d'autres astuces dans ce genre qui permettent de gagner du temps en général ?

C'était juste pour te donner un exemple. En pratique, on ne fait pas ce test car la probabilité qu'il soit vrai est extrêmement faible. Des astuces de programmation, tu peux en imaginer plein : c'est difficile de les lister.

je pensais à |dx/x|<=cond.|db/b|

D'accord. Ce résultat ne dit pas que c'est le conditionnement qui détermine la précision. Il dit que cette erreur est majoré par le produit du conditionnement et d'une "backward error" qui correspond ici à la norme relative du résidu calculé. Il est important de remarquer que le conditionnement ne dépend que des données tandis que la backward error dépend de l'algorithme de résolution.

j'avais vu ceci sur un lien internet (peut etre pas super correct....)
c'etait pour le gradient conjugué preconditionné

Alors ce devait être pour une famille bien particulière de préconditionneurs et/ou de matrices.

Oui, je comprends.
j'avais lu une fois un document qui disais que les directions peuvent ne plus être conjugués à cause des erreurs d'arrondi et du coup il faut de temps en temps (toute les x itérations) réinitialiser la direction de recherche où un truc dans ce genre....

Oui, on appelle ça des cas de "breakdown". Cela dit, pour le cas particulier des gradients conjugués appliquées au systèmes linéaires, on ne fait généralement pas de redémarrage.

aie, j'ai un peu de mal là car pour moi les deux sont liés...
par exemple si je prends des equadiff (EDO) la stabilité du schéma dépend du solveur mais aussi de l'EDO :
- si j'ai une EDO pas du tout raide il se peut que mon schéma soit stable
- mais si je prends le meme solveur mais une EDO raide alors il se peut que ça ne soit pas stable

Dans le cas présent, ton EDO est réduite à un système linéaire qui à mon avis dans ton esprit est inversible et n'admet qu'un seul point d'équilibre trivial. La théorie de Lyapunov n'est pas utile ici et n'a rien à voir avec la stabilité numérique de ton solveur (qui n'est pas un schéma).

c'est vrai pour la précision je n'y avais pas pensé....
en fait ce n'est pas si critique que ceci le conditionnement :

si j'ai des variables avec une precision de 1e-14 (double) et que j'ai conditionnement énorme 1e6 alors j'aurais quand même une précision relative de 1e-8 ce qui est super.

ps: j'ai pas fais une erreur de raisonnement là ?

Aucune erreur de raisonnement à ceci près que l'erreur relative sera au pire de 1e-8 mais elle peut très bien être à 1e-16 en réalité : c'est juste qu'on n'en sait rien. C'est là qu'intervient la précision machine. Tu considères qu'un conditionnement de 1e6 est énorme mais ramené à la précision machine ce n'est pas forcément évident. Si ta précision machine est 1e-8 c'est très mal conditionné mais si elle est à 1e-16 c'est en fait plutôt bien conditionné.

[membreComplexe12]

Bonjour Aleph69 et merci beaucoup pour ton aide

C'était juste pour te donner un exemple. En pratique, on ne fait pas ce test car la probabilité qu'il soit vrai est extrêmement faible. Des astuces de programmation, tu peux en imaginer plein : c'est difficile de les lister.

OK, je comprends

D'accord. Ce résultat ne dit pas que c'est le conditionnement qui détermine la précision. Il dit que cette erreur est majoré par le produit du conditionnement et d'une "backward error" qui correspond ici à la norme relative du résidu calculé. Il est important de remarquer que le conditionnement ne dépend que des données tandis que la backward error dépend de l'algorithme de résolution.

OKKKKKK j'ai enfin compris la nuance
merci beaucoup pour cette explication, c'est très clair à présent
=> juste une dernière question : on a moyen d'estimer cette "backward error" et son evolution au cours des itérations ?

Alors ce devait être pour une famille bien particulière de préconditionneurs et/ou de matrices.

peut être, je ne sais plus désolé

Oui, on appelle ça des cas de "breakdown". Cela dit, pour le cas particulier des gradients conjugués appliquées au systèmes linéaires, on ne fait généralement pas de redémarrage.

OK, c'est noté

Dans le cas présent, ton EDO est réduite à un système linéaire qui à mon avis dans ton esprit est inversible et n'admet qu'un seul point d'équilibre trivial. La théorie de Lyapunov n'est pas utile ici et n'a rien à voir avec la stabilité numérique de ton solveur (qui n'est pas un schéma).

Je ne comprends pas vraiment ce que tu veux dire car un système d'EDO n'est pas un système linéaire pour moi (et Lyapunov connais pas du tout)....
mais ce n'est pas grave je me pencherai sur ce genre de chose une autre fois.

Aucune erreur de raisonnement à ceci près que l'erreur relative sera au pire de 1e-8 mais elle peut très bien être à 1e-16 en réalité : c'est juste qu'on n'en sait rien. C'est là qu'intervient la précision machine. Tu considères qu'un conditionnement de 1e6 est énorme mais ramené à la précision machine ce n'est pas forcément évident. Si ta précision machine est 1e-8 c'est très mal conditionné mais si elle est à 1e-16 c'est en fait plutôt bien conditionné.

Ouf, je suis bien content, si mon raisonnement aurait été faux ça voulais dire que j'avais rien compris.
Je comprends tout à fait la remarque que tu as fais ici, je suis tout à fait d'accord!

MERCI POUR TOUT

[Aleph69]

Bonjour,

=> juste une dernière question : on a moyen d'estimer cette "backward error" et son evolution au cours des itérations ?

Oui. Tu peux lire cet article qui rappelle les résultats élémentaires de la théorie et fournit des références historiques pour les systèmes linéaires. Pour d'autres algorithmes, se référer au livre du même auteur "Accuracy and stability of numerical algorithms".

[membreComplexe12]

merci beaucoup pour tout !!!!