[Matlab] Problema Fattorizzazione LU con pivoting parziale

[HeroGian]

Dovrei implementare su Matlab una funzione che calcoli la fattorizzazione LU di una matrice con pivoting parziale.
Solo che c'è qualche errore che non riesco ad individuare nel codice, perchè dopo aver fattorizzato, scrivendo su Matlab inv(P)*L*U non mi torna fuori la matrice A di partenza..


function [L, U, P] = GaussPP(A, b)
n = length(A);
L = eye(n);
P = eye(n);

for k=1 : n-1
[~, idx] = max(abs(A(k:n, k)));
idx = idx+k-1;

if idx > k
s1 = A(k, :);
A(k, :) = A(idx, :);
A(idx, :) = s1;

s1 = P(k, :);
P(k, :) = P(idx, :);
P(idx, :) = s1;

s2 = b(k);
b(k) = b(idx);
b(idx) = s2;
end

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


L'errore è sicuramente nella parte del pivoting perchè senza quello funziona..