Matlab Décomposition LU pivotante partielle

Question

Je suis en train de travailler avec ma décomposition lu en grande partie basée sur LU decomposition with partial pivoting Matlab

function [L,U,P] = lup(A) 
n = length(A); 
L = eye(n); 
U = zeros(n); 
P = eye(n); 

for k=1:n-1 
% find the entry in the left column with the largest abs value (pivot) 
[~,r] = max(abs(A(k:end,k))); 
r = n-(n-k+1)+r; 

A([k r],:) = A([r k],:); 
P([k r],:) = P([r k],:); 
L([k r],:) = L([r k],:); 

% from the pivot down divide by the pivot 
L(k+1:n,k) = A(k+1:n,k)/A(k,k); 

U(k,1:n) = A(k,1:n); 
A(k+1:n,1:n) = A(k+1:n,1:n) - L(k+1:n,k)*A(k,1:n); 

end 
U(:,end) = A(:,end); 

end 

Il semble fonctionner pour la plupart des matrices (ce qui correspond à la fonction Matlab lu), mais la matrice suivante semble pour produire des résultats différents:

A = [ 
3 -7 -2  2 
-3  5  1  0 
6 -4  0 -5 
-9  5 -5 12 
]; 

Je n'arrive pas à comprendre ce qui ne va pas. Il semble bien fonctionner sur les matrices mentionnées dans le poste lié

Answer 1

vous étiez assez proche. J'ai changé trois lignes au total

for k=1:n-1 est devenu for k=1:n nous ne faisons pas -1 parce que nous voulons aussi L(n,n)=u(n,n)/u(n,n)=1 avec votre méthode que nous quittions ceci

L(k+1:n,k) = A(k+1:n,k)/A(k,k); est devenu L(k:n,k) = A(k:n,k)/A(k,k); parce que vous partiez à L(k,k)=A(k,k)/A(k,k)=1

parce que le changement k+1 nous ne devons commencer par une matrice d'identité l puisque nous reproduisons maintenant les 1 sur les diagonales si L=eyes(n); est devenu L=zeros(n);

et le code complété

function [L,U,P] = lup(A) 
% lup factorization with partial pivoting 
% [L,U,P] = lup(A) returns unit lower triangular matrix L, upper 
% triangular matrix U, and permutation matrix P so that P*A = L*U. 
    n = length(A); 
    L = zeros(n); 
    U = zeros(n); 
    P = eye(n); 


    for k=1:n 
     % find the entry in the left column with the largest abs value (pivot) 
     [~,r] = max(abs(A(k:end,k))); 
     r = n-(n-k+1)+r;  

     A([k r],:) = A([r k],:); 
     P([k r],:) = P([r k],:); 
     L([k r],:) = L([r k],:); 

     % from the pivot down divide by the pivot 
     L(k:n,k) = A(k:n,k)/A(k,k); 

     U(k,1:n) = A(k,1:n); 
     A(k+1:n,1:n) = A(k+1:n,1:n) - L(k+1:n,k)*A(k,1:n); 

    end 
    U(:,end) = A(:,end); 

end