Polynôme caractéristique

[pascal16]

Wello,
est-ce que quelqu'un pourrait vérifier le polynôme caractéristique de cette matrice :


c'est pour vérifier un algo

normalement, c'est :
x⁶ + 23 x⁵ + 2 x⁴ - 3473 x³ - 45444 x² - 141872 x¹ + 374788

[Pisigma]

Bonjour,

sauf erreur d'encodage de ma part, c'est la réponse fournie par https://www.dcode.fr/polynome-caracteristique-matrice

mais tu as peut-être utilisé le même logiciel

[aviateur]

Oui c'est bon, j'ai vérifié. Est ce que ton code peut calculer le polynôme car. d'une matrice
formelle de taille 14? C'est à dire que par exemple tu remplaces au moins un nombre par une lettre.

[pascal16]

Merci à tous les deux
Cool pour dcode, j'avais pas vu le petit encart qui affiche la taille de la matrice, je pensais qu'il était limité à 3.

j'ai juste fait un algo "det(xId-M)", j'ai transformé la matrice réelle en matrice de polynômes de degré 1 et développé par rapport à la première colonne de façon récursive pour constituer le polynôme.
La classe faite n'a pas de limite, 14*14 sans problème, mais sans variable, c'est pas du calcul formel. Ca serait quoi le but ? Ca serait assez chronophage de repenser toute la structure.

Au passage avec la transformée LU->U*L->...->LU, j'ai bien retrouvé 4 valeurs propres réelles.

On peut aller bien plus haut et bien plus vite en utilisant le calcul par bloc 14*14 revient 2 det 7*7 simples et 2 det avec un x en 7*7

[Ben314]

Salut,

On peut aller bien plus haut et bien plus vite en utilisant le calcul par bloc 14*14 revient 2 det 7*7 simples et 2 det avec un x en 7*7

Tu peut détailler un peu ton truc ?
Parce que si tu compte utiliser un tel "calcul par bloc" sur une matrice quelconque, j'ai un peu des doutes...

[aviateur]

j'ai juste fait un algo "det(xId-M)", j'ai transformé la matrice réelle en matrice de polynômes de degré 1 et développé par rapport à la première colonne de façon récursive pour constituer le polynôme.
La classe faite n'a pas de limite, 14*14 sans problème, mais sans variable, c'est pas du calcul formel. Ca serait quoi le but ? Ca serait assez chronophage de repenser toute la structure

Rebonjour, le but de ma question était en fait caché.
Parce que calculer le polynôme caractéristique d'une matrice carrée c'est quelque part calculer le déterminant
d'une matrice formelle puisque tu introduis la variable -x.
Mais alors quand on écrit un algorithme parfois certains oublient le temps de calcul.
Ce temps de calcul est souvent une fonction croissante de n ( ici taille de la matrice).
Il se peut que ce temps augmente trop vite de sorte que pour des valeurs de n à peine grande le temps
le temps dépasse l'entendement.
Justement est-ce qua tu as tenu compte de cela. J'ai pris comme exemple n=14. Mais n=20?
Autrement si tu donnes à manger à ton algorithme une matrice de taille de cette ordre (n,14,8,20.?) , c'est quoi la réponse.

[pascal16]

