diff --git a/src/addition.jl b/src/addition.jl index 94f8797..ab4cbc0 100644 --- a/src/addition.jl +++ b/src/addition.jl @@ -2,7 +2,7 @@ # between pdmat and pdmat +(a::PDMat, b::AbstractPDMat) = PDMat(a.mat + Matrix(b)) -+(a::PDiagMat, b::AbstractPDMat) = PDMat(_adddiag!(Matrix(b), a.diag)) ++(a::PDiagMat, b::AbstractPDMat) = PDMat(_adddiag!(Matrix(b), a.diag, true)) +(a::ScalMat, b::AbstractPDMat) = PDMat(_adddiag!(Matrix(b), a.value)) if HAVE_CHOLMOD +(a::PDSparseMat, b::AbstractPDMat) = PDMat(a.mat + Matrix(b)) diff --git a/src/utils.jl b/src/utils.jl index 3e1453f..c85d938 100644 --- a/src/utils.jl +++ b/src/utils.jl @@ -30,7 +30,7 @@ function _adddiag!(a::Union{Matrix, SparseMatrixCSC}, v::Real) return a end -function _adddiag!(a::Union{Matrix, SparseMatrixCSC}, v::Vector, c::Real) +function _adddiag!(a::Union{Matrix, SparseMatrixCSC}, v::AbstractVector, c::Real) @check_argdims eachindex(v) == axes(a, 1) == axes(a, 2) if c == one(c) for i in eachindex(v) @@ -45,8 +45,8 @@ function _adddiag!(a::Union{Matrix, SparseMatrixCSC}, v::Vector, c::Real) end _adddiag(a::Union{Matrix, SparseMatrixCSC}, v::Real) = _adddiag!(copy(a), v) -_adddiag(a::Union{Matrix, SparseMatrixCSC}, v::Vector, c::Real) = _adddiag!(copy(a), v, c) -_adddiag(a::Union{Matrix, SparseMatrixCSC}, v::Vector{T}) where {T<:Real} = _adddiag!(copy(a), v, one(T)) +_adddiag(a::Union{Matrix, SparseMatrixCSC}, v::AbstractVector, c::Real) = _adddiag!(copy(a), v, c) +_adddiag(a::Union{Matrix, SparseMatrixCSC}, v::AbstractVector{T}) where {T<:Real} = _adddiag!(copy(a), v, one(T)) function wsumsq(w::AbstractVector, a::AbstractVector) diff --git a/test/addition.jl b/test/addition.jl index 9de57f2..32b48e1 100644 --- a/test/addition.jl +++ b/test/addition.jl @@ -2,6 +2,18 @@ using PDMats + +# New AbstractPDMat type for the tests below +# Supports only functions needed in the tests below +struct ScalMat3D{T<:Real} <: AbstractPDMat{T} + value::T +end +Base.Matrix(a::ScalMat3D) = Matrix(Diagonal(fill(a.value, 3))) +Base.size(::ScalMat3D) = (3, 3) +# Not generally correct +Base.:*(a::ScalMat3D, c::Real) = ScalMat3D(a.value * c) +Base.getindex(a::ScalMat3D, i::Int, j::Int) = i == j ? a.value : zero(a.value) + @testset "addition" begin for T in (Float64, Float32) printstyled("Testing addition with eltype = $T\n"; color=:blue) @@ -11,26 +23,35 @@ using PDMats pm1 = PDMat(M) pm2 = PDiagMat(V) - pm3 = ScalMat(3, X) - pm4 = X * I + pm3 = PDiagMat(sparse(V)) + pm4 = ScalMat(3, X) pm5 = PDSparseMat(sparse(M)) + pm6 = ScalMat3D(X) - pmats = Any[pm1, pm2, pm3] #, pm5] + pmats = Any[pm1, pm2, pm3, pm4, pm5, pm6] for p1 in pmats, p2 in pmats pr = p1 + p2 @test size(pr) == size(p1) @test Matrix(pr) ≈ Matrix(p1) + Matrix(p2) - pr = pdadd(p1, p2, convert(T, 1.5)) - @test size(pr) == size(p1) - @test Matrix(pr) ≈ Matrix(p1) + Matrix(p2) * convert(T, 1.5) + if p1 isa ScalMat3D + @test_broken pdadd(p1, p2, convert(T, 1.5)) + else + pr = pdadd(p1, p2, convert(T, 1.5)) + @test size(pr) == size(p1) + @test Matrix(pr) ≈ Matrix(p1) + Matrix(p2) * convert(T, 1.5) + end end for p1 in pmats - pr = p1 + pm4 - @test size(pr) == size(p1) - @test Matrix(pr) ≈ Matrix(p1) + pm4 + if p1 isa ScalMat3D + @test_broken p1 + X * I + else + pr = p1 + X * I + @test size(pr) == size(p1) + @test Matrix(pr) ≈ Matrix(p1) + X * I + end end end end