Refactor code to avoid out of bounds operation after deflation

src/svd.jl

@@ -167,10 +167,7 @@ function __svd!(
  e = B.ev
  iteration = 0
 
- # See LAWN3 page 6 and 22
- σ⁻ = estimate_σ⁻!(d, e, n1, n2, tol)
- fudge = n
- thresh = tol * σ⁻
+ thresh = tol * zero(first(d))
 
  if B.uplo === 'U'
  while true
@@ -199,17 +196,21 @@ function __svd!(
  break
  end
  end
+ @debug "Active submatrix" iteration n1 n2 d[n1] d[n2] e[n1] e[n2 - 1]
 
  # Deflation check. See LAWN3 p21
  # The Demmel-Kahan iteration moves the zero to the end and produces a
  # zero super diagonal element as well
  for i = n1:n2
- if d[i] == 0
+ if B.dv[i] == 0
  @debug "Deflation! Exact zero in diagonal element" i
  svdDemmelKahan!(B, n1, n2, U, Vᴴ)
 
  # We have now moved a zero to the end so the problem is one smaller
  n2 -= 1
+
+ # We'll start over to find the relevant submatrix after deflation
+ @goto top
  end
  end
 
@@ -225,7 +226,7 @@ function __svd!(
  fudge = n2 - n1 + 1
  thresh = tol * σ⁻
 
- @debug "__svd!" iteration n1 n2 d[n1] d[n2] e[n1] e[n2-1] thresh
+ @debug "estimated quantities" σ⁻ σ⁺ fudge thresh
 
  if fudge * tol * σ⁻ <= eps(σ⁺)
  svdDemmelKahan!(B, n1, n2, U, Vᴴ)

test/svd.jl

@@ -248,4 +248,9 @@ using Test, GenericLinearAlgebra, LinearAlgebra, Quaternions, DoubleFloats
  @test F.U == I
  @test F.Vt == I
  end
+
+ @testset "Issue 121" begin
+ @test svdvals(BigFloat[0 0; 1 -1]) ≈ [sqrt(2), 0]
+ @test svdvals(BigFloat[1 0 0; 0 0 0; 0 1 -1]) ≈ [sqrt(2), 1, 0]
+ end
 end