% Valeurs propres et vecteurs propres d'une matrice
clear all   % supprimer toutes les variables A

format short

disp('-----------------------------------------------');
fprintf(' \n')
disp('Si vous n''avez pas accès au Symbolic Math Toolbox, voici comment')
disp('faires les calculs nécessaires en MATLAB pour reproduire nos étapes')
disp('à la main qui donnent la forme de Jordan d''une matrice.')

fprintf(' \n')

disp('Le premier cas:  aucun erreure d''arrondissement.')

disp('>> A= [3 1 -1 ; 0 1 4; 0 -1 5]  ')

A =[     3     1    -1
     0     1     4
     0    -1     5]

disp('----------------------------------------------- (taper une touche)---');
pause
fprintf(' \n')

disp('>> [V,D]=eig(A)')

[V,D]=eig(A)

disp('Puisque V a des colonne repétées, on reconnait le cas d''une matrice défectueuse.')
disp('La seule valeur propre est 3.')
disp('----------------------------------------------- (taper une touche)---');
pause
fprintf(' \n')

disp('Selon la méthode de Jordan, il faut résoudre (A-3I)w = v pour un vecteur')
disp('propre v.  Ici, on n''a que le seul, donc on fait:')

disp('>> v=V(:,1); [A - 3*eye(3) v]')

v=V(:,1); [A - 3*eye(3) v]

disp('>> R= rref([A-3*eye(3) v])')

R= rref([A-3*eye(3) v])

disp('Or qu''on peut lire le résultat du système.')

disp('-----------------------------------------------');
pause
fprintf(' \n')

disp('Mais c''est un désastre en générale.  Par exemple:')

disp('>> B=    [ 0     0     1;     -7     2     3;   -10    -1     7]')

B=    [ 0     0     1;     -7     2     3;   -10    -1     7]

disp('C''est une matrice semblable à la précédente.')
disp('-----------------------------------------------');
pause
fprintf(' \n')

disp('>> [V,D]=eig(B); V    --- pour voir les vecteurs propres')
[V,D]=eig(B); V

disp('>> v=V(:,1); rref([B-3*eye(3) v])')

v=V(:,1); rref([B-3*eye(3) v])

disp('Qui est un système incohérent, à cause des erreures d''arrondissement.')
disp('-----------------------------------------------');
pause
fprintf(' \n')

disp('Par contre, on sait qu''une solution à (B-3I)w=v existe.  Utiliser ')
disp('le pseudoinverse pour arrondir à une solution proche:')

disp('>> w= pinv(B-3*eye(3))*v')
w= pinv(B-3*eye(3))*v

disp('>> (B-3*eye(3))*w-v    --- vérification')

(B-3*eye(3))*w-v

disp('-----------------------------------------------');
pause
fprintf(' \n')

disp('Une autre option:')

fprintf('\n')
disp('>> T=[ B-3*eye(3) eye(3) v]  --- afin de forcer MATLAB de ne pas toucher col 7 ')
format short
T=[ B-3*eye(3) eye(3) v]

disp('-----------------------------------------------');
pause
fprintf(' \n')

disp('>> R= rref(T)')

R=rref(T)

disp('et on voit que la dernière colonne représente maintenant un système cohérent')
disp('(à arrondissement près)!  C''est qu''on n''a pas permis à MATLAB de chercher ')
disp('un pivot dans cette dernière colonne.  Alors c''est comme si on l''aurait réduit')
disp('à la main.')
disp('-----------------------------------------------');
pause
fprintf(' \n')

disp('Mais parfois, la difficulté c''est qu''il faut résoudre (A-3I)w=v')
disp('pour un vecteur qui est une combinaison linéaire inconnue d''une base')
disp('de vecteurs propres.  Remarquons:')

disp('C =1     4    -2;1     1     1;4    -8     7]')
C =[1     4    -2;1     1     1;4    -8     7]
     
disp('>> [V,D]=eig(C); V    --- voir les vecteurs propres (toute val prp=3)')
  [V,D]=eig(C); V   

  disp('-----------------------------------------------');
pause
fprintf(' \n')

disp('>> Z=[ C-3*eye(3) eye(3) V ]; R= rref(Z)')

Z=[ C-3*eye(3) eye(3) V ]; R= rref(Z)

disp('Effacer les trois colonnes qui proviennent du terme eye(3); le reste ')
disp('représente ce qu''on aurait trouvé, jusqu''à arrondissement près, si ')
disp('on l''aurait réduit à la main.')

disp('-----------------------------------------------');
pause
fprintf(' \n')

disp('Donc ici on voit que c''est V(:,1) qui donne un système cohérent.')
disp('En générale, il faudra choisir une combinaison linéaire des dernières')
disp('colonnes qui donnera un système cohérent.')


fprintf(' \n')
disp('(Les vrais algorithmes de MATLAB sont formidables.) ')

disp('Fin de Jordanalamain') 


