Skip to content

Commit 42a435f

Browse files
committed
Tolerance increase, test matrix change, and system test tolerance
1 parent 750b7ff commit 42a435f

File tree

2 files changed

+60
-66
lines changed

2 files changed

+60
-66
lines changed

src/lobpcg.jl

Lines changed: 14 additions & 15 deletions
Original file line numberDiff line numberDiff line change
@@ -743,6 +743,8 @@ function (iterator::LOBPCGIterator{Generalized})(residualTolerance, log) where {
743743
end
744744
end
745745

746+
default_tolerance(::Type{T}) where {T<:Number} = eps(real(T))^(real(T)(3)/10)
747+
746748
"""
747749
The Locally Optimal Block Preconditioned Conjugate Gradient Method (LOBPCG)
748750
@@ -771,7 +773,7 @@ Finds the `nev` extremal eigenvalues and their corresponding eigenvectors satisf
771773
772774
- `maxiter`: maximum number of iterations; default is 200;
773775
774-
- `tol::Number`: tolerance to which residual vector norms must be under.
776+
- `tol::Real`: tolerance to which residual vector norms must be under.
775777
776778
# Output
777779
@@ -808,7 +810,7 @@ end
808810
809811
- `maxiter`: maximum number of iterations; default is 200;
810812
811-
- `tol::Number`: tolerance to which residual vector norms must be under.
813+
- `tol::Real`: tolerance to which residual vector norms must be under.
812814
813815
# Output
814816
@@ -818,10 +820,9 @@ function lobpcg(A, largest::Bool, X0; kwargs...)
818820
lobpcg(A, nothing, largest, X0; kwargs...)
819821
end
820822
function lobpcg(A, B, largest, X0;
821-
not_zeros=false, log=false, P=nothing,
822-
C=nothing, tol=nothing, maxiter=200)
823+
not_zeros=false, log=false, P=nothing, maxiter=200,
824+
C=nothing, tol::Real=default_tolerance(eltype(X0)))
823825
X = copy(X0)
824-
T = eltype(X)
825826
n = size(X, 1)
826827
sizeX = size(X, 2)
827828
sizeX > n && throw("X column dimension exceeds the row dimension")
@@ -848,15 +849,15 @@ end
848849
849850
- `maxiter`: maximum number of iterations; default is 200;
850851
851-
- `tol::Number`: tolerance to which residual vector norms must be under.
852+
- `tol::Real`: tolerance to which residual vector norms must be under.
852853
853854
# Output
854855
855856
- `results`: a `LOBPCGResults` struct. `r.λ` and `r.X` store the eigenvalues and eigenvectors.
856857
857858
"""
858-
function lobpcg!(iterator::LOBPCGIterator; log=false, tol=nothing, maxiter=200, not_zeros=false)
859-
T = eltype(iterator.XBlocks.block)
859+
function lobpcg!(iterator::LOBPCGIterator; log=false, maxiter=200, not_zeros=false,
860+
tol::Real=default_tolerance(eltype(iterator.XBlocks.block)))
860861
X = iterator.XBlocks.block
861862
iterator.constr!(iterator.XBlocks.block, iterator.tempXBlocks.block)
862863
if !not_zeros
@@ -869,10 +870,9 @@ function lobpcg!(iterator::LOBPCGIterator; log=false, tol=nothing, maxiter=200,
869870
end
870871
n = size(X, 1)
871872
sizeX = size(X, 2)
872-
residualTolerance = (tol isa Void) ? (eps(real(T)))^(real(T)(4)/10) : real(tol)
873873
iterator.iteration[] = 1
874874
while iterator.iteration[] <= maxiter
875-
state = iterator(residualTolerance, log)
875+
state = iterator(tol, log)
876876
if log
877877
push!(iterator.trace, state)
878878
end
@@ -881,7 +881,7 @@ function lobpcg!(iterator::LOBPCGIterator; log=false, tol=nothing, maxiter=200,
881881
end
882882
@inbounds iterator.λ .= view(iterator.ritz_values, 1:sizeX)
883883

884-
results = LOBPCGResults(iterator.λ, X, residualTolerance, iterator.residuals, iterator.iteration[], maxiter, all((x)->(norm(x)<=residualTolerance), view(iterator.residuals, 1:sizeX)), iterator.trace)
884+
results = LOBPCGResults(iterator.λ, X, tol, iterator.residuals, iterator.iteration[], maxiter, all((x)->(norm(x)<=tol), view(iterator.residuals, 1:sizeX)), iterator.trace)
885885

886886
return results
887887
end
@@ -910,7 +910,7 @@ lobpcg(A, [B,] largest, X0, nev; kwargs...) -> results
910910
911911
- `maxiter`: maximum number of iterations; default is 200;
912912
913-
- `tol::Number`: tolerance to which residual vector norms must be under.
913+
- `tol::Real`: tolerance to which residual vector norms must be under.
914914
915915
# Output
916916
@@ -920,9 +920,8 @@ function lobpcg(A, largest::Bool, X0, nev::Int; kwargs...)
920920
lobpcg(A, nothing, largest, X0, nev; kwargs...)
921921
end
922922
function lobpcg(A, B, largest::Bool, X0, nev::Int;
923-
not_zeros=false, log=false, P=nothing,
924-
C=nothing, tol=nothing, maxiter=200)
925-
T = eltype(X0)
923+
not_zeros=false, log=false, P=nothing, maxiter=200,
924+
C=nothing, tol::Real=default_tolerance(eltype(X0)))
926925
n = size(X0, 1)
927926
sizeX = size(X0, 2)
928927
nev > n && throw("Number of eigenvectors desired exceeds the row dimension.")

test/lobpcg.jl

Lines changed: 46 additions & 51 deletions
Original file line numberDiff line numberDiff line change
@@ -33,10 +33,9 @@ end
3333
@testset "Matrix{$T}" for T in (Float32, Float64, Complex64, Complex128)
3434
@testset "largest = $largest" for largest in (true, false)
3535
A = rand(T, n, n)
36-
A = A' * A + I
36+
A = A' + A + 20I
3737
b = rand(T, n, 1)
38-
tol = eps(real(T))
39-
38+
tol = IterativeSolvers.default_tolerance(T)
4039
r = lobpcg(A, largest, b; tol=tol, maxiter=Inf, log=false)
4140
λ, X = r.λ, r.X
4241
@test norm(A*X - X*λ) tol
@@ -51,12 +50,11 @@ end
5150
@testset "Matrix{$T}" for T in (Float32, Float64, Complex64, Complex128)
5251
@testset "largest = $largest" for largest in (true, false)
5352
A = rand(T, n, n)
54-
A = A' * A + I
53+
A = A' + A + 20I
5554
B = rand(T, n, n)
56-
B = B' * B + I
55+
B = B' + B + 20I
5756
b = rand(T, n, 1)
58-
tol = eps(real(T))
59-
57+
tol = IterativeSolvers.default_tolerance(T)
6058
r = lobpcg(A, B, largest, b; tol=tol, maxiter=Inf, log=true)
6159
λ, X = r.λ, r.X
6260
@test max_err(A*X - B*X*λ) tol
@@ -72,8 +70,7 @@ end
7270
A = laplace_matrix(Float64, 20, 2)
7371
rhs = randn(size(A, 2), 1)
7472
scale!(rhs, inv(norm(rhs)))
75-
tol = 1e-5
76-
73+
tol = IterativeSolvers.default_tolerance(Float64)
7774
@testset "Matrix" begin
7875
@testset "largest = $largest" for largest in (true, false)
7976
r = lobpcg(A, largest, rhs; tol=tol, maxiter=Inf)
@@ -87,10 +84,9 @@ end
8784
@testset "Matrix{$T}" for T in (Float32, Float64, Complex64, Complex128)
8885
@testset "largest = $largest" for largest in (true, false)
8986
A = rand(T, n, n)
90-
A = A' * A + I
87+
A = A' + A + 20I
9188
b = zeros(T, n, 1)
92-
tol = eps(real(T))
93-
89+
tol = IterativeSolvers.default_tolerance(T)
9490
r = lobpcg(A, largest, b; tol=tol, maxiter=Inf, log=false)
9591
λ, X = r.λ, r.X
9692
@test norm(A*X - X*λ) tol
@@ -101,11 +97,11 @@ end
10197
@testset "Matrix{$T}" for T in (Float32, Float64, Complex64, Complex128)
10298
@testset "largest = $largest" for largest in (true, false)
10399
A = rand(T, n, n)
104-
A = A' * A + I
100+
A = A' + A + 20I
105101
B = rand(T, n, n)
106-
B = B' * B + I
102+
B = B' + B + 20I
107103
b = zeros(T, n, 1)
108-
tol = eps(real(T))
104+
tol = IterativeSolvers.default_tolerance(T)
109105

110106
r = lobpcg(A, B, largest, b; tol=tol, maxiter=Inf, log=true)
111107
λ, X = r.λ, r.X
@@ -119,8 +115,8 @@ end
119115
@testset "Matrix{$T}" for T in (Float32, Float64, Complex64, Complex128)
120116
@testset "largest = $largest" for largest in (true, false)
121117
A = rand(T, n, n)
122-
A = A' * A + I
123-
tol = eps(real(T))
118+
A = A' + A + 20I
119+
tol = IterativeSolvers.default_tolerance(T)
124120

125121
r = lobpcg(A, largest, 1; tol=tol, maxiter=Inf, log=false)
126122
λ, X = r.λ, r.X
@@ -132,10 +128,10 @@ end
132128
@testset "Matrix{$T}" for T in (Float32, Float64, Complex64, Complex128)
133129
@testset "largest = $largest" for largest in (true, false)
134130
A = rand(T, n, n)
135-
A = A' * A + I
131+
A = A' + A + 20I
136132
B = rand(T, n, n)
137-
B = B' * B + I
138-
tol = eps(real(T))
133+
B = B' + B + 20I
134+
tol = IterativeSolvers.default_tolerance(T)
139135

140136
r = lobpcg(A, B, largest, 1; tol=tol, maxiter=Inf, log=true)
141137
λ, X = r.λ, r.X
@@ -149,8 +145,8 @@ end
149145
@testset "Matrix{$T}" for T in (Float32, Float64, Complex64, Complex128)
150146
@testset "largest = $largest" for largest in (true, false)
151147
A = rand(T, n, n)
152-
A = A' * A + I
153-
tol = eps(real(T))
148+
A = A' + A + 20I
149+
tol = IterativeSolvers.default_tolerance(T)
154150
b = rand(T, n, 1)
155151
itr = LOBPCGIterator(A, largest, b)
156152

@@ -164,11 +160,11 @@ end
164160
@testset "Matrix{$T}" for T in (Float32, Float64, Complex64, Complex128)
165161
@testset "largest = $largest" for largest in (true, false)
166162
A = rand(T, n, n)
167-
A = A' * A + I
163+
A = A' + A + 20I
168164
B = rand(T, n, n)
169-
B = B' * B + I
165+
B = B' + B + 20I
170166
b = rand(T, n, 1)
171-
tol = eps(real(T))
167+
tol = IterativeSolvers.default_tolerance(T)
172168
itr = LOBPCGIterator(A, B, largest, b)
173169

174170
r = lobpcg!(itr; tol=tol, maxiter=Inf, log=true)
@@ -183,8 +179,8 @@ end
183179
@testset "Matrix{$T}" for T in (Float32, Float64, Complex64, Complex128)
184180
@testset "largest = $largest" for largest in (true, false)
185181
A = rand(T, n, n)
186-
A = A' * A + I
187-
tol = eps(real(T))
182+
A = A' + A + 20I
183+
tol = IterativeSolvers.default_tolerance(T)
188184
P = JacobiPrec(diag(A))
189185
r = lobpcg(A, largest, 1; P=P, tol=tol, maxiter=Inf, log=false)
190186
λ, X = r.λ, r.X
@@ -196,11 +192,11 @@ end
196192
@testset "Matrix{$T}" for T in (Float32, Float64, Complex64, Complex128)
197193
@testset "largest = $largest" for largest in (true, false)
198194
A = rand(T, n, n)
199-
A = A' * A + I
195+
A = A' + A + 20I
200196
P = JacobiPrec(diag(A))
201197
B = rand(T, n, n)
202-
B = B' * B + I
203-
tol = eps(real(T))
198+
B = B' + B + 20I
199+
tol = IterativeSolvers.default_tolerance(T)
204200

205201
r = lobpcg(A, B, largest, 1; P=P, tol=tol, maxiter=Inf, log=true)
206202
λ, X = r.λ, r.X
@@ -214,8 +210,8 @@ end
214210
@testset "Matrix{$T}" for T in (Float32, Float64, Complex64, Complex128)
215211
@testset "largest = $largest" for largest in (true, false)
216212
A = rand(T, n, n)
217-
A = A' * A + I
218-
tol = eps(real(T))
213+
A = A' + A + 20I
214+
tol = IterativeSolvers.default_tolerance(T)
219215
r = lobpcg(A, largest, 1; tol=tol, maxiter=Inf, log=false)
220216
λ1, X1 = r.λ, r.X
221217
r = lobpcg(A, largest, 1; C=copy(r.X), tol=tol, maxiter=Inf, log=false)
@@ -229,10 +225,10 @@ end
229225
@testset "Matrix{$T}" for T in (Float32, Float64, Complex64, Complex128)
230226
@testset "largest = $largest" for largest in (true, false)
231227
A = rand(T, n, n)
232-
A = A' * A + I
228+
A = A' + A + 20I
233229
B = rand(T, n, n)
234-
B = B' * B + I
235-
tol = eps(real(T))^0.4
230+
B = B' + B + 20I
231+
tol = IterativeSolvers.default_tolerance(T)
236232
r = lobpcg(A, B, largest, 1; tol=tol, maxiter=Inf, log=false)
237233
λ1, X1 = r.λ, r.X
238234
r = lobpcg(A, B, largest, 1; C=copy(r.X), tol=tol, maxiter=Inf, log=false)
@@ -251,9 +247,9 @@ end
251247
@testset "Matrix{$T}" for T in (Float32, Float64, Complex64, Complex128)
252248
@testset "largest = $largest" for largest in (true, false)
253249
A = rand(T, n, n)
254-
A = A' * A + I
250+
A = A' + A + 20I
255251
b = rand(T, n, 2)
256-
tol = eps(real(T))
252+
tol = IterativeSolvers.default_tolerance(T)
257253

258254
r = lobpcg(A, largest, b; tol=tol, maxiter=Inf, log=false)
259255
λ, X = r.λ, r.X
@@ -269,11 +265,11 @@ end
269265
@testset "Matrix{$T}" for T in (Float32, Float64, Complex64, Complex128)
270266
@testset "largest = $largest" for largest in (true, false)
271267
A = rand(T, n, n)
272-
A = A' * A + I
268+
A = A' + A + 20I
273269
B = rand(T, n, n)
274-
B = B' * B + I
270+
B = B' + B + 20I
275271
b = rand(T, n, 2)
276-
tol = eps(real(T))^(real(T)(4/10))
272+
tol = IterativeSolvers.default_tolerance(T)
277273
r = lobpcg(A, B, largest, b; tol=tol, maxiter=Inf, log=true)
278274
λ, X = r.λ, r.X
279275
@test max_err(A*X - B*X*diagm(λ)) tol
@@ -292,8 +288,8 @@ end
292288
@testset "Matrix{$T}" for T in (Float32, Float64, Complex64, Complex128)
293289
@testset "largest = $largest" for largest in (true, false)
294290
A = rand(T, n, n)
295-
A = A' * A + I
296-
tol = eps(real(T))^0.4
291+
A = A' + A + 20I
292+
tol = IterativeSolvers.default_tolerance(T)
297293
X0 = rand(T, n, block_size)
298294
r = lobpcg(A, largest, X0, 3, tol=tol, maxiter=Inf, log=true)
299295
λ, X = r.λ, r.X
@@ -306,11 +302,10 @@ end
306302
@testset "Matrix{$T}" for T in (Float32, Float64, Complex64, Complex128)
307303
@testset "largest = $largest" for largest in (true, false)
308304
A = rand(T, n, n)
309-
A = A' * A + I
305+
A = A' + A + 20I
310306
B = rand(T, n, n)
311-
B = B' * B + I
312-
tol = eps(real(T))^0.4
313-
307+
B = B' + B + 20I
308+
tol = IterativeSolvers.default_tolerance(T)
314309
X0 = rand(T, n, block_size)
315310
r = lobpcg(A, B, largest, X0, 3, tol=tol, maxiter=Inf, log=true)
316311
λ, X = r.λ, r.X
@@ -324,8 +319,8 @@ end
324319
@testset "Matrix{$T}" for T in (Float32, Float64, Complex64, Complex128)
325320
@testset "largest = $largest" for largest in (true, false)
326321
A = rand(T, n, n)
327-
A = A' * A + I
328-
tol = eps(real(T))
322+
A = A' + A + 20I
323+
tol = IterativeSolvers.default_tolerance(T)
329324
r = lobpcg(A, largest, 1; tol=tol, maxiter=Inf, log=false)
330325
λ1, X1 = r.λ, r.X
331326

@@ -342,10 +337,10 @@ end
342337
@testset "Matrix{$T}" for T in (Float32, Float64, Complex64, Complex128)
343338
@testset "largest = $largest" for largest in (true, false)
344339
A = rand(T, n, n)
345-
A = A' * A + 2I
340+
A = A' + A + 20I
346341
B = rand(T, n, n)
347-
B = B' * B + 2I
348-
tol = eps(real(T))^0.4
342+
B = B' + B + 20I
343+
tol = IterativeSolvers.default_tolerance(T)
349344
r = lobpcg(A, B, largest, 1; tol=tol, maxiter=Inf, log=false)
350345
λ1, X1 = r.λ, r.X
351346

0 commit comments

Comments
 (0)