diff --git a/Project.toml b/Project.toml index 2c53cb3..8e60927 100644 --- a/Project.toml +++ b/Project.toml @@ -1,6 +1,6 @@ name = "ArrayLayouts" uuid = "4c555306-a7a7-4459-81d9-ec55ddd5c99a" -version = "1.15" +version = "1.15.1" authors = ["Sheehan Olver "] [deps] diff --git a/ext/ArrayLayoutsSparseArraysExt.jl b/ext/ArrayLayoutsSparseArraysExt.jl index 52e9c03..32b543e 100644 --- a/ext/ArrayLayoutsSparseArraysExt.jl +++ b/ext/ArrayLayoutsSparseArraysExt.jl @@ -6,6 +6,7 @@ using SparseArrays # Specifying the full namespace is necessary because of https://github.com/JuliaLang/julia/issues/48533 # See https://github.com/JuliaStats/LogExpFunctions.jl/pull/63 import ArrayLayouts.LinearAlgebra +using ArrayLayouts.LinearAlgebra: Adjoint, Transpose, Symmetric, Hermitian, Diagonal, UpperTriangular, UnitUpperTriangular, LowerTriangular, UnitLowerTriangular import Base: copyto!, == @@ -23,6 +24,22 @@ copyto!(dest::SubArray{<:Any,2,<:LayoutMatrix}, src::SparseArrays.AbstractSparse ArrayLayouts.dot(a,b) +# ambiguity between \ for sparse matrices and the layout overloads of \ (see @layoutldiv). +# SparseArrays' solvers need strided right-hand sides, so convert to Array. +# Symmetric/Hermitian are included since CHOLMOD does not support ldiv! with other arrays. +const LayoutLdivRHS = Any[:LayoutVector, :LayoutMatrix] +for Tri in (:UpperTriangular, :UnitUpperTriangular, :LowerTriangular, :UnitLowerTriangular), + Par in (:(LayoutMatrix{T}), :(SubArray{T,2,<:LayoutMatrix{T}}), :(Adjoint{T,<:LayoutMatrix{T}}), :(Transpose{T,<:LayoutMatrix{T}})) + push!(LayoutLdivRHS, :($Tri{T,<:$Par} where T)) +end +for Sp in (:(SparseArrays.AbstractSparseMatrixCSC), + :(Adjoint{<:Any,<:SparseArrays.AbstractSparseMatrixCSC}), :(Transpose{<:Any,<:SparseArrays.AbstractSparseMatrixCSC}), + :(Symmetric{<:Any,<:SparseArrays.AbstractSparseMatrixCSC}), :(Hermitian{<:Any,<:SparseArrays.AbstractSparseMatrixCSC})), + Typ in LayoutLdivRHS + @eval Base.:\(A::$Sp, B::$Typ) = A \ convert(Array, B) + Typ === :LayoutVector && @eval Base.:\(A::$Sp, B::Diagonal{<:Any,<:LayoutVector}) = A \ Diagonal(convert(Array, B.diag)) +end + # ambiguity from SparseArrays.ReadOnly if isdefined(SparseArrays, :ReadOnly) ==(A::SparseArrays.ReadOnly, B::LayoutVector) = equals_layout(A, B) diff --git a/test/test_layoutarray.jl b/test/test_layoutarray.jl index dac1dc2..38d58b4 100644 --- a/test/test_layoutarray.jl +++ b/test/test_layoutarray.jl @@ -619,6 +619,18 @@ BroadcastStyle(::Type{<:MyVector2{T}}) where {T} = MyBroadcastStyle{1}() copyto!(view(M, :, :), S) @test S == M2 end + @testset "ldiv" begin + S = sparse([4.0 1 0; 1 4 1; 0 1 4]) + b = [1.0, 2, 3] + X = [1.0 2 0; 2 3 1; 3 4 2] + for A in (S, S', transpose(S), Symmetric(S), Hermitian(S), Symmetric(S, :L)) + for B in (MyVector(b), MyMatrix(X), Diagonal(MyVector(b)), + UpperTriangular(MyMatrix(X)), UnitLowerTriangular(view(MyMatrix(X), :, :)), + LowerTriangular(MyMatrix(X)'), UnitUpperTriangular(transpose(MyMatrix(X)))) + @test A \ B ≈ Matrix(A) \ Array(B) + end + end + end end @testset "mul! with subarrays" begin