function [S,D,P,info,idx] = EGPP(A)
//
// Funzione EGPP: cerca una matrice di permutazione P ed una fattorizzazione
//                LR della matrice PA.
//
// A: matrice N x N;
// S, D: fattorizzazione LR di PA (S triangolare inferiore con S(k,k) = 1,
//       D triangolare superiore, SD = PA);
//
// info = 0: pivot nullo al passo idx;
//      = 1: la procedura ha trovato P e determinato la fattorizzazione.
//
// Se info = 0 si pone S = D = P = [].
//
[NR, NC] = size(A);
N = NR;
info = 1; // resta tale solo se la procedura termina correttamente
//
P = eye(N,N); // inizialmente P = I
for k = 1:N-1, // colonne da 1 a N-1
    [mx,idx] = max(abs(A(k:N,k))); idx = idx + (k-1);
    if idx ~= k then
       riga = P(k,:); P(k,:) = P(idx,:); P(idx,:) = riga;
       riga = A(k,:); A(k,:) = A(idx,:); A(idx,:) = riga;
    end;
    if A(k,k) == 0 then
            info = 0; idx = k;
            printf("\n *** Errore: pivot nullo al passo %d\n\n",k);
            S = []; D = []; P = [];
            return;
       else for j = k+1:N, // righe da k+1 a N
                A(j,k) = A(j,k) / A(k,k);
                A(j,k+1:N) = A(j,k+1:N) - A(j,k)*A(k,k+1:NC);
            end;
       end;
end;
if info == 1 then
   D = triu(A);
   S = eye(A) + tril(A,-1);
end;
endfunction
