clear all;clc;
format long;
n = 60;
m = 80;
A = randn(n,m);
A = A * A';
b = randn(n,1);
tol=1e-4;
  
D=diag(diag(A));
L=tril(A)-D;
U=triu(A)-D;

GS=max(abs(eig(-inv(D+L)*U))) 
if (GS<1)
    disp('Gauss Seidel method converges')
else disp('Gauss Seidel method diverges')
end

%% Gauss Seidel Method

itrG=0;
x=ones(n,1);
erG(1)=norm(A*x-b);
normValG=10;

while (normValG>tol && itrG<300)
    x_old=x;
    
    for i=1:n
        
        sigma=0;
        
        for j=1:i-1
            sigma=sigma+A(i,j)*x(j);
        end
        
        for j=i+1:n
            sigma=sigma+A(i,j)*x_old(j);
        end
        
        x(i)=(b(i)-sigma)/A(i,i);
    end
    
    normValG=norm(x_old-x);
    
    itrG=itrG+1;
end

xGS=x
itrG

%% CG

[xCG iterCG] = conjgrad(A,b,tol)