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.

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

LinAlg.fillslots! -> LinAlg.fillstored! #25030

Merged
merged 4 commits into from
Dec 14, 2017

Conversation

fredrikekre
Copy link
Member

@fredrikekre fredrikekre commented Dec 11, 2017

See #17670

@fredrikekre fredrikekre added deprecation This change introduces or involves a deprecation linear algebra Linear algebra labels Dec 11, 2017
# 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)
(xT == zero(eltype(A)) || _small_enough(A)) && return fillstored!(A, xT)
Copy link
Member

Choose a reason for hiding this comment

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

iszero(xT)? :)

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
Copy link
Member

Choose a reason for hiding this comment

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

Perhaps nix the A? The user will not be aware of the array's name in the method signature :).

@@ -309,18 +309,18 @@ end
@test promote(C,A) isa Tuple{Tridiagonal, Tridiagonal}
end

import Base.LinAlg: fillslots!, UnitLowerTriangular
@testset "fill! and fillslots!" begin
import Base.LinAlg: fillstored!, UnitLowerTriangular
Copy link
Member

Choose a reason for hiding this comment

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

Perhaps import -> using while touching this code?

import Base.LinAlg: fillslots!, UnitLowerTriangular
@testset "fill! and fillslots!" begin
import Base.LinAlg: fillstored!, UnitLowerTriangular
@testset "fill! and fillstored!" begin
let #fill!
Copy link
Member

Choose a reason for hiding this comment

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

This enclosing let seems vestigial; perhaps remove and de-indent while touching this code?

@@ -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)
Copy link
Member

Choose a reason for hiding this comment

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

I wonder whether we should constrain the mutated part of A.data in case something else aliases A.data?

Copy link
Member Author

Choose a reason for hiding this comment

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

It would be very weird to expect the following:

julia> A = fill(2., 2, 2);

julia> T = UpperTriangular(A);

julia> LinAlg.fillstored!(T, 0);

julia> A
2×2 Array{Float64,2}:
 0.0  0.0
 2.0  0.0

I think those wrappers should have the license to do whatever they like with the ignored part. (Also, the ignored part is still stored so the function name implies that we will fill that too :trollface: )

Copy link
Member

Choose a reason for hiding this comment

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

The behavior you highlight does not strike me as weird, but rather as desirable :). Consider, for example, that packing the factors from a decomposition into a single matrix is common. Were you using <:AbstractTriangular wrappers to work with / manipulate those wrappers externally (which we do in LinAlg IIRC), the behavior above is correct, whereas the implemented behavior silently corrupts other objects. Another angle from which to consider this point: Should fillstored!(view(A, 1:2, 1:2), 0) corrupt all other data in A? :)

Copy link
Member

@KristofferC KristofferC Dec 12, 2017

Choose a reason for hiding this comment

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

A view clearly does not own it's wrapped Array while I would argue a type like UpperTriangular does.
We have no way in Julia to express these semantics but giving an Array to something like UpperTriangular and keep reusing that Array while having some expectations on how this type will modify it seems very suspicious and error prone. Ideally, there would be a way where the original data buffers could be "moved" to another binding while making it illegal to use the prevous one (std::movestyle).

However, since we don't have that, it is imo fair game for the type to do whatever you want with the passed in Array and no consideration has to be taken to any future use of this Array.

Copy link
Member

Choose a reason for hiding this comment

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

A view clearly does not own it's wrapped Array

<:AbstractTriangular are also views, and just as a SubArray does not own its parent's storage beyond the SubArray's purview, a <:AbstractTriangular should not either. Neither SubArrays nor <:AbstractTriangular allow mutation of the parent beyond the view's purview via setindex! and, insofar as I am aware, no other operations over such views allow/cause such mutation either. Implementations of operations over such views that do perform such mutation break the abstraction, and in this case unnecessarily so. Moreover, as mentioned above there are good reasons not to perform such mutation: Packed storage is common in numerical linear algebra, as in e.g. most dense triangular or orthogonal decompositions. See, e.g., base/linalg/lu.jl for operations over triangular views of subsections of packed decompositions that such mutation could jeopardize. Best!

Copy link
Member

@Sacha0 Sacha0 left a comment

Choose a reason for hiding this comment

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

Generally looks great! Thanks @fredrikekre! :)

(The only non-negligible review comment concerns the fillstored!(A::AbstractTriangular, x) definition. Not certain what's best.)

@fredrikekre fredrikekre force-pushed the fe/fillstored branch 2 times, most recently from 68865e2 to 60f23df Compare December 13, 2017 17:17
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."))
Copy link
Member

@Sacha0 Sacha0 Dec 14, 2017

Choose a reason for hiding this comment

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

Perhaps nix the x=, given that what x refers to is not user-visible?

@testset "fill! and fillstored!" begin
let # fillstored!
A = Tridiagonal(randn(2), randn(3), randn(2))
@test fillstored!(A, 3) == Tridiagonal([3, 3.], [3, 3, 3.], [3, 3.])
Copy link
Member

Choose a reason for hiding this comment

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

No need for the trailing .s given that == does not check type equality?

]
for A in exotic_arrays
fill!(A, 0)
for a in A
Copy link
Member

@Sacha0 Sacha0 Dec 14, 2017

Choose a reason for hiding this comment

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

Perhaps @test iszero(fill!(A, 0))?

Copy link
Member

@Sacha0 Sacha0 left a comment

Choose a reason for hiding this comment

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

Looks great! :) I will give this another pass tomorrow (with benefit of sleep).

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
Copy link
Member

Choose a reason for hiding this comment

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

Perhaps simplify to the following?

@test fill!(Bidiagonal(fill(0, 1, 1), :U), 2) == fill(2, 1, 1)
@test fill!(SymTridiagonal(fill(0, 1, 1)), 2) == fill(2, 1, 1)

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])
Copy link
Member

Choose a reason for hiding this comment

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

Perhaps the other triangular types are worth testing as well?

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)
Copy link
Member

Choose a reason for hiding this comment

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

Perhaps simplify to @test fill!(x, 0) == ...?

Copy link
Member

@Sacha0 Sacha0 left a comment

Choose a reason for hiding this comment

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

Looks great! Thanks @fredrikekre! Mergeworthy insofar as I am concerned :).

@fredrikekre fredrikekre merged commit 1f33cdd into JuliaLang:master Dec 14, 2017
@fredrikekre fredrikekre deleted the fe/fillstored branch December 14, 2017 23:41
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
deprecation This change introduces or involves a deprecation linear algebra Linear algebra
Projects
None yet
Development

Successfully merging this pull request may close these issues.

3 participants