\documentclass{article}

\textheight 23cm
%\textwidth 16.5cm \columnsep 10pt \columnseprule 0pt

%\usepackage[pdftex]{hyperref}
\usepackage[latin1]{inputenc}
\usepackage[T1]{fontenc}
\usepackage[francais]{babel}
\begin{document}
\newtheorem{rem}{Remarque}
\newcommand{\on}{\framebox{ON}}
\newcommand{\tr}{\mbox{trace }}
\newcommand{\enter}{\framebox{ENTER}}
\newcommand{\haut}{\framebox{$\uparrow$}}
\newcommand{\bas}{\framebox{$\downarrow$}}
\newcommand{\droit}{\framebox{$\rightarrow$}}
\newcommand{\gauche}{\framebox{$\leftarrow$}}
\def\fb{\framebox}
\def\ls{shift gauche }
\def\rs{shift droit }
\newcommand{\atan}{\mbox{atan\,}}

\title{TP 3: Alg\`ebre lin\'eaire.}
\maketitle

L'alg\`ebre lin\'eaire se pr\^ete bien au calcul sur machine : on aborde
ici les m\'ethodes de r\'eduction des endomorphismes (diagonalisation).

\section{Rappels et d\'efinitions.} 
\label{sec:scolaire}
Pour diagonaliser un endomorphisme \`a la main,
on factorise le polyn\^ome caract\'eristique det$(\lambda I-A)$
et on r\'esoud $(\lambda I-A)v=0$ pour les racines $\lambda $ trouv\'ees.
On va pr\'esenter d'autres m\'ethodes plus efficaces \`a la machine ainsi 
qu'une m\'ethode num\'erique.\\
D\'efinitions :\\
{\bf matrice caract\'eristique de  $A$} c'est $\lambda I-A$\\ 
{\bf polyn\^ome caract\'eristique de $A$} c'est det$(\lambda I-A)$=$P(\lambda)$ \\
{\bf matrice adjointe de $A$} est la matrice des cofacteurs  de la transpos\'ee
 de $A$.\\
 On a donc si $B(\lambda)$ = matrice adjointe de  $\lambda I-A$ :\\
$B(\lambda)*(\lambda I-A)=(\lambda I-A)*B(\lambda)=P(\lambda)*I $ 


\section{Logiciels} \label{sec:expe}

\subsection{MuPAD/maple} \label{sec:mupad}
Ne pas oublier de taper \verb|export(linalg)| (\verb|MuPAD|)
ou {\tt with(linalg)} (maple) dans votre session.

Quelques instructions d'alg\`ebre lin\'eaire:
\begin{enumerate}
\item {\tt charpoly(A,x)}: polyn\^ome 
caract\'eristique
\item {\tt det(A)}: d\'eterminant d'une matrice $A$,
\item {\tt scalarProduct(A,B)}/{\tt dotprod(A,B)}: 
produit scalaire des deux vecteurs $A$ et $B$,
\item {\tt eigenvalues(A)}/{\tt eigenvals(A)}: 
valeurs propres de la matrice $A$. Si vous voulez
des valeurs num\'eriques (et si $A$ est \`a valeurs num\'eriques) 
\'ecrivez au moins un des coefficients de $A$ comme nombre r\'eel ou complexe.
\item {\tt eigenvectors(A)}/{\tt eigenvects(A)}: vecteurs propres de $A$,
\item {\tt matlinsolve(A,b)}/{\tt linsolve(A,b)}: 
r\'esout l'\'equation lin\'eaire $Ax=b$,
\item {\tt f:=(i,j)->(1/(i+j-1));matrix(2,2,f)}:
pour cr\'eer une matrice (ici de Hilbert).
\item \verb|(i,j) -> if i=j then 1; else 0; end_if;|: utile pour cr\'eer
la matrice identit\'e.
\item \verb|ncols(A)|/\verb|coldim(A)|: nombre de colonnes d'une matrice,
\verb|trace(A)|: trace de $A$
\item {\tt transpose(A)}: transpos\'ee de $A$
%\item {\tt ludecomp(A)}: donne la d\'ecomposition $LU$ de $A$.
\item \verb|?linalg| pour la liste des commandes d'alg\`ebre lin\'eaire
\end{enumerate}

On rappelle \'egalement l'instruction {\tt gcdex} pour l'algorithme
de B\'ezout pour les polyn\^omes.

\subsection{Calculatrices} \label{sec:calc}
Les calculatrices formelles Casio ne disposent pas de commande
d'algebre lin\'eaire symbolique. On pr\'esente ici les commandes pour
les HP49 et TI92/89.

\subsubsection{Commandes de la HP49} \label{sec:erable}
Pour entrer une matrices $A$ vous pouvez utiliser le MatrixWriter 
(\verb|MTRW| sur le clavier). On peut aussi cr\'eer une matrice
dont l'\'el\'ement $m_{i,j}$ est donn\'e par une fonction
$(i,j) \rightarrow m_{i,j}$ \`a l'aide de \verb|LC2M|.

Pour les commandes utiles se r\'eferrer \`a ''La HP49G en r\'esum\'e''.

\subsubsection{TI92/89} \label{sec:ti}
Commandes utiles:
\begin{itemize}
\item \verb|identity(n)| et \verb|randmat(n,m)| pour cr\'eer 
des matrices identit\'e et al\'eatoire
\item  \verb|rref| : r\'eduction de Gau\ss\
\item  \verb|egv| (TI89 et 92+ seulement) : valeurs propres/vecteurs propres. Attention, le calcul est num\'erique uniquement.
Sur tous les mod\`eles on peut bien sur calculer le polynome
caract\'eristique avec \verb|det(A-xI)| ou $I$ d\'esigne la matrice
identit\'e, puis factoriser avec \verb|factor| et calculer les vecteurs
propres \`a la main en r\'esolvant le syst\`eme.
\end{itemize}


\section{R\'eduction exacte des endomorphismes par le polynome minimal} 
\label{sec:diag}


\subsection{Calcul du polyn\^ome minimal} \label{sec:calcul}
Soit $A$ une matrice $n\times n$. Nous allons d\'eterminer le polyn\^ome
minimal de $A$, c'est-\`a-dire un polyn\^ome non nul de degr\'e minimal
qui annule $A$. On consid\`ere pour cela $A$ comme un {\em vecteur}\/ 
de l'espace des matrices ayant $n^2$ coordonn\'ees (par exemple la matrice
identité en dimension 2 a pour coordonnées (1,0,0,1))
et on va chercher
par la m\'ethode du pivot de Gau\ss\ une relation non triviale
entre les puissances successives de $A$ (on peut se limiter \`a $I$,
$A$, .., $A^n$ car le polyn\^ome caract\'eristique de $A$ annule $A$).\\

{\bf Exercice 1} (à rendre à la fin de la 1ère séance de ce TP)\\
Soit la matrice~:
\[A=\left(\begin{array}{rrrr}
 0 & 1 & 1  \\
 1 & 0 & 1  \\
 1 & 1 & 0  \\
\end{array}\right)\] 
Calculer \`a la calculatrice $A^2,\ A^3$ et appliquer le pivot de Gau\ss\ 
sur les 4 lignes constitu\'ees par $I,\ A,\ A^2,\ A^3$ 
(si une colonne a des 
z\'eros sur et sous la diagonale, on passe \`a la colonne suivante en 
cherchant le pivot \`a partir de la m\^eme ligne).\\
%mettant des z\'eros sous $\alpha_{i,i+1}$).\\
 N'oubliez pas de garder une trace des op\'erations effectu\'ees en notant 
en bout de lignes l'expression des lignes ($I, A...A^3$).\\ 
Trouver le polyn\^ome minimal de $A$.\\
Trouver par la m\^eme m\'ethode mais en faisant les calculs avec un logiciel
(HP49 instruction \verb|PMINI| ou avec le programme MuPAD ci-dessous), 
le polyn\^ome minimal des matrices d\'efinies par :
\[B= \left(\begin{array}{cccc}
-63 & 70 & -12 & -19 \\
-46 & 52 & -8 & -14 \\
45 & -48 & 10 & 15 \\
25 & -26 & 4 & 13
\end{array}\right) , \quad
C= \left(\begin{array}{rrrr}
-8 &  8 & 8 & 6 \\
 5 &-2 &  0 & 5 \\
 0 & 1 & -4 & 6 \\
 -6 & -8 & -6& -6
\end{array}\right)\]

\subsection{Programme pour le polyn\^ome minimal.} 
\label{sec:mupad2}
Ce programme est écrit en {\tt MuPAD}, pour les utilisateurs de maple,
il y a quelques remplacements à faire indiqués entre parenthèses.
Attention, MuPAD utilise \verb|//| pour un commentaire (comme en C++),
en maple utilisez le signe \verb|#|.
\begin{verbatim}
// programme polynome minimal a sauver sous le nom minimal_poly.mu
// utiliser read("minimal_poly.mu") a l'interieur de MuPAD

// fonction iequalj utile pour definir la matrice identite
iequalj:=(l,k)->if (l=k) then 1; else 0; end_if; ; // (maple) fi

// definition de la procedure
// en argument, on place la matrice A dont on cherche le 
// polynome minimal, la procedure retourne la matrice reduite 
// ou on lit le polynome minimal
minimal_poly:=proc(A)
local n,Id,B,Aj;            // variables locales
begin  // (maple) pas de begin
 n:=ncols(A);   // dimension de la matrice A, (maple) coldim
 Id:=matrix(n,n,iequalj);   // matrice identite
 B:=matrix(n+1,n^2+1,0);    // cree la matrice resultat
 Aj:=Id;

 // on remplit la 1ere ligne avec les coordonnees de l'identite
 for i from 1 to n^2 do 
   B[1,i]:=Id[((i-1)div n)+1,((i-1) mod n)+1]: 
 end_for: // (maple) od

 // on remplit les lignes 2 a n+1 avec 
 // les coordonnees des puissances de A
 for j from 1 to n do
  Aj:=Aj*A;
  for i from 1 to n^2 do 
    B[j+1,i]:=(Aj)[((i-1)div n)+1,((i-1) mod n)+1]: // (maple) iquo (pour div)
  end_for; // (maple) od
 end_for; // (maple) od

// on remplit la derniere colonne avec les puissances de A
 for i from 1 to n+1 do B[i,n^2+1]:=a^(i-1): end_for; // (maple) od

// on reduit la matrice
 gaussElim(B); // (maple) gausselim(B)
end_proc: // (maple) end
\end{verbatim}

Pour utiliser cette procedure, tapez par exemple les commandes:
\begin{verbatim}
export(linalg); // (maple) with(linalg);
read("minimal_poly.mu");
A:=matrix(2,2,[[1,2],[3,4]]);
minimal_poly(A);
\end{verbatim}


\subsection{Recherche de vecteurs propres} \label{sec:propres}
On suppose que le polyn\^ome minimal $P_A$ est sans facteur carr\'e
({\em square-free}\/) ce qui se v\'erifie en cherchant le pgcd
de $P_A$ avec sa d\'eriv\'ee.
\`A l'aide du polyn\^ome minimal, il est facile de d\'eterminer les
vecteurs propres de $A$ correspondant aux racines de $P_A$ que l'on
sait d\'eterminer car si $P_A(X)=(X-\lambda_1)\times Q(X)$ on a~:
\[ \mbox{Im }Q(A) = \mbox{Ker }(A-\lambda_1 I) \]
En effet puisque  $(A-\lambda_1I)\times Q(A)=P_A(A)=0$ on a~ :
\[ \mbox{Im }Q(A) \subset \mbox{Ker }(A-\lambda_1 I) \]
$(X-\lambda_1)$ et $Q(X)$ sont premiers entre eux car $P_A$ est sans facteur carr\'e donc d'apr\`es B\'ezout,
 il existe deux polyn\^omes $U$ et $V$ tels que~:
\[ U(X)\times (X-\lambda_1)+Q(X)\times V(X)=1\]
on a donc~:
\[ U(A)\times (A-\lambda_1I)+Q(A)\times V(A)=I\] 
donc :
\[ \mbox{Im }Q(A) \supset \mbox{Ker }(A-\lambda_1 I) \]
%Exemple :\\
%$P_A(X)=(X-2)*Q(X)$\\
%Puisque $P_A$ est sans facteur carr\'e $X-2$ et $Q(X)$ sont premiers entre eux
%donc il existe $U(X)$ et $V(X)$  v\'erifiant :\\
%$(X-2).U(X)+Q(X).V(X)=1$\\
%On a donc :\\
%$(A-2I).U(A)+Q(A).V(A)=I$\\
%donc tous les vecteurs $w$ se d\'ecomposent en $w=w_1+w_2$ :\\
%$w_1=(A-2I).U(A).w\ \in \mbox{Ker }Q(A)$\\
%$w_2=Q(A).V(A).w \ \in \mbox{Ker }(A-2I)$\\ 
%et $\mbox{Ker }Q(A) \cap \mbox{Ker }(A-2I)=\{0\}$ (d'apr\`es B`'ezout)\\
%donc si $w \ \in \mbox{Ker }(A-2I)$ on a $w_1=0$ et $w=w_2$ (l'application
%$w->w_2$ est donc surjective sur le sous-espace propre associ\'e \`a 2 ce qui permet d'avoir une base du sous-espace propre associ\'e \`a 2.\\
Lorsque le polyn\^ome minimal n'est pas {\em square-free}\/, l'identit\'e
de B\'ezout donne les projecteurs sur les {\em espaces caract\'eristiques}
($Q(A)$ envoie $K^n$ sur Ker$(A-\lambda_1 I)^{p_1}$ lorsque $\lambda_1$ est 
d'ordre $p_1$).
%En fait, en appliquant l'identit\'e de B\'ezout, on peut m\^eme
%construire les projecteurs sur l'espace propre correspondant \`a $\lambda $
%et la somme directe des autres espaces propres.\\

{\bf Exercice 2} (à rendre au plus tard le jour de l'examen)\\
Montrer que $B$ est diagonalisable et
trouver une base de vecteurs propres de la matrice $B$ de l'exercice 1
en utilisant la m\'ethode ci-dessus.

\section{M\'ethode itératives}
\subsection{Méthode de la puissance.} 
\label{sec:puissance}
La m\'ethode de la puissance est une m\'ethode num\'erique qui permet de 
d\'eterminer la valeur propre de module maximal d'une matrice \`a coefficients
r\'eels (en supposant que $A$ poss\`ede
une seule valeur propre de module maximal qui est alors r\'eelle). 
On prend un vecteur 
colonne $v$ au hasard et on calcule la suite r\'ecurrente:
\[ v_0=v\ ,\ v_{n+1}= Av_n/||Av_n|| \]
Si la composante de $v_0$ sur l'espace propre correspondant \`a la valeur
propre de plus grand module n'est pas nulle, $\pm v_n$ tend vers un vecteur
(norm\'e) de cet espace propre.

D\'emonstration~:\\
si la base propre est $w_1,...,w_n$ (avec $Aw_i=\lambda_i w_i$ avec
$|\lambda_1|>|\lambda_2|>...>|\lambda_n|$) on a :
\[ v_0=a_1w_1+a_2w_2+....+a_nw_n \] 
donc~:
\[ Av_0=a_1\lambda_1w_1+a_2\lambda_2w_2+....+a_n\lambda_nw_n \]
et on pose $||Av_0||^{-1}=\mu_1$.
Comme $v_1=\mu_1 Av_0$, on a~: 
\[ Av_1=\mu_1 A^2v_0=
\mu_1 (a_1\lambda_1^2w_1+a_2\lambda_2^2w_2+....+a_n\lambda_n^2w_n) \]
et on pose  $||Av_1||^{-1}=\mu_2$ donc~:
\[ Av_2=
\mu_2\mu_1(a_1\lambda_1^2w_1+a_2\lambda_2^2w_2+....+a_n\lambda_n^2w_n)\]
.....\\
si on pose $||Av_{k-1}||^{-1}=\mu_k$ et $C_k=\mu_k..\mu_1$ on a :\\
\[ v_k=\mu_kAv_{k-1} \]
donc~:
\begin{eqnarray*} 
v_k &=&\mu_k..\mu_1  A^kv_0=\mu_k..\mu_1 
(a_1\lambda_1^kw_1+a_2\lambda_2^kw_2+....+a_n\lambda_n^kw_n) \\
&=&C_k\lambda_1^k(a_1w_1+a_2(\frac{\lambda_2}{\lambda_1})^kw_2+...
a_n(\frac{\lambda_n}{\lambda_1})^kw_n)\\
&\simeq &C_k\lambda_1^ka_1w_1 
\end{eqnarray*}
On a bien pour $k$ assez grand, $v_k$ colin\'eaire au vecteur propre $w_1$ 
donc $Av_k\simeq\lambda_1v_k $  et puisque $v_k$ est de norme 1,
$|\lambda_1|\simeq||Av_k||$. 
Si la premi\`ere coordonn\'ee $(v_k)_1$ de $v_k$ est non nulle, on a~:
\[ \lambda_1 \simeq \frac{(Av_k)_1}{(v_k)_1} \]
d'autre part~:
\[ v_{k+1} \simeq \frac{\lambda_1}{|\lambda_1|}v_k \]
donc $\pm v_k$ converge vers un vecteur propre, en pratique on teste
si $||v_k-v_{k+1}||$ ou $||v_k+v_{k+1}||$ est inf\'erieur \`a $\varepsilon$,
auquel cas on arr\^ete le calcul des $v_k$.
 
{\bf Exercice 3} (à rendre à la fin de la 2ème séance de ce TP)\\
\'Ecrire un programme MuPAD ou calculatrice mettant en oeuvre
cet algorithme. Utilisez ce programme pour trouver une valeur approch\'ee 
de la valeur propre de norme maximale par la m\'ethode
de la puissance de la matrice $C$ de l'exercice 1 puis
d'une matrice al\'eatoire.\\
Pour g\'en\'erer une matrice al\'eatoire:
\begin{itemize}
\item Avec \verb|MuPAD| utilisez {\tt linalg::randomMatrix(m, n,Dom::Integer)}
\item Avec {\tt maple} utilisez {\tt linalg::randmatrix}
\item HP49: \verb|RAND({4,4})|, TI89/92: \verb|RandMat(4,4)|.
\end{itemize}
Attention, votre matrice al\'eatoire peut avoir un couple de complexe
conjugu\'es comme valeur propre de plus grand module (en pratique
le test d'arr\^et du programme ne sera jamais atteint). Dans ce cas
r\'eessayez avec une autre matrice.

\subsection{Méthode des itérations inverses.}
La  m\'ethode des it\'erations inverses consiste lorsque l'endomorphisme $f$ 
est inversible \`a appliquer l'algorithme de la puissance \`a $f^{-1}$ 
en effet:\\
si $\lambda$ est une valeur propre de $f $ et si $w$ est un vecteur propre
associ\'e \`a $\lambda$ on a :\\
$\displaystyle f(w)=\lambda w \ \Leftrightarrow \ f^{-1}w=\frac{1}{\lambda}w$\\
La  m\'ethode des it\'erations inverses permet donc de trouver la valeur 
propre de plus petit module (\`a condition d'avoir invers\'e la matrice 
$A$ associ\'ee \`a $f$).\\
La  m\'ethode des it\'erations inverses est utile
lorsqu'on connait une valeur approch\'ee $\tau$ d'une valeur propre $\lambda$.\\ On
peut alors am\'eliorer $\tau$ en utilisant des it\'erations inverses, 
puisqu'alors la matice $B=A-\tau I$ est inversible et poss\`ede $\lambda-\tau$
(qui est tr\`es petit) comme valeur propre.\\
On cherche l'inverse de $B=A-\tau I $,
% pour cela on r\'esoud \[ (A-\tau I)y=b \]
puis on pose :\\
$y_0=B^{-1}b$ o\`u $b$ est un vecteur al\'eatoire  ($y_0$ est alors
 proche d'un vecteur propre correspondant \`a $\lambda \approx \tau$.\\
 On it\`ere ensuite la proc\'edure.

{\bf Exercice 4} (à rendre au plus tard le jour de l'examen)\\
Trouvez la plus petite valeur propre des matrices de l'exercice 3
(changez de matrice al\'eatoire si n\'ecessaire).

Pour trouver les autres valeurs propres/vecteurs propres,
il faut pouvoir \'eliminer la valeur propre trouv\'ee.\\
On sait en particulier le faire quand $A$ est sym\'etrique,
car il suffit de remplacer $A$ par $A'=A-\lambda_1 w_1 \ ^t w_1$.\\
En effet, on prend une base orthonormale $(w_1,...,w_n)$ de vecteurs
propres de $A$. $A'$ a les m\^emes vecteurs propres que $A$ et pour
valeurs propres correspondantes $0$ et les $\lambda_k$ ($k>1$) car pour
$k>1$~:
\[ A' w_k=A w_k +0 = \lambda_k w_k \]
puisque $w_k$ est orthogonal \`a $w_1$

{\bf Exercice 5} (à rendre au plus tard le jour de l'examen)\\
Trouvez toutes les valeurs propres de la matrice~:
\[D= \left(\begin{array}{rrr}
 9 & -1 & 2  \\
 -1 & 5 & -2  \\
 2 & -2 & -7 
\end{array}\right) \]

\subsection{Cas o\`u il y un couple de complexe conjugu\'e de module maximal.}
Soit $(\lambda,\overline{\lambda})$ le couple de valeurs propres
correspondant au couple de vecteurs propres $(w,\overline{w})$.
Montrer que $v_k$ s'approche de
l'espace vectoriel engendr\'e par $w$ et $\overline{w}$. En d\'eduire que
$A^2v_{k}$ est approximativement combinaison lin\'eaire \`a coefficients
r\'eels de $Av_k$ et $v_{k}$ lorsque $k$ est grand~:
\[ v_{k+2}=av_{k+1}+bv_k , \mbox{ avec } \lambda^2=a \lambda +b \]

{\bf Exercice 6} (à rendre au plus tard le jour de l'examen)\\
Trouver le couple de valeurs propres de plus grand module de la matrice~:
\[E= \left(\begin{array}{rrrr}
 -3 & -8 & -3 & 4  \\
 7 & 8 & 6 & 5  \\
 -1 & 4 & -4 & 1 \\
7  & -4 & -5 & 7 
\end{array}\right) \]
Faites de m\^eme pour le couple de plus petit module.

\end{document}

\section{Bonus: diagonalisation exacte.} 
\label{sec:souriau}

\subsection{Calcul du polyn\^ome caract\'eristique et de la matrice adjointe}
\label{sec:adj}
\subsubsection{Les formules de Newton et m\'ethode de Leverrier}
On note $P$ le polyn\^ome caract\'eristique de $A$ d\'efini par :\\
$P(\lambda)=\det(\lambda I-A)=(\lambda-\lambda_1)(\lambda-\lambda_2)...(\lambda-\lambda_n)=\lambda^n-p_1\lambda^{n-1}-p_2\lambda^{n-2}...-p_n$\\
On pose :\\
$s_1=\mbox{tr}(A)=\Sigma \lambda_i=\Sigma _{i=1}^n a_{i,i}$\\
$s_k=\mbox{tr}(A^k)=\Sigma \lambda_i^k \ (k=1..n)$\\
On montre alors les formules de Newton :\\
$p_1=s_1$\\
$2.p_2=s_2-p_1.s_1$\\
$k.p_k=s_k-p_1.s_{k-1}-...-p_{k-1}.s_1 \ (k=1..n)$\\ 
Ces formules permettent de calculer les  coefficients $p_k \ (k=1..n)$ du  
polyn\^ome caract\'eristique de $A$ lorsqu'on a calcul\'e les traces des 
matrices $A\ A^2..A^n$.\\
La matrice adjointe $B(\lambda)$ est un polyn\^ome en $\lambda$ de degr\'e 
$n-1$ \`a coefficients matriciels :\\
$B(\lambda)=\lambda^{n-1}I+\lambda^{n-2}B_1+...B_{n-1}$\\
puisque $B(\lambda) \times (\lambda I-A)=P(\lambda) I$ on a la relation de  r\'ecurrence :\\
$B_0=I$, $B_k=AB_{k-1}-p_kI$ \\
donc
 $B_k=A^k-p_1A^{k-1}-p_2A^{k-2}...-p_kI$
\subsubsection{La m\'ethode de Faddeev}
La m\'ethode de Faddeev permet le calcul simultan\'e des coefficients $p_i \ (i=1..n)$ du polyn\^ome caract\'eristique $P(\lambda)=\det(\lambda I-A)$ 
et des coefficients matriciels
$B_i \ (i=1..n-1)$ du polyn\^ome en $\lambda$ donnant la matrice adjointe.\\
Au lieu de calculer les traces des puissaces de $A$, Faddeev calcule les traces
 des matrices $A_i \ (i=1..n)$ d\'efinies par :\\
$A_1=A,\ p_1=\mbox{tr}(A),\ B_1=A_1-p_1I$\\  
$A_2=AB_1,\ p_2=\frac{1}{2}\mbox{tr}(A_2),\ B_2=A_2-p_2I$\\ 
$A_k=AB_{k-1},\ p_k=\frac{1}{k}\mbox{tr}(A_k),\ B_k=A_k-p_kI$\\  
On montre facilement que :\\
$B_n=A_n-p_nI=0$\\
les nombres $p_1,p_2...p_n$ et les matrices $B_1,B_2,...B_n$ ainsi d\'efinis 
sont bien les coefficients de $P(\lambda)$ et de $B(\lambda)$ puisque on
retrouve les formules de Newton en prenant les traces des \'egalit\'es :
$A_k=A^k-p_1A^{k-1}-...-p_{k-1}A$\\
$B_k=A^k-p_1A^{k-1}-...-p_{k-1}A-p_kI$
\subsubsection{D\'emonstration de la  m\'ethode de Faddeev}
Th\'eor\`eme : $P'(\lambda)$, la d\'eriv\'ee  du polyn\^ome caract\'eristique, 
est \'egal \`a la trace de  $B(\lambda)$ ($B(\lambda)$ est la matrice adjointe
de $\lambda I-A$ et on note $b_{i,j}(\lambda)$ ses coefficients).\\
D\'emonstration :\\
si on note $V_1(\lambda),...V_n(\lambda)$ les vecteurs colonnes de 
$\lambda I-A$, on a :
$\det(\lambda I-A)-\det(\lambda_0 I-A)=\det(V_1(\lambda)-V_1(\lambda_0),V_2(\lambda),,...,V_n(\lambda))+$\\$\det(V_1(\lambda_0),V_2(\lambda)-V_2(\lambda_0),...,V_n(\lambda))+...+\det(V_1(\lambda_0),V_2(\lambda_0),...,V_n(\lambda)-V_n(\lambda_0))$\\
On a donc :\\
$\displaystyle  P'(\lambda_0)=\lim_{\lambda \rightarrow \lambda_0}\frac{P(\lambda)-P(\lambda_0)}{\lambda-\lambda_0}=\det(V'_1(\lambda_0),V_2(\lambda_0),...,V_n(\lambda_0) )+$\\
$\det(V_1(\lambda_0),V'_2(\lambda_0),...,V_n(\lambda_0) )+...+\det(V_1(\lambda_0),V_2(\lambda_0),...,V'_n(\lambda_0) )$ \\
Pour finir la d\'emonstration, il suffit de remarquer que :\\
$V'_i(\lambda_0)=e_i$ \\
 $\det(V_1(\lambda_0),V_2(\lambda_0),...,V'_i(\lambda_0),...,V_n(\lambda_0) )=b_{i,i}(\lambda_0)$ \\
et donc que  $P'(\lambda_0)=\sum_{i=1}^n b_{i,i}(\lambda_0)=\tr(B(\lambda_0)$\\

Les relations de r\'ecurrence : \\
On note :\\
$P(\lambda)=\lambda^n-p_1\lambda^{n-1}-p_2\lambda^{n-2}...-p_n$\\
$B(\lambda)=\lambda^{n-1}I+\lambda^{n-2}B_1+...+B_{n-1}$\\
Donc  :\\
$P'(\lambda)=\tr(B(\lambda))=\lambda^{n-1}\tr(I)+\lambda^{n-2}\tr(B_1)+...+\tr(B_{n-1})$\\
ou encore \\
$n \lambda^{n-1}-p_1(n-1)\lambda^{n-2}-p_2(n-2)\lambda^{n-3}...-p_{n-1}=\lambda^{n-1}\tr(I)+\lambda^{n-2}\tr(B_1)+...+\tr(B_{n-1})$\\
ou encore $\tr(B_i)=-p_i(n-i)$\\
Comme $\ B(\lambda) \times (\lambda I-A)=P(\lambda) I \ $ on a les relations de  r\'ecurrence :\\
$B_0=I$,.., $B_i=AB_{i-1}-p_iI$ \\
donc
 $\tr(B_i)=-p_i(n-i)=\tr(AB_{i-1})-np_i$\\
$p_i=\frac{\tr(AB_{i-1})}{i}$

\pagebreak

On calcule donc :\\
$p_1=\tr(AB_0)$ et $B_1=AB_0-p_1I$\\
$p_2=\frac{\tr(AB_1)}{2}$ et $B_2=AB_1-p_2I$\\
....\\
$p_i=\frac{\tr(AB_{i-1})}{i}$ et $B_i=AB_{i-1}-p_iI$ \\
\subsection{Exercice 6} 
Soit la matrice :\\
\[ A=\left(\begin{array}{cccc}
 2& -1& 1&2\\
 0 & 1 & 1&0\\
 -1 & 1 & 1&1\\
 1 & 1 & 1&0\\
 \end{array}\right),\]
D\'eterminer $p_1,p_2,p_3p_4$ et $B_1,B_2,B_3,B_4$ :\\
- par la  m\'ethode de Leverrier\\
- par la  m\'ethode de Faddeev\\
En d\'eduire le polyn\^ome caract\'eristique de $A$ et calculer $A^{-1}$.\\
V\'erifier votre r\'esultat \`a l'aide des commandes {\tt MuPAD}.


\section{Programme pour la m\'ethode de Fadeev.} \label{sec:Fadeev}
\begin{verbatim}
// programme a sauver sous le nom fadeev.mu
// fonction iequalj utile pour definir la matrice identite
iequalj:=(k,l)->(if (k=l) then 1; else 0; end_if;); 

fadeev:=proc(A)
local Aj,AAj,Id,coef,n,pcara,lmat;
begin
 n:=ncols(A);
 Id:=Dom::Matrix()(n,n,iequalj);  // matrice identite
 Aj:=Id;
 lmat:=[];                        // liste vide
 pcara:=[1];                       // coefficient de plus grand degre
 for j from 1 to n do
  lmat:=append(lmat,Aj):          // rajoute Aj a la liste de matrices
  AAj:=Aj*A:
  coef:=-tr(AAj)/j:
  pcara:=append(pcara,coef):        // rajoute coef au polynome caracteristique
  Aj:=AAj+coef*Id:
 end_for;
 lmat,pcara;                       // resultat
end_proc;
\end{verbatim}
Exemple d'utilisation:
\begin{verbatim}
export(linalg);
read("fadeev.mu");
A:=Dom::Matrix()(2,2,[[1,2],[3,4]]);
res:=fadeev(A);
(res[1])[2]*A;
\end{verbatim}


\end{document}

\subsection{Avec la HP49G mode RPN}
\begin{verbatim}
<<
  DUP IDN DUP           @ A Id A_j (initialise a Id pour j=0)
  DUP SIZE 1 GET        @ A Id A_j n
  { } { 1 }             @ A Id A_j n lmat pcara (liste des matrices et poly car)
  -> A ID AJ N LMAT PCARA
  <<
     1 N FOR J          @ (boucle J de 1 a N)
       LMAT AJ +        @ lmat
       'LMAT' STO       @ (rajoute A_j a la liste des matrices)
       AJ A *           @ A A*A_j
       DUP TRACE        @ A*A_j trace
       J NEG /          @ A*A_j -trace/j
       PCARA OVER +      @ A*A_j -trace/j pcara
       'PCARA' STO       @ A*A_j -trace/j (ajoute le coefficient a pcara)
       ID * +           @ A_j+1
       'AJ' STO         @ (stocke le resultat)
     NEXT
     LMAT PCARA          @ rappelle lmat et pcara sur la pile
  >>
>> 'FADEEV' STO
\end{verbatim}
Pour obtenir le polynome caracteristique sur la pile sous forme symbolique
tapez \verb|X PEVAL|.

\section{Quelques r\'ef\'erences} \label{sec:ref}

\begin{itemize}
\item 
Press et al.: \\
Numerical Recipies in Pascal (exite aussi en C)


\item Davenport, Siret, Tournier:\\
Calcul formel: Syst\`emes et algorithmes de manipulations alg\'ebriques


\item Gantmacher:\\
Th\'eorie des matrices

\end{itemize}

\section{M\'ethode de la puissance.} \label{sec:puissance2}
Programme pour effectuer une it\'eration:
\begin{verbatim}
export(linalg);
M:=Dom::Matrix()(2,2,[[1,2],[3,4]]);
v:=Dom::Matrix()(2,1,[1.2,2.3]);
iter:=proc()
begin
 v:=normalize(M*v);
end_proc;
\end{verbatim}
Modifier $M$ et $v$ puis tapez \verb|iter();|.

\begin{verbatim}
<<
  M SWAP *
  DUP ABS /
>> 'ITER' STO
\end{verbatim}
Sur la \verb|49G mode RPN| , on stocke la matrice dans une variable (\verb|'M' STO|)
puis on tape un vecteur, puis \verb|ITER| plusieurs fois.

(\verb|[ [ 1 2 ] [ 3 4] ]| ou shift droit \verb|ENTER| pour 
utiliser matrixwriter). Stockez la dans la variable $A$ en
tapant \verb|A STO|.
L'instruction \verb|ARRAY->| permet de d\'ecomposer une matrice 
sur la pile, \verb|EVAL * ->LIST| la recompose alors en 
un vecteur liste.
Il suffit ensuite de rajouter la puissance de \verb|a| correspondante
au vecteur liste par exemple \verb|'a^2' +|.
Lorsque tous les vecteurs listes sont cr\'e\'es sur la pile, entrez
le nombre de vecteurs listes puis tapez \verb|->LIST| pour cr\'eer
la matrice, puis \verb|REF| pour la r\'eduire.

On peut \'egalement programmer cet algorithme en utilisant le langage
de la HP48. Le principe general est que tous les arguments d'une
fonction sont pris sur la pile et tous les r\'esultats sont renvoy\'es
sur la pile:
\begin{itemize}
\item pour cr\'eer un programme, utiliser les guillemets fran\c{c}ais 
\verb|<< >>|, taper le texte du programme puis stockez le r\'esultat
(\verb|'PROG' STO|), pour l'ex\'ecuter, taper le(s) argument(s) puis
le nom du programme.
\item pour d\'eclarer des variables locales:\verb|-> A B << ... >>|
ici le niveau 2 de la pile est assign\'e \`a \verb|A| et le niveau 1
\`a \verb|B| (et ces 2 niveaux sont supprim\'es). On peut utiliser
\verb|A| et \verb|B| entre les guillemets qui suivent leur d\'eclaration.
\item test: \verb|IF .. THEN .. ELSE .. END|: la clause \verb|THEN|
est ex\'ecut\'ee lorsque le niveau 1 de la pile est non nul (comme en
C, une valeur nulle est consid\'er\'ee comme {\em false}\/ et une
valeur non nulle {\em true}\/).
\item boucle ind\'efinie: \verb|1 10 FOR I .. NEXT|, 
\item boucle d\'efinie: \verb|WHILE .. REPEAT .. END|
(ou \verb|DO .. UNTIL .. END|). La boucle est ex\'ecut\'ee tant que
le niveau 1 de la pile n'est pas nul lorsque \verb|REPEAT| s'ex\'ecute.
\end{itemize}









