From 858ccac1f93ab74aecd6f6c217e2cb553e7866ae Mon Sep 17 00:00:00 2001 From: Tim Holy Date: Tue, 28 Apr 2026 05:45:28 -0500 Subject: [PATCH 1/2] Drop `getindex` support, Af' returns FactoredMatrix Formerly `getindex` was supported for `Adjoint{FactoredMatrix}` and `Transpose{FactoredMatrix}`, but not `FactoredMatrix` itself. Since `getindex` should not be supported by `Factorization` objects, just delete these methods. Also makes the lowercase `adjoint(A)` and `transpose(A)` return `FactoredMatrix` objects. This was made possible by #9 which generalized the types for `U` and `V`. Note the uppercase `Adjoint` and `Transpose` type-wrappers still return the expected types. --- src/FactoredMatrices.jl | 8 +--- test/runtests.jl | 92 ++++++++++++++++++++++------------------- 2 files changed, 51 insertions(+), 49 deletions(-) diff --git a/src/FactoredMatrices.jl b/src/FactoredMatrices.jl index 0fb1a6f..946df3f 100644 --- a/src/FactoredMatrices.jl +++ b/src/FactoredMatrices.jl @@ -79,12 +79,8 @@ length(A::FactoredMatrix) = size(A.U, 1) * size(A.V, 2) Array(A::FactoredMatrix) = A.U * A.V -adjoint(A::FactoredMatrix) = Adjoint(A) -transpose(A::FactoredMatrix) = Transpose(A) -getindex(A::Adjoint{T,FactoredMatrix{T}}, i, j) where T = - conj(sum(A.parent.U[j,:] .* A.parent.V[:,i])) -getindex(A::Transpose{T,FactoredMatrix{T}}, i, j) where T = - sum(A.parent.U[j,:] .* A.parent.V[:,i]) +adjoint(A::FactoredMatrix{T}) where T = FactoredMatrix{T}(adjoint(A.V), adjoint(A.U)) +transpose(A::FactoredMatrix{T}) where T = FactoredMatrix{T}(transpose(A.V), transpose(A.U)) # --- mul! --- diff --git a/test/runtests.jl b/test/runtests.jl index aa38414..7edae70 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -15,15 +15,15 @@ function check_vec_alloc(Mf, b, ws, c) end function check_allocs(Mf, ws, C, D, E, r15x5, r10x5, r5x10, r5x15) - @test @allocated(mul!(r15x5, Mf, C; cache=ws)) == 0 - @test @allocated(mul!(r10x5, Mf', D; cache=ws)) == 0 - @test @allocated(mul!(r10x5, transpose(Mf), D; cache=ws)) == 0 - @test @allocated(mul!(r5x10, D', Mf; cache=ws)) == 0 - @test @allocated(mul!(r5x10, transpose(D), Mf; cache=ws)) == 0 - @test @allocated(mul!(r15x5, Mf, E'; cache=ws)) == 0 - @test @allocated(mul!(r15x5, Mf, transpose(E); cache=ws)) == 0 - @test @allocated(mul!(r5x15, E, Mf'; cache=ws)) == 0 - @test @allocated(mul!(r5x15, E, transpose(Mf); cache=ws)) == 0 + @test @allocated(mul!(r15x5, Mf, C; cache=ws)) == 0 + @test @allocated(mul!(r10x5, Adjoint(Mf), D; cache=ws)) == 0 + @test @allocated(mul!(r10x5, Transpose(Mf), D; cache=ws)) == 0 + @test @allocated(mul!(r5x10, D', Mf; cache=ws)) == 0 + @test @allocated(mul!(r5x10, transpose(D), Mf; cache=ws)) == 0 + @test @allocated(mul!(r15x5, Mf, E'; cache=ws)) == 0 + @test @allocated(mul!(r15x5, Mf, transpose(E); cache=ws)) == 0 + @test @allocated(mul!(r5x15, E, Adjoint(Mf); cache=ws)) == 0 + @test @allocated(mul!(r5x15, E, Transpose(Mf); cache=ws)) == 0 end @testset "FactoredMatrices" begin @@ -34,18 +34,24 @@ end Mf = FactoredMatrix(U, V) M = Array(Mf) @test M == U*V + # lowercase adjoint(Mf) and transpose(Mf) re-wrap `U` and `V` + @test Mf' isa FactoredMatrix{T} + @test transpose(Mf) isa FactoredMatrix{T} + # uppercase return the expected types + @test Adjoint(Mf) isa Adjoint{T, <:FactoredMatrix{T}} + @test Transpose(Mf) isa Transpose{T, <:FactoredMatrix{T}} C = myrand(T, 10, 5) # for Mf*C (15×10)*(10×5) → 15×5 D = myrand(T, 15, 5) # for Mf'*D (10×15)*(15×5) → 10×5 E = myrand(T, 5, 10) # for E*Mf' (5×10)*(10×15) → 5×15 @test M*C ≈ Mf*C - @test M'*D ≈ Mf'*D - @test transpose(M)*D ≈ transpose(Mf)*D + @test M'*D ≈ Mf'*D ≈ Adjoint(Mf)*D + @test transpose(M)*D ≈ transpose(Mf)*D ≈ Transpose(Mf)*D @test D'*M ≈ D'*Mf @test transpose(D)*M ≈ transpose(D)*Mf @test M*adjoint(E) ≈ Mf*E' @test M*transpose(E) ≈ Mf*transpose(E) @test E*adjoint(M) ≈ E*Mf' - @test E*transpose(M) ≈ E*transpose(Mf) + @test E*transpose(M) ≈ E*transpose(Mf) ≈ E*Transpose(Mf) @test M'*D ≈ Mf'*D @test D'*M ≈ D'*Mf @test M*E' ≈ Mf*E' @@ -92,23 +98,23 @@ end r15x8 = Matrix{T}(undef, 15, 8) r10x8 = Matrix{T}(undef, 10, 8) @test Mf * Mf_A ≈ M * MA - @test Mf' * Mf_B ≈ M' * MB + @test Adjoint(Mf) * Mf_B ≈ M' * MB @test transpose(Mf) * Mf_B ≈ transpose(M) * MB - @test Mf * Mf_C' ≈ M * MC' - @test Mf * transpose(Mf_C) ≈ M * transpose(MC) - @test Mf' * Mf_D' ≈ M' * MD' - @test Mf' * transpose(Mf_D) ≈ M' * transpose(MD) - @test transpose(Mf) * Mf_D' ≈ transpose(M) * MD' - @test transpose(Mf) * transpose(Mf_D) ≈ transpose(M) * transpose(MD) + @test Mf * Adjoint(Mf_C) ≈ M * MC' + @test Mf * Transpose(Mf_C) ≈ M * transpose(MC) + @test Adjoint(Mf) * Mf_D' ≈ M' * MD' + @test Adjoint(Mf) * Transpose(Mf_D) ≈ M' * transpose(MD) + @test Transpose(Mf) * Mf_D' ≈ transpose(M) * MD' + @test Transpose(Mf) * Transpose(Mf_D) ≈ transpose(M) * transpose(MD) mul!(r15x8, Mf, Mf_A); @test r15x8 ≈ M * MA - mul!(r10x8, Mf', Mf_B); @test r10x8 ≈ M' * MB + mul!(r10x8, Adjoint(Mf), Mf_B); @test r10x8 ≈ M' * MB mul!(r10x8, transpose(Mf), Mf_B); @test r10x8 ≈ transpose(M) * MB - mul!(r15x8, Mf, Mf_C'); @test r15x8 ≈ M * MC' - mul!(r15x8, Mf, transpose(Mf_C)); @test r15x8 ≈ M * transpose(MC) - mul!(r10x8, Mf', Mf_D'); @test r10x8 ≈ M' * MD' - mul!(r10x8, Mf', transpose(Mf_D)); @test r10x8 ≈ M' * transpose(MD) - mul!(r10x8, transpose(Mf), Mf_D'); @test r10x8 ≈ transpose(M) * MD' - mul!(r10x8, transpose(Mf), transpose(Mf_D)); @test r10x8 ≈ transpose(M) * transpose(MD) + mul!(r15x8, Mf, Adjoint(Mf_C)); @test r15x8 ≈ M * MC' + mul!(r15x8, Mf, Transpose(Mf_C)); @test r15x8 ≈ M * transpose(MC) + mul!(r10x8, Adjoint(Mf), Mf_D'); @test r10x8 ≈ M' * MD' + mul!(r10x8, Adjoint(Mf), Transpose(Mf_D)); @test r10x8 ≈ M' * transpose(MD) + mul!(r10x8, Transpose(Mf), Mf_D'); @test r10x8 ≈ transpose(M) * MD' + mul!(r10x8, Transpose(Mf), Transpose(Mf_D)); @test r10x8 ≈ transpose(M) * transpose(MD) end # mul! FM×FM→AM with cache @@ -119,8 +125,8 @@ end r10x8 = Matrix{T}(undef, 10, 8) ws3x3 = FactoredMatrices.Workspace{T}(3, 3) mul!(r15x8, Mf, Mf_A; cache=ws3x3); @test r15x8 ≈ M * MA - mul!(r10x8, Mf', Mf_B; cache=ws3x3); @test r10x8 ≈ M' * MB - mul!(r10x8, transpose(Mf), Mf_B; cache=ws3x3); @test r10x8 ≈ transpose(M) * MB + mul!(r10x8, Adjoint(Mf), Mf_B; cache=ws3x3); @test r10x8 ≈ M' * MB + mul!(r10x8, Transpose(Mf), Mf_B; cache=ws3x3); @test r10x8 ≈ transpose(M) * MB end # mul! with FactoredMatrix output: FM×AM, AM×FM, AM×AM, adj/transpose variants @@ -131,25 +137,25 @@ end # FM × AM → FM Cf = FactoredMatrix(Matrix{T}(undef, 15, 3), Matrix{T}(undef, 3, 5)) - mul!(Cf, Mf, C_mat); @test Array(Cf) ≈ M * C_mat + mul!(Cf, Mf, C_mat); @test Array(Cf) ≈ M * C_mat # adj(FM) × AM → FM Cf2 = FactoredMatrix(Matrix{T}(undef, 10, 3), Matrix{T}(undef, 3, 5)) - mul!(Cf2, Mf', D_mat); @test Array(Cf2) ≈ M' * D_mat - mul!(Cf2, transpose(Mf), D_mat); @test Array(Cf2) ≈ transpose(M) * D_mat + mul!(Cf2, Adjoint(Mf), D_mat); @test Array(Cf2) ≈ M' * D_mat + mul!(Cf2, Transpose(Mf), D_mat); @test Array(Cf2) ≈ transpose(M) * D_mat # AM × FM → FM A_mat = myrand(T, 5, 15) Cf3 = FactoredMatrix(Matrix{T}(undef, 5, 3), Matrix{T}(undef, 3, 10)) - mul!(Cf3, A_mat, Mf); @test Array(Cf3) ≈ A_mat * M + mul!(Cf3, A_mat, Mf); @test Array(Cf3) ≈ A_mat * M # AM × adj(FM) → FM Cf4 = FactoredMatrix(Matrix{T}(undef, 5, 3), Matrix{T}(undef, 3, 15)) - mul!(Cf4, E_mat, Mf'); @test Array(Cf4) ≈ E_mat * M' - mul!(Cf4, E_mat, transpose(Mf)); @test Array(Cf4) ≈ E_mat * transpose(M) + mul!(Cf4, E_mat, Adjoint(Mf)); @test Array(Cf4) ≈ E_mat * M' + mul!(Cf4, E_mat, Transpose(Mf)); @test Array(Cf4) ≈ E_mat * transpose(M) # AM × AM → FM A2 = myrand(T, 5, 3); B2 = myrand(T, 3, 7) Cf5 = FactoredMatrix(Matrix{T}(undef, 5, 3), Matrix{T}(undef, 3, 7)) - mul!(Cf5, A2, B2); @test Array(Cf5) ≈ A2 * B2 + mul!(Cf5, A2, B2); @test Array(Cf5) ≈ A2 * B2 end # mul! FM×FM→FM (without and with cache, plain and adj/transpose) @@ -162,15 +168,15 @@ end Cf10x8 = FactoredMatrix(Matrix{T}(undef, 10, 3), Matrix{T}(undef, 3, 8)) mul!(Cf15x8, Mf, Mf_A); @test Array(Cf15x8) ≈ M * MA - mul!(Cf10x8, Mf', Mf_B); @test Array(Cf10x8) ≈ M' * MB - mul!(Cf10x8, transpose(Mf), Mf_B); @test Array(Cf10x8) ≈ transpose(M) * MB - mul!(Cf15x8, Mf, Mf_C'); @test Array(Cf15x8) ≈ M * MC' - mul!(Cf15x8, Mf, transpose(Mf_C)); @test Array(Cf15x8) ≈ M * transpose(MC) + mul!(Cf10x8, Adjoint(Mf), Mf_B); @test Array(Cf10x8) ≈ M' * MB + mul!(Cf10x8, Transpose(Mf), Mf_B); @test Array(Cf10x8) ≈ transpose(M) * MB + mul!(Cf15x8, Mf, Adjoint(Mf_C)); @test Array(Cf15x8) ≈ M * MC' + mul!(Cf15x8, Mf, Transpose(Mf_C)); @test Array(Cf15x8) ≈ M * transpose(MC) mul!(Cf15x8, Mf, Mf_A; cache=ws3x3); @test Array(Cf15x8) ≈ M * MA - mul!(Cf10x8, Mf', Mf_B; cache=ws3x3); @test Array(Cf10x8) ≈ M' * MB - mul!(Cf10x8, transpose(Mf), Mf_B; cache=ws3x3); @test Array(Cf10x8) ≈ transpose(M) * MB - mul!(Cf15x8, Mf, Mf_C'; cache=ws3x3); @test Array(Cf15x8) ≈ M * MC' + mul!(Cf10x8, Adjoint(Mf), Mf_B; cache=ws3x3); @test Array(Cf10x8) ≈ M' * MB + mul!(Cf10x8, Transpose(Mf), Mf_B; cache=ws3x3); @test Array(Cf10x8) ≈ transpose(M) * MB + mul!(Cf15x8, Mf, Adjoint(Mf_C); cache=ws3x3); @test Array(Cf15x8) ≈ M * MC' check_fmfm_fm_alloc(Mf, Mf_A, ws3x3, Cf15x8) end @@ -178,7 +184,7 @@ end let b = vec(myrand(T, 10, 1)) c = Vector{T}(undef, 15) - mul!(c, Mf, b); @test c ≈ M * b + mul!(c, Mf, b); @test c ≈ M * b ws_vec = FactoredMatrices.Workspace(Mf, 1) mul!(c, Mf, b; cache=ws_vec); @test c ≈ M * b check_vec_alloc(Mf, b, ws_vec, c) From b393b2cddfdb74029b09ccfb49dbce46ae3fe529 Mon Sep 17 00:00:00 2001 From: Tim Holy Date: Tue, 28 Apr 2026 05:58:58 -0500 Subject: [PATCH 2/2] Bump to v1 --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index 6489d7d..d43d2a0 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "FactoredMatrices" uuid = "13d07ad6-c0d4-11e8-35b1-91a532e61cc8" authors = ["Tim Holy "] -version = "0.1.3" +version = "1.0.0" [deps] LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"