Skip to content

Commit a58422f

Browse files
araujomsdkarrasch
authored andcommitted
fix infinite loop in cond (#682)
Closes #680 For `n==2` this code was trying to find three non-parallel +-1 vectors, which is mathematically impossible.
1 parent 4d4a98d commit a58422f

2 files changed

Lines changed: 15 additions & 8 deletions

File tree

src/linalg.jl

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -1745,7 +1745,7 @@ function opnormestinv(A::AbstractSparseMatrixCSC{T}, t::Integer = min(2,maximum(
17451745
repeated = true
17461746
end
17471747
end
1748-
if !repeated
1748+
if !repeated && 2^(n-1) 2t #we need enough non-parallel ±1 vectors
17491749
saux2 = S[1:n,j]' * S_old[1:n,1:t]
17501750
if _any_abs_eq(saux2,n)
17511751
repeated = true

test/linalg_solvers.jl

Lines changed: 14 additions & 7 deletions
Original file line numberDiff line numberDiff line change
@@ -129,20 +129,24 @@ end
129129

130130
@testset "sparse matrix cond" begin
131131
local A = sparse(reshape([1.0], 1, 1))
132-
Ac = sprandn(20, 20,.5) + im*sprandn(20, 20,.5)
133-
Ar = sprandn(20, 20,.5) + eps()*I
134132
@test cond(A, 1) == 1.0
135-
# For a discussion of the tolerance, see #14778
136-
@test 0.99 <= cond(Ar, 1) \ opnorm(Ar, 1) * opnorm(inv(Array(Ar)), 1) < 3
137-
@test 0.99 <= cond(Ac, 1) \ opnorm(Ac, 1) * opnorm(inv(Array(Ac)), 1) < 3
138-
@test 0.99 <= cond(Ar, Inf) \ opnorm(Ar, Inf) * opnorm(inv(Array(Ar)), Inf) < 3
139-
@test 0.99 <= cond(Ac, Inf) \ opnorm(Ac, Inf) * opnorm(inv(Array(Ac)), Inf) < 3
140133
@test_throws ArgumentError cond(A,2)
141134
@test_throws ArgumentError cond(A,3)
142135
Arect = spzeros(10, 6)
143136
@test_throws DimensionMismatch cond(Arect, 1)
144137
@test_throws ArgumentError cond(Arect,2)
145138
@test_throws DimensionMismatch cond(Arect, Inf)
139+
Ac = sprandn(20, 20,.5) + im*sprandn(20, 20,.5)
140+
Ar = sprandn(20, 20,.5) + eps()*I
141+
# For a discussion of the tolerance, see #14778
142+
@test 0.99 <= cond(Ar, 1) \ opnorm(Ar, 1) * opnorm(inv(Array(Ar)), 1) < 3
143+
@test 0.99 <= cond(Ac, 1) \ opnorm(Ac, 1) * opnorm(inv(Array(Ac)), 1) < 3
144+
@test 0.99 <= cond(Ar, Inf) \ opnorm(Ar, Inf) * opnorm(inv(Array(Ar)), Inf) < 3
145+
@test 0.99 <= cond(Ac, Inf) \ opnorm(Ac, Inf) * opnorm(inv(Array(Ac)), Inf) < 3
146+
#issue 680
147+
A22 = sparse(randn(2,2))
148+
@test 0.99 cond(Array(A22), 1) / cond(A22, 1) < 3
149+
@test 0.99 cond(Array(A22), Inf) / cond(A22, Inf) < 3
146150
end
147151

148152
@testset "sparse matrix opnormestinv" begin
@@ -158,6 +162,9 @@ end
158162
@test_throws ArgumentError SparseArrays.opnormestinv(Ac,0)
159163
@test_throws ArgumentError SparseArrays.opnormestinv(Ac,21)
160164
@test_throws DimensionMismatch SparseArrays.opnormestinv(sprand(3,5,.9))
165+
#issue 680
166+
A33 = sparse(randn(3,3))
167+
@test SparseArrays.opnormestinv(A33,3) opnorm(inv(Array(A33)),1) atol=1e-4
161168
end
162169

163170
@testset "factorization" begin

0 commit comments

Comments
 (0)