function [x,iter,relresid] = jacobi(A,b,x0,tol,maxiter)
%JACOBI Jacobi iterative method for solving Ax = b.
%
%   [x,iter,relresid] = jacobi(A,b,x0,tol,maxiter) computes
%   the solution x of Ax = b using the Jacobi iterative method.
%   x0 is the initial approximation to the solution,
%   tol is the error tolerance, and maxiter is the number of
%   iterations to execute. If the Jacobi method converges, iter
%   contains the number of iterations required to meet the
%   error tolerance. If the Jacobi iteration does not converge,
%   iter = -1. relresid is the relative residual obtained by
%   the iteration.

[m,n] = size(A);
if m~=n
	disp('matrix A is not square')  ;
	return;
end;

xnew = zeros(n,1);
x = x0;
for k = 1:maxiter
	for i = 1:n
		sum = 0;
		for j = 1:n
			if (j ~= i)
				sum = sum + A(i,j)*x(j);
			end
		end
		xnew(i) = (b(i)-sum)/A(i,i);
   end
   relresid = norm(b-A*xnew)/norm(b);
	if relresid < tol
		iter = k;
		x = xnew;
		return;
	end
	x = xnew;
end
x = xnew;
iter = -1;