Fix bug where number type in stopping criterion was derived from `tol…

…` instead of real(eltype(H)) (#137)

Previously ArnoldiMethod.jl would accidentially pick `eps(typeof(tol))` instead `eps(real(typeof(residual)))`

So, with bigfloat and tol = 1e-30 you would actually just get eps(Float64) accuracy, whereas `tol = big(1e-30)` would have given the requested tolerance.

This fixes that.

src/run.jl

@@ -127,17 +127,19 @@ be scale invariant: the matrix `B = αA` for some constant `α` has the same
 eigenvectors with eigenvalue λα, so this scaling with `α` cancels in the 
 inequality.
 """
-struct IsConverged{RV<:RitzValues,T}
-    ritz::RV
-    tol::T
-    H_frob_norm::RefValue{T}
+struct IsConverged{Tv,Tr,Tt}
+    ritz::RitzValues{Tv,Tr}
+    tol::Tt
+    H_frob_norm::RefValue{Tr}
 
-    IsConverged(ritz::R, tol::T) where {R<:RitzValues,T} =
-        new{R,T}(ritz, tol, RefValue(zero(T)))
+    IsConverged(ritz::RitzValues{Tv,Tr}, tol::Tt) where {Tv,Tr,Tt} =
+        new{Tv,Tr,Tt}(ritz, tol, RefValue(zero(Tr)))
+end
+
+function (r::IsConverged{Tv,Tr})(i::Integer) where {Tv,Tr}
+    @inbounds r.ritz.rs[i] <= max(eps(Tr) * r.H_frob_norm[], r.tol * abs(r.ritz.λs[i]))
 end
 
-(r::IsConverged{RV,T})(i::Integer) where {RV,T} = @inbounds return r.ritz.rs[i] <=
-                 max(eps(T) * r.H_frob_norm[], r.tol * abs(r.ritz.λs[i]))
 
 """
     History(mvproducts, nconverged, converged, nev)

test/partial_schur.jl

@@ -26,6 +26,18 @@ using LinearAlgebra
     end
 end
 
+@testset "Stopping criterion specified in different number type is fine" begin
+    A = spdiagm(
+        -1 => fill(big(-1.0), 99),
+        0 => fill(big(2.0), 100),
+        1 => fill(big(-1.0), 99),
+    )
+    tol = 1e-30 # specified as Float64
+    schur, history = partialschur(A, tol = tol, maxdim = 30, nev = 2)
+    @test history.converged
+    @test norm(A * schur.Q - schur.Q * schur.R) < size(A, 1) * tol
+end
+
 @testset "Right number type" begin
     A = [rand(Bool) ? 1 : 0 for i = 1:10, j = 1:10]
     @inferred partialschur(A, nev = 2, mindim = 3, maxdim = 8)