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

transition linalg to Adjoint/Transpose externally #25083

Merged
merged 21 commits into from
Dec 16, 2017

Conversation

Sacha0
Copy link
Member

@Sacha0 Sacha0 commented Dec 14, 2017

#24969 transitioned linalg to Adjoint/Transpose internally; this pull request transitions linalg to Adjoint/Transpose externally, completing item seven in #24969's OP and constituting the next major step towards JuliaLang/LinearAlgebra.jl#57.

In greater detail, this pull request contains commits that : (1) give Adjoint/Transpose all of ConjRowVector/RowVector's behaviors, making the former a drop-in replacement for the latter; (2) shield ConjRowVector/RowVector from downstream changes to adjoint(...), transpose(...), .', ', and A[ct]_(mul|ldiv|rdiv)_B[ct]; and (3) subsequently transition all forms of construction of ConjRowVector/RowVector (apart from their constructors) to Adjoint/Transpose and clean up the fallout. All intermediate commits pass the linalg, sparse, and stdlib test suites locally.

Ref. JuliaLang/LinearAlgebra.jl#57. Best!

@Sacha0
Copy link
Member Author

Sacha0 commented Dec 15, 2017

AV i686 failure is #25096, fixed by #25103. Otherwise CI appears happy 🎉!

Copy link
Member

@timholy timholy left a comment

Choose a reason for hiding this comment

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

This was such a heroic effort, can't thank you enough for tackling it @Sacha0!

I haven't had the time to follow this carefully, but I noticed one concern. Upon (very) superficial reading the rest of this seems fine and definitely should be part of 1.0.

# instead of simply mapping the caller's operation over the set of unwrapped vectors,
# we map Adjoint/Transpose composed with the caller's operationt over the set of unwrapped
# vectors. then re-wrapping the result vector yields a wrapped vector with the correct entries.
map(f, avs::AdjointAbsVec...) = Adjoint(map(Adjoint∘f, vectorfyall(avs...)...))
Copy link
Member

Choose a reason for hiding this comment

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

I don't understand why getindex doesn't just solve this for you. Are you trying to return an Adjoint-wrapped object? If so, why?

Copy link
Member Author

Choose a reason for hiding this comment

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

Indeed, returning an Adjoint/Transpose-wrapped vector from higher order operations over Adjoint/Transpose-wrapped vectors is these methods' purpose :).

Why? The superficial reason is preservation of ConjRowVector/RowVector's behavior in Adjoint/Transpose-wrapped vectors. The reason for ConjRowVector/RowVector's behavior: Without re-wrapping, the result of these operations is a single-row Matrix, which has different semantics from ConjRowVector/RowVector. Consequently, for example, foo.(anadjointvector) * avector would produce a one-by-one matrix rather than a scalar.

In other words, re-wrapping preserves 4774 semantics through higher order operations insofar as possible.

Does that clarify? Thoughts? :)

Copy link
Member

Choose a reason for hiding this comment

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

(Sidenote: typo in the comment operationt)

Copy link
Member Author

Choose a reason for hiding this comment

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

Good catch! Thanks! :)

Copy link
Member

Choose a reason for hiding this comment

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

What about Adjoint(map(Adjoint∘f∘Adjoint, parent.(avs)...))?

Copy link
Member

Choose a reason for hiding this comment

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

(so you don't have to create a temporary)

Copy link
Member Author

Choose a reason for hiding this comment

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

Avoiding the temporaries would be lovely! Deciphering the above expression is more than I have mental bandwidth for right at the moment though 😄. Perhaps a post-freeze optimization? :)

Copy link
Member Author

Choose a reason for hiding this comment

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

@timholy, revisiting this rewrite briefly, given Adjoint is unary I imagine you instead meant e.g. Adjoint(map((xs...) -> Adjoint(f(Adjoint.(xs)...)), parent.(avs)...))? Thanks again! :)

Copy link
Member

Choose a reason for hiding this comment

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

Yes, that looks right.

@Sacha0 Sacha0 changed the title [WIP] transition linalg to Adjoint/Transpose externally transition linalg to Adjoint/Transpose externally Dec 15, 2017
@Sacha0
Copy link
Member Author

Sacha0 commented Dec 15, 2017

Much thanks for the kind words @timholy! :)

A tangential observation: The lazification of adjoint and transpose creates substantial future optimization opportunity. A simple and obvious example, efficient mul!(Adjoint(C), A, B) (which the A_mul_B system could not express) becomes possible. Less obviously, operations like efficient foo(A::Adjoint{<:Any,<:AbstractMatrix}) for foo some matrix operation (e.g. log, sqrt, pow, lufact, qrfact, and so on) become possible as well (instead of foo(adjoint(A)) with adjoint(A) eager).

@Sacha0 Sacha0 added this to the 1.0 milestone Dec 16, 2017
…ose(...), .', ', and A[ct]_(mul|ldiv|rdiv)_B[ct][!].
…ition associated tests to Adjoint/Transpose.
@Sacha0
Copy link
Member Author

Sacha0 commented Dec 16, 2017

Rebased and settled a couple todos. Assuming CI remains happy, I plan to merge this pull request tomorrow morning and carry on with downstream steps towards JuliaLang/LinearAlgebra.jl#57. Thanks all! :)

@Sacha0
Copy link
Member Author

Sacha0 commented Dec 16, 2017

The AV i686 failure is unrelated, the FreeBSD logs are not accessible (but FreeBSD passed formerly IIRC), and (given CircleCI did not run over night) continuing to wait for CircleCI does not seem worthwhile. Thanks all! :)

@martinholters
Copy link
Member

I haven't bisected, but most probably this PR or one of its cousins is responsible for

julia> v = SparseVector([1.0])'
1×1 Adjoint{Float64,SparseVector{Float64,Int64}}:
 1.0

julia> v - v
1×1 Adjoint{Any,SparseVector{Any,Int64}}:
Error showing value of type Adjoint{Any,SparseVector{Any,Int64}}:
ERROR: MethodError: no method matching zero(::Type{Any})

Note that the error by itself is only a consequence of the element type becoming Any, which the the real problem. (This is also a consequence of JuliaSparse/SparseArrays.jl#59, but we'd want to have the element type be inferable in this easy case anyway.)

@Sacha0 do you have anything in the pipeline that will fix this, or should I open a dedicated issue?

@martinholters
Copy link
Member

This seems to fix it:

--- a/base/linalg/adjtrans.jl
+++ b/base/linalg/adjtrans.jl
@@ -158,8 +158,8 @@ map(f, tvs::TransposeAbsVec...) = Transpose(map(Transpose∘f, vectorfyall(tvs..
 
 # broadcast over collections of Adjoint/Transpose-wrapped vectors and numbers
 # similar explanation for these definitions as for map above
-broadcast(f, avs::Union{Number,AdjointAbsVec}...) = Adjoint(broadcast(Adjoint∘f, vectorfyall(avs...)...))
-broadcast(f, tvs::Union{Number,TransposeAbsVec}...) = Transpose(broadcast(Transpose∘f, vectorfyall(tvs...) ...))
+broadcast(f, avs::Union{Number,AdjointAbsVec}...) = Adjoint(broadcast((x->Adjoint(x))∘f, vectorfyall(avs...)...))
+broadcast(f, tvs::Union{Number,TransposeAbsVec}...) = Transpose(broadcast((x->Transpose(x))∘f, vectorfyall(tvs...) ...))

I don't have the time to prepare a PR right now, so if anyone would be so kind to pick this up...

@Sacha0
Copy link
Member Author

Sacha0 commented Dec 20, 2017

Great catch @martinholters! Perhaps JuliaLang/LinearAlgebra.jl#492 will resolve this issue as a side effect. I will have a look later today. Thanks! :)

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
linear algebra Linear algebra
Projects
None yet
Development

Successfully merging this pull request may close these issues.

4 participants