gauss elimination

function x = gauss_elimination(A, B)
    [n, m] = size(A);

    if n ~= m
        error('Not a square matrix');
    end

    for i = 1 : n - 1
        for j = i + 1 : n
            for k = i + 1 : n
                A(j, k) = A(j, k) - (A(j, i) / A(i, i)) * A(i, k);
            end

            B(j, 1) = B(j, 1) - (A(j, i) / A(i, i)) * B(i, 1);
        end
    end

    x = zeros(1, n);
    x(1, n) = B(n) / A(n, n);

    for i = n - 1 : -1 : 1
        x(1, i) = B(i, 1);
        for j = i + 1 : n
            x(1, i) = x(1, i) - A(i, j) * x(1, j);
        end
        x(1, i) = x(1, i) / A(i, i);
    end