From dbefc57abc8f21e025e05909216d81cc87792f5f Mon Sep 17 00:00:00 2001 From: Jishnu Bhattacharya Date: Wed, 26 Feb 2025 12:52:24 +0530 Subject: [PATCH 1/3] Forward `Diagonal` `rmul!`/`lmul!` for adj/trans to parent --- src/diagonal.jl | 14 ++++++++++++++ test/diagonal.jl | 14 ++++++++++++++ 2 files changed, 28 insertions(+) diff --git a/src/diagonal.jl b/src/diagonal.jl index 90572d7b..b5642800 100644 --- a/src/diagonal.jl +++ b/src/diagonal.jl @@ -360,6 +360,13 @@ function rmul!(A::AbstractMatrix, D::Diagonal) end return A end +# A' = A' * D => A = D' * A +# This uses the fact that D' is a Diagonal +function rmul!(A::AdjOrTransAbsMat, D::Diagonal) + f = wrapperop(A) + lmul!(f(D), f(A)) + A +end # T .= T * D function rmul!(T::Tridiagonal, D::Diagonal) matmul_size_check(size(T), size(D)) @@ -395,6 +402,13 @@ function lmul!(D::Diagonal, T::Tridiagonal) end return T end +# A' = D * A' => A = A * D' +# This uses the fact that D' is a Diagonal +function lmul!(D::Diagonal, A::AdjOrTransAbsMat) + f = wrapperop(A) + rmul!(f(A), f(D)) + A +end @inline function __muldiag_nonzeroalpha!(out, D::Diagonal, B, alpha::Number, beta::Number) @inbounds for j in axes(B, 2) diff --git a/test/diagonal.jl b/test/diagonal.jl index 3dd334aa..bceb2c4b 100644 --- a/test/diagonal.jl +++ b/test/diagonal.jl @@ -1257,6 +1257,20 @@ end end end +@testset "rmul!/lmul! for adj/trans" begin + A = rand(5,4); B = similar(A) + for f in (adjoint, transpose) + D = Diagonal(rand(size(A,1))) + B .= A + rmul!(f(B), D) + @test f(B) == f(A) * D + D = Diagonal(rand(size(A,2))) + B .= A + lmul!(D, f(B)) + @test f(B) == D * f(A) + end +end + struct SMatrix1{T} <: AbstractArray{T,2} elt::T end From fda040168e47737b106e25034b7810d754d2fe08 Mon Sep 17 00:00:00 2001 From: Jishnu Bhattacharya Date: Wed, 26 Feb 2025 12:59:45 +0530 Subject: [PATCH 2/3] Add tests for complex arrays --- test/diagonal.jl | 22 ++++++++++++---------- 1 file changed, 12 insertions(+), 10 deletions(-) diff --git a/test/diagonal.jl b/test/diagonal.jl index bceb2c4b..6f6bbabb 100644 --- a/test/diagonal.jl +++ b/test/diagonal.jl @@ -1258,16 +1258,18 @@ end end @testset "rmul!/lmul! for adj/trans" begin - A = rand(5,4); B = similar(A) - for f in (adjoint, transpose) - D = Diagonal(rand(size(A,1))) - B .= A - rmul!(f(B), D) - @test f(B) == f(A) * D - D = Diagonal(rand(size(A,2))) - B .= A - lmul!(D, f(B)) - @test f(B) == D * f(A) + for T in (Float64, ComplexF64) + A = rand(T,5,4); B = similar(A) + for f in (adjoint, transpose) + D = Diagonal(rand(size(A,1))) + B .= A + rmul!(f(B), D) + @test f(B) == f(A) * D + D = Diagonal(rand(size(A,2))) + B .= A + lmul!(D, f(B)) + @test f(B) == D * f(A) + end end end From a9ee0f59feabd19b57e970b2b712db86d49647c1 Mon Sep 17 00:00:00 2001 From: Jishnu Bhattacharya Date: Wed, 26 Feb 2025 14:46:17 +0530 Subject: [PATCH 3/3] Test for complex Diagonal --- test/diagonal.jl | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/test/diagonal.jl b/test/diagonal.jl index 6f6bbabb..4675e6f2 100644 --- a/test/diagonal.jl +++ b/test/diagonal.jl @@ -1261,11 +1261,11 @@ end for T in (Float64, ComplexF64) A = rand(T,5,4); B = similar(A) for f in (adjoint, transpose) - D = Diagonal(rand(size(A,1))) + D = Diagonal(rand(T, size(A,1))) B .= A rmul!(f(B), D) @test f(B) == f(A) * D - D = Diagonal(rand(size(A,2))) + D = Diagonal(rand(T, size(A,2))) B .= A lmul!(D, f(B)) @test f(B) == D * f(A)