function [x,n,nr,reshist]=it_jacobi(b,A,N,tol)
% rescaling
D=diag(diag(A));
norm(eye(size(A,1))-inv(D)*A,Inf)
A1=inv(D)*A;
b1=D\b;
nb=norm(b1);
% use 0 as initial guess
x=zeros(size(b));
for i=1:N
    r=b1-A1*x;
    nr=norm(r);
    reshist(i)=nr/nb;
    if nr<tol*nb
        break;
    end
    x=x+r;
end
nr=norm(r)/nb;
n=i;
