给我一个完整的LU分解法的matlab代码

function [L,U] = LU_decomposition(A)
% LU decomposition using Gaussian elimination with partial pivoting
% A: input matrix
% L: lower triangular matrix
% U: upper triangular matrix

[m,n] = size(A);
if m ~= n
error('Input matrix must be square');
end

L = eye(n); % initialize L as identity matrix
U = A; % initialize U as A

for k = 1:n-1 % loop through columns
[~,p] = max(abs(U(k:n,k))); % find pivot row
p = p+k-1;
if U(p,k) == 0
error('LU decomposition failed: matrix is singular');
end
if p ~= k
% swap rows p and k in U
temp = U(k,:);
U(k,:) = U(p,:);
U(p,:) = temp;
% swap rows p and k in L
temp = L(k,1:k-1);
L(k,1:k-1) = L(p,1:k-1);
L(p,1:k-1) = temp;
end
% eliminate entries below pivot in U
for i = k+1:n
L(i,k) = U(i,k)/U(k,k);
U(i,:) = U(i,:) - L(i,k)*U(k,:);
end
end

end