From e5b3073017f00fc866b1a4873aa416b88410520c Mon Sep 17 00:00:00 2001 From: Katharine Hyatt Date: Mon, 22 Dec 2025 13:22:11 -0500 Subject: [PATCH 1/3] Use TestSuite for projections tests --- test/amd/projections.jl | 104 ------------------------ test/cuda/projections.jl | 104 ------------------------ test/projections.jl | 113 +++++--------------------- test/runtests.jl | 12 +-- test/testsuite/TestSuite.jl | 1 + test/testsuite/projections.jl | 144 ++++++++++++++++++++++++++++++++++ 6 files changed, 169 insertions(+), 309 deletions(-) delete mode 100644 test/amd/projections.jl delete mode 100644 test/cuda/projections.jl create mode 100644 test/testsuite/projections.jl diff --git a/test/amd/projections.jl b/test/amd/projections.jl deleted file mode 100644 index a06b152c..00000000 --- a/test/amd/projections.jl +++ /dev/null @@ -1,104 +0,0 @@ -using MatrixAlgebraKit -using Test -using TestExtras -using StableRNGs -using LinearAlgebra: LinearAlgebra, Diagonal, norm -using AMDGPU - -const BLASFloats = (Float32, Float64, ComplexF32, ComplexF64) - -@testset "project_(anti)hermitian! for T = $T" for T in BLASFloats - rng = StableRNG(123) - m = 54 - noisefactor = eps(real(T))^(3 / 4) - for alg in (NativeBlocked(blocksize = 16), NativeBlocked(blocksize = 32), NativeBlocked(blocksize = 64)) - for A in (ROCArray(randn(rng, T, m, m)), Diagonal(ROCArray(randn(rng, T, m)))) - Ah = (A + A') / 2 - Aa = (A - A') / 2 - Ac = copy(A) - - Bh = project_hermitian(A, alg) - @test ishermitian(Bh) - @test Bh ≈ Ah - @test A == Ac - Bh_approx = Bh + noisefactor * Aa - # this is still hermitian for real Diagonal: |A - A'| == 0 - @test !ishermitian(Bh_approx) || norm(Aa) == 0 - @test ishermitian(Bh_approx; rtol = 10 * noisefactor) - - Ba = project_antihermitian(A, alg) - @test isantihermitian(Ba) - @test Ba ≈ Aa - @test A == Ac - Ba_approx = Ba + noisefactor * Ah - @test !isantihermitian(Ba_approx) - # this is never anti-hermitian for real Diagonal: |A - A'| == 0 - @test isantihermitian(Ba_approx; rtol = 10 * noisefactor) || norm(Aa) == 0 - - Bh = project_hermitian!(Ac, alg) - @test Bh === Ac - @test ishermitian(Bh) - @test Bh ≈ Ah - - copy!(Ac, A) - Ba = project_antihermitian!(Ac, alg) - @test Ba === Ac - @test isantihermitian(Ba) - @test Ba ≈ Aa - end - end -end - -@testset "project_isometric! for T = $T" for T in BLASFloats - rng = StableRNG(123) - m = 54 - @testset "size ($m, $n)" for n in (37, m) - k = min(m, n) - svdalgs = (ROCSOLVER_QRIteration(), ROCSOLVER_Jacobi()) - algs = (PolarViaSVD.(svdalgs)...,) # PolarNewton()) # TODO - @testset "algorithm $alg" for alg in algs - A = ROCArray(randn(rng, T, m, n)) - W = project_isometric(A, alg) - @test isisometric(W) - W2 = project_isometric(W, alg) - @test W2 ≈ W # stability of the projection - @test W * (W' * A) ≈ A - - Ac = similar(A) - W2 = @constinferred project_isometric!(copy!(Ac, A), W, alg) - @test W2 === W - @test isisometric(W) - - # test that W is closer to A then any other isometry - for k in 1:10 - δA = ROCArray(randn(rng, T, size(A)...)) - W = project_isometric(A, alg) - W2 = project_isometric(A + δA / 100, alg) - @test norm(A - W2) >= norm(A - W) - end - end - - m == n && @testset "DiagonalAlgorithm" begin - A = Diagonal(ROCArray(randn(rng, T, m))) - alg = PolarViaSVD(DiagonalAlgorithm()) - W = project_isometric(A, alg) - @test isisometric(W) - W2 = project_isometric(W, alg) - @test W2 ≈ W # stability of the projection - @test W * (W' * A) ≈ A - - Ac = similar(A) - W2 = @constinferred project_isometric!(copy!(Ac, A), W, alg) - @test W2 === W - @test isisometric(W) - - # test that W is closer to A then any other isometry - for k in 1:10 - δA = Diagonal(ROCArray(randn(rng, T, m))) - W = project_isometric(A, alg) - W2 = project_isometric(A + δA / 100, alg) - @test norm(A - W2) >= norm(A - W) - end - end - end -end diff --git a/test/cuda/projections.jl b/test/cuda/projections.jl deleted file mode 100644 index 677ef520..00000000 --- a/test/cuda/projections.jl +++ /dev/null @@ -1,104 +0,0 @@ -using MatrixAlgebraKit -using Test -using TestExtras -using StableRNGs -using LinearAlgebra: LinearAlgebra, Diagonal, norm -using CUDA - -const BLASFloats = (Float32, Float64, ComplexF32, ComplexF64) - -@testset "project_(anti)hermitian! for T = $T" for T in BLASFloats - rng = StableRNG(123) - m = 54 - noisefactor = eps(real(T))^(3 / 4) - for alg in (NativeBlocked(blocksize = 16), NativeBlocked(blocksize = 32), NativeBlocked(blocksize = 64)) - for A in (CuArray(randn(rng, T, m, m)), Diagonal(CuArray(randn(rng, T, m)))) - Ah = (A + A') / 2 - Aa = (A - A') / 2 - Ac = copy(A) - - Bh = project_hermitian(A, alg) - @test ishermitian(Bh) - @test Bh ≈ Ah - @test A == Ac - Bh_approx = Bh + noisefactor * Aa - # this is still hermitian for real Diagonal: |A - A'| == 0 - @test !ishermitian(Bh_approx) || norm(Aa) == 0 - @test ishermitian(Bh_approx; rtol = 10 * noisefactor) - - Ba = project_antihermitian(A, alg) - @test isantihermitian(Ba) - @test Ba ≈ Aa - @test A == Ac - Ba_approx = Ba + noisefactor * Ah - @test !isantihermitian(Ba_approx) - # this is never anti-hermitian for real Diagonal: |A - A'| == 0 - @test isantihermitian(Ba_approx; rtol = 10 * noisefactor) || norm(Aa) == 0 - - Bh = project_hermitian!(Ac, alg) - @test Bh === Ac - @test ishermitian(Bh) - @test Bh ≈ Ah - - copy!(Ac, A) - Ba = project_antihermitian!(Ac, alg) - @test Ba === Ac - @test isantihermitian(Ba) - @test Ba ≈ Aa - end - end -end - -@testset "project_isometric! for T = $T" for T in BLASFloats - rng = StableRNG(123) - m = 54 - @testset "size ($m, $n)" for n in (37, m) - k = min(m, n) - svdalgs = (CUSOLVER_SVDPolar(), CUSOLVER_QRIteration(), CUSOLVER_Jacobi()) - algs = (PolarViaSVD.(svdalgs)...,) # PolarNewton()) # TODO - @testset "algorithm $alg" for alg in algs - A = CuArray(randn(rng, T, m, n)) - W = project_isometric(A, alg) - @test isisometric(W) - W2 = project_isometric(W, alg) - @test W2 ≈ W # stability of the projection - @test W * (W' * A) ≈ A - - Ac = similar(A) - W2 = @constinferred project_isometric!(copy!(Ac, A), W, alg) - @test W2 === W - @test isisometric(W) - - # test that W is closer to A then any other isometry - for k in 1:10 - δA = CuArray(randn(rng, T, size(A)...)) - W = project_isometric(A, alg) - W2 = project_isometric(A + δA / 100, alg) - @test norm(A - W2) >= norm(A - W) - end - end - - m == n && @testset "DiagonalAlgorithm" begin - A = Diagonal(CuArray(randn(rng, T, m))) - alg = PolarViaSVD(DiagonalAlgorithm()) - W = project_isometric(A, alg) - @test isisometric(W) - W2 = project_isometric(W, alg) - @test W2 ≈ W # stability of the projection - @test W * (W' * A) ≈ A - - Ac = similar(A) - W2 = @constinferred project_isometric!(copy!(Ac, A), W, alg) - @test W2 === W - @test isisometric(W) - - # test that W is closer to A then any other isometry - for k in 1:10 - δA = Diagonal(CuArray(randn(rng, T, m))) - W = project_isometric(A, alg) - W2 = project_isometric(A + δA / 100, alg) - @test norm(A - W2) >= norm(A - W) - end - end - end -end diff --git a/test/projections.jl b/test/projections.jl index 9a6644cf..16044f63 100644 --- a/test/projections.jl +++ b/test/projections.jl @@ -1,106 +1,35 @@ using MatrixAlgebraKit -using MatrixAlgebraKit: check_hermitian, default_hermitian_tol using Test using TestExtras using StableRNGs using LinearAlgebra: LinearAlgebra, Diagonal, norm, normalize! +using CUDA, AMDGPU -const BLASFloats = (Float32, Float64, ComplexF32, ComplexF64) +BLASFloats = (Float32, Float64, ComplexF32, ComplexF64) +GenericFloats = (BigFloat, Complex{BigFloat}) -@testset "project_(anti)hermitian! for T = $T" for T in BLASFloats - rng = StableRNG(123) - m = 54 - noisefactor = eps(real(T))^(3 / 4) +@isdefined(TestSuite) || include("testsuite/TestSuite.jl") +using .TestSuite - mat0 = zeros(T, (1, 1)) - @test ishermitian(mat0) - @test ishermitian(mat0; atol = default_hermitian_tol(mat0)) - @test isnothing(check_hermitian(mat0)) +is_buildkite = get(ENV, "BUILDKITE", "false") == "true" - for alg in (NativeBlocked(blocksize = 16), NativeBlocked(blocksize = 32), NativeBlocked(blocksize = 64)) - for A in (randn(rng, T, m, m), Diagonal(randn(rng, T, m))) - Ah = (A + A') / 2 - Aa = (A - A') / 2 - Ac = copy(A) - Bh = project_hermitian(A, alg) - @test ishermitian(Bh) - @test Bh ≈ Ah - @test A == Ac - Bh_approx = Bh + noisefactor * Aa - # this is still hermitian for real Diagonal: |A - A'| == 0 - @test !ishermitian(Bh_approx) || norm(Aa) == 0 - @test ishermitian(Bh_approx; rtol = 10 * noisefactor) - - Ba = project_antihermitian(A, alg) - @test isantihermitian(Ba) - @test Ba ≈ Aa - @test A == Ac - Ba_approx = Ba + noisefactor * Ah - @test !isantihermitian(Ba_approx) - # this is never anti-hermitian for real Diagonal: |A - A'| == 0 - @test isantihermitian(Ba_approx; rtol = 10 * noisefactor) || norm(Aa) == 0 - - Bh = project_hermitian!(Ac, alg) - @test Bh === Ac - @test ishermitian(Bh) - @test Bh ≈ Ah - - copy!(Ac, A) - Ba = project_antihermitian!(Ac, alg) - @test Ba === Ac - @test isantihermitian(Ba) - @test Ba ≈ Aa +m = 54 +for T in (BLASFloats..., GenericFloats...) + TestSuite.seed_rng!(123) + if T ∈ BLASFloats + if CUDA.functional() + TestSuite.test_projections(CuMatrix{T}, (m, m); test_blocksize = false) + TestSuite.test_projections(Diagonal{T, CuVector{T}}, m; test_blocksize = false) end - end - - # test approximate error calculation - A = normalize!(randn(rng, T, m, m)) - Ah = project_hermitian(A) - Aa = project_antihermitian(A) - - Ah_approx = Ah + noisefactor * Aa - ϵ = norm(project_antihermitian(Ah_approx)) - @test !ishermitian(Ah_approx; atol = (999 // 1000) * ϵ) - @test ishermitian(Ah_approx; atol = (1001 // 1000) * ϵ) - - Aa_approx = Aa + noisefactor * Ah - ϵ = norm(project_hermitian(Aa_approx)) - @test !isantihermitian(Aa_approx; atol = (999 // 1000) * ϵ) - @test isantihermitian(Aa_approx; atol = (1001 // 1000) * ϵ) -end - -@testset "project_isometric! for T = $T" for T in BLASFloats - rng = StableRNG(123) - m = 54 - @testset "size ($m, $n)" for n in (37, m) - k = min(m, n) - if LinearAlgebra.LAPACK.version() < v"3.12.0" - svdalgs = (LAPACK_DivideAndConquer(), LAPACK_QRIteration(), LAPACK_Bisection()) - else - svdalgs = (LAPACK_DivideAndConquer(), LAPACK_QRIteration(), LAPACK_Bisection(), LAPACK_Jacobi()) - end - algs = (PolarViaSVD.(svdalgs)..., PolarNewton()) - @testset "algorithm $alg" for alg in algs - A = randn(rng, T, m, n) - W = project_isometric(A, alg) - @test isisometric(W) - W2 = project_isometric(W, alg) - @test W2 ≈ W # stability of the projection - @test W * (W' * A) ≈ A - - Ac = similar(A) - W2 = @constinferred project_isometric!(copy!(Ac, A), W, alg) - @test W2 === W - @test isisometric(W) - - # test that W is closer to A then any other isometry - for k in 1:10 - δA = randn(rng, T, m, n) - W = project_isometric(A, alg) - W2 = project_isometric(A + δA / 100, alg) - @test norm(A - W2) > norm(A - W) - end + if AMDGPU.functional() + TestSuite.test_projections(ROCMatrix{T}, (m, m); test_blocksize = false) + TestSuite.test_projections(Diagonal{T, ROCVector{T}}, m; test_blocksize = false) end end + if !is_buildkite + TestSuite.test_projections(T, (m, m)) + AT = Diagonal{T, Vector{T}} + TestSuite.test_projections(AT, m) + end end diff --git a/test/runtests.jl b/test/runtests.jl index a32f26d6..a2613823 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -7,9 +7,6 @@ if !is_buildkite @safetestset "Algorithms" begin include("algorithms.jl") end - @safetestset "Projections" begin - include("projections.jl") - end @safetestset "Truncate" begin include("truncate.jl") end @@ -71,12 +68,12 @@ end @safetestset "Polar Decomposition" begin include("polar.jl") end +@safetestset "Projections" begin + include("projections.jl") +end using CUDA if CUDA.functional() - @safetestset "CUDA Projections" begin - include("cuda/projections.jl") - end @safetestset "CUDA SVD" begin include("cuda/svd.jl") end @@ -93,9 +90,6 @@ end using AMDGPU if AMDGPU.functional() - @safetestset "AMDGPU Projections" begin - include("amd/projections.jl") - end @safetestset "AMDGPU SVD" begin include("amd/svd.jl") end diff --git a/test/testsuite/TestSuite.jl b/test/testsuite/TestSuite.jl index 737c8c68..f2f55bbc 100644 --- a/test/testsuite/TestSuite.jl +++ b/test/testsuite/TestSuite.jl @@ -72,5 +72,6 @@ is_pivoted(alg::MatrixAlgebraKit.LQViaTransposedQR) = is_pivoted(alg.qr_alg) include("qr.jl") include("lq.jl") include("polar.jl") +include("projections.jl") end diff --git a/test/testsuite/projections.jl b/test/testsuite/projections.jl new file mode 100644 index 00000000..b9f47df9 --- /dev/null +++ b/test/testsuite/projections.jl @@ -0,0 +1,144 @@ +using TestExtras +using MatrixAlgebraKit: ishermitian +using LinearAlgebra: Diagonal, normalize +using GenericLinearAlgebra + +function test_projections(T::Type, sz; kwargs...) + summary_str = testargs_summary(T, sz) + return @testset "projections $summary_str" begin + test_project_antihermitian(T, sz; kwargs...) + test_project_hermitian(T, sz; kwargs...) + test_project_isometric(T, sz; kwargs...) + end +end + +function test_project_antihermitian( + T::Type, sz; + test_blocksize = true, + atol::Real = 0, rtol::Real = precision(T), + kwargs... + ) + summary_str = testargs_summary(T, sz) + return @testset "project_antihermitian! $summary_str" begin + noisefactor = eps(real(eltype(T)))^(3 / 4) + A = instantiate_matrix(T, sz) + Ac = deepcopy(A) + Ah = (A + A') / 2 + Aa = (A - A') / 2 + + # do we support `blocksize = Int`? + if test_blocksize + Ba = project_antihermitian(A; blocksize = 16) + @test isantihermitian(Ba) + @test Ba ≈ Aa + @test A == Ac + Ba_approx = Ba + noisefactor * Ah + @test !isantihermitian(Ba_approx) + # this is never anti-hermitian for real Diagonal: |A - A'| == 0 + @test isantihermitian(Ba_approx; rtol = 10 * noisefactor) || norm(Aa) == 0 + + copy!(Ac, A) + Ba = project_antihermitian!(Ac; blocksize = 16) + @test Ba === Ac + @test isantihermitian(Ba) + @test Ba ≈ Aa + end + + # test approximate error calculation + A = normalize(A) + Ah = project_hermitian(A) + Aa = project_antihermitian(A) + + Ah_approx = Ah + noisefactor * Aa + ϵ = norm(project_antihermitian(Ah_approx)) + # this is never off-hermitian for real Diagonal: |A - A'| == 0 + @test !ishermitian(Ah_approx; atol = (999 // 1000) * ϵ) || norm(Aa) == 0 + @test ishermitian(Ah_approx; atol = (1001 // 1000) * ϵ) + + Aa_approx = Aa + noisefactor * Ah + ϵ = norm(project_hermitian(Aa_approx)) + @test !isantihermitian(Aa_approx; atol = (999 // 1000) * ϵ) + @test isantihermitian(Aa_approx; atol = (1001 // 1000) * ϵ) + end +end + +function test_project_hermitian( + T::Type, sz; + test_blocksize = true, + atol::Real = 0, rtol::Real = precision(T), + kwargs... + ) + summary_str = testargs_summary(T, sz) + return @testset "project_hermitian! $summary_str" begin + noisefactor = eps(real(eltype(T)))^(3 / 4) + A = instantiate_matrix(T, sz) + Ac = deepcopy(A) + + # do we support `blocksize = Int`? + if test_blocksize + Ah = (A + A') / 2 + Aa = (A - A') / 2 + + Bh = project_hermitian(A; blocksize = 16) + @test ishermitian(Bh) + @test Bh ≈ Ah + @test A == Ac + Bh_approx = Bh + noisefactor * Aa + # this is still hermitian for real Diagonal: |A - A'| == 0 + @test !ishermitian(Bh_approx) || norm(Aa) == 0 + @test ishermitian(Bh_approx; rtol = 10 * noisefactor) + + Bh = project_hermitian!(Ac; blocksize = 16) + @test Bh === Ac + @test ishermitian(Bh) + @test Bh ≈ Ah + end + + # test approximate error calculation + A = normalize(A) + Ah = project_hermitian(A) + Aa = project_antihermitian(A) + + Ah_approx = Ah + noisefactor * Aa + ϵ = norm(project_antihermitian(Ah_approx)) + # this is still hermitian for real Diagonal: |A - A'| == 0 + @test !ishermitian(Ah_approx; atol = (999 // 1000) * ϵ) || norm(Aa) == 0 + @test ishermitian(Ah_approx; atol = (1001 // 1000) * ϵ) + + Aa_approx = Aa + noisefactor * Ah + ϵ = norm(project_hermitian(Aa_approx)) + @test !isantihermitian(Aa_approx; atol = (999 // 1000) * ϵ) + @test isantihermitian(Aa_approx; atol = (1001 // 1000) * ϵ) + end +end + +function test_project_isometric( + T::Type, sz; + atol::Real = 0, rtol::Real = precision(T), + kwargs... + ) + summary_str = testargs_summary(T, sz) + return @testset "project_isometric! $summary_str" begin + A = instantiate_matrix(T, sz) + Ac = deepcopy(A) + k = min(size(A)...) + W = project_isometric(A) + @test isisometric(W) + W2 = project_isometric(W) + @test W2 ≈ W # stability of the projection + @test W * (W' * A) ≈ A + + W2 = @testinferred project_isometric!(Ac, W) + @test W2 === W + @test isisometric(W) + + # test that W is closer to A then any other isometry + for k in 1:10 + δA = instantiate_matrix(T, sz) + W = project_isometric(A) + W2 = project_isometric(A + δA / 100) + # must be ≥ for real Diagonal case + @test norm(A - W2) ≥ norm(A - W) + end + end +end From 9d19d33f9a445c0c4c281b715f8a1e209d590506 Mon Sep 17 00:00:00 2001 From: Katharine Hyatt Date: Tue, 23 Dec 2025 07:31:46 +0100 Subject: [PATCH 2/3] Go back to inplace normalize --- test/testsuite/projections.jl | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/test/testsuite/projections.jl b/test/testsuite/projections.jl index b9f47df9..7d80a8a7 100644 --- a/test/testsuite/projections.jl +++ b/test/testsuite/projections.jl @@ -1,6 +1,6 @@ using TestExtras using MatrixAlgebraKit: ishermitian -using LinearAlgebra: Diagonal, normalize +using LinearAlgebra: Diagonal, normalize! using GenericLinearAlgebra function test_projections(T::Type, sz; kwargs...) @@ -45,7 +45,7 @@ function test_project_antihermitian( end # test approximate error calculation - A = normalize(A) + A = normalize!(A) Ah = project_hermitian(A) Aa = project_antihermitian(A) @@ -95,7 +95,7 @@ function test_project_hermitian( end # test approximate error calculation - A = normalize(A) + A = normalize!(A) Ah = project_hermitian(A) Aa = project_antihermitian(A) From 66201279b9220448a85033410343249505d0ff6b Mon Sep 17 00:00:00 2001 From: Katharine Hyatt Date: Tue, 23 Dec 2025 07:48:19 +0100 Subject: [PATCH 3/3] Bypass normalize for now --- test/testsuite/projections.jl | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/test/testsuite/projections.jl b/test/testsuite/projections.jl index 7d80a8a7..3147084a 100644 --- a/test/testsuite/projections.jl +++ b/test/testsuite/projections.jl @@ -1,6 +1,5 @@ using TestExtras using MatrixAlgebraKit: ishermitian -using LinearAlgebra: Diagonal, normalize! using GenericLinearAlgebra function test_projections(T::Type, sz; kwargs...) @@ -37,7 +36,6 @@ function test_project_antihermitian( # this is never anti-hermitian for real Diagonal: |A - A'| == 0 @test isantihermitian(Ba_approx; rtol = 10 * noisefactor) || norm(Aa) == 0 - copy!(Ac, A) Ba = project_antihermitian!(Ac; blocksize = 16) @test Ba === Ac @test isantihermitian(Ba) @@ -45,7 +43,8 @@ function test_project_antihermitian( end # test approximate error calculation - A = normalize!(A) + normA = norm(A) + A ./= normA Ah = project_hermitian(A) Aa = project_antihermitian(A) @@ -95,7 +94,8 @@ function test_project_hermitian( end # test approximate error calculation - A = normalize!(A) + normA = norm(A) + A ./= normA Ah = project_hermitian(A) Aa = project_antihermitian(A)