clear;

%--------------------------------------------------------
% Metodo QR per il calcolo di tutti gli autovalori di una
% matrice data.
%--------------------------------------------------------

A=[30 1 2 3;4 15 -4 -2; -1 0 3 5; -3  5 0 -1];
kmax=100; tol=1.e-3;
[Q,R]=qr(A);
T=Q*R;
% QQ è la matrice contenente il prodotto comulativo
% delle matrici ortogonali Q_k. Ovvero
% QQ=\prod_{k=1}^M Q_k. 
QQ=Q;
k=1;
while convergenzaQR(T,tol)==0 & k <= kmax,
    T1=R*Q;
    [Q,R]=qr(T1);
    QQ=QQ*Q;
    T=Q*R;
   % pause
    k=k+1;
end
disp('Numero iterazioni ' )
k
% 
% Al termine del programma,
% T rappresenta la matrice triangolare superiore 
% avente gli autovalori di A lungo la diagonale.
% Inoltre QQ'*A*QQ=T.

% Verifica
disp('Calcolo del residuo diag(T)-eig(A) ' )
norm(diag(T),2)-norm(eig(A),2)

function nn=convergenzaQR(T,tol)
%---------------------------------------------------------
% Controlla che gli elementi extradiagonali della matrice
% T sono sotto la tolleranza
%----------------------------------------------------------
[n1,n2]=size(T);
if n1 ~= n2,
    error('Attenzione .... matrice non quadrata');
    return
end
if sum(abs(diag(T,-1))) <= tol,
    nn=1;
else
    nn=0;
end
return






%
% ESERCIZIO 2, del 06/02/06
%
clear
%---------------------------------------------------
% Metodo QR con shift per il calcolo di tutti gli autovalori di 
% una data matrice quadrata A
%---------------------------------------------------
tol=1.e-10;
A=[30 1 2 3;4 15 -4 -2; -1 0 3 5; -3  5 0 -1];

%alfa=input('Dammi alfa = ')
%A=[alfa 2 3 10; 5 12 10 7; 9 7 6 13; 4 16 18 0];
%------------------------------------------------
% -> Matrice testo Quarteroni e Saleri, 
%                                  pag. 149 <--
%A=[alfa 2 3 13; 5 11 10 8; 9 7 6 12; 4 14 15 1];
%------------------------------------------------
% -> Matrice testo Quarteroni e Saleri,
%                                   pag. 155 <--
% A=[30 1 2 3; 4 15 -4 -2; -1 0 3 5; -3 5 0 -1];
%------------------------------------------------
%------------------------------------------------
% -> Matrice testo Quarteroni, Sacco e Saleri, function nn=convergenzaQRShift(T,tol,k)
%---------------------------------------------------------
%                                   pag. 178 <--
% A=[17 24 1 8 15; 23 5 7 14 16; 4 6 13 20 22; 10 12 19 21 3; 11 18 25 2 9];
%------------------------------------------------
iter_max=1000;
n=max(size(A));
[Q,R]=qr(A);
T=Q*R;
iter=1;
for k=n:-1:2,
    I=eye(k);
    while convergenzaQRShift(T,tol,k)==0 & iter <= iter_max,
     mu=T(k,k);
     [Q,R]=qr(T(1:k,1:k)-mu*I);
     T(1:k,1:k)=R*Q+mu*I;
     iter=iter+1;
    end
T(k,k-1)=0;
end

% Alla fine T rappresenta la matrice triangolare superiore 
% avente gli autovalori di A lungo la diagonale.
iter

% Verifica
disp('Calcolo del residuo diag(T)-eig(A) ' )
norm(diag(T),2)-norm(eig(A),2)



