From e2d474ea3a405e0d1b42331ff1475af812eb2a51 Mon Sep 17 00:00:00 2001
From: Romeo Valentin <romeov@stanford.edu>
Date: Fri, 20 Oct 2023 17:12:45 -0700
Subject: [PATCH 1/2] Add retrying normalization if necessary.

See #25.
---
 Project.toml |  2 ++
 src/svd.jl   | 16 ++++++++++++++++
 2 files changed, 18 insertions(+)

diff --git a/Project.toml b/Project.toml
index c051292..bd26cb4 100644
--- a/Project.toml
+++ b/Project.toml
@@ -4,7 +4,9 @@ version = "0.4.3"
 
 [deps]
 Adapt = "79e6a3ab-5dfb-504d-930d-738a2a938a0e"
+KrylovKit = "0b1a1467-8014-51b9-945f-bf0ae24f4b77"
 LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
+Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c"
 
 [compat]
 Adapt = "2.1, 3.0"
diff --git a/src/svd.jl b/src/svd.jl
index 881add7..3544dc9 100644
--- a/src/svd.jl
+++ b/src/svd.jl
@@ -1,3 +1,5 @@
+using KrylovKit
+using Random: rand!
 function _biLanczosIterations!(A, stepsize, αs, βs, U, V, μs, νs, maxνs, maxμs, τ, reorth_in, tolreorth, debug)
 
     m, n = size(A)
@@ -50,6 +52,20 @@ function _biLanczosIterations!(A, stepsize, αs, βs, U, V, μs, νs, maxνs, ma
                 nReorthVecs += 1
             end
             α = norm(v)
+
+            if α < τ * maximum(size(A)) * eps()  # orthogonalization failed, see https://github.com/poulson/PROPACK/blob/2465f89d5b1fba56de71c3e69e27d017c3dc2295/double/dlanbpro.F#L384
+                @warn "Restart orthogonalization"
+                b = V[1:(j-1)]
+                B = KrylovKit.OrthonormalBasis(b ./ norm.(b))
+                for _ in 1:3
+                    v = rand!(v) * 2 .- 1
+                    KrylovKit.orthonormalize!(v, B, KrylovKit.ModifiedGramSchmidt())
+                    α = norm(v)
+                    if !(α < τ * maximum(size(A)) * eps())
+                        break
+                    end
+                end
+            end
         end
 
         ## update the result vectors

From d3cc117f0675bcb089a14a744c5f80223a510b77 Mon Sep 17 00:00:00 2001
From: Romeo Valentin <romeov@stanford.edu>
Date: Fri, 3 Nov 2023 02:56:34 -0700
Subject: [PATCH 2/2] Remove @warn, instead print when debug is true.

---
 src/svd.jl | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/src/svd.jl b/src/svd.jl
index 3544dc9..44abcd7 100644
--- a/src/svd.jl
+++ b/src/svd.jl
@@ -54,7 +54,7 @@ function _biLanczosIterations!(A, stepsize, αs, βs, U, V, μs, νs, maxνs, ma
             α = norm(v)
 
             if α < τ * maximum(size(A)) * eps()  # orthogonalization failed, see https://github.com/poulson/PROPACK/blob/2465f89d5b1fba56de71c3e69e27d017c3dc2295/double/dlanbpro.F#L384
-                @warn "Restart orthogonalization"
+                debug && println("Restart orthogonalization")
                 b = V[1:(j-1)]
                 B = KrylovKit.OrthonormalBasis(b ./ norm.(b))
                 for _ in 1:3