使用 LDL 分解更新 cholrank1

Question

我有对称正定 (SPD) 矩阵的 cholrank1 更新程序 (wikipedia)。

function [L] = cholupdate(L,x)
    p = length(x);
    for k=1:p
        r = sqrt(L(k,k)^2 + x(k)^2);
        c = r / L(k, k);
        s = x(k) / L(k, k);
        L(k, k) = r;
        L(k+1:p,k) = (L(k+1:p,k) + s*x(k+1:p)) / c;
        x(k+1:p) = c*x(k+1:p) - s*L(k+1:p,k);
    end
end

它适用于 LL 分解。我尝试像这样修复程序以使用 LDL 分解（即不调用 sqrt）：

function [L] = cholupdate_ldl(L,x)
    p = length(x);
    for k=1:p
        r = L(k,k) + x(k)^2;
        c = r / L(k, k);
        s = x(k) / L(k, k);
        L(k, k) = r;
        L(k+1:p,k) = (L(k+1:p,k) + s*x(k+1:p)) / c;
        x(k+1:p) = sqrt(c)*(x(k+1:p) - x(k)*L(k+1:p,k));
    end
end

它工作正常，但我被迫使用 sqrt。如何在完全不使用 sqrt 的情况下更新 LDL 分解？

Answer 1

有很多方法。参见 Gill、Golub、Murray 和 Saunders (1974)：计算数学中的Methods for Modifying Matrix Factorizations。为了正式总结您的问题，我引用论文：

终于到了算法：

这是我在 MATLAB 中的实现：

function [L1,D1] = ldlt_update(L0,D0,z)

n = size(L0,1) ;
D1 = zeros(n,n) ;
L1 = zeros(n,n) ;

a = 1 ;
w = z ;
for jj = 1:n
    p = w(jj) ;
    D1(jj,jj) = D0(jj,jj) + a*p^2 ;
    b = p*a/D1(jj,jj) ;
    a = D0(jj,jj)*a/D1(jj,jj) ;
    for r = jj+1:n
        w(r) = w(r) - p*L0(r,jj) ;
        L1(r,jj) = L0(r,jj) + b*w(r) ;
    end
end
end

上面引用的论文和 Gill、Murray 和 Wright (1982)：Practical Optimization. Brian Borchers has a complete set of MATLAB code for working with real symmetric positive definite LDLT factorizations as defined in Golub and Van Loan (2013): Matrix Computations and on his web site 中有一个替代算法。

使用 LDL 分解更新 cholrank1

cholrank1 update with LDL decomposition

matlab

matrix

linear-algebra