LinearAlgebra 包中的 LDLt 函数是否仅适用于 SymTriDiagonal 矩阵?

Does LDLt function from LinearAlgebra package only work with SymTriDiagonal matrices?

我正在研究对称矩阵的 LDL^t 分解。我自己的代码工作正常,但是,当我想使用 LinearAlgebra 包中的 LDLt 函数时,以下代码不起作用

using LinearAlgebra

A = [5.0 1 0 0 0; 1 5 1 0 0; 0 1 5 1 0; 0 0 1 5 1; 0 0 0 1 5]

ldltS = LDLt(A)

display(ldltS)

Julia 报告错误“错误:LoadError:类型数组没有字段 dv”。

如果我用 SymTriDiagonal 构建 A 而不是上面代码中的 A,则上面的代码有效。

你的代码有效,试试:

ldltS = LDLt(A);

然而 LDLt 的文档显示:

Matrix factorization type of the LDLt factorization of a real SymTridiagonal matrix S

如果你想强制它使用常规矩阵,问题是 show 显示假设它是 SymTridiagonal 矩阵的对象。参见LinearAlgebra\src\ldlt.jl中的getproperty方法:

function getproperty(F::LDLt, d::Symbol)      
    Fdata = getfield(F, :data)                
    if d === :d                               
        return Fdata.dv                       
    elseif d === :D                           
        return Diagonal(Fdata.dv)             
    elseif d === :L                           
        return UnitLowerTriangular(Fdata)     
    elseif d === :Lt                          
        return UnitUpperTriangular(Fdata)     
    else                                      
        return getfield(F, d)                 
    end                                       
end

行:

return Diagonal(Fdata.dv)             

可以替换为

return Diagonal(Fdata)             

你可以覆盖这个函数:

julia> function Base.getproperty(F::LDLt, d::Symbol)
           Fdata = getfield(F, :data)
           if d === :d
               return Fdata.dv
           elseif d === :D
               return Diagonal(Fdata)
           elseif d === :L
               return UnitLowerTriangular(Fdata)
           elseif d === :Lt
               return UnitUpperTriangular(Fdata)
           else
               return getfield(F, d)
           end
       end

现在可以使用了(关于 LDLt 中的代数假设这是一个 SymTridiagonal 矩阵):

julia> ldltS
LDLt{Float64, Matrix{Float64}}
L factor:
5×5 UnitLowerTriangular{Float64, Matrix{Float64}}:
 1.0   ⋅    ⋅    ⋅    ⋅
 1.0  1.0   ⋅    ⋅    ⋅
 0.0  1.0  1.0   ⋅    ⋅
 0.0  0.0  1.0  1.0   ⋅
 0.0  0.0  0.0  1.0  1.0
D factor:
5×5 Diagonal{Float64, Vector{Float64}}:
 5.0   ⋅    ⋅    ⋅    ⋅
  ⋅   5.0   ⋅    ⋅    ⋅
  ⋅    ⋅   5.0   ⋅    ⋅
  ⋅    ⋅    ⋅   5.0   ⋅
  ⋅    ⋅    ⋅    ⋅   5.0

您可以使用 LDLFactorizations.ldl:

julia> using LDLFactorizations

julia> M = [55 225 730; 225 979 3162; 730 3162 10216.0]
3×3 Matrix{Float64}:
  55.0   225.0    730.0
 225.0   979.0   3162.0
 730.0  3162.0  10216.0

julia> LDLT = ldl(M)
LDLFactorizations.LDLFactorization{Float64, Int64, Int64, Int64}(true, true, false, 3, [2, 3, -1], [2, 1, 0], [3, 3, 3], [1, 2, 3], [1, 2, 3], [1, 3, 4, 4], Int64[], Int64[], [2, 3, 3], [4.090909090909091, 13.272727272727273, 2.9999999999999942], [55.0, 58.54545454545462, 5.684341886080801e-13], [0.0, 0.0, 0.0], [1, 1, 2], 0.0, 0.0, 0.0, 3)

julia> LDLT.P # identity, so we can ignore
3-element Vector{Int64}:
 1
 2
 3

julia> (LDLT.L+I) * LDLT.D * LinearAlgebra.transpose(LDLT.L+I)
3×3 SparseArrays.SparseMatrixCSC{Float64, Int64} with 9 stored entries:
  55.0   225.0    730.0
 225.0   979.0   3162.0
 730.0  3162.0  10216.0