Problems with posdef symmetric matrix and LinearSolve

I am trying to move some FEM code to use the LinearSolve package. Consider the code

using LinearAlgebra, SparseArrays, LinearSolve

# Create a posdef symmetric matrix
A = sprand(100,100,0.01); A = A + A' + 100*I;

# rhs
b=rand(100);

# Set the problem
prob = LinearProblem(A,b)
sol = solve(prob)

# Check algorithm
sol.alg 

# ir returns KLUFactorization(true, true)

# Check solution...OK
norm(A * sol.u .- b)

# Enforce symmetry to use Cholesky, since A is symmetric and posdef
prob2 = LinearProblem(Symmetric(A),b)
sol2 = solve(prob2)

# Gives the error 
# = 
 ERROR: MethodError: no method matching ldiv!(::SparseArrays.CHOLMOD.Factor{Float64}, # ::Vector{Float64})

Closest candidates are:
  ldiv!(::Any, ::ChainRulesCore.AbstractThunk)
   @ ChainRulesCore ~/.julia/packages/ChainRulesCore/C73ay/src/tangent_types/thunks.jl:95
  ldiv!(::Any, ::LinearSolve.InvPreconditioner, ::Any)
   @ LinearSolve ~/.julia/packages/LinearSolve/HbBUp/src/preconditioners.jl:30
  ldiv!(::Any, ::LinearSolve.ComposePreconditioner, ::Any)
   @ LinearSolve ~/.julia/packages/LinearSolve/HbBUp/src/preconditioners.jl:17
  ...
=#

Am I doing something wrong with LinearSolve or is it a bug ?

Thank you very much,
Eduardo.

It was a bug. Fixed in patch v1.32.2 which should be released within an hour.

3 Likes

Wow! That was fast!

Thank you very much for your help and for all the packages.

Released