function [x,flag,k]=GS(A,b,toll)
P=tril(A);
N=P-A;
xold=b;
flag=0;
for k=1:1000
    x=P\(b+N*xold);  
    inc=norm(x-xold);
    if inc<toll
        return
    end
    xold=x;
end
flag=1;