%jacobi iteration for approximating the solution to Ax = b

A = [7,1,-1,2;1,8,0,-2;-1,0,4,-1;2,-2,-1,6];
b = [3;-5;4;-3];
n = length(b);

x_0 = [1;1;1;1];        %initial guess of solution
I = eye(n);           %identity matrix
N = 10;               %number of iterations

%Jacobi method

Q = zeros(n,n);
for i=1:n
    Q(i,i) = A(i,i);  %Q = diagonal matrix with same diagonal entries as A
end

x = x_0;
for i=1:N
    y = (Q-A)*x + b;
    for j=1:n
        x(j) = y(j)/Q(j,j);   %solve Qx = y
    end
end

disp('jacobi iteration estimate is...')
x

disp('exact solution is...')
x_exact = linsolve(A,b)
