Skip to content

Commit

Permalink
Refactor tests II
Browse files Browse the repository at this point in the history
  • Loading branch information
@lkdvos
lkdvos committed Mar 7, 2024
1 parent d31b115 commit 7a27157
Show file tree
Hide file tree
Showing 9 changed files with 413 additions and 1,058 deletions.
285 changes: 85 additions & 200 deletions test/eigsolve.jl
View file

Large diffs are not rendered by default.

196 changes: 54 additions & 142 deletions test/expintegrator.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -12,17 +12,24 @@ function ϕ(A, v, p)
return exp(A′)[1:m, end]
end

@testset "Lanczos - expintegrator full" begin
@testset for T in (Float32, Float64, ComplexF32, ComplexF64)
@testset for orth in (cgs2, mgs2, cgsr, mgsr)
@testset "Lanczos - expintegrator full ($mode)" for mode in (:vector, :inplace, :outplace)
scalartypes = mode === :vector ? (Float32, Float64, ComplexF32, ComplexF64) :
(ComplexF64,)
orths = mode === :vector ? (cgs2, mgs2, cgsr, mgsr) : (mgsr,)
@testset for T in scalartypes
@testset for orth in orths
A = rand(T, (n, n)) .- one(T) / 2
A = (A + A') / 2
V = one(A)
W = zero(A)
alg = Lanczos(; orth=orth, krylovdim=n, maxiter=2, tol=tolerance(T),
verbosity=2)
for k in 1:n
W[:, k] = first(@constinferred exponentiate(A, 1, view(V, :, k), alg))
W[:, k] = unwrapvec(first(@constinferred exponentiate(wrapop(A, Val(mode)),
1,
wrapvec(view(V, :, k),
Val(mode)),
alg)))
end
@test W ≈ exp(A)

Expand All @@ -33,60 +40,37 @@ end
randn(real(T)) + im * randn(real(T)))
for p in 1:pmax
u = ntuple(i -> rand(T, n), p + 1)
w, info = @constinferred expintegrator(A, t, u, alg)
w, info = @constinferred expintegrator(wrapop(A, Val(mode)), t,
wrapvec.(u, Ref(Val(mode))), alg)
w2 = exp(t * A) * u[1]
for j in 1:p
w2 .+= t^j * ϕ(t * A, u[j + 1], j)
end
@test info.converged > 0
@test w2 ≈ w
@test w2 ≈ unwrapvec(w)
end
end
end
end

@testset "MinimalVec{$IP}" for IP in (true, false)
T = ComplexF64
A = rand(T, (n, n)) .- one(T) / 2
A = (A + A') / 2
V = one(A)
W = zero(A)
alg = Lanczos(; krylovdim=n, maxiter=2, tol=tolerance(T), verbosity=2)
for k in 1:n
W[:, k] = unwrap(first(@constinferred exponentiate(wrapop(A), 1,
MinimalVec{IP}(view(V, :, k)),
alg)))
end
@test W ≈ exp(A)

pmax = 5
alg = Lanczos(; krylovdim=n, maxiter=2, tol=tolerance(T), verbosity=1)
for t in (rand(real(T)), -rand(real(T)), im * randn(real(T)),
randn(real(T)) + im * randn(real(T)))
for p in 1:pmax
u = ntuple(i -> MinimalVec{IP}(rand(T, n)), p + 1)
w, info = @constinferred expintegrator(wrapop(A), t, u, alg)
w2 = exp(t * A) * unwrap(u[1])
for j in 1:p
w2 .+= t^j * ϕ(t * A, unwrap(u[j + 1]), j)
end
@test info.converged > 0
@test w2 ≈ unwrap(w)
end
end
end
end

@testset "Arnoldi - expintegrator full" begin
@testset for T in (Float32, Float64, ComplexF32, ComplexF64)
@testset for orth in (cgs2, mgs2, cgsr, mgsr)
@testset "Arnoldi - expintegrator full ($mode)" for mode in (:vector, :inplace, :outplace)
scalartypes = mode === :vector ? (Float32, Float64, ComplexF32, ComplexF64) :
(ComplexF64,)
orths = mode === :vector ? (cgs2, mgs2, cgsr, mgsr) : (mgsr,)
@testset for T in scalartypes
@testset for orth in orths
A = rand(T, (n, n)) .- one(T) / 2
V = one(A)
W = zero(A)
alg = Arnoldi(; orth=orth, krylovdim=n, maxiter=2, tol=tolerance(T),
verbosity=2)
for k in 1:n
W[:, k] = first(@constinferred exponentiate(A, 1, view(V, :, k), alg))
W[:, k] = unwrapvec(first(@constinferred exponentiate(wrapop(A, Val(mode)),
1,
wrapvec(view(V, :, k),
Val(mode)),
alg)))
end
@test W ≈ exp(A)

Expand All @@ -97,52 +81,27 @@ end
randn(real(T)) + im * randn(real(T)))
for p in 1:pmax
u = ntuple(i -> rand(T, n), p + 1)
w, info = @constinferred expintegrator(A, t, u, alg)
w, info = @constinferred expintegrator(wrapop(A, Val(mode)), t,
wrapvec.(u, Ref(Val(mode))), alg)
w2 = exp(t * A) * u[1]
for j in 1:p
w2 .+= t^j * ϕ(t * A, u[j + 1], j)
end
@test info.converged > 0
@test w2 ≈ w
end
end
end
end

