From 0b47749f295ed1d1fdd4a6593d6bfde8223154a7 Mon Sep 17 00:00:00 2001 From: "Viral B. Shah" Date: Mon, 26 Oct 2026 13:25:04 +0000 Subject: [PATCH] Return `Inf` from sparse `cond` of a singular matrix instead of throwing `cond(A, 1)` and `cond(A, Inf)` caught no `SingularException` from the factorization or the triangular solve in `opnormestinv`; dense `cond` returns `Inf` for a singular matrix, and now so does sparse. Co-Authored-By: Claude Opus 5.5 --- src/linalg.jl | 18 ++++++++++-------- test/solvers/issues.jl | 6 ++++++ 2 files changed, 16 insertions(+), 8 deletions(-) diff --git a/src/linalg.jl b/src/linalg.jl index 4d59b8da..446f071d 100644 --- a/src/linalg.jl +++ b/src/linalg.jl @@ -990,14 +990,16 @@ small enough to store dense. function cond(A::AbstractSparseMatrixCSC, p::Real=2) # as for a dense matrix; the estimator of the inverse's norm needs at least one column (p == 1 || p == Inf) && checksquare(A) == 0 && return zero(real(float(eltype(A)))) - if p == 1 - normAinv = opnormestinv(A) - normA = opnorm(A, 1) - return normA * normAinv - elseif p == Inf - normAinv = opnormestinv(copy(A')) - normA = opnorm(A, Inf) - return normA * normAinv + if p == 1 || p == Inf + # a singular matrix has an infinite condition number, as for a dense matrix; a + # triangular or diagonal `factorize(A)` reports it only in the solve + normAinv = try + p == 1 ? opnormestinv(A) : opnormestinv(copy(A')) + catch e + e isa SingularException || rethrow() + return convert(float(real(eltype(A))), Inf) + end + return opnorm(A, p) * normAinv elseif p == 2 throw(ArgumentError("only 1- and Inf-norm condition numbers are implemented for sparse matrices, for 2-norm try cond(Array(A), 2) instead")) else diff --git a/test/solvers/issues.jl b/test/solvers/issues.jl index 51ca3b60..60e094f9 100644 --- a/test/solvers/issues.jl +++ b/test/solvers/issues.jl @@ -254,4 +254,10 @@ end @test factorize(sparse([0.0 1; 1 0])) \ b ≈ [0.0 1; 1 0] \ b end +@testset "cond of a singular matrix is Inf" begin + @test cond(sparse([1.0 2; 3 6]), 1) === Inf + # the triangular solve, not the factorization, finds this one singular + @test cond(sparse([1.0 0; 2 0]), Inf) === Inf +end + end # module