Define functions for `Cholesky`

[devmotion]

This PR defines

• dim
• whiten, whiten!
• unwhiten, unwhiten!
• quad, quad!
• invquad, invquad!
• X_A_Xt, Xt_A_X, X_invA_Xt, Xt_invA_X

for Cholesky.

Ultimately, I think we should forward (some of) the generic fallbacks for AbstractMatrix to these functions (I noticed they could also be cleaned a bit since e.g. https://github.com/JuliaStats/PDMats.jl/blob/8d88c15f7c9931038ab4b741fea92f7fc990da4d/src/generics.jl#L124 is never hit for x::AbstractMatrix due to the definition in https://github.com/JuliaStats/PDMats.jl/blob/8d88c15f7c9931038ab4b741fea92f7fc990da4d/src/generics.jl#L125; additionally, it is problematic e.g. for StaticArrays that out-of-place methods fall back to in-place methods). However, I assumed it would be cleaner to not touch the dispatch pipeline in this PR but only add and test definitions for Cholesky.

This PR also provides a solution to https://github.com/JuliaStats/PDMats.jl/issues/140 @st--: If one does not want to construct the full matrix, it is possible to work with Cholesky directly (at least when using the PDMats-specific API).

[codecov-commenter]

Codecov Report

Attention: 8 lines in your changes are missing coverage. Please review.

Files	Coverage Δ
src/chol.jl	80.48% <75.75%> (-19.52%)	:arrow_down:

:loudspeaker: Thoughts on this report? Let us know!.

[btmit]

It might be worth extending addition as well. Below uses a full rank update to add two matrices directly using their Cholesky decompositions. Obviously slower than just adding the matrices if you have them, but if you're starting and ending with the Cholesky this seems the way to go. IIRC this is more numerically stable as well. Unfortunately, this doesn't work for StaticArrays #https://github.com/JuliaArrays/StaticArrays.jl/issues/930

using LinearAlgebra, PDMats, BenchmarkTools

function fullrankupdate(c1::Cholesky, c2::Cholesky)
    cout = copy(c1)
    return fullrankupdate!(cout, c2)
end

function fullrankupdate!(c1::Cholesky, c2::Cholesky)
    for c in eachcol(PDMats.chol_lower(c2))
        lowrankupdate!(c1, c)
    end
    return c1
end

c1 = cholesky([1.0 0.5; 0.5 2.0])
c2 = cholesky([2.0 0.2; 0.2 1.0])

@btime PDMat($c1) + PDMat($c2);

@btime fullrankupdate($c1, $c2);

@assert PDMats.chol_lower(cholesky(PDMat(c1) + PDMat(c2))) ≈ PDMats.chol_lower(fullrankupdate(c1, c2))

458.122 ns (12 allocations: 688 bytes)
70.585 ns (5 allocations: 240 bytes)
true

[devmotion]

I don't think this should be done in this PR. In general, I'm not sure if it's a good idea at all to define addition for Cholesky at all but if it is desired it should probably be done in base.