% Jacobi.m - Jacobi method for iterative solution of Ax=b
% 
% INPUT:
%   A       - coefficient matrix
%   b       - right-hand side
%   x0      - initial guess
%   tol     - tolerance for estimated relative error in infinity norm
%   maxiter - maximum number of Jacobi iterations to perform
%
% OUTPUT:
%   x       - final iteration
%   relerr  - estimated relative error in infinity norm
%   niter   - number of iterations executed
%
% Notes:
%  - If A has many zero components, pass A in sparse format.
%
% Written by Douglas Meade
% Created 19 September 2005
% Revised 27 September 2005

function [x,relerr,niter] = jacobi(A,b,x0,tol,maxiter)

relerr = inf;
niter = 1;

S = diag( diag(A) );
T = S - A;

while relerr >= tol & niter < maxiter,
   x=S \ (b+T*x0);
   relerr = norm( x-x0,inf )/norm( x,inf );
   x0 = x;
   niter = niter+1;
end