Skip to content

Commit

Permalink
Merge pull request #24372 from Sacha0/constructI
Browse files Browse the repository at this point in the history
constructors for Matrix and SparseMatrixCSC from UniformScaling
  • Loading branch information
Sacha0 authored Oct 31, 2017
2 parents 1536aab + f5bdacb commit d2d995d
Show file tree
Hide file tree
Showing 4 changed files with 60 additions and 15 deletions.
10 changes: 10 additions & 0 deletions base/linalg/uniformscaling.jl
Original file line number Diff line number Diff line change
Expand Up @@ -371,3 +371,13 @@ UniformScaling{Float64}
```
"""
chol(J::UniformScaling, args...) = ((C, info) = _chol!(J, nothing); @assertposdef C info)


## Matrix construction from UniformScaling
Matrix{T}(s::UniformScaling, dims::Dims{2}) where {T} = setindex!(zeros(T, dims), T(s.λ), diagind(dims...))
Matrix{T}(s::UniformScaling, m::Integer, n::Integer) where {T} = Matrix{T}(s, Dims((m, n)))
Matrix(s::UniformScaling, m::Integer, n::Integer) = Matrix(s, Dims((m, n)))
Matrix(s::UniformScaling, dims::Dims{2}) = Matrix{eltype(s)}(s, dims)
# convenience variations that accept a single integer to specify dims
Matrix{T}(s::UniformScaling, m::Integer) where {T} = Matrix{T}(s, m, m)
Matrix(s::UniformScaling, m::Integer) = Matrix(s, m, m)
38 changes: 23 additions & 15 deletions base/sparse/sparsematrix.jl
Original file line number Diff line number Diff line change
Expand Up @@ -1529,29 +1529,37 @@ if not specified.
`sparse(α*I, m, n)` can be used to efficiently create a sparse
multiple `α` of the identity matrix.
"""
speye(::Type{T}, m::Integer, n::Integer) where {T} = speye_scaled(T, oneunit(T), m, n)
speye(::Type{T}, m::Integer, n::Integer) where {T} = SparseMatrixCSC{T}(UniformScaling(one(T)), Dims((m, n)))
sparse(s::UniformScaling, m::Integer, n::Integer=m) = SparseMatrixCSC(s, Dims((m, n)))

function one(S::SparseMatrixCSC{T}) where T
m,n = size(S)
if m != n; throw(DimensionMismatch("multiplicative identity only defined for square matrices")); end
speye(T, m)
end

speye_scaled(diag, m::Integer, n::Integer) = speye_scaled(typeof(diag), diag, m, n)
## SparseMatrixCSC construction from UniformScaling
SparseMatrixCSC{Tv,Ti}(s::UniformScaling, m::Integer, n::Integer) where {Tv,Ti} = SparseMatrixCSC{Tv,Ti}(s, Dims((m, n)))
SparseMatrixCSC{Tv}(s::UniformScaling, m::Integer, n::Integer) where {Tv} = SparseMatrixCSC{Tv}(s, Dims((m, n)))
SparseMatrixCSC(s::UniformScaling, m::Integer, n::Integer) = SparseMatrixCSC(s, Dims((m, n)))
SparseMatrixCSC{Tv}(s::UniformScaling, dims::Dims{2}) where {Tv} = SparseMatrixCSC{Tv,Int}(s, dims)
SparseMatrixCSC(s::UniformScaling, dims::Dims{2}) = SparseMatrixCSC{eltype(s)}(s, dims)
function SparseMatrixCSC{Tv,Ti}(s::UniformScaling, dims::Dims{2}) where {Tv,Ti}
@boundscheck first(dims) < 0 && throw(ArgumentError("first dimension invalid ($(first(dims)) < 0)"))
@boundscheck last(dims) < 0 && throw(ArgumentError("second dimension invalid ($(last(dims)) < 0)"))
iszero(s.λ) && return spzeros(Tv, Ti, dims...)
m, n, k = dims..., min(dims...)
nzval = fill!(Vector{Tv}(k), Tv(s.λ))
rowval = copy!(Vector{Ti}(k), 1:k)
colptr = copy!(Vector{Ti}(n + 1), 1:(k + 1))
for i in (k + 2):(n + 1) colptr[i] = (k + 1) end
SparseMatrixCSC{Tv,Ti}(dims..., colptr, rowval, nzval)
end
# convenience variations that accept a single integer to specify dims
SparseMatrixCSC{Tv,Ti}(s::UniformScaling, m::Integer) where {Tv,Ti} = SparseMatrixCSC{Tv,Ti}(s, m, m)
SparseMatrixCSC{Tv}(s::UniformScaling, m::Integer) where {Tv} = SparseMatrixCSC{Tv}(s, m, m)
SparseMatrixCSC(s::UniformScaling, m::Integer) = SparseMatrixCSC(s, m, m)

function speye_scaled(::Type{T}, diag, m::Integer, n::Integer) where T
((m < 0) || (n < 0)) && throw(ArgumentError("invalid array dimensions"))
if iszero(diag)
return SparseMatrixCSC(m, n, ones(Int, n+1), Vector{Int}(0), Vector{T}(0))
end
nnz = min(m,n)
colptr = Vector{Int}(1+n)
colptr[1:nnz+1] = 1:nnz+1
colptr[nnz+2:end] = nnz+1
SparseMatrixCSC(Int(m), Int(n), colptr, Vector{Int}(1:nnz), fill!(Vector{T}(nnz), diag))
end

sparse(S::UniformScaling, m::Integer, n::Integer=m) = speye_scaled(S.λ, m, n)

Base.iszero(A::SparseMatrixCSC) = iszero(view(A.nzval, 1:(A.colptr[size(A, 2) + 1] - 1)))

Expand Down
11 changes: 11 additions & 0 deletions test/linalg/uniformscaling.jl
Original file line number Diff line number Diff line change
Expand Up @@ -190,6 +190,17 @@ end
end
end

@testset "Matrix construction from UniformScaling" begin
@test Matrix(2I, 3)::Matrix{Int} == 2*eye(3)
@test Matrix(2I, 3, 3)::Matrix{Int} == 2*eye(3)
@test Matrix(2I, 3, 4)::Matrix{Int} == 2*eye(3, 4)
@test Matrix(2I, 4, 3)::Matrix{Int} == 2*eye(4, 3)
@test Matrix(2.0I, 3, 3)::Matrix{Float64} == 2*eye(3)
@test Matrix{Real}(2I, 3)::Matrix{Real} == 2*eye(3)
@test Matrix{Real}(2I, 3, 3)::Matrix{Real} == 2*eye(3)
@test Matrix{Float64}(2I, 3, 3)::Matrix{Float64} == 2*eye(3)
end

@testset "equality comparison of matrices with UniformScaling" begin
# AbstractMatrix methods
diagI = Diagonal(fill(1, 3))
Expand Down
16 changes: 16 additions & 0 deletions test/sparse/sparse.jl
Original file line number Diff line number Diff line change
Expand Up @@ -36,6 +36,22 @@ end
@test sparse(Any[1,2,3], Any[1,2,3], Any[1,1,1], 5, 4) == sparse([1,2,3], [1,2,3], [1,1,1], 5, 4)
end

@testset "SparseMatrixCSC construction from UniformScaling" begin
@test_throws ArgumentError SparseMatrixCSC(I, -1)
@test_throws ArgumentError SparseMatrixCSC(I, -1, 3)
@test_throws ArgumentError SparseMatrixCSC(I, 3, -1)
@test SparseMatrixCSC(2I, 3)::SparseMatrixCSC{Int,Int} == 2*eye(3)
@test SparseMatrixCSC(2I, 3, 3)::SparseMatrixCSC{Int,Int} == 2*eye(3)
@test SparseMatrixCSC(2I, 3, 4)::SparseMatrixCSC{Int,Int} == 2*eye(3, 4)
@test SparseMatrixCSC(2I, 4, 3)::SparseMatrixCSC{Int,Int} == 2*eye(4, 3)
@test SparseMatrixCSC(2.0I, 3, 3)::SparseMatrixCSC{Float64,Int} == 2*eye(3)
@test SparseMatrixCSC{Real}(2I, 3)::SparseMatrixCSC{Real,Int} == 2*eye(3)
@test SparseMatrixCSC{Real}(2I, 3, 3)::SparseMatrixCSC{Real,Int} == 2*eye(3)
@test SparseMatrixCSC{Float64}(2I, 3, 3)::SparseMatrixCSC{Float64,Int} == 2*eye(3)
@test SparseMatrixCSC{Float64,Int32}(2I, 3, 3)::SparseMatrixCSC{Float64,Int32} == 2*eye(3)
@test SparseMatrixCSC{Float64,Int32}(0I, 3, 3)::SparseMatrixCSC{Float64,Int32} == spzeros(Float64, Int32, 3, 3)
end

se33 = speye(3)
do33 = ones(3)

Expand Down

0 comments on commit d2d995d

Please sign in to comment.