% Jacobi iteration scheme, applied component by component

A=[7 2 -3 1; -4 8 1 -2; 3 -2 9 2; 1 -1 1 4];
b=[-1;16;-3;5];

x=[1;1;1;1]; % initial guess

itermax=50;
tol=1e-8;

for k=1:itermax
    % assign most recent x to xprev
    xprev=x;
    % calculate new x
    x(1)=( b(1) - (A(1,2)*xprev(2)+A(1,3)*xprev(3)+A(1,4)*xprev(4)) )/A(1,1);
    x(2)=( b(2) - (A(2,1)*xprev(1)+A(2,3)*xprev(3)+A(2,4)*xprev(4)) )/A(2,2);
    x(3)=( b(3) - (A(3,1)*xprev(1)+A(3,2)*xprev(2)+A(3,4)*xprev(4)) )/A(3,3);
    x(4)=( b(4) - (A(4,1)*xprev(1)+A(4,2)*xprev(2)+A(4,3)*xprev(3)) )/A(4,4);
    % print iteration number and the x just computed
    fprintf('%2d  %12.8f  %12.8f  %12.8f  %12.8f\n',k,x(1),x(2),x(3),x(4))
    % check convergence
    if norm( x - xprev ) <= tol
        break
    end
end