Skip to content

Commit

Permalink
rewrite fillstored!(::SpecialMatrix, x) without generated function
Browse files Browse the repository at this point in the history
move implementations to more suitable places
  • Loading branch information
fredrikekre committed Dec 12, 2017
1 parent 8148227 commit bdae975
Show file tree
Hide file tree
Showing 4 changed files with 67 additions and 83 deletions.
34 changes: 0 additions & 34 deletions base/linalg/bidiag.jl
Original file line number Diff line number Diff line change
Expand Up @@ -601,37 +601,3 @@ function eigvecs(M::Bidiagonal{T}) where T
Q #Actually Triangular
end
eigfact(M::Bidiagonal) = Eigen(eigvals(M), eigvecs(M))

# fill! methods
_valuefields(::Type{<:Diagonal}) = [:diag]
_valuefields(::Type{<:Bidiagonal}) = [:dv, :ev]
_valuefields(::Type{<:Tridiagonal}) = [:dl, :d, :du]
_valuefields(::Type{<:SymTridiagonal}) = [:dv, :ev]
_valuefields(::Type{<:AbstractTriangular}) = [:data]

const SpecialArrays = Union{Diagonal,Bidiagonal,Tridiagonal,SymTridiagonal,AbstractTriangular}

function fillstored!(A::SpecialArrays, x)
xT = convert(eltype(A), x)
if @generated
quote
$([ :(fill!(A.$field, xT)) for field in _valuefields(A) ]...)
end
else
for field in _valuefields(A)
fill!(getfield(A, field), xT)
end
end
return A
end

_small_enough(A::Bidiagonal) = size(A, 1) <= 1
_small_enough(A::Tridiagonal) = size(A, 1) <= 2
_small_enough(A::SymTridiagonal) = size(A, 1) <= 2

