function [x,errs,iter]=sor(A,x,b,w,max_iter,tol)
% Metodo SOR per la risoluzione di sistemi lineari
% con parametro di rilassamento w
%
% Input:
%            A          matrice n x n non singolare
%            x          approssimazione iniziale della soluzione 
%            b          vettore noto
%            w          parametro di rilassamento (0 < w < 2)
%            max_iter   numer massimo di iterazioni
%            tol        tolleranza del metodo
%Output:  
%            x          approssimazione della soluzione
%            errs:      errore relativo norm(x_new-x_old)/norm(x_new) di uno step 
%            iter:      numero di iterazioni del metodo

D = diag(diag(A));

if det(D) == 0
    error('Metodo SOR non applicabile');
end

E = -tril( A, -1 );
F = -triu( A,  1 );
M =  D - w*E;
N = (1-w)*D+w*F;

for iter=1:max_iter                       
    x_old = x;
    x     = M\(N*x_old+w*b);                % nuova iterazione
    errs(iter) = norm(x-x_old)/norm(x);   % stima errore relativo
    if (errs(iter) <= tol)
        return
    end  
end
end