Polynôme caractéristique

  1. Créer une matrice A de taille n:=10 à coefficients dans \([-20,20]\). Elle servira pour toute la suite. On pourra faire varier n ou considérer la matrice Ap = A%1009 pour changer de cadre.

    solution ⇲ ⇱ solution
    n := 10;
    A := matrix(n,n,(j,k)->rand(41)-20);
    

Méthode de Lagrange

  1. Calculer la liste des déterminants des matrices \(kI_n-A\) pour \(0\leq k\leq n\).

    solution ⇲ ⇱ solution
    ly := seq( det(k-A), k, 0, n );
    
  2. En déduire le polynôme caractéristique de \(A\) en utilisant la fonction lagrange.

    solution ⇲ ⇱ solution
    lx := seq( k, k, 0, n );
    normal(lagrange(lx,ly));
    

Méthode de Le Verrier (sommes de Newton)

Fonctions symétriques élémentaires et sommes de Newton

Soient \(x_1,\dots x_n\) des indéterminées.

Il existe plusieurs familles d’expressions invariantes par permutation des \(x_i\), en particulier

  • les fonction symétriques élémentaires

    \[\sigma_k(x_1,\dots x_n) = \sum_{1\leq i_1< i_2<\dots i_k\leq n} x_{i_1}x_{i_2}\dots x_{i_k}\]

    qui tirent leur nom d”«élémentaires» du fait que tout polynôme \(f(x_1,\dots x_n)\) invariant par permutation des \(x_i\) s’exprime en fonction des \(\sigma_k\) ;

  • les sommes de Newton

    \[S_k(x_1,\dots x_n) = \sum_{i=1}^n x_i^k.\]

Ces familles ont une expression élégante en terme de série génératrice : si l’on considère le polynôme

\[P(X) = \prod_{i=1}^n (X-x_i)\]

les fonction symétriques élémentaires sont (à un signe près) les coefficients du polynôme développé

\[P(X) = \sum_{k=0}^n (-1)^k \sigma_k(x_1,\dots x_n) x^{n-k}\]

tandis que les sommes de Newton sont données par le développement en série