function fill!(A::Union{Bidiagonal,Tridiagonal,SymTridiagonal}, x)
xT = convert(eltype(A), x)
(xT == zero(eltype(A)) || _small_enough(A)) && return fillstored!(A, xT)
throw(ArgumentError("array A of type $(typeof(A)) and size $(size(A)) can
not be filled with x=$x, since some of its entries are constrained."))
end
18 changes: 18 additions & 0 deletions base/linalg/special.jl
Original file line number Diff line number Diff line change
Expand Up @@ -124,3 +124,21 @@ end

A_mul_Bc!(A::AbstractTriangular, B::Union{QRCompactWYQ,QRPackedQ}) = A_mul_Bc!(full!(A), B)
A_mul_Bc(A::AbstractTriangular, B::Union{QRCompactWYQ,QRPackedQ}) = A_mul_Bc(copy!(similar(parent(A)), A), B)

# fill[stored]! methods
fillstored!(A::Diagonal, x) = (fill!(A.diag, x); A)
fillstored!(A::Bidiagonal, x) = (fill!(A.dv, x); fill!(A.ev, x); A)
fillstored!(A::Tridiagonal, x) = (fill!(A.dl, x); fill!(A.d, x); fill!(A.du, x); A)
fillstored!(A::SymTridiagonal, x) = (fill!(A.dv, x); fill!(A.ev, x); A)

_small_enough(A::Bidiagonal) = size(A, 1) <= 1
_small_enough(A::Tridiagonal) = size(A, 1) <= 2
_small_enough(A::SymTridiagonal) = size(A, 1) <= 2

# TODO: Add Diagonal to this method when 0.7 deprecations are removed
function fill!(A::Union{Bidiagonal,Tridiagonal,SymTridiagonal}, x)
xT = convert(eltype(A), x)
(iszero(xT) || _small_enough(A)) && return fillstored!(A, xT)
throw(ArgumentError("array of type $(typeof(A)) and size $(size(A)) can
not be filled with x=$x, since some of its entries are constrained."))
end
2 changes: 2 additions & 0 deletions base/linalg/triangular.jl
Original file line number Diff line number Diff line change
Expand Up @@ -418,6 +418,8 @@ end
scale!(A::Union{UpperTriangular,LowerTriangular}, c::Number) = scale!(A,A,c)
scale!(c::Number, A::Union{UpperTriangular,LowerTriangular}) = scale!(A,c)

fillstored!(A::AbstractTriangular, x) = (fill!(A.data, x); A)

# Binary operations
+(A::UpperTriangular, B::UpperTriangular) = UpperTriangular(A.data + B.data)
+(A::LowerTriangular, B::LowerTriangular) = LowerTriangular(A.data + B.data)
Expand Down
96 changes: 47 additions & 49 deletions test/linalg/bidiag.jl
Original file line number Diff line number Diff line change
Expand Up @@ -309,58 +309,56 @@ end
@test promote(C,A) isa Tuple{Tridiagonal, Tridiagonal}
end

import Base.LinAlg: fillstored!, UnitLowerTriangular
using Base.LinAlg: fillstored!, UnitLowerTriangular
@testset "fill! and fillstored!" begin
let #fill!
let # fillstored!
A = Tridiagonal(randn(2), randn(3), randn(2))
@test fillstored!(A, 3) == Tridiagonal([3, 3.], [3, 3, 3.], [3, 3.])
B = Bidiagonal(randn(3), randn(2), :U)
@test fillstored!(B, 2) == Bidiagonal([2.,2,2], [2,2.], :U)
S = SymTridiagonal(randn(3), randn(2))
@test fillstored!(S, 1) == SymTridiagonal([1,1,1.], [1,1.])
Ult = UnitLowerTriangular(randn(3,3))
@test fillstored!(Ult, 3) == UnitLowerTriangular([1 0 0; 3 1 0; 3 3 1])
end
let # fill!(exotic, 0)
exotic_arrays = Any[Tridiagonal(randn(3), randn(4), randn(3)),
Bidiagonal(randn(3), randn(2), rand([:U,:L])),
SymTridiagonal(randn(3), randn(2)),
sparse(randn(3,4)),
# Diagonal(randn(5)), # Diagonal fill! deprecated, see below
sparse(rand(3)),
# LowerTriangular(randn(3,3)), # AbstractTriangular fill! deprecated, see below
# UpperTriangular(randn(3,3)) # AbstractTriangular fill! deprecated, see below
]
for A in exotic_arrays
fill!(A, 0)
for a in A
@test a == 0
end
let # fillstored!
A = Tridiagonal(randn(2), randn(3), randn(2))
@test fillstored!(A, 3) == Tridiagonal([3, 3.], [3, 3, 3.], [3, 3.])
B = Bidiagonal(randn(3), randn(2), :U)
@test fillstored!(B, 2) == Bidiagonal([2.,2,2], [2,2.], :U)
S = SymTridiagonal(randn(3), randn(2))
@test fillstored!(S, 1) == SymTridiagonal([1,1,1.], [1,1.])
Ult = UnitLowerTriangular(randn(3,3))
@test fillstored!(Ult, 3) == UnitLowerTriangular([1 0 0; 3 1 0; 3 3 1])
end
let # fill!(exotic, 0)
exotic_arrays = Any[Tridiagonal(randn(3), randn(4), randn(3)),
Bidiagonal(randn(3), randn(2), rand([:U,:L])),
SymTridiagonal(randn(3), randn(2)),
sparse(randn(3,4)),
# Diagonal(randn(5)), # Diagonal fill! deprecated, see below
sparse(rand(3)),
# LowerTriangular(randn(3,3)), # AbstractTriangular fill! deprecated, see below
# UpperTriangular(randn(3,3)) # AbstractTriangular fill! deprecated, see below
]
for A in exotic_arrays
fill!(A, 0)
for a in A
@test a == 0
end
# Diagonal and AbstractTriangular fill! were defined as fillstored!,
# not matching the general behavior of fill!, and so have been deprecated.
# In a future dev cycle, these fill! methods should probably be reintroduced
# with behavior matching that of fill! for other structured matrix types.
# In the interm, equivalently test fillstored! below
@test iszero(fillstored!(Diagonal(fill(1, 3)), 0))
@test iszero(fillstored!(LowerTriangular(fill(1, 3, 3)), 0))
@test iszero(fillstored!(UpperTriangular(fill(1, 3, 3)), 0))
end
let # fill!(small, x)
val = randn()
b = Bidiagonal(randn(1,1), :U)
st = SymTridiagonal(randn(1,1))
for x in (b, st)
@test Array(fill!(x, val)) == fill!(Array(x), val)
end
b = Bidiagonal(randn(2,2), :U)
st = SymTridiagonal(randn(3), randn(2))
t = Tridiagonal(randn(3,3))
for x in (b, t, st)
@test_throws ArgumentError fill!(x, val)
@test Array(fill!(x, 0)) == fill!(Array(x), 0)
end
# Diagonal and AbstractTriangular fill! were defined as fillstored!,
# not matching the general behavior of fill!, and so have been deprecated.
# In a future dev cycle, these fill! methods should probably be reintroduced
# with behavior matching that of fill! for other structured matrix types.
# In the interm, equivalently test fillstored! below
@test iszero(fillstored!(Diagonal(fill(1, 3)), 0))
@test iszero(fillstored!(LowerTriangular(fill(1, 3, 3)), 0))
@test iszero(fillstored!(UpperTriangular(fill(1, 3, 3)), 0))
end
let # fill!(small, x)
val = randn()
b = Bidiagonal(randn(1,1), :U)
st = SymTridiagonal(randn(1,1))
for x in (b, st)
@test Array(fill!(x, val)) == fill!(Array(x), val)
end
b = Bidiagonal(randn(2,2), :U)
st = SymTridiagonal(randn(3), randn(2))
t = Tridiagonal(randn(3,3))
for x in (b, t, st)
@test_throws ArgumentError fill!(x, val)
@test Array(fill!(x, 0)) == fill!(Array(x), 0)
end
end
end
Expand Down

0 comments on commit bdae975

Please sign in to comment.