numerical-methods - 带有部分旋转的 LU 分解如何工作？

Question

我正在使用最初将 L 的主对角线上的元素设置为 1 的方法（认为这是 Doolittle 的方法，但不确定，因为我看到它的名称不同）。我知道有大量的文档、论文和书籍，但我找不到不使用高斯消元法寻找 L 的简单示例。

Answer 1

与完全透视相比，部分透视仅使用行交换，而完全透视也对列进行透视。如下图和代码所示部分旋转的主要目的是交换行以找到那里的最大值 u 以避免在 for 循环中除以一个非常小的值，这会导致很大的条件数。

如果您尝试 LU 分解和一些病态矩阵（如任意对角占优矩阵）的简单实现，它应该会爆炸。

function [L,U,P] = my_lu_piv(A)
n = size(A,1);
I = eye(n);
O = zeros(n);
L = I;
U = O;
P = I;
    function change_rows(k,p)
        x = P(k,:); P(k,:) = P(p,:); P(p,:) = x;
        x = A(k,:); A(k,:) = A(p,:); A(p,:) = x;
        x = v(k); v(k) = v(p); v(p) = x;
    end

    function change_L(k,p)
        x = L(k,1:k-1); L(k,1:k-1) = L(p,1:k-1);
        L(p,1:k-1) = x;
    end
for k = 1:n
    if k == 1, v(k:n) = A(k:n,k);
    else
        z = L(1:k-1,1:k -1)\ A(1:k-1,k);
        U(1:k-1,k) = z;
        v(k:n) = A(k:n,k)-L(k:n,1:k-1)*z;
    end
    if k<n
        x = v(k:n); p = (k-1)+find(abs(x) == max(abs(x))); % find index p
        change_rows(k,p);
        L(k+1:n,k) = v(k+1:n)/v(k);
        if k > 1, change_L(k,p); end
    end
    U(k,k) = v(k);
end
end