\[-\frac{\tilde P'(X)}{\tilde P(X)} = \sum_i \frac{x_i}{1-x_iX} = \sum_{k\geq 0} S_{k+1} X^k\]

\(\tilde P(x)=x^nP(\frac1x)\) est le polynôme réciproque de \(P\).

  1. Vérifier les identités ci-dessus sur l’exemple du polynôme P:=(x-a)(x-b)(x-c)(x-d).

  2. On pose \(P(x)=x^5-6x^4-7x^3+10x-9\). Déterminer la liste de ses racines complexes.

    solution ⇲ ⇱ solution
    P:=x^5-6x^4-7x^3+10x-9
    R:=col(roots(P),0)
    
  3. Calculer les 15 premières sommes de Newton de \(P\) à l’aide des racines complexes. Que constate-t-on sur la précision, pourquoi ?

    solution ⇲ ⇱ solution
    SP:=[ seq( sum(R[j]^k,j,0,4)  ,k,0,15) ]
    

    On voit que la précision absolue diminue au fur et à mesure : au début on est proche de la valeur entière, mais les sommes sont de plus en plus grandes et à la fin on ne peut en connaître qu’une valeur approchée. Cela vient du fait que \(P\) a une «grande» racine qui absorbe le comportement des autres.

  4. Déterminer le polynôme réciproque \(Q(X)=\tilde P(X)\), puis calculer les quinze premières sommes de Newton de \(P\) à l’aide d’un développement en série. Qu’en est-il de la stabilité numérique ?

    solution ⇲ ⇱ solution
     := normal(x^5*subst(P,x,1/x))
    series(diff(Q,x)/Q,x,0,15)
    

    À présent les sommes sont calculées comme des entiers exacts.

La série génératrice des sommes de Newton sert aussi dans l’autre sens pour retrouver \(P\) : en effet l’équation différentielle se résout en

\[ \begin{align}\begin{aligned}\frac{\tilde P(x)'}{\tilde P(x)} = -\sum_{k\geq 0} S_{k+1} x^k\\\Leftrightarrow\\\tilde P(x) = \exp( -\sum_{k\geq 1} S_{k}\frac{x^k}k )\end{aligned}\end{align} \]
  1. Retrouver \(P\) à l’aide de ses 6 premières sommes de Newton:

    SP := 6+50x+342x^2+2362x^3+16551x^4+115394x^5;
    

    On utilisera int( ,x) pour intégrer des polynômes, series( ,x) pour obtenir un développement en serie, et convert(  ,polynom) pour passer d’une série à un polynôme.

    solution ⇲ ⇱ solution
    ISP := int(SP,x)
    SQ := series(exp(-ISP),x)
    Q := convert(SQ,polynom)
    

Méthode de Le Verrier

Soit \(A\in M_n(\k)\), et \(\chi_A(x)\) son polynôme caractéristique.

  1. Montrer que les sommes de Newton de \(\chi_A\) sont données, pour tout \(k\geq 0\), par

    \[S_k(\chi_A) = \trace(A^k)\]
  2. Calculer les dix premières sommes de Newton de \(\chi_A\).

  3. En déduire un algorithme de calcul du polynôme caractéristique. Quelle est sa complexité, lors de quelle étape est-elle atteinte ?

    solution ⇲ ⇱ solution

    On peut écrire un algorithme en deux étapes :

    • calculer \(n\) sommes de Newton, pour une complexité de \(n\) multiplications matricielles, soit \(O(n^{1+2.81})\) (avec Strassen).
    • reconstruire le polynôme caractéristique, en complexité \(n^2\)

    de sorte que c’est l’étape de calcul des sommes qui est la plus lourde et donne la complexité de l’algorithme.

    charpoly_newton(A) = {
      my(d = matsize(A)[1]);
      my(B = matid(d));
      my(S = x*Ser(vector(d,i,B*=A;-trace(B)/i)));
      polrecip(Pol(exp(S)));
    };
    

    Remarque: il est possible de calculer les traces des \(A^k\) plus rapidement à l’aide d’une stratégie «pas de bébé-pas de géant» : en effet, si l’on décompose \(A^k = A^{aq+r}\), où \(a=\floor{\sqrt{n}}\), la trace de \(A^k\) nécessite \(O(n^2)\) opérations (seuls les coefficients de la diagonale nous intéressent).

  4. Programmer cet algorithme, vérifier en pratique sa complexité.

  5. Comparer avec l’algorithme ci-dessus pour des matrices à coefficients dans \(\Z/p\Z\) (la matrice Ap).

Méthode de Faddeev

Soit \(A\in M_n(\k)\), on définit \(M=x I_n -A\in M_n(\k[x])\) de sorte que \(\det(M)=\chi_A(x)\) soit le polynôme caractéristique de \(A\).

On dispose de deux équations :

  • d’une part l’expression du produit matrice-comatrice :

    \[M \times \Com(M)^T = \det(M).I_n\]
  • de l’autre la dérivée du déterminant, obtenue en utilisant la multilinéarité (on somme les dérivées partielles selon chaque colonne) :

    \[\frac{\d}{\d x}\det(M) = \trace( \Com(M) )\]

Si l’on décompose ces deux équations dans \(M_n(\k)[x]\) en introduisant les coefficients \(c_k\in\k\) et \(A_k\in M_n(\k)\)

\[ \begin{align}\begin{aligned}\begin{split}\det(M) &= \chi_A(x) = \sum_{k=0}^n c_k x^{n-k} \\\end{split}\\\Com(M) &= \sum_{k=0}^{n-1} C_k x^{n-1-k}\end{aligned}\end{align} \]

on obtient les formules de récurrence suivantes

\[ \begin{align}\begin{aligned}c_k &= \frac{\trace(C_k)}{n-k}\\C_k &= AC_{k-1} + c_k I_n\end{aligned}\end{align} \]

en partant de \(C_0=I_n\).

  1. Utiliser cette méthode pour calculer \(\chi_A(x)\).

    solution ⇲ ⇱ solution

    Puisque les expressions sont imbriquées (\(C_k\) dépend de \(c_k\) qui dépend de \(C_k\)), on résout la première ligne en \(c_k=-\trace(AC_{k-1})/k\) ce qui invite à décomposer l’itération de la manière suivante

    Faddev(A,x) := {
      local n, c, C, P, k;
      n := size(A);
      C := identity(n);
      P := x^n;
      for(k:=1;k<n;k++) {
        C := A*C;
        c := -trace(C)/k;
        C := C+c;
        P += c*x^(n-k);
      }
      return(P);
    }
    

Note

exercice

La technique de Faddeev et celle de Le Verrier sont en réalité rigoureusement identiques si l’on retrouve itérativement les fonctions symétriques à partir des sommes de Newton.

La vision «Le Verrier» des choses est plus globale et se prête mieux à l’adoption d’algorithmes rapides pour la calcul simultané des traces de \(A^k\) et l’obtention de \(\chi_A\) à partir de ses sommes de Newton.

De son côté, la présentation «Faddeev» fournit gratuitement l’inverse de \(A-XI_n\).