From 1e8906315619e9f3d8896580c2fdcd91bcaf5659 Mon Sep 17 00:00:00 2001 From: raghav9-97 Date: Mon, 17 Dec 2018 15:50:06 +0530 Subject: [PATCH] Fixed sparse cholesky to return 1D-Vector --- stdlib/SuiteSparse/src/cholmod.jl | 10 +++++++++- stdlib/SuiteSparse/test/cholmod.jl | 5 +++++ 2 files changed, 14 insertions(+), 1 deletion(-) diff --git a/stdlib/SuiteSparse/src/cholmod.jl b/stdlib/SuiteSparse/src/cholmod.jl index 9260b648f83a2..c8dc6e83532cf 100644 --- a/stdlib/SuiteSparse/src/cholmod.jl +++ b/stdlib/SuiteSparse/src/cholmod.jl @@ -1686,10 +1686,18 @@ end (\)(L::Factor, B::SparseVecOrMat) = sparse(spsolve(CHOLMOD_A, L, Sparse(B, 0))) \(adjL::Adjoint{<:Any,<:Factor}, B::Dense) = (L = adjL.parent; solve(CHOLMOD_A, L, B)) -\(adjL::Adjoint{<:Any,<:Factor}, B::VecOrMat) = (L = adjL.parent; Matrix(solve(CHOLMOD_A, L, Dense(B)))) \(adjL::Adjoint{<:Any,<:Factor}, B::Sparse) = (L = adjL.parent; spsolve(CHOLMOD_A, L, B)) \(adjL::Adjoint{<:Any,<:Factor}, B::SparseVecOrMat) = (L = adjL.parent; \(adjoint(L), Sparse(B))) +function \(adjL::Adjoint{<:Any,<:Factor}, b::StridedVector) + L = adjL.parent + return Vector(solve(CHOLMOD_A, L, Dense(b))) +end +function \(adjL::Adjoint{<:Any,<:Factor}, B::StridedMatrix) + L = adjL.parent + return Matrix(solve(CHOLMOD_A, L, Dense(B))) +end + const RealHermSymComplexHermF64SSL = Union{ Symmetric{Float64,SparseMatrixCSC{Float64,SuiteSparse_long}}, Hermitian{Float64,SparseMatrixCSC{Float64,SuiteSparse_long}}, diff --git a/stdlib/SuiteSparse/test/cholmod.jl b/stdlib/SuiteSparse/test/cholmod.jl index f7e46ec3d4106..7c220267a1d59 100644 --- a/stdlib/SuiteSparse/test/cholmod.jl +++ b/stdlib/SuiteSparse/test/cholmod.jl @@ -676,6 +676,11 @@ end @test_throws ArgumentError logdet(Fnew) end +@testset "Issue #28985" begin + @test typeof(cholesky(sparse(I, 4, 4))'\rand(4)) == Array{Float64, 1} + @test typeof(cholesky(sparse(I, 4, 4))'\rand(4,1)) == Array{Float64, 2} +end + @testset "Issue with promotion during conversion to CHOLMOD.Dense" begin @test CHOLMOD.Dense(fill(1, 5)) == fill(1, 5, 1) @test CHOLMOD.Dense(fill(1f0, 5)) == fill(1, 5, 1)