Matlab code for Permuted LU decomposition

Shubham Shakya

2020 9 19

0 个回答

7 次查看（30 天）

0 个投票

This is a permuted LU decomposition

but its showing error as not enough input arguments

function [P, L, U, sign] = splu(A)
% splu  Square PA=LU factorization *with row exchanges*.
%
% [P, L, U] = splu(A), for a square, invertible matrix A,
% uses Gaussian elimination to compute a permutation
% matrix P, a lower triangular matrix L and 
% an upper triangular matrix U so that P*A = L*U.
% P, L and U are the same size as A.
% sign = det(P); it is 1 or -1.
%
% See also slu, lu, rref, partic, nulbasis, colbasis.
[m, n] = size(A);
if m ~= n
    error('Matrix must be square.')
end
P = eye(n, n);
L = eye(n, n);
U = zeros(n, n);
tol = sqrt(eps);
sign = 1;
for k = 1:n
    if abs(A(k, k)) < tol
        for r = k:n
            if abs(A(r, k)) >= tol
                break
            end
            if r == n
                if nargout == 4
                    sign = 0;
                    return
                else
                    disp('A is singular within tolerance.')
                    error(['No pivot in column ' int2str(k) '.'])
                end
            end
        end
        A([r k], 1:n) = A([k r], 1:n);
        if k > 1, L([r k], 1:k-1) = L([k r], 1:k-1); end
        P([r k], 1:n) = P([k r], 1:n);
        sign = -sign;
    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
    end
    for j = k:n
        U(k, j) = A(k, j) * (abs(A(k, j)) >= tol);
    end
end
if nargout < 4
    roworder = P*(1:n)';
    disp('Pivots in rows:'), disp(roworder'); end
end