diff --git a/Project.toml b/Project.toml index e5c83b2..f729a83 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "ArrayLayouts" uuid = "4c555306-a7a7-4459-81d9-ec55ddd5c99a" authors = ["Sheehan Olver "] -version = "1.11.0" +version = "1.11.1" [deps] FillArrays = "1a297f60-69ca-5386-bcde-b61e274b549b" diff --git a/src/diagonal.jl b/src/diagonal.jl index 8a8b2be..a2c83ea 100644 --- a/src/diagonal.jl +++ b/src/diagonal.jl @@ -63,28 +63,47 @@ _similar(A::SymTridiagonal) = similar(Tridiagonal(A.ev, A.dv, A.ev)) _copy_diag(M::T, ::T) where {T<:Rmul} = copyto!(_similar(M.A), M) _copy_diag(M::T, ::T) where {T<:Lmul} = copyto!(_similar(M.B), M) _copy_diag(M, _) = copy(M) +_bidiagonal(A::Bidiagonal) = A +function _bidiagonal(A) + # we assume that the matrix is indeed bidiagonal, + # so that the conversion is lossless + if iszero(view(A, diagind(A, -1))) + uplo = :U + else + uplo = :L + end + Bidiagonal(A, uplo) +end function copy(M::Rmul{<:BidiagonalLayout,<:DiagonalLayout}) - A = convert(Bidiagonal, M.A) + A = _bidiagonal(M.A) _copy_diag(Rmul(A, M.B), M) end function copy(M::Lmul{<:DiagonalLayout,<:BidiagonalLayout}) - B = convert(Bidiagonal, M.B) + B = _bidiagonal(M.B) _copy_diag(Lmul(M.A, B), M) end +# we assume that the matrix is indeed tridiagonal, +# so that the conversion is lossless +_tridiagonal(A::Tridiagonal) = A +_tridiagonal(A) = Tridiagonal(A) function copy(M::Rmul{<:TridiagonalLayout,<:DiagonalLayout}) - A = convert(Tridiagonal, M.A) + A = _tridiagonal(M.A) _copy_diag(Rmul(A, M.B), M) end function copy(M::Lmul{<:DiagonalLayout,<:TridiagonalLayout}) - B = convert(Tridiagonal, M.B) + B = _tridiagonal(M.B) _copy_diag(Lmul(M.A, B), M) end +# we assume that the matrix is indeed symmetric tridiagonal, +# so that the conversion is lossless +_symtridiagonal(A::SymTridiagonal) = A +_symtridiagonal(A) = SymTridiagonal(A) function copy(M::Rmul{<:SymTridiagonalLayout,<:DiagonalLayout}) - A = convert(SymTridiagonal, M.A) + A = _symtridiagonal(M.A) _copy_diag(Rmul(A, M.B), M) end function copy(M::Lmul{<:DiagonalLayout,<:SymTridiagonalLayout}) - B = convert(SymTridiagonal, M.B) + B = _symtridiagonal(M.B) _copy_diag(Lmul(M.A, B), M) end diff --git a/test/test_layoutarray.jl b/test/test_layoutarray.jl index 3655aa0..60a664a 100644 --- a/test/test_layoutarray.jl +++ b/test/test_layoutarray.jl @@ -338,18 +338,24 @@ Base.copy(A::MyVector) = MyVector(copy(A.A)) @testset "Diagonal * Bidiagonal/Tridiagonal with structured diags" begin n = size(D,1) - B = Bidiagonal(map(MyVector, (rand(n), rand(n-1)))..., :U) - MB = MyMatrix(B) + BU = Bidiagonal(map(MyVector, (rand(n), rand(n-1)))..., :U) + MBU = MyMatrix(BU) + BL = Bidiagonal(map(MyVector, (rand(n), rand(n-1)))..., :L) + MBL = MyMatrix(BL) S = SymTridiagonal(map(MyVector, (rand(n), rand(n-1)))...) MS = MyMatrix(S) T = Tridiagonal(map(MyVector, (rand(n-1), rand(n), rand(n-1)))...) MT = MyMatrix(T) - DA, BA, SA, TA = map(Array, (D, B, S, T)) + DA, BUA, BLA, SA, TA = map(Array, (D, BU, BL, S, T)) if VERSION >= v"1.11" - @test D * B ≈ DA * BA - @test B * D ≈ BA * DA - @test D * MB ≈ DA * BA - @test MB * D ≈ BA * DA + @test D * BU ≈ DA * BUA + @test BU * D ≈ BUA * DA + @test D * MBU ≈ DA * BUA + @test MBU * D ≈ BUA * DA + @test D * BL ≈ DA * BLA + @test BL * D ≈ BLA * DA + @test D * MBL ≈ DA * BLA + @test MBL * D ≈ BLA * DA end if VERSION >= v"1.12.0-DEV.824" @test D * S ≈ DA * SA