Grand principe fondamental de la programmation :
-> bien, faire le choses pour qu'elles marchent
-> PUIS optimiser (le code optimisé n'est pas forcément le plus rapide pour les cas simples ni le plus facile à débogger).
Et on peut alors tester ses optimisation avec un base qui donne des valeurs.

Le but étant de manier des classes, et comme l'interface s'arrête à 8*8, je ne suis pas aller plus loin, et comme il est fait, il est lent,(instantané pour 8*8, 3 à 5 sec pour 10*10..) plus lent qu'un gauss (instantané pour n=30).

Comme beaucoup de gens ici, j'aime bien trouver l'algo qui va bien, comme pour du Mandelbrot que j'ai fait il y a 6 mois.

[pascal16]

pour les matrice de plus grande taille, une optimisation :
on extrait des lignes car l'accès est très lent, le découpage en sous-matrice trop chronophage.
on transforme la matrice par pivot Gauss en
AB
0C
on passe de une matrice de taille 2n à deux de taille n.
pour un calcul de det, quand n est impaire on peut rajouter une ligne et une colonne de 0 avec un 1 sur sur la diagonale

[aviateur]

Bonjour
Quand tu transformes ta matrice en pivot de Gauss sous la forme que tu dis.
Si je comprends bien tu fais un calcul par bloc: Ta matrice est de taille 2n par exemple,
et elle s'écrit par blocs:
M={{m1,m2},{m3,m4}} et tu arrives à un blocs P={{A,B},{0,C}}
Admettons que tu n'as pas eu de problème du genre pivot nul.
Mais alors ton algo te permet de calculer le polynôme caractéristique de P,
mais comment fais-tu pour retrouver celui de M?

[pascal16]

Gauss, je l'utilise de façon classique seulement pour l'ordre de la matrice.
La méthode de mettre les grosses matrices en 4 sous-matrices dont 1 nulle peut servir pour les déterminants classiques.
Si tu veux retrouver M, l'algo LU convient mieux, mais à priori, on peut générer pour chaque élément un polynome de degré n.

pour le polynôme caractéristique, j'utilise le développement par la première colonne
la matrice est transformée en matrice de polynômes d'ordre 1 (X.Id-M)
je ne multiplie donc les polynomes en retour de l'ordre n-1 que par une constante ou un terme ax+b.
je connais à l'avance le degré max du polynome de l'étape n.

[Edit] j'ai tenté de modifier mon pivot de Gauss pour s'adapter au polynome caractéristique, mais il faudrait que je réécrive toute la classe et les méthodes pour passer en matrice de polynome (dont on connais le degré max, c'est faisable dans la journée).

test de vitesse pour le polynôme caractéristique ( l'interface graphique fait de 8*8, c'est instantané)
n=10, développement itératif par rapport à la première colonne -> 2 sec
n=11 -> 24sec
n=11, quelques optimisations et un parallel.for -> 11 secondes
n=20, par Gauss, buggé, -> 2secondes, soit vers 5 sec avec une réécriture complète de la classe

[pascal16]

Polynôme caractéristique par algo de Le Verrier
instantané pour matrice 20*20, 1 seconde pour matrice 50*50.
Programme fait en quelques dizaines de minutes car la classe matrice était déjà faite.

Code: Tout sélectionner

/// <summary> Polynome caractéristique, Algo de Le Verrier </summary>
      public double[] Pol_Car_LeVerrier()
      {
         //Pour A ∈ Mn(K), on pose : A0 = A et ∀k ∈ { 1, . . . , n}          Ak = A × (Ak−1 − 1/k Tr(Ak−1)·In)
         //χA(X) = X ^ n − somme de 1 à n(1 / k Tr(Ak−1) · X ^ (n−k)
         int n = GetnLignes();
         double[] resultat = new double[n + 1]; // indices de 0 à n
         Matrice A = new Matrice(this);
         Matrice A0 = new Matrice(this);
         Matrice I = new Matrice(Identite(n));
         resultat[n] = 1;
         for (int k = 1; k <= n; k++)
         {
            resultat[n - k] = -(1 / (double)k) * A.Trace();
            A = Multiplier(A0, Soustraire(A, I.Multiplier(A.Trace() / (double)k)));
         }
         return resultat;
      }

[aviateur]

Ok, alors j'ai une une question. Si tu as une matrice à coefficients entiers, alors les coefficients du pol. caract. sont entiers.
Et je vois que dans ton algorithme tu fais au moins une division (1/k Trace(A^{k-1}..). Il me semble alors que la machine fait une petite erreur et à la fin est-ce que tu retrouves les coefficients entiers exactement?
Même question si tu as une matrice de rationnels, est-ce que à la fin le résultat est exact ou bien approché?

[Ben314]

Salut,
Je sais pas ce que c'est comme langage, mais le type double, en général, c'est un type flottant, donc forcément c'est du calcul approximatif.
Et ça pose parfois problème dans les cas de long calculs itératifs, par exemple (et justement) dans le cas du calcul de déterminant. Par exemple, l'algo. du Pivot de Gauss est connu pour être "peu stable" (i.e. gros risque d'une erreur assez importante sur le résultat par cumul des erreurs faites à chaque opération).
Concernant l'algo. de Le Verrier, je sais pas si c'est "relativement stable" ou pas.
De mémoire, il me semble qu'en passant par la décomposition LU, c'est "assez stable", mais ça reste à vérifier.

[aviateur]

Justement c'était un peu le but de ma question. Car il me semble bien que la méthode de Leverrier est très instable.

[Ben314]

Je viens de trouver ça sur le Net :
https://www-fourier.ujf-grenoble.fr/~pa ... c/algo.pdf
où ils parlent des méthodes numériques pour calculer déterminant et polynôme minimal/caractéristique dans les pages 292 et suivantes et ils parlent de la méthode de Leverrier-Faddeev-Souriau page 304 et suivantes mais sans indiquer si c'est "relativement stable" ou pas...

[aviateur]

Bonjour
C'est instable. Je l'ai vérifié sur un exemple avec n=15.
La matrice A est une matrice de rationnels crée aléatoirement:
1. Voici les coefficients calculés exactement puis arrondi avec une précision de 8 chiffres après la virgule.
{1., -0.402778, -2.18306, 5.61796, -14.2378, 23.5361, -30.9791, \
80.3896, 313.056, -580.252, 543.516, -536.719, 324.124, -785.548, \
651.943, -278.624}
2. voici le même algo mais à chaque étape le coefficient est arrondi à 8 chiffres après la virgule
{1, -0.40277778, -2.1830633, 5.617959, -14.23782, 23.53609, -30.9791, \
80.390, 313.06, -580.3, 544., -5.4*10^2, 3.*10^2, -0.*10^2, 0.*10^2,
0.*10^3}
à partir du 7ème coefficient on voit bien l'instabilité.

Pour info voici la première ligne de A.

[pascal16]

c'est du C#
double = réel double précision

Le pivot de Gauss : le choix du pivot le plus grand possible 'stabilise' le calcul.

[aviateur]

Et bien, si je passe en double précision (-i.e 10^(-16)) pour n=15 cela tient la route. Mais pour n=20, cela ne va plus.
L'instabilité ne permet pas de considérer de grandes matrices.

[pascal16]

on peut créer une classe de calcul de rationnels et faire une matrice de rationnels au résultat exact.
le calcul par les entiers peut donner aussi des résultats astronomiques sur de grandes matrices
cf ma matrice qui en est déjà à 300 000 et dons les arrondis en double se sont tous très bien passés, et pour une matrice 50*50, avec es entiers entre -10 et 10, je monte à 1E54 .