From 97e8fe0fbd997f8954578b2ff983e3bf4da1a79d Mon Sep 17 00:00:00 2001 From: Katharine Hyatt Date: Tue, 23 Dec 2025 12:15:05 +0100 Subject: [PATCH 1/4] Use TestSuite for eig --- ext/MatrixAlgebraKitGenericSchurExt.jl | 13 +- test/cuda/eig.jl | 108 --------------- test/eig.jl | 156 +++++---------------- test/runtests.jl | 10 +- test/testsuite/TestSuite.jl | 1 + test/testsuite/eig.jl | 184 +++++++++++++++++++++++++ 6 files changed, 237 insertions(+), 235 deletions(-) delete mode 100644 test/cuda/eig.jl create mode 100644 test/testsuite/eig.jl diff --git a/ext/MatrixAlgebraKitGenericSchurExt.jl b/ext/MatrixAlgebraKitGenericSchurExt.jl index 686ca87b..7ee0939b 100644 --- a/ext/MatrixAlgebraKitGenericSchurExt.jl +++ b/ext/MatrixAlgebraKitGenericSchurExt.jl @@ -13,7 +13,18 @@ MatrixAlgebraKit.initialize_output(::typeof(eig_full!), A::AbstractMatrix, ::GS_ MatrixAlgebraKit.initialize_output(::typeof(eig_vals!), A::AbstractMatrix, ::GS_QRIteration) = nothing function MatrixAlgebraKit.eig_full!(A::AbstractMatrix, DV, ::GS_QRIteration) - D, V = GenericSchur.eigen!(A) + D_, V_ = GenericSchur.eigen!(A) + D, V = DV + if !isnothing(D) + copyto!(D, Diagonal(D_)) + else + D = D_ + end + if !isnothing(V) + copyto!(V, V_) + else + V = V_ + end return Diagonal(D), V end diff --git a/test/cuda/eig.jl b/test/cuda/eig.jl deleted file mode 100644 index e40ee9b8..00000000 --- a/test/cuda/eig.jl +++ /dev/null @@ -1,108 +0,0 @@ -using MatrixAlgebraKit -using LinearAlgebra: Diagonal -using Test -using TestExtras -using StableRNGs -using CUDA -using MatrixAlgebraKit: TruncatedAlgorithm, diagview, norm - -include(joinpath("..", "utilities.jl")) - -BLASFloats = (Float32, Float64, ComplexF32, ComplexF64) - -@testset "eig_full! for T = $T" for T in BLASFloats - rng = StableRNG(123) - m = 54 - for alg in (CUSOLVER_Simple(), :CUSOLVER_Simple, CUSOLVER_Simple) - A = CuArray(randn(rng, T, m, m)) - Tc = complex(T) - - D, V = @constinferred eig_full(A; alg = ($alg)) - @test eltype(D) == eltype(V) == Tc - @test A * V ≈ V * D - - alg′ = @constinferred MatrixAlgebraKit.select_algorithm(eig_full!, A, $alg) - - Ac = similar(A) - D2, V2 = @constinferred eig_full!(copy!(Ac, A), (D, V), alg′) - @test D2 === D - @test V2 === V - @test A * V ≈ V * D - - Dc = @constinferred eig_vals(A, alg′) - @test eltype(Dc) == Tc - @test parent(D) ≈ Dc - end -end - -#= -@testset "eig_trunc! for T = $T" for T in (Float32, Float64, ComplexF32, ComplexF64) - rng = StableRNG(123) - m = 54 - for alg in (CUSOLVER_Simple(),) - A = CuArray(randn(rng, T, m, m)) - A *= A' # TODO: deal with eigenvalue ordering etc - # eigenvalues are sorted by ascending real component... - D₀ = sort!(eig_vals(A); by=abs, rev=true) - rmin = findfirst(i -> abs(D₀[end - i]) != abs(D₀[end - i - 1]), 1:(m - 2)) - r = length(D₀) - rmin - - D1, V1, ϵ1 = @constinferred eig_trunc(A; alg, trunc=truncrank(r)) - @test length(D1.diag) == r - @test A * V1 ≈ V1 * D1 - - s = 1 + sqrt(eps(real(T))) - trunc = trunctol(; atol=s * abs(D₀[r + 1])) - D2, V2, ϵ2 = @constinferred eig_trunc(A; alg, trunc) - @test length(diagview(D2)) == r - @test A * V2 ≈ V2 * D2 - - # trunctol keeps order, truncrank might not - # test for same subspace - @test V1 * ((V1' * V1) \ (V1' * V2)) ≈ V2 - @test V2 * ((V2' * V2) \ (V2' * V1)) ≈ V1 - end -end - -@testset "eig_trunc! specify truncation algorithm T = $T" for T in BLASFloats - rng = StableRNG(123) - m = 4 - atol = sqrt(eps(real(T))) - V = randn(rng, T, m, m) - D = Diagonal(real(T)[0.9, 0.3, 0.1, 0.01]) - A = V * D * inv(V) - alg = TruncatedAlgorithm(LAPACK_Simple(), truncrank(2)) - D2, V2, ϵ2 = @constinferred eig_trunc(A; alg) - @test diagview(D2) ≈ diagview(D)[1:2] - @test ϵ2 ≈ norm(diagview(D)[3:4]) atol = atol - @test_throws ArgumentError eig_trunc(A; alg, trunc = (; maxrank = 2)) - - alg = TruncatedAlgorithm(LAPACK_Simple(), truncerror(; atol = 0.2, p = 1)) - D3, V3, ϵ3 = @constinferred eig_trunc(A; alg) - @test diagview(D3) ≈ diagview(D)[1:2] - @test ϵ3 ≈ norm(diagview(D)[3:4]) atol = atol -end -=# - -@testset "eig for Diagonal{$T}" for T in BLASFloats - rng = StableRNG(123) - m = 54 - Ad = CuArray(randn(rng, T, m)) - A = Diagonal(Ad) - atol = sqrt(eps(real(T))) - - D, V = @constinferred eig_full(A) - @test D isa Diagonal{T} && size(D) == size(A) - @test V isa Diagonal{T} && size(V) == size(A) - @test A * V ≈ V * D - - D2 = @constinferred eig_vals(A) - @test D2 isa AbstractVector{T} && length(D2) == m - @test diagview(D) ≈ D2 - - #=A2 = Diagonal(T[0.9, 0.3, 0.1, 0.01]) - alg = TruncatedAlgorithm(DiagonalAlgorithm(), truncrank(2)) - D2, V2, ϵ2 = @constinferred eig_trunc(A2; alg) - @test diagview(D2) ≈ diagview(A2)[1:2] - @test ϵ2 ≈ norm(diagview(A2)[3:4]) atol = atol=# -end diff --git a/test/eig.jl b/test/eig.jl index d709bd20..01a3bd94 100644 --- a/test/eig.jl +++ b/test/eig.jl @@ -4,127 +4,45 @@ using TestExtras using StableRNGs using LinearAlgebra: Diagonal using MatrixAlgebraKit: TruncatedAlgorithm, diagview, norm +using CUDA, AMDGPU BLASFloats = (Float32, Float64, ComplexF32, ComplexF64) -GenericFloats = (Float16, BigFloat, Complex{BigFloat}) - -@testset "eig_full! for T = $T" for T in BLASFloats - rng = StableRNG(123) - m = 54 - for alg in (LAPACK_Simple(), LAPACK_Expert(), :LAPACK_Simple, LAPACK_Simple) - A = randn(rng, T, m, m) - Tc = complex(T) - - D, V = @constinferred eig_full(A; alg = ($alg)) - @test eltype(D) == eltype(V) == Tc - @test A * V ≈ V * D - - alg′ = @constinferred MatrixAlgebraKit.select_algorithm(eig_full!, A, $alg) - - Ac = similar(A) - D2, V2 = @constinferred eig_full!(copy!(Ac, A), (D, V), alg′) - @test D2 === D - @test V2 === V - @test A * V ≈ V * D - - Dc = @constinferred eig_vals(A, alg′) - @test eltype(Dc) == Tc - @test D ≈ Diagonal(Dc) +GenericFloats = (BigFloat, Complex{BigFloat}) + +@isdefined(TestSuite) || include("testsuite/TestSuite.jl") +using .TestSuite + +is_buildkite = get(ENV, "BUILDKITE", "false") == "true" + +m = 54 +for T in (BLASFloats..., GenericFloats...) + TestSuite.seed_rng!(123) + if T ∈ BLASFloats + if CUDA.functional() + TestSuite.test_eig(CuMatrix{T}, (m, m); test_trunc = false) + TestSuite.test_eig_algs(CuMatrix{T}, (m, m), (CUSOLVER_Simple(),)) + TestSuite.test_eig(Diagonal{T, CuVector{T}}, m) + TestSuite.test_eig_algs(Diagonal{T, CuVector{T}}, m, (DiagonalAlgorithm(),)) + end + #= not yet supported + if AMDGPU.functional() + TestSuite.test_eig(ROCMatrix{T}, (m, m); test_blocksize = false) + TestSuite.test_eig_algs(ROCMatrix{T}, (m, m), (ROCSOLVER_Simple(),)) + TestSuite.test_eig(Diagonal{T, ROCVector{T}}, m; test_blocksize = false) + TestSuite.test_eig_algs(Diagonal{T, ROCVector{T}}, m, (DiagonalAlgorithm(),)) + end=# end -end - -@testset "eig_trunc! for T = $T" for T in BLASFloats - rng = StableRNG(123) - m = 54 - for alg in (LAPACK_Simple(), LAPACK_Expert()) - A = randn(rng, T, m, m) - A *= A' # TODO: deal with eigenvalue ordering etc - # eigenvalues are sorted by ascending real component... - D₀ = sort!(eig_vals(A); by = abs, rev = true) - rmin = findfirst(i -> abs(D₀[end - i]) != abs(D₀[end - i - 1]), 1:(m - 2)) - r = length(D₀) - rmin - atol = sqrt(eps(real(T))) - - D1, V1, ϵ1 = @constinferred eig_trunc(A; alg, trunc = truncrank(r)) - @test length(diagview(D1)) == r - @test A * V1 ≈ V1 * D1 - @test ϵ1 ≈ norm(view(D₀, (r + 1):m)) atol = atol - - s = 1 + sqrt(eps(real(T))) - trunc = trunctol(; atol = s * abs(D₀[r + 1])) - D2, V2, ϵ2 = @constinferred eig_trunc(A; alg, trunc) - @test length(diagview(D2)) == r - @test A * V2 ≈ V2 * D2 - @test ϵ2 ≈ norm(view(D₀, (r + 1):m)) atol = atol - - s = 1 - sqrt(eps(real(T))) - trunc = truncerror(; atol = s * norm(@view(D₀[r:end]), 1), p = 1) - D3, V3, ϵ3 = @constinferred eig_trunc(A; alg, trunc) - @test length(diagview(D3)) == r - @test A * V3 ≈ V3 * D3 - @test ϵ3 ≈ norm(view(D₀, (r + 1):m)) atol = atol - - s = 1 - sqrt(eps(real(T))) - trunc = truncerror(; atol = s * norm(@view(D₀[r:end]), 1), p = 1) - D4, V4 = @constinferred eig_trunc_no_error(A; alg, trunc) - @test length(diagview(D4)) == r - @test A * V4 ≈ V4 * D4 - # trunctol keeps order, truncrank might not - # test for same subspace - @test V1 * ((V1' * V1) \ (V1' * V2)) ≈ V2 - @test V2 * ((V2' * V2) \ (V2' * V1)) ≈ V1 - @test V1 * ((V1' * V1) \ (V1' * V3)) ≈ V3 - @test V3 * ((V3' * V3) \ (V3' * V1)) ≈ V1 + if !is_buildkite + TestSuite.test_eig(T, (m, m)) + if T ∈ BLASFloats + LAPACK_EIG_ALGS = (LAPACK_Simple(), LAPACK_Expert()) + TestSuite.test_eig_algs(T, (m, m), LAPACK_EIG_ALGS) + elseif T ∈ GenericFloats + GS_EIG_ALGS = (GS_QRIteration(),) + TestSuite.test_eig_algs(T, (m, m), GS_EIG_ALGS) + end + AT = Diagonal{T, Vector{T}} + TestSuite.test_eig(AT, m) + TestSuite.test_eig_algs(AT, m, (DiagonalAlgorithm(),)) end end - -@testset "eig_trunc! specify truncation algorithm T = $T" for T in BLASFloats - rng = StableRNG(123) - m = 4 - atol = sqrt(eps(real(T))) - V = randn(rng, T, m, m) - D = Diagonal(real(T)[0.9, 0.3, 0.1, 0.01]) - A = V * D * inv(V) - alg = TruncatedAlgorithm(LAPACK_Simple(), truncrank(2)) - D2, V2, ϵ2 = @constinferred eig_trunc(A; alg) - @test diagview(D2) ≈ diagview(D)[1:2] - @test ϵ2 ≈ norm(diagview(D)[3:4]) atol = atol - @test_throws ArgumentError eig_trunc(A; alg, trunc = (; maxrank = 2)) - - alg = TruncatedAlgorithm(LAPACK_Simple(), truncerror(; atol = 0.2, p = 1)) - D3, V3, ϵ3 = @constinferred eig_trunc(A; alg) - @test diagview(D3) ≈ diagview(D)[1:2] - @test ϵ3 ≈ norm(diagview(D)[3:4]) atol = atol - - alg = TruncatedAlgorithm(LAPACK_Simple(), truncerror(; atol = 0.2, p = 1)) - D4, V4 = @constinferred eig_trunc_no_error(A; alg) - @test diagview(D4) ≈ diagview(D)[1:2] -end - -@testset "eig for Diagonal{$T}" for T in (BLASFloats..., GenericFloats...) - rng = StableRNG(123) - m = 54 - Ad = randn(rng, T, m) - A = Diagonal(Ad) - atol = sqrt(eps(real(T))) - - D, V = @constinferred eig_full(A) - @test D isa Diagonal{T} && size(D) == size(A) - @test V isa Diagonal{T} && size(V) == size(A) - @test A * V ≈ V * D - - D2 = @constinferred eig_vals(A) - @test D2 isa AbstractVector{T} && length(D2) == m - @test diagview(D) ≈ D2 - - A2 = Diagonal(T[0.9, 0.3, 0.1, 0.01]) - alg = TruncatedAlgorithm(DiagonalAlgorithm(), truncrank(2)) - D2, V2, ϵ2 = @constinferred eig_trunc(A2; alg) - @test diagview(D2) ≈ diagview(A2)[1:2] - @test ϵ2 ≈ norm(diagview(A2)[3:4]) atol = atol - - A3 = Diagonal(T[0.9, 0.3, 0.1, 0.01]) - alg = TruncatedAlgorithm(DiagonalAlgorithm(), truncrank(2)) - D3, V3 = @constinferred eig_trunc_no_error(A3; alg) - @test diagview(D3) ≈ diagview(A3)[1:2] -end diff --git a/test/runtests.jl b/test/runtests.jl index 7325d410..062801c0 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -16,9 +16,6 @@ if !is_buildkite @safetestset "Hermitian Eigenvalue Decomposition" begin include("eigh.jl") end - @safetestset "General Eigenvalue Decomposition" begin - include("eig.jl") - end @safetestset "Generalized Eigenvalue Decomposition" begin include("gen_eig.jl") end @@ -51,7 +48,6 @@ if !is_buildkite @safetestset "Hermitian Eigenvalue Decomposition" begin include("genericlinearalgebra/eigh.jl") end - end @safetestset "QR / LQ Decomposition" begin @@ -67,15 +63,15 @@ end @safetestset "Schur Decomposition" begin include("schur.jl") end +@safetestset "General Eigenvalue Decomposition" begin + include("eig.jl") +end using CUDA if CUDA.functional() @safetestset "CUDA SVD" begin include("cuda/svd.jl") end - @safetestset "CUDA General Eigenvalue Decomposition" begin - include("cuda/eig.jl") - end @safetestset "CUDA Hermitian Eigenvalue Decomposition" begin include("cuda/eigh.jl") end diff --git a/test/testsuite/TestSuite.jl b/test/testsuite/TestSuite.jl index deee4e12..6e13567c 100644 --- a/test/testsuite/TestSuite.jl +++ b/test/testsuite/TestSuite.jl @@ -74,5 +74,6 @@ include("lq.jl") include("polar.jl") include("projections.jl") include("schur.jl") +include("eig.jl") end diff --git a/test/testsuite/eig.jl b/test/testsuite/eig.jl new file mode 100644 index 00000000..feafb769 --- /dev/null +++ b/test/testsuite/eig.jl @@ -0,0 +1,184 @@ +using TestExtras +using MatrixAlgebraKit: TruncatedAlgorithm +using LinearAlgebra: I +using GenericSchur + +function test_eig(T::Type, sz; test_trunc = true, kwargs...) + summary_str = testargs_summary(T, sz) + return @testset "eig $summary_str" begin + test_eig_full(T, sz; kwargs...) + test_trunc && test_eig_trunc(T, sz; kwargs...) + end +end + +function test_eig_algs(T::Type, sz, algs; test_trunc = true, kwargs...) + summary_str = testargs_summary(T, sz) + return @testset "eig algorithms $summary_str" begin + test_eig_full_algs(T, sz, algs; kwargs...) + test_trunc && test_eig_trunc_algs(T, sz, algs; kwargs...) + end +end + +function test_eig_full( + T::Type, sz; + atol::Real = 0, rtol::Real = precision(T), + kwargs... + ) + summary_str = testargs_summary(T, sz) + return @testset "eig_full! $summary_str" begin + A = instantiate_matrix(T, sz) + Ac = deepcopy(A) + Tc = isa(A, Diagonal) ? eltype(T) : complex(eltype(T)) + D, V = @testinferred eig_full(A) + @test eltype(D) == eltype(V) == Tc + @test A * V ≈ V * D + + D2, V2 = @testinferred eig_full!(Ac, (D, V)) + @test D2 === D + @test V2 === V + @test A * V ≈ V * D + + Dc = @testinferred eig_vals(A) + @test eltype(Dc) == Tc + @test D ≈ Diagonal(Dc) + end +end + +function test_eig_full_algs( + T::Type, sz, algs; + atol::Real = 0, rtol::Real = precision(T), + kwargs... + ) + summary_str = testargs_summary(T, sz) + return @testset "eig_full! algorithm $alg $summary_str" for alg in algs + A = instantiate_matrix(T, sz) + Ac = deepcopy(A) + Tc = isa(A, Diagonal) ? eltype(T) : complex(eltype(T)) + D, V = @testinferred eig_full(A; alg) + @test eltype(D) == eltype(V) == Tc + @test A * V ≈ V * D + + D2, V2 = @testinferred eig_full!(Ac, (D, V); alg) + @test D2 === D + @test V2 === V + @test A * V ≈ V * D + + Dc = @testinferred eig_vals(A; alg) + @test eltype(Dc) == Tc + @test D ≈ Diagonal(Dc) + end +end + +function test_eig_trunc( + T::Type, sz; + test_blocksize = true, + atol::Real = 0, rtol::Real = precision(T), + kwargs... + ) + summary_str = testargs_summary(T, sz) + return @testset "eig_trunc! $summary_str" begin + A = instantiate_matrix(T, sz) + A *= A' # TODO: deal with eigenvalue ordering etc + Ac = deepcopy(A) + Tc = complex(eltype(T)) + # eigenvalues are sorted by ascending real component... + D₀ = sort!(eig_vals(A); by = abs, rev = true) + m = size(A, 1) + rmin = findfirst(i -> abs(D₀[end - i]) != abs(D₀[end - i - 1]), 1:(m - 2)) + r = length(D₀) - rmin + atol = sqrt(eps(real(eltype(T)))) + + D1, V1, ϵ1 = @testinferred eig_trunc(A; trunc = truncrank(r)) + @test length(diagview(D1)) == r + @test A * V1 ≈ V1 * D1 + @test ϵ1 ≈ norm(view(D₀, (r + 1):m)) atol = atol + + s = 1 + sqrt(eps(real(eltype(T)))) + trunc = trunctol(; atol = s * abs(D₀[r + 1])) + D2, V2, ϵ2 = @testinferred eig_trunc(A; trunc) + @test length(diagview(D2)) == r + @test A * V2 ≈ V2 * D2 + @test ϵ2 ≈ norm(view(D₀, (r + 1):m)) atol = atol + + s = 1 - sqrt(eps(real(eltype(T)))) + trunc = truncerror(; atol = s * norm(@view(D₀[r:end]), 1), p = 1) + D3, V3, ϵ3 = @testinferred eig_trunc(A; trunc) + @test length(diagview(D3)) == r + @test A * V3 ≈ V3 * D3 + @test ϵ3 ≈ norm(view(D₀, (r + 1):m)) atol = atol + + D4, V4 = @testinferred eig_trunc_no_error(A; trunc) + @test length(diagview(D4)) == r + @test A * V4 ≈ V4 * D4 + + # trunctol keeps order, truncrank might not + # test for same subspace + @test V1 * ((V1' * V1) \ (V1' * V2)) ≈ V2 + @test V2 * ((V2' * V2) \ (V2' * V1)) ≈ V1 + @test V1 * ((V1' * V1) \ (V1' * V3)) ≈ V3 + @test V3 * ((V3' * V3) \ (V3' * V1)) ≈ V1 + @test V4 * ((V4' * V4) \ (V4' * V1)) ≈ V1 + + m4 = 4 + atol = sqrt(eps(real(eltype(T)))) + V = T <: Diagonal ? I : randn!(similar(A, eltype(T), m4, m4)) + Ddiag = similar(A, real(eltype(T)), m4) + copyto!(Ddiag, real(eltype(T))[0.9, 0.3, 0.1, 0.01]) + D = Diagonal(Ddiag) + A = T <: Diagonal ? Diagonal(V * D * inv(V)) : V * D * inv(V) + alg = TruncatedAlgorithm(MatrixAlgebraKit.default_eig_algorithm(A), truncrank(2)) + D2, V2, ϵ2 = @testinferred eig_trunc(A; alg) + @test diagview(D2) ≈ diagview(D)[1:2] + @test ϵ2 ≈ norm(diagview(D)[3:4]) atol = atol + @test_throws ArgumentError eig_trunc(A; alg, trunc = (; maxrank = 2)) + + alg = TruncatedAlgorithm(MatrixAlgebraKit.default_eig_algorithm(A), truncerror(; atol = 0.2, p = 1)) + D3, V3, ϵ3 = @testinferred eig_trunc(A; alg) + @test diagview(D3) ≈ diagview(D)[1:2] + @test ϵ3 ≈ norm(diagview(D)[3:4]) atol = atol + + D4, V4 = @testinferred eig_trunc_no_error(A; alg) + @test diagview(D4) ≈ diagview(D)[1:2] + end +end + +function test_eig_trunc_algs( + T::Type, sz, algs; + atol::Real = 0, rtol::Real = precision(T), + kwargs... + ) + summary_str = testargs_summary(T, sz) + return @testset "eig_trunc! algorithm $alg $summary_str" for alg in algs + A = instantiate_matrix(T, sz) + A *= A' # TODO: deal with eigenvalue ordering etc + Ac = deepcopy(A) + Tc = complex(eltype(T)) + # eigenvalues are sorted by ascending real component... + D₀ = sort!(eig_vals(A; alg); by = abs, rev = true) + m = size(A, 1) + rmin = findfirst(i -> abs(D₀[end - i]) != abs(D₀[end - i - 1]), 1:(m - 2)) + r = length(D₀) - rmin + atol = sqrt(eps(real(eltype(T)))) + + m4 = 4 + atol = sqrt(eps(real(eltype(T)))) + V = T <: Diagonal ? I : randn!(similar(A, eltype(T), m4, m4)) + Ddiag = similar(A, real(eltype(T)), m4) + copyto!(Ddiag, real(eltype(T))[0.9, 0.3, 0.1, 0.01]) + D = Diagonal(Ddiag) + A = T <: Diagonal ? Diagonal(V * D * inv(V)) : V * D * inv(V) + truncalg = TruncatedAlgorithm(alg, truncrank(2)) + D2, V2, ϵ2 = @testinferred eig_trunc(A; alg = truncalg) + @test diagview(D2) ≈ diagview(D)[1:2] + @test ϵ2 ≈ norm(diagview(D)[3:4]) atol = atol + @test_throws ArgumentError eig_trunc(A; alg = truncalg, trunc = (; maxrank = 2)) + + truncalg = TruncatedAlgorithm(alg, truncerror(; atol = 0.2, p = 1)) + D3, V3, ϵ3 = @testinferred eig_trunc(A; alg = truncalg) + @test diagview(D3) ≈ diagview(D)[1:2] + @test ϵ3 ≈ norm(diagview(D)[3:4]) atol = atol + + D4, V4 = @testinferred eig_trunc_no_error(A; alg = truncalg) + @test diagview(D4) ≈ diagview(D)[1:2] + end +end From 94dcdaed165764aad86720ce892bdcd87a3a4237 Mon Sep 17 00:00:00 2001 From: Katharine Hyatt Date: Tue, 23 Dec 2025 14:42:26 +0100 Subject: [PATCH 2/4] Don't test trunc for CUDA algs --- test/eig.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/test/eig.jl b/test/eig.jl index 01a3bd94..13d703e8 100644 --- a/test/eig.jl +++ b/test/eig.jl @@ -20,7 +20,7 @@ for T in (BLASFloats..., GenericFloats...) if T ∈ BLASFloats if CUDA.functional() TestSuite.test_eig(CuMatrix{T}, (m, m); test_trunc = false) - TestSuite.test_eig_algs(CuMatrix{T}, (m, m), (CUSOLVER_Simple(),)) + TestSuite.test_eig_algs(CuMatrix{T}, (m, m), (CUSOLVER_Simple(),); test_trunc = false) TestSuite.test_eig(Diagonal{T, CuVector{T}}, m) TestSuite.test_eig_algs(Diagonal{T, CuVector{T}}, m, (DiagonalAlgorithm(),)) end From 71f61da89e0b0099c9c0e84f781e45821137eead Mon Sep 17 00:00:00 2001 From: Katharine Hyatt Date: Tue, 23 Dec 2025 14:48:25 +0100 Subject: [PATCH 3/4] No trunc for any CUDA --- test/eig.jl | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/test/eig.jl b/test/eig.jl index 13d703e8..7cc54c5d 100644 --- a/test/eig.jl +++ b/test/eig.jl @@ -21,8 +21,8 @@ for T in (BLASFloats..., GenericFloats...) if CUDA.functional() TestSuite.test_eig(CuMatrix{T}, (m, m); test_trunc = false) TestSuite.test_eig_algs(CuMatrix{T}, (m, m), (CUSOLVER_Simple(),); test_trunc = false) - TestSuite.test_eig(Diagonal{T, CuVector{T}}, m) - TestSuite.test_eig_algs(Diagonal{T, CuVector{T}}, m, (DiagonalAlgorithm(),)) + TestSuite.test_eig(Diagonal{T, CuVector{T}}, m; test_trunc = false) + TestSuite.test_eig_algs(Diagonal{T, CuVector{T}}, m, (DiagonalAlgorithm(),); test_trunc = false) end #= not yet supported if AMDGPU.functional() From ea4a72386737f8c38c359a79bbf8abb7310ee1bb Mon Sep 17 00:00:00 2001 From: Katharine Hyatt Date: Tue, 23 Dec 2025 13:15:13 -0500 Subject: [PATCH 4/4] Comments --- ext/MatrixAlgebraKitGenericSchurExt.jl | 13 +------------ test/testsuite/eig.jl | 8 ++------ 2 files changed, 3 insertions(+), 18 deletions(-) diff --git a/ext/MatrixAlgebraKitGenericSchurExt.jl b/ext/MatrixAlgebraKitGenericSchurExt.jl index 7ee0939b..686ca87b 100644 --- a/ext/MatrixAlgebraKitGenericSchurExt.jl +++ b/ext/MatrixAlgebraKitGenericSchurExt.jl @@ -13,18 +13,7 @@ MatrixAlgebraKit.initialize_output(::typeof(eig_full!), A::AbstractMatrix, ::GS_ MatrixAlgebraKit.initialize_output(::typeof(eig_vals!), A::AbstractMatrix, ::GS_QRIteration) = nothing function MatrixAlgebraKit.eig_full!(A::AbstractMatrix, DV, ::GS_QRIteration) - D_, V_ = GenericSchur.eigen!(A) - D, V = DV - if !isnothing(D) - copyto!(D, Diagonal(D_)) - else - D = D_ - end - if !isnothing(V) - copyto!(V, V_) - else - V = V_ - end + D, V = GenericSchur.eigen!(A) return Diagonal(D), V end diff --git a/test/testsuite/eig.jl b/test/testsuite/eig.jl index feafb769..2dbea8b9 100644 --- a/test/testsuite/eig.jl +++ b/test/testsuite/eig.jl @@ -34,9 +34,7 @@ function test_eig_full( @test A * V ≈ V * D D2, V2 = @testinferred eig_full!(Ac, (D, V)) - @test D2 === D - @test V2 === V - @test A * V ≈ V * D + @test A * V2 ≈ V2 * D2 Dc = @testinferred eig_vals(A) @test eltype(Dc) == Tc @@ -59,9 +57,7 @@ function test_eig_full_algs( @test A * V ≈ V * D D2, V2 = @testinferred eig_full!(Ac, (D, V); alg) - @test D2 === D - @test V2 === V - @test A * V ≈ V * D + @test A * V2 ≈ V2 * D2 Dc = @testinferred eig_vals(A; alg) @test eltype(Dc) == Tc