%====================================================================
% Script per determinare la soluzione di un sistema lineare
% con il metodo di Jacobi
% Necessita della function jacobi
%====================================================================

clear all

% Definizione dati 
n = 50;
A = diag(4*ones(n,1)) + diag(-2*ones(n-1,1),1) + diag(ones(n-1,1),-1);
sol = ones(n,1);
b = A*sol;
toll = 1e-5;
x0 = zeros(n,1);
kmax = 50;

% Controlli generali
[mA,nA] = size(A);
if (mA~=nA)
    error('La matrice A non e'' quadrata')
end

if length(b) ~= nA
    error('Le dimensioni di A e b non sono compatibili')
end

% Risoluzione

% Metodo di Jacobi
[x,k,flag] = jacobi (A,b,x0,toll,kmax);

% Controlli
if flag ~= 0  
    error('Non si puo` applicare il metodo di Jacobi: D singolare!') 
end

disp('Numero di iterate')
disp(k)
if k == kmax
    disp('Massimo numero di iterazioni raggiunto')
end

% Uscita risultati
disp('Dimensione del sistema')
disp(mA)
disp('Tolleranza')
disp(toll)

% Grafico

for j = 1:k+1
    err(j) = norm(x(:,j)-sol);
end
xk = 0:k;
semilogy(xk,err,'r');
title('Metodo di Jacobi')
xlabel('iterate')
ylabel('norma dell''errore')