function d = CompDet(At,l)
% Compute the determinant using the output from the forward elimination
% procedure Gauss using Gaussian elimination with scaled partial pivoting
%
% INPUT: At (Matrix resulting from forward elimination procedure using
%           scaled partial pivoting. Matrix At is upper triangular if the 
%           rows are considered in the order given in l.)
%        l (Index vector. Element i of the vector gives the row of At used
%           as the i-th pivot row.)
% OUTPUT d (Determinant of the original coefficient matrix A)
%

%% Initialize stuff
d = At(l(1),1); % Set determinant equal to first "diagonal" element of At
count = 0; % Counter of row interchanges
ltemp = l; % Copy of index vector used to count row interchanges.
n = length(At); % Find number of rows and columns in A

%% Compute determinant
for i=1:n-1
    
        % Count if number of row swaps are odd or even
        if ltemp(i) ~= i
            ltemp(ltemp == i) = ltemp(i);
            ltemp(i) = i;
            count = count + 1;
        end
        
    % Multiply current determinant value by next "diagonal" element of At
    d = d*At(l(i+1),i+1);
end

% Change the sign if number of row swaps is odd.
d = d*(-1)^count;