`cholesky(::Symmetric{<:Any,<:SMatrix})` doesn't respect `uplo`

[danielwe]

StaticArrays' implementation of Cholesky decomposition assumes that uplo == 'U' in Symmetric and Hermitian wrappers. This gives incorrect results when uplo == 'L' and the underlying array is not actually symmetric/hermitian.

MWE:

julia> using LinearAlgebra, StaticArrays

julia> A = Symmetric(SA[1.5 -1; 1 2], :L)
2×2 Symmetric{Float64, SMatrix{2, 2, Float64, 4}}:
 1.5  1.0
 1.0  2.0

julia> C = cholesky(A);

julia> C.L * C.U
2×2 SMatrix{2, 2, Float64, 4} with indices SOneTo(2)×SOneTo(2):
  1.5  -1.0
 -1.0   2.0

julia> C.L * C.U ≈ A
false

I suppose the easiest solution is to replace the following line:

https://github.com/JuliaArrays/StaticArrays.jl/blob/0feac146bb71ff48ca0a8465739f6b442395292c/src/cholesky.jl#L9

with something like

    if A.uplo == 'L'
        return _cholesky(Size(A), A.data', check)
    else
        return _cholesky(Size(A), A.data, check)
    end

However, it would be even nicer if uplo were preserved by the factorization, i.e., if cholesky(A).uplo == A.uplo, the way it currently works in Base. That would require a change slightly deeper in the code, but could use the same basic trick: if uplo == 'L', take the conjugate transpose, carry out the factorization as if uplo == 'U', and finally take the conjugate transpose of the output and use it to instantiate a Cholesky instance with uplo == 'L'.