function xnew = jacobi(A, b, x0, m)
[p, n] = size(A);
if not(p==n)
    error('matrice non quadrata');
end
if not(length(b)==n) || not(length(x0)==n)
    error('lunghezza vettori sbagliata');
end
xnew = zeros(n, 1);
xold = x0;
for k = 1:m
  for i = 1:n
    somma = 0;
    for j = 1:n
      if not(j==i)
       somma = somma + A(i,j)*xold(j);
      end
    end
    xnew(i) = (b(i) - somma) / A(i,i);
  end
  % qui ho sia xnew che xold a disposizione
  % volendo, posso controllare se norm(xold-xnew)
  % e' piccola e fermare l'iterazione
  xold = xnew;
end