Skip to content

Commit

Permalink
SparseArrays: Allow division by diagonal matrices
Browse files Browse the repository at this point in the history
  • Loading branch information
eschnett committed Apr 23, 2020
1 parent a8275e5 commit 5b740f5
Show file tree
Hide file tree
Showing 3 changed files with 12 additions and 2 deletions.
4 changes: 2 additions & 2 deletions stdlib/LinearAlgebra/src/diagonal.jl
Original file line number Diff line number Diff line change
Expand Up @@ -483,8 +483,8 @@ rdiv!(A::AbstractMatrix{T}, adjD::Adjoint{<:Any,<:Diagonal{T}}) where {T} =
rdiv!(A::AbstractMatrix{T}, transD::Transpose{<:Any,<:Diagonal{T}}) where {T} =
(D = transD.parent; rdiv!(A, D))

(/)(A::Union{StridedMatrix, AbstractTriangular}, D::Diagonal) =
rdiv!((typeof(oneunit(eltype(D))/oneunit(eltype(A)))).(A), D)
(/)(A::Union{StridedMatrix, AbstractTriangular}, D::Diagonal) =
rdiv!((typeof(oneunit(eltype(D))/oneunit(eltype(A)))).(A), D)

(\)(F::Factorization, D::Diagonal) =
ldiv!(F, Matrix{typeof(oneunit(eltype(D))/oneunit(eltype(F)))}(D))
Expand Down
4 changes: 4 additions & 0 deletions stdlib/SparseArrays/src/linalg.jl
Original file line number Diff line number Diff line change
Expand Up @@ -174,6 +174,10 @@ function (*)(A::AbstractSparseMatrixCSC, D::Diagonal)
T = Base.promote_op(*, eltype(D), eltype(A))
mul!(LinearAlgebra.copy_oftype(A, T), A, D)
end
function (/)(A::AbstractSparseMatrixCSC, D::Diagonal)
T = typeof(oneunit(eltype(A))/oneunit(eltype(D)))
rdiv!(LinearAlgebra.copy_oftype(A, T), D)
end

# Sparse matrix multiplication as described in [Gustavson, 1978]:
# http://dl.acm.org/citation.cfm?id=355796
Expand Down
6 changes: 6 additions & 0 deletions stdlib/SparseArrays/test/sparse.jl
Original file line number Diff line number Diff line change
Expand Up @@ -519,6 +519,9 @@ dA = Array(sA)
@test lmul!(Diagonal(bi), copy(dA)) ldiv!(Diagonal(b), copy(sA))
@test lmul!(Diagonal(bi), copy(dA)) ldiv!(transpose(Diagonal(b)), copy(sA))
@test lmul!(Diagonal(conj(bi)), copy(dA)) ldiv!(adjoint(Diagonal(b)), copy(sA))
Aob = Diagonal(b) \ sA
@test Aob == ldiv!(Diagonal(b), copy(sA))
@test issparse(Aob)
@test_throws DimensionMismatch ldiv!(Diagonal(fill(1., length(b)+1)), copy(sA))
@test_throws LinearAlgebra.SingularException ldiv!(Diagonal(zeros(length(b))), copy(sA))

Expand All @@ -527,6 +530,9 @@ dA = Array(sA)
@test rmul!(copy(dAt), Diagonal(bi)) rdiv!(copy(sAt), Diagonal(b))
@test rmul!(copy(dAt), Diagonal(bi)) rdiv!(copy(sAt), transpose(Diagonal(b)))
@test rmul!(copy(dAt), Diagonal(conj(bi))) rdiv!(copy(sAt), adjoint(Diagonal(b)))
Atob = sAt / Diagonal(b)
@test Atob == rdiv!(copy(dAt), Diagonal(b))
@test issparse(Atob)
@test_throws DimensionMismatch rdiv!(copy(sAt), Diagonal(fill(1., length(b)+1)))
@test_throws LinearAlgebra.SingularException rdiv!(copy(sAt), Diagonal(zeros(length(b))))
end
Expand Down

0 comments on commit 5b740f5

Please sign in to comment.