Cholesky decomposition of low-rank positive-semidefinite matrix

General Usage
,
3

For example, here is a rank-3 5x5 positive semidefinite matrix:

julia> B = rand(3,5); A = Hermitian(B'B);

cholesky(A) throws PosDefException, and cholesky(A, Val(true)) throws RankDeficientException. However, passing check=false forces the factorization to proceed even if it is rank-deficient:

julia> F = cholesky(A, Val(true); check=false)
CholeskyPivoted{Float64, Matrix{Float64}}
U factor with rank 4:
5×5 UpperTriangular{Float64, Matrix{Float64}}:
 1.26174  0.830057  1.04195    0.862304     0.981338
  ⋅       0.708177  0.242843  -0.34719     -0.0424352
  ⋅        ⋅        0.573781   0.122706     0.375544
  ⋅        ⋅         ⋅         1.05367e-8  -1.05367e-8
  ⋅        ⋅         ⋅          ⋅          -2.22045e-16
permutation:
5-element Vector{Int64}:
 3
 1
 4
 5
 2

Notice the slightly negative pivot due to roundoff errors, as well as the small 4th pivot (on the order of the square root of the precision — recall that Cholesky takes square roots of pivots). This is still a factorization of A (up to roundoff), taking into account the permutation F.p:

julia> F.L * F.L' - A[F.p, F.p]
5×5 Matrix{Float64}:
  0.0          0.0  -2.22045e-16  0.0  0.0
  0.0          0.0   0.0          0.0  0.0
 -2.22045e-16  0.0   0.0          0.0  0.0
  0.0          0.0   0.0          0.0  0.0
  0.0          0.0   0.0          0.0  2.22045e-16

However, in general the pivoted-Cholesky algorithm only factorizes a submatrix of A[F.p, F.p], up to the point where the algorithm terminates.

5 Likes