@testset "MinimalVec{$IP}" for IP in (true, false)
T = ComplexF64
A = rand(T, (n, n)) .- one(T) / 2
V = one(A)
W = zero(A)
alg = Arnoldi(; krylovdim=n, maxiter=2, tol=tolerance(T), verbosity=2)
for k in 1:n
W[:, k] = unwrap(first(@constinferred exponentiate(wrapop(A), 1,
MinimalVec{IP}(view(V, :, k)),
alg)))
end
@test W ≈ exp(A)

pmax = 5
alg = Arnoldi(; krylovdim=n, maxiter=2, tol=tolerance(T), verbosity=1)
for t in (rand(real(T)), -rand(real(T)), im * randn(real(T)),
randn(real(T)) + im * randn(real(T)))
for p in 1:pmax
u = ntuple(i -> MinimalVec{IP}(rand(T, n)), p + 1)
w, info = @constinferred expintegrator(wrapop(A), t, u, alg)
w2 = exp(t * A) * unwrap(u[1])
for j in 1:p
w2 .+= t^j * ϕ(t * A, unwrap(u[j + 1]), j)
@test w2 ≈ unwrapvec(w)
end
@test info.converged > 0
@test w2 ≈ unwrap(w)
end
end
end
end

@testset "Lanczos - expintegrator iteratively" begin
@testset for T in (Float32, Float64, ComplexF32, ComplexF64)
@testset for orth in (cgs2, mgs2, cgsr, mgsr)
@testset "Lanczos - expintegrator iteratively ($mode)" for mode in
(:vector, :inplace, :outplace)
scalartypes = mode === :vector ? (Float32, Float64, ComplexF32, ComplexF64) :
(ComplexF64,)
orths = mode === :vector ? (cgs2, mgs2, cgsr, mgsr) : (mgsr,)
@testset for T in scalartypes
@testset for orth in orths
A = rand(T, (N, N)) .- one(T) / 2
A = (A + A') / 2
s = norm(eigvals(A), 1)
Expand All @@ -152,56 +111,34 @@ end
randn(real(T)) + im * randn(real(T)))
for p in 1:pmax
u = ntuple(i -> rand(T, N), p + 1)
w1, info = @constinferred expintegrator(A, t, u...;
w1, info = @constinferred expintegrator(wrapop(A, Val(mode)), t,
wrapvec.(u, Ref(Val(mode)))...;
maxiter=100, krylovdim=n,
eager=true)
@assert info.converged > 0
w2 = exp(t * A) * u[1]
for j in 1:p
w2 .+= t^j * ϕ(t * A, u[j + 1], j)
end
@test w2 ≈ w1
w1, info = @constinferred expintegrator(A, t, u...;
@test w2 ≈ unwrapvec(w1)
w1, info = @constinferred expintegrator(wrapop(A, Val(mode)), t,
wrapvec.(u, Ref(Val(mode)))...;
maxiter=100, krylovdim=n,
tol=1e-3, eager=true)
@test w1 ≈ w2 atol = 1e-2 * abs(t)
end
end
end
end

@testset "MinimalVec{$IP}" for IP in (true, false)
T = ComplexF64
A = rand(T, (N, N)) .- one(T) / 2
A = (A + A') / 2
s = norm(eigvals(A), 1)
rmul!(A, 1 / (10 * s))
pmax = 5
for t in (rand(real(T)), -rand(real(T)), im * randn(real(T)),
randn(real(T)) + im * randn(real(T)))
for p in 1:pmax
u = ntuple(i -> MinimalVec{IP}(rand(T, N)), p + 1)
w1, info = @constinferred expintegrator(wrapop(A), t, u...;
maxiter=100, krylovdim=n,
eager=true)
@assert info.converged > 0
w2 = exp(t * A) * unwrap(u[1])
for j in 1:p
w2 .+= t^j * ϕ(t * A, unwrap(u[j + 1]), j)
@test unwrapvec(w1) ≈ w2 atol = 1e-2 * abs(t)
end
@test w2 ≈ unwrap(w1)
w1, info = @constinferred expintegrator(wrapop(A), t, u...;
maxiter=100, krylovdim=n,
tol=1e-3, eager=true)
@test unwrap(w1) ≈ w2 atol = 1e-2 * abs(t)
end
end
end
end

@testset "Arnoldi - expintegrator iteratively" begin
@testset for T in (Float32, Float64, ComplexF32, ComplexF64)
@testset for orth in (cgs2, mgs2, cgsr, mgsr)
@testset "Arnoldi - expintegrator iteratively ($mode)" for mode in
(:vector, :inplace, :outplace)
scalartypes = mode === :vector ? (Float32, Float64, ComplexF32, ComplexF64) :
(ComplexF64,)
orths = mode === :vector ? (cgs2, mgs2, cgsr, mgsr) : (mgsr,)
@testset for T in scalartypes
@testset for orth in orths
A = rand(T, (N, N)) .- one(T) / 2
s = norm(eigvals(A), 1)
rmul!(A, 1 / (10 * s))
Expand All @@ -210,47 +147,22 @@ end
randn(real(T)) + im * randn(real(T)))
for p in 1:pmax
u = ntuple(i -> rand(T, N), p + 1)
w1, info = @constinferred expintegrator(A, t, u...;
w1, info = @constinferred expintegrator(wrapop(A, Val(mode)), t,
wrapvec.(u, Ref(Val(mode)))...;
maxiter=100, krylovdim=n,
eager=true)
@test info.converged > 0
w2 = exp(t * A) * u[1]
for j in 1:p
w2 .+= t^j * ϕ(t * A, u[j + 1], j)
end
@test w2 ≈ w1
w1, info = @constinferred expintegrator(A, t, u...;
@test w2 ≈ unwrapvec(w1)
w1, info = @constinferred expintegrator(wrapop(A, Val(mode)), t,
wrapvec.(u, Ref(Val(mode)))...;
maxiter=100, krylovdim=n,
tol=1e-3, eager=true)
@test w1 ≈ w2 atol = 1e-2 * abs(t)
end
end
end
end

@testset "MinimalVec{$IP}" for IP in (true, false)
T = ComplexF64
A = rand(T, (N, N)) .- one(T) / 2
s = norm(eigvals(A), 1)
rmul!(A, 1 / (10 * s))
pmax = 5
for t in (rand(real(T)), -rand(real(T)), im * randn(real(T)),
randn(real(T)) + im * randn(real(T)))
for p in 1:pmax
u = ntuple(i -> MinimalVec{IP}(rand(T, N)), p + 1)
w1, info = @constinferred expintegrator(wrapop(A), t, u...;
maxiter=100, krylovdim=n,
eager=true)
@test info.converged > 0
w2 = exp(t * A) * unwrap(u[1])
for j in 1:p
w2 .+= t^j * ϕ(t * A, unwrap(u[j + 1]), j)
@test unwrapvec(w1) ≈ w2 atol = 1e-2 * abs(t)
end
@test w2 ≈ unwrap(w1)
w1, info = @constinferred expintegrator(wrapop(A), t, u...;
maxiter=100, krylovdim=n,
tol=1e-3, eager=true)
@test unwrap(w1) ≈ w2 atol = 1e-2 * abs(t)
end
end
end
Expand Down
Loading

0 comments on commit 7a27157

Please sign in to comment.