%SOR 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
w = 1.1;              %acceleration factor

%SOR (successive overrelaxation) method

x = x_0;
for k=1:N
    for i=1:n
        x(i) = x(i)*A(i,i)*(1-w)/w;
        x(i) = x(i)+b(i);
        for j=1:i-1
            x(i) = x(i)-A(i,j)*x(j);
        end
        for j=i+1:n
            x(i) = x(i)-A(i,j)*x(j);
        end
        x(i) = w*x(i)/A(i,i);
    end
end
            
disp('SOR iteration estimate is...')
x

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