Skip to content
Draft
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 2 additions & 0 deletions src/bidiag.jl
Original file line number Diff line number Diff line change
Expand Up @@ -435,6 +435,8 @@ function diag(M::Bidiagonal, n::Integer=0)
return v
end

isstoredband(A::Bidiagonal, k::Integer) = k == 0 || k == _offdiagind(A.uplo)

function +(A::Bidiagonal, B::Bidiagonal)
if A.uplo == B.uplo || length(A.dv) == 0
Bidiagonal(A.dv+B.dv, A.ev+B.ev, A.uplo)
Expand Down
12 changes: 12 additions & 0 deletions src/dense.jl
Original file line number Diff line number Diff line change
Expand Up @@ -237,6 +237,18 @@ end

fillstored!(A::AbstractMatrix, v) = fill!(A, v)

"""
isstoredband(A::AbstractMatrix, k::Integer)

Return whether the `k`-th band of `A` is stored as a vector, as opposed to
being generated during indexing.
For example, only the principal diagonal would be stored for a `Diagonal`.

!!! note
This is a conservative check that may have false negatives but should not have false positives.
"""
isstoredband(A::AbstractMatrix, k::Integer) = false

diagind(m::Integer, n::Integer, k::Integer=0) = diagind(IndexLinear(), m, n, k)
diagind(::IndexLinear, m::Integer, n::Integer, k::Integer=0) =
k <= 0 ? range(1-k, step=m+1, length=min(m+k, n)) : range(k*m+1, step=m+1, length=min(m, n-k))
Expand Down
3 changes: 3 additions & 0 deletions src/diagonal.jl
Original file line number Diff line number Diff line change
Expand Up @@ -940,6 +940,9 @@ function diag(D::Diagonal, k::Integer=0)
end
return v
end

isstoredband(::Diagonal, k::Integer) = k == 0

tr(D::Diagonal) = sum(tr, D.diag)
det(D::Diagonal) = prod(det, D.diag)
function logdet(D::Diagonal{<:Complex}) # make sure branch cut is correct
Expand Down
2 changes: 2 additions & 0 deletions src/hessenberg.jl
Original file line number Diff line number Diff line change
Expand Up @@ -667,3 +667,5 @@ function logdet(F::Hessenberg)
d,s = logabsdet(F)
return d + log(s)
end

isstoredband(U::UpperHessenberg, k::Integer) = k >= -1 && isstoredband(parent(U), k)
5 changes: 5 additions & 0 deletions src/triangular.jl
Original file line number Diff line number Diff line change
Expand Up @@ -238,6 +238,11 @@ Base.isassigned(A::UpperOrLowerTriangular, i::Int, j::Int) =
Base.isstored(A::UpperOrLowerTriangular, i::Int, j::Int) =
_shouldforwardindex(A, i, j) ? Base.isstored(A.data, i, j) : false

isstoredband(U::UpperTriangular, k::Integer) = k >= 0 && isstoredband(parent(U), k)
isstoredband(L::LowerTriangular, k::Integer) = k <= 0 && isstoredband(parent(L), k)
isstoredband(U::UnitUpperTriangular, k::Integer) = k > 0 && isstoredband(parent(U), k)
isstoredband(L::UnitLowerTriangular, k::Integer) = k < 0 && isstoredband(parent(L), k)

@propagate_inbounds function getindex(A::Union{UnitLowerTriangular{T}, UnitUpperTriangular{T}}, i::Int, j::Int) where {T}
if _shouldforwardindex(A, i, j)
A.data[i,j]
Expand Down
2 changes: 2 additions & 0 deletions src/tridiag.jl
Original file line number Diff line number Diff line change
Expand Up @@ -1230,3 +1230,5 @@ function fillband!(T::SymTridiagonal, x, l, u)
end
return T
end

isstoredband(::Union{SymTridiagonal, Tridiagonal}, k::Integer) = -1 <= k <= 1
27 changes: 27 additions & 0 deletions test/hessenberg.jl
Original file line number Diff line number Diff line change
Expand Up @@ -320,4 +320,31 @@ end
@test U == U2
end

@testset "isstoredband" begin
U = UpperHessenberg(Diagonal(1:4))
@test LinearAlgebra.isstoredband(U, 0)
@test !LinearAlgebra.isstoredband(U, 1)
@test !LinearAlgebra.isstoredband(U, -1)

U = UpperHessenberg(Bidiagonal(1:4, 1:3, :U))
@test LinearAlgebra.isstoredband(U, 0)
@test LinearAlgebra.isstoredband(U, 1)
@test !LinearAlgebra.isstoredband(U, 2)
@test !LinearAlgebra.isstoredband(U, -1)

U = UpperHessenberg(Bidiagonal(1:4, 1:3, :L))
@test LinearAlgebra.isstoredband(U, 0)
@test LinearAlgebra.isstoredband(U, -1)
@test !LinearAlgebra.isstoredband(U, 1)
@test !LinearAlgebra.isstoredband(U, 2)

for Tri in (Tridiagonal(1:3, 1:4, 1:3), SymTridiagonal(1:4, 1:3))
U = UpperHessenberg(Tri)
@test LinearAlgebra.isstoredband(U, 0)
@test LinearAlgebra.isstoredband(U, 1)
@test LinearAlgebra.isstoredband(U, -1)
@test !LinearAlgebra.isstoredband(U, 2)
end
end

end # module TestHessenberg
43 changes: 43 additions & 0 deletions test/triangular.jl
Original file line number Diff line number Diff line change
Expand Up @@ -1110,4 +1110,47 @@ end
end
end

@testset "isstoredband" begin
U = UpperTriangular(Diagonal(1:4))
@test LinearAlgebra.isstoredband(U, 0)
@test !LinearAlgebra.isstoredband(U, 1)
@test !LinearAlgebra.isstoredband(U, -1)
L = LowerTriangular(Diagonal(1:4))
@test LinearAlgebra.isstoredband(L, 0)
@test !LinearAlgebra.isstoredband(L, 1)
@test !LinearAlgebra.isstoredband(L, -1)
for T in (UnitUpperTriangular, UnitLowerTriangular)
U = T(Diagonal(1:4))
@test !LinearAlgebra.isstoredband(U, 0)
@test !LinearAlgebra.isstoredband(U, 1)
@test !LinearAlgebra.isstoredband(U, -1)
end

U = UpperTriangular(Bidiagonal(1:4, 1:3, :U))
@test LinearAlgebra.isstoredband(U, 0)
@test LinearAlgebra.isstoredband(U, 1)
@test !LinearAlgebra.isstoredband(U, 2)
@test !LinearAlgebra.isstoredband(U, -1)

U = UpperTriangular(Bidiagonal(1:4, 1:3, :L))
@test LinearAlgebra.isstoredband(U, 0)
@test !LinearAlgebra.isstoredband(U, 1)
@test !LinearAlgebra.isstoredband(U, 2)
@test !LinearAlgebra.isstoredband(U, -1)

for Tri in (Tridiagonal(1:3, 1:4, 1:3), SymTridiagonal(1:4, 1:3))
U = UpperTriangular(Tri)
@test LinearAlgebra.isstoredband(U, 0)
@test LinearAlgebra.isstoredband(U, 1)
@test !LinearAlgebra.isstoredband(U, 2)
@test !LinearAlgebra.isstoredband(U, -1)

L = LowerTriangular(Tri)
@test LinearAlgebra.isstoredband(L, 0)
@test LinearAlgebra.isstoredband(L, -1)
@test !LinearAlgebra.isstoredband(L, -2)
@test !LinearAlgebra.isstoredband(L, 1)
end
end

end # module TestTriangular