@@ -2,8 +2,16 @@ using IterativeSolvers
2
2
using LinearMaps
3
3
using Base. Test
4
4
5
+ import Base. A_ldiv_B!
6
+
5
7
include (" laplace_matrix.jl" )
6
8
9
+ struct JacobiPrec{TD}
10
+ diagonal:: TD
11
+ end
12
+
13
+ A_ldiv_B! (y, P:: JacobiPrec , x) = y .= x ./ P. diagonal
14
+
7
15
function max_err (R)
8
16
r = zeros (real (eltype (R)), size (R, 2 ))
9
17
for j in 1 : length (r)
18
26
@testset " Locally Optimal Block Preconditioned Conjugate Gradient" begin
19
27
srand (1234321 )
20
28
@testset " Single eigenvalue" begin
29
+ n = 10
21
30
@testset " Small full system" begin
22
- n = 10
23
31
@testset " Simple eigenvalue problem" begin
24
32
@testset " Matrix{$T }" for T in (Float32, Float64, Complex64, Complex128)
25
33
@testset " largest = $largest " for largest in (true , false )
59
67
end
60
68
end
61
69
end
62
-
63
70
@testset " Sparse Laplacian" begin
64
71
A = laplace_matrix (Float64, 20 , 2 )
65
72
rhs = randn (size (A, 2 ), 1 )
74
81
end
75
82
end
76
83
end
84
+ @testset " Zero initial solution" begin
85
+ @testset " Simple eigenvalue problem" begin
86
+ @testset " Matrix{$T }" for T in (Float32, Float64, Complex64, Complex128)
87
+ @testset " largest = $largest " for largest in (true , false )
88
+ A = rand (T, n, n)
89
+ A = A' * A + I
90
+ b = zeros (T, n, 1 )
91
+ tol = √ eps (real (T))
92
+
93
+ r = lobpcg (A, largest, b; tol= tol, maxiter= Inf , log= false )
94
+ λ, X = r. λ, r. X
95
+ @test norm (A* X - X* λ) ≤ tol
96
+ end
97
+ end
98
+ end
99
+ @testset " Generalized eigenvalue problem" begin
100
+ @testset " Matrix{$T }" for T in (Float32, Float64, Complex64, Complex128)
101
+ @testset " largest = $largest " for largest in (true , false )
102
+ A = rand (T, n, n)
103
+ A = A' * A + I
104
+ B = rand (T, n, n)
105
+ B = B' * B + I
106
+ b = zeros (T, n, 1 )
107
+ tol = √ eps (real (T))
108
+
109
+ r = lobpcg (A, B, largest, b; tol= tol, maxiter= Inf , log= true )
110
+ λ, X = r. λ, r. X
111
+ @test max_err (A* X - B* X* λ) ≤ tol
112
+ end
113
+ end
114
+ end
115
+ end
116
+ @testset " No initial solution" begin
117
+ @testset " Simple eigenvalue problem" begin
118
+ @testset " Matrix{$T }" for T in (Float32, Float64, Complex64, Complex128)
119
+ @testset " largest = $largest " for largest in (true , false )
120
+ A = rand (T, n, n)
121
+ A = A' * A + I
122
+ tol = √ eps (real (T))
123
+
124
+ r = lobpcg (A, largest, 1 ; tol= tol, maxiter= Inf , log= false )
125
+ λ, X = r. λ, r. X
126
+ @test norm (A* X - X* λ) ≤ tol
127
+ end
128
+ end
129
+ end
130
+ @testset " Generalized eigenvalue problem" begin
131
+ @testset " Matrix{$T }" for T in (Float32, Float64, Complex64, Complex128)
132
+ @testset " largest = $largest " for largest in (true , false )
133
+ A = rand (T, n, n)
134
+ A = A' * A + I
135
+ B = rand (T, n, n)
136
+ B = B' * B + I
137
+ tol = √ eps (real (T))
138
+
139
+ r = lobpcg (A, B, largest, 1 ; tol= tol, maxiter= Inf , log= true )
140
+ λ, X = r. λ, r. X
141
+ @test max_err (A* X - B* X* λ) ≤ tol
142
+ end
143
+ end
144
+ end
145
+ end
146
+ @testset " Inplace" begin
147
+ @testset " Simple eigenvalue problem" begin
148
+ @testset " Matrix{$T }" for T in (Float32, Float64, Complex64, Complex128)
149
+ @testset " largest = $largest " for largest in (true , false )
150
+ A = rand (T, n, n)
151
+ A = A' * A + I
152
+ tol = √ eps (real (T))
153
+ b = rand (T, n, 1 )
154
+ itr = LOBPCGIterator (A, b, largest)
155
+
156
+ r = lobpcg! (itr; tol= tol, maxiter= Inf , log= false )
157
+ λ, X = r. λ, r. X
158
+ @test norm (A* X - X* λ) ≤ tol
159
+ end
160
+ end
161
+ end
162
+ @testset " Generalized eigenvalue problem" begin
163
+ @testset " Matrix{$T }" for T in (Float32, Float64, Complex64, Complex128)
164
+ @testset " largest = $largest " for largest in (true , false )
165
+ A = rand (T, n, n)
166
+ A = A' * A + I
167
+ B = rand (T, n, n)
168
+ B = B' * B + I
169
+ b = rand (T, n, 1 )
170
+ tol = √ eps (real (T))
171
+ itr = LOBPCGIterator (A, B, b, largest)
172
+
173
+ r = lobpcg! (itr; tol= tol, maxiter= Inf , log= true )
174
+ λ, X = r. λ, r. X
175
+ @test max_err (A* X - B* X* λ) ≤ tol
176
+ end
177
+ end
178
+ end
179
+ end
180
+ @testset " Jacobi preconditioner" begin
181
+ @testset " Simple eigenvalue problem" begin
182
+ @testset " Matrix{$T }" for T in (Float32, Float64, Complex64, Complex128)
183
+ @testset " largest = $largest " for largest in (true , false )
184
+ A = rand (T, n, n)
185
+ A = A' * A + I
186
+ tol = √ eps (real (T))
187
+ P = JacobiPrec (diag (A))
188
+ r = lobpcg (A, largest, 1 ; P= P, tol= tol, maxiter= Inf , log= false )
189
+ λ, X = r. λ, r. X
190
+ @test norm (A* X - X* λ) ≤ tol
191
+ end
192
+ end
193
+ end
194
+ @testset " Generalized eigenvalue problem" begin
195
+ @testset " Matrix{$T }" for T in (Float32, Float64, Complex64, Complex128)
196
+ @testset " largest = $largest " for largest in (true , false )
197
+ A = rand (T, n, n)
198
+ A = A' * A + I
199
+ P = JacobiPrec (diag (A))
200
+ B = rand (T, n, n)
201
+ B = B' * B + I
202
+ tol = √ eps (real (T))
203
+
204
+ r = lobpcg (A, B, largest, 1 ; P= P, tol= tol, maxiter= Inf , log= true )
205
+ λ, X = r. λ, r. X
206
+ @test max_err (A* X - B* X* λ) ≤ tol
207
+ end
208
+ end
209
+ end
210
+ end
211
+
212
+ @testset " Constraint" begin
213
+ @testset " Simple eigenvalue problem" begin
214
+ @testset " Matrix{$T }" for T in (Float32, Float64, Complex64, Complex128)
215
+ @testset " largest = $largest " for largest in (true , false )
216
+ A = rand (T, n, n)
217
+ A = A' * A + I
218
+ tol = √ eps (real (T))
219
+ r = lobpcg (A, largest, 1 ; tol= tol, maxiter= Inf , log= false )
220
+ λ1, X1 = r. λ, r. X
221
+ r = lobpcg (A, largest, 1 ; C= copy (r. X), tol= tol, maxiter= Inf , log= false )
222
+ λ2, X2 = r. λ, r. X
223
+ @test norm (A* X2 - X2* λ2) ≤ tol
224
+ @test isapprox (real (Ac_mul_B (X1, X2)[1 ,1 ]), 0 , atol= n* tol)
225
+ end
226
+ end
227
+ end
228
+ @testset " Generalized eigenvalue problem" begin
229
+ @testset " Matrix{$T }" for T in (Float32, Float64, Complex64, Complex128)
230
+ @testset " largest = $largest " for largest in (true , false )
231
+ A = rand (T, n, n)
232
+ A = A' * A + I
233
+ B = rand (T, n, n)
234
+ B = B' * B + I
235
+ tol = √ eps (real (T))
236
+ r = lobpcg (A, B, largest, 1 ; tol= tol, maxiter= Inf , log= false )
237
+ λ1, X1 = r. λ, r. X
238
+ r = lobpcg (A, B, largest, 1 ; C= copy (r. X), tol= tol, maxiter= Inf , log= false )
239
+ λ2, X2 = r. λ, r. X
240
+ @test norm (A* X2 - B* X2* λ2) ≤ tol
241
+ @test isapprox (real (Ac_mul_B (X1, B* X2)[1 ,1 ]), 0 , atol= n* tol)
242
+ end
243
+ end
244
+ end
245
+ end
77
246
end
78
247
@testset " Two eigenvalues" begin
79
248
@testset " Small full system" begin
104
273
B = rand (T, n, n)
105
274
B = B' * B + I
106
275
b = rand (T, n, 2 )
107
- tol = √ eps (real (T))
108
-
276
+ tol = eps (real (T))^ (real (T)(4 / 10 ))
109
277
r = lobpcg (A, B, largest, b; tol= tol, maxiter= Inf , log= true )
110
278
λ, X = r. λ, r. X
111
279
@test max_err (A* X - B* X* diagm (λ)) ≤ tol
0 commit comments