diff --git a/Project.toml b/Project.toml
index 51168952a..d0895a1a6 100644
--- a/Project.toml
+++ b/Project.toml
@@ -81,7 +81,7 @@ HYPRE = "1.4.0"
 InteractiveUtils = "1.10"
 IterativeSolvers = "0.9.3"
 JET = "0.8.28"
-KLU = "0.5"
+KLU = "0.6"
 KernelAbstractions = "0.9.16"
 Krylov = "0.9"
 KrylovKit = "0.6"
diff --git a/src/factorization.jl b/src/factorization.jl
index fcdd80bbc..98c431d9e 100644
--- a/src/factorization.jl
+++ b/src/factorization.jl
@@ -855,23 +855,17 @@ function SciMLBase.solve!(cache::LinearCache, alg::KLUFactorization; kwargs...)
     A = convert(AbstractMatrix, A)
     if cache.isfresh
         cacheval = @get_cacheval(cache, :KLUFactorization)
-        if cacheval !== nothing && alg.reuse_symbolic
+        if alg.reuse_symbolic
             if alg.check_pattern && !(SparseArrays.decrement(SparseArrays.getcolptr(A)) ==
                  cacheval.colptr &&
-                 SparseArrays.decrement(SparseArrays.getrowval(A)) == cacheval.rowval)
-                fact = KLU.klu(SparseMatrixCSC(size(A)..., getcolptr(A), rowvals(A),
-                    nonzeros(A)))
-            else
-                # If we have a cacheval already, run umfpack_symbolic to ensure the symbolic factorization exists
-                # This won't recompute if it does.
-                KLU.klu_analyze!(cacheval)
-                copyto!(cacheval.nzval, nonzeros(A))
-                if cacheval._numeric === C_NULL # We MUST have a numeric factorization for reuse, unlike UMFPACK.
-                    KLU.klu_factor!(cacheval)
-                end
-                fact = KLU.klu!(cacheval,
+                 SparseArrays.decrement(SparseArrays.getrowval(A)) ==
+                 cacheval.rowval)
+                fact = KLU.klu(
                     SparseMatrixCSC(size(A)..., getcolptr(A), rowvals(A),
-                        nonzeros(A)))
+                        nonzeros(A)),
+                    check = false)
+            else
+                fact = KLU.klu!(cacheval, nonzeros(A), check = false)
             end
         else
             # New fact each time since the sparsity pattern can change
@@ -882,9 +876,14 @@ function SciMLBase.solve!(cache::LinearCache, alg::KLUFactorization; kwargs...)
         cache.cacheval = fact
         cache.isfresh = false
     end
-
-    y = ldiv!(cache.u, @get_cacheval(cache, :KLUFactorization), cache.b)
-    SciMLBase.build_linear_solution(alg, y, nothing, cache)
+    F = @get_cacheval(cache, :KLUFactorization)
+    if F.common.status == KLU.KLU_OK
+        y = ldiv!(cache.u, F, cache.b)
+        SciMLBase.build_linear_solution(alg, y, nothing, cache)
+    else
+        SciMLBase.build_linear_solution(
+            alg, cache.u, nothing, cache; retcode = ReturnCode.Infeasible)
+    end
 end
 
 ## CHOLMODFactorization
diff --git a/test/default_algs.jl b/test/default_algs.jl
index 5093953ca..8d3ac261e 100644
--- a/test/default_algs.jl
+++ b/test/default_algs.jl
@@ -39,7 +39,7 @@ A = spzeros(2, 2)
 # test that solving a singular problem doesn't error
 prob = LinearProblem(A, ones(2))
 @test solve(prob, UMFPACKFactorization()).retcode == ReturnCode.Infeasible
-@test_broken solve(prob, KLUFactorization()).retcode == ReturnCode.Infeasible
+@test solve(prob, KLUFactorization()).retcode == ReturnCode.Infeasible
 
 @test LinearSolve.defaultalg(sprand(10^4, 10^4, 1e-5) + I, zeros(1000)).alg ===
       LinearSolve.DefaultAlgorithmChoice.KLUFactorization