% Gauss-Seidel iterative solver for systems of linear equations
% Input parameters, consider system: Ax=b
%   A: 
%   x: current solution values, can start with zero vector
%   b: right hand side, b as in Ax=b
%   tolerance: iterations limit as residual norm
%   maxit: after this amount of iterations we will stop
% Return values:
%   solution 
%   iterations 

function [result, iterations] = gseidel2(A, x, b, tolerance, maxit)
    [ilim, jlim] = size(A);
    if (ilim ~= jlim) 
        disp('Provide square matrix');
        return;
    end
    N=ilim;
    
    %using L2 norm sqrt(1/N * sum( residual(j) ^ 2)) for each j=row
    tolerance_limit = tolerance^2 * jlim;
    %iterations_limit = tolerance;
    iterations = 0;
    while(1)
        iterations=iterations+1;
        for i = 1 : N
            if (A(i,i)~=0) 
                x(i) = (b(i) + A(i,i)*x(i) - dot(A(i,:), x)) / A(i,i);
            end
        end

        %calculating L2 norm
        norm=0;
        for i = 1 : N 
           tmp = b(i) - dot(A(i,:), x);
           norm = norm + tmp^2;
        end
        
        %checking tolerance with L2 norm
        if (norm<tolerance_limit) 
            break;
        end

        %checking iterations limit
        if (iterations>maxit)
            break;
        end
    end	
    result = x;
end