clear all
clc

A=[1 2 -2; 1 1 1; 2 2 1];
b=sum(A,2);


[x,e,k]=JACOBI(A,b,100,eps);
iter=linspace(1,k,k);
plot(e)
hold on
plot(iter,e(iter), '*')

function [x,e,k] = JACOBI(A,b,Nmax,tol)

%metodo di iterazione di Jacobi
%A=D+E, dove D è diagonale con in diagonale la diagonale di A
%B=-inv(D)*E;
%Jacobi converge se il raggio spettrale di B è <1
n=length(A);
x0=zeros(n,1);
D=diag(diag(A));
E=A-D;
B=-D\E;
h=D\b;
x=B*x0+h;
e(1)=norm(x-x0)/norm(x);
RHO=max(abs(eig(B)));

if RHO<1
    disp('Il metodo converge!');
    k=1;
    while k<Nmax && e(k)>eps
        x0=x;
        x=B*x + h;
        k=k+1;
        e(k)=norm(x-x0)/norm(x);
    end
else
    disp('Il metodo di Jacobi non converge!');
end
end