Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Make convert(::Type{(Lower|Upper)Triangular}, A::Bidiagonal) preserve storage structure. Add equivalent methods for Unit(Lower|Upper)Triangular. #17656

Closed
wants to merge 1 commit into from
Closed

Make convert(::Type{(Lower|Upper)Triangular}, A::Bidiagonal) preserve storage structure. Add equivalent methods for Unit(Lower|Upper)Triangular. #17656

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
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
34 changes: 32 additions & 2 deletions base/linalg/special.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -6,8 +6,38 @@
convert{T}(::Type{Bidiagonal}, A::Diagonal{T})=Bidiagonal(A.diag, zeros(T, size(A.diag,1)-1), true)
convert{T}(::Type{SymTridiagonal}, A::Diagonal{T})=SymTridiagonal(A.diag, zeros(T, size(A.diag,1)-1))
convert{T}(::Type{Tridiagonal}, A::Diagonal{T})=Tridiagonal(zeros(T, size(A.diag,1)-1), A.diag, zeros(T, size(A.diag,1)-1))
convert(::Type{LowerTriangular}, A::Bidiagonal) = !A.isupper ? LowerTriangular(full(A)) : throw(ArgumentError("Bidiagonal matrix must have lower off diagonal to be converted to LowerTriangular"))
convert(::Type{UpperTriangular}, A::Bidiagonal) = A.isupper ? UpperTriangular(full(A)) : throw(ArgumentError("Bidiagonal matrix must have upper off diagonal to be converted to UpperTriangular"))

# methods for conversion from Bidiagonal to [Unit](Upper|Lower)Triangular
function convert(::Type{LowerTriangular}, A::Bidiagonal)
if A.isupper
throw(ArgumentError("upper Bidiagonal matrices cannot be converted to LowerTriangular"))
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I guess it could if all the off-diagonals are zero?

Copy link
Member Author

@Sacha0 Sacha0 Jul 27, 2016 •
edited
Loading

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It could, which segues to a related issue.

On the one hand, existing {constructors for annotation types} that accept an A<:AbstractMatrix as the primary argument do not require that A 'conform' to the annotation. For example, LowerTriangular(Bidiagonal(rand(5), rand(4), true)) and Symmetric(rand(5,5)) are always valid.

On the other hand, convert(::Type{$(AnnotationType)}, A::$(SpecialMatrixType))-like methods inconsistently require that A 'conform' to the annotation type either in storage structure (e.g. as in the code you identified) or in value (e.g. as in your suggestion). And equivalent convert methods from other storage types (e.g. Matrix, SparseMatrixCSC) to annotation types do not appear to exist. For example, convert(LowerTriangular, rand(4, 4)) throws a MethodError.

Do convert methods from special matrix types to annotation types make sense generally? If so, how should they behave? Also if so, shouldn't equivalent convert methods from other storage types to annotation types exist? And otherwise, why do these particular methods exist?

Thoughts? Best!

Edit: Maybe I should open a dedicated issue?
Edit: grep -rn 'convert(LowerTriangular base/` yields nothing, FWIW.

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

You're very much getting into #8001 territory here. If you think of the "annotation wrapper type" as a view into whatever subset of the underlying array happens to be important for that type of annotation, then I guess the structure and the contents of the "inactive" part of the wrapper-type view don't actually matter and could be anything.

Is wrapping around an object always the meaning we should give for convert here? Maybe not as that could get confusing. But then since constructors and convert are pretty closely tied together at the moment, I'm not sure what syntax you would use as an alternative to the constructors.

end
LowerTriangular(A)
end
function convert(::Type{UpperTriangular}, A::Bidiagonal)
if !A.isupper
throw(ArgumentError("lower Bidiagonal matrices cannot be converted to UpperTriangular"))
end
UpperTriangular(A)
end
function convert(::Type{UnitLowerTriangular}, A::Bidiagonal)
if A.isupper
throw(ArgumentError("upper Bidiagonal matrices cannot be converted to UnitLowerTriangular"))
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

and here if all the off-diagonals are zero and all the diagonals are 1?

elseif !all(A.dv .== one(eltype(A)))
throw(ArgumentError(string("Bidiagonal matrices with non-one entries on the ",
"diagonal cannot be converted to UnitLowerTriangular")))
end
UnitLowerTriangular(A)
end
function convert(::Type{UnitUpperTriangular}, A::Bidiagonal)
if !A.isupper
throw(ArgumentError("lower Bidiagonal matrices cannot be converted to UnitUpperTriangular"))
elseif !all(A.dv .== one(eltype(A)))
throw(ArgumentError(string("Bidiagonal matrices with non-one entries on the ",
"diagonal cannot be converted to UnitUpperTriangular")))
end
UnitUpperTriangular(A)
end

function convert(::Type{UnitUpperTriangular}, A::Diagonal)
if !all(A.diag .== one(eltype(A)))
Expand Down
25 changes: 25 additions & 0 deletions test/linalg/special.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -128,3 +128,28 @@ for typ in [UpperTriangular,LowerTriangular,Base.LinAlg.UnitUpperTriangular,Base
@test Base.LinAlg.A_mul_Bc(atri,qrb[:Q]) ≈ full(atri) * qrb[:Q]'
@test Base.LinAlg.A_mul_Bc!(copy(atri),qrb[:Q]) ≈ full(atri) * qrb[:Q]'
end

# Test conversion from Bidiagonal to <:AbstractTriangular
let
lowerbidiagmat = Bidiagonal(rand(3), rand(2), false)
upperbidiagmat = Bidiagonal(rand(3), rand(2), true)
unitlowerbidiagmat = Bidiagonal(ones(3), rand(2), false)
unitupperbidiagmat = Bidiagonal(ones(3), rand(2), true)
# test that conversion from upper(lower) bidiagonal to lower(upper) [unit]triangular throws
@test_throws ArgumentError convert(UpperTriangular, unitlowerbidiagmat)
@test_throws ArgumentError convert(LowerTriangular, unitupperbidiagmat)
@test_throws ArgumentError convert(Base.LinAlg.UnitUpperTriangular, unitlowerbidiagmat)
@test_throws ArgumentError convert(Base.LinAlg.UnitLowerTriangular, unitupperbidiagmat)
# test that conversion from non-unit bidiag to unit triangular throws
@test_throws ArgumentError convert(Base.LinAlg.UnitUpperTriangular, upperbidiagmat)
@test_throws ArgumentError convert(Base.LinAlg.UnitLowerTriangular, lowerbidiagmat)
# test that conversion from bidiagonal to triangular preserves bidiagonal storage structure
@test typeof(convert(UpperTriangular, upperbidiagmat)) ==
UpperTriangular{eltype(upperbidiagmat),Bidiagonal{eltype(upperbidiagmat)}}
@test typeof(convert(LowerTriangular, lowerbidiagmat)) ==
LowerTriangular{eltype(lowerbidiagmat),Bidiagonal{eltype(lowerbidiagmat)}}
@test typeof(convert(Base.LinAlg.UnitUpperTriangular, unitupperbidiagmat)) ==
Base.LinAlg.UnitUpperTriangular{eltype(unitupperbidiagmat),Bidiagonal{eltype(unitupperbidiagmat)}}
@test typeof(convert(Base.LinAlg.UnitLowerTriangular, unitlowerbidiagmat)) ==
Base.LinAlg.UnitLowerTriangular{eltype(unitlowerbidiagmat),Bidiagonal{eltype(unitlowerbidiagmat)}}
end