Jeff Bezanson PhD

[stevengj]

One of the core Julia authors, @JeffBezanson, has become a problematic developer. He needs to graduate from MIT, ideally by January 2015. Dependencies:

• [x] Thesis proposal.
• [x] Thesis outline.
• [x] Meeting with thesis committee.
• [x] Rough draft sent to advisor(s).
• [x] Schedule defense.
• [x] Final draft sent to committee.
• [x] Defense.
• [x] Alcohol.

This is a priority issue, to ensure that arms are not broken and to guarantee long-term viability of the Julia project.

cc: @alanedelman, @jiahao, @StefanKarpinski, @ViralBShah, @samanamarasinghe, @gjs

Edit (VS): This issue is closed with the following thesis. I am putting it up here, since many people will be interested in finding it. https://github.com/JeffBezanson/phdthesis/blob/master/main.pdf

[johnmyleswhite]

cc @fperez who also is interested in this outstanding issue

[jiahao]

Supporting information attached.

[JeffBezanson]

As a procedural issue, closing this might require work in a separate repository. The thesis should perhaps be included in Base to make it easier for others to contribute.

Also the ordering of the last task is misleading; it will actually recur frequently throughout the process.

[stevengj]

+1 for including it in Base, or at least in julia/doc/thesis. Or maybe theses to allow for future needs.

(Please go ahead and open up a thesis branch, Jeff.)

[jiahao]

Also the ordering of the last task is misleading; it will actually recur frequently throughout the process.

*has also already recurred

[stevengj]

I'm looking forward to being present at the Close issue ceremony.

[timholy]

I can't reproduce this issue locally; is it MIT-specific?

[jakebolewski]

One of the core Julia authors, @JeffBezanson, has become a problematic developer academic

[tknopp]

Is this the github version of a PHD thesis? Jeff has to open a PR with his proposal and the committee will decide whether to merge or not...

[ViralBShah]

+Inf for speedy resolution of this one!

[Carreau]

I had the same problem on IPython repo a few month ago, hopefully it was fixed 32 days ago. I'm pretty sure it involved coffee, annoying paperwork and last minute change of plan because of jackhammers.

Good luck !

[stevengj]

Updated: Jeff met with his thesis committee and gave us a rough outline.

[timholy]

Glad to hear progress is being made!

But git log --date=short --all --since=22.days.ago --author="Jeff Bezanson" still makes one wonder how he has time for writing a thesis. Either that, or he's a superhero. Actually, scratch that: we all know he is a superhero, so never mind.

[jiahao]

The commits involving juliatypes.jl record our attempts to describe Julia's type system, which is directly relevant thesis work.

[johnmyleswhite]

The type system work seems to already be hitting some nerves: https://twitter.com/plt_hulk/status/535045242920378369

[StefanKarpinski]

I kind of doubt that's directly in response to Jeff's work, although I could be wrong. Hilarious tweet either way though.

[timholy]

@jiahao, my comment was mostly tongue-in-cheek---I kinda wondered about that very thing. I, for one, tend to have a lot of commits when I'm polishing something up for presentation.

[jiahao]

@timholy humor noted. :)

It would be remiss to not mention the lovely Magritte homage our local theory collaborator @jeanqasaur made and posted on twitter:

"The Treachery of Types" has a nice ring to it, no?

[tonyhffong]

That's quite funny.

[timholy]

Love it!

[stevengj]

Call for help: Jeff is looking for nice examples that show off multiple dispatch (and maybe staged functions), things which would be much harder/slower in languages without those features.

[tonyhffong]

er, show?

[timholy]

(and maybe staged functions)

subarray.jl and subarray2.jl should serve rather nicely. Design document is at http://docs.julialang.org/en/latest/devdocs/subarrays/

[johnmyleswhite]

I think the Distributions package really makes multiple dispatch seem useful. Having things like rand(Gamma(1, 1), 5, 5) vs rand(Normal(0, 1), 3) is a huge gain in expressivity at no performance cost because of multiple dispatch.

[jakebolewski]

I don't see how that is the best example as it really is showing off single-dispatch. How is it different than Gamma(1,1).rand(5,5) which you would do in a more traditional OO language like Python or Java?

On Fri, Dec 19, 2014 at 1:39 PM, John Myles White notifications@github.com wrote:

I think the Distributions package really makes multiple dispatch seem useful. Having things like rand(Gamma(1, 1), 5, 5) vs rand(Normal(0, 1), 3) is a huge gain in expressivity at no performance cost because of multiple dispatch.

— Reply to this email directly or view it on GitHub https://github.com/JuliaLang/julia/issues/8839#issuecomment-67678367.

[johnmyleswhite]

Ok. Replace that with examples of computing KL-divergences using analytic results: kl(Normal(0, 1), Normal(0, 1)) vs kl(Normal(0, 1), Gamma(1, 1)).

[timholy]

I should have also added, there were some potentially-useful statistics on what life would have been like without staged functions in my initial post in #8235. The take-home message: generating all methods through dimension 8 resulted in > 5000 separate methods, and required over 4 minutes of parsing&lowering time (i.e., a 4-minute delay while compiling julia). By comparison, the stagedfunction implementation loads in a snap, and of course can go even beyond 8 dimensions.

[vtjnash]

Vs single-dispatch, it still demonstrates the unification of what other OO languages describe as: functions vs methods. You could constrast to python's sorted(a) vs a.sort(). In comparison to "traditional" OO languages, it dramatically changes what it means for a function to be "associated with" a class.

You could point out how it replaces the need for static vs instance variables, since you can dispatch on the instance or the type. I might have some more ideas from some recent IRC conversations, when I can get to a computer.

[eschnett]

I have an implementation of fmap. This traverses several arrays and applies a function to each set of elements. This implementation is actually very slow since the number of arrays can be arbitrary. To make this useful, I have manually created specialization of this for various numbers of arguments yss. I've always wanted to write a staged function for this, but haven't done so yet.

The staged function would need to evaluate in particular the map call that produces the arguments to f.

-erik

function fmap{T,D}(f::Function, xs::Array{T,D}, yss::Array...;
                   R::Type=eltype(f))
    [@assert size(ys) == size(xs) for ys in yss]
    rs = similar(xs, R)
    @simd for i in 1:length(xs)
        @inbounds rs[i] = f(xs[i], map(ys->ys[i], yss)...)
    end
    rs::Array{R,D}
end

On Fri, Dec 19, 2014 at 12:46 PM, Steven G. Johnson < notifications@github.com> wrote:

Call for help: Jeff is looking for nice examples that show off multiple dispatch (and maybe staged functions), things which would be much harder/slower in languages without those features.

Reply to this email directly or view it on GitHub https://github.com/JuliaLang/julia/issues/8839#issuecomment-67671331.

Erik Schnetter schnetter@cct.lsu.edu http://www.perimeterinstitute.ca/personal/eschnetter/

[eschnett]

Sorry, ignore the call to eltype(f) in the function signature in my code, that's non-standard.

-erik

On Fri, Dec 19, 2014 at 3:08 PM, Erik Schnetter schnetter@cct.lsu.edu wrote:

I have an implementation of fmap. This traverses several arrays and applies a function to each set of elements. This implementation is actually very slow since the number of arrays can be arbitrary. To make this useful, I have manually created specialization of this for various numbers of arguments yss. I've always wanted to write a staged function for this, but haven't done so yet.

The staged function would need to evaluate in particular the map call that produces the arguments to f.

-erik

function fmap{T,D}(f::Function, xs::Array{T,D}, yss::Array...;
                   R::Type=eltype(f))
    [@assert size(ys) == size(xs) for ys in yss]
    rs = similar(xs, R)
    @simd for i in 1:length(xs)
        @inbounds rs[i] = f(xs[i], map(ys->ys[i], yss)...)
    end
    rs::Array{R,D}
end

On Fri, Dec 19, 2014 at 12:46 PM, Steven G. Johnson < notifications@github.com> wrote:

Call for help: Jeff is looking for nice examples that show off multiple dispatch (and maybe staged functions), things which would be much harder/slower in languages without those features.

Reply to this email directly or view it on GitHub https://github.com/JuliaLang/julia/issues/8839#issuecomment-67671331.

Erik Schnetter schnetter@cct.lsu.edu http://www.perimeterinstitute.ca/personal/eschnetter/

Erik Schnetter schnetter@cct.lsu.edu http://www.perimeterinstitute.ca/personal/eschnetter/

[stevengj]

@timholy, since Matlab and NumPy have subarrays/slices too, why can we argue that multiple dispatch is essential here?

[ivarne]

Ease of implementation? As far as I can tell you can simulate multiple dispatch in any language, so it isn't essential for anything.

Maybe not good to sugest something we haven't yet decided we want. In https://github.com/JuliaLang/julia/issues/9297, there is a proposal that allows us to have both efficient UTF-8 buffer positions and convenient indexing without holes so that you can do convenient arithmetic when you want that. Regex and search will return the wrapped internal index, but s[2] can give the second character regardless of how many bytes that were used to encode the first.

[timholy]

Can they make efficient subarrays/slices of AbstractArrays, or does their implementation only work for contiguous blocks of memory? This isn't a terribly hard problem to solve if you can assume your parent array has contiguous memory; it gets more interesting when you don't make that assumption.

[stevengj]

Yes, that's the key feature we are looking for: not just that something can be done nicely with multiple dispatch and/or staged functions, but that the lack of these features in other languages made implementation of the feature much harder (ideally, so much harder that no one has even attempted it).

@timholy, a NumPy array is characterized by a fixed stride for each dimension, not necessarily contiguity (essentially equivalent to our DenseArray). This property is preserved under slicing, so slices themselves can be sliced etcetera.

[vtjnash]

aberrant had some good questions on IRC along these lines. I've tried to extract just the relevant bits of comments (from among the unrelated conversations and notifications) below:

2014-12-10 (EST)
11:41 aberrant: “Organizing methods into function objects rather than having named bags of methods “inside” each object ends up being a highly beneficial aspect of the language design.”
11:41 aberrant: why?
12:20 Travisty: aberrant: I can’t speak for them, but I imagine that the argument is that it’s a nice separation of concerns. I have data (which I will represent with types) and routines for operating on that data (which I will represent as functions), and rather than having some routines belong to specific types, they are kept separate
12:21 aberrant: Travisty: I sort of understand the argument, but I’m not sure I agree with it :)
12:22 Travisty: aberrant: Yeah, sure. This is the sort of thing that may be hard to argue about from first principles, and it may be useful to look at examples. I think one place where this design simplified things was in impementing the standard mathematical functions on all of the numeric types, and dealing with conversion
12:22 Travisty: I’m not sure, but I think the solution in julia is quite elegant because of this design and it would be a bit trickier to do it in the traditional OOP setting
12:23 aberrant: Travisty: perhaps. I need to think about it some more. I really like pure OO, and this is a bit of a change that I need to wrap my head around.
...
12:54 vtjnash: julia has a convention that a method name will end in a ! to signify that the method will mutate one of it's arguments
12:56 aberrant: that’s one thing I sorely miss in python. foo.sort() vs foo.sorted() always confused me.
12:57 vtjnash: except that in python, isn't it sort(foo) vs sorted(foo) ?
12:57 aberrant: it might be :)
12:58 aberrant: no
12:58 aberrant: it’s foo.sort vs sorted(foo)
12:58 vtjnash: ah
12:58 aberrant: foo.sort modifies foo in place.
12:58 aberrant: see?
12:58 aberrant: that’s what I mean.
12:58 vtjnash: well, anyways, that's an unintentional example of why . oriented programming is a pain
12:58 aberrant: sort(foo) vs sort!(foo) makes much more sense.
12:59 vtjnash: python made a reasonable choice there
12:59 vtjnash: and tries to help you remember
12:59 vtjnash: but it still was forced to make some decision

2014-12-14 (EST)
15:25 aberrant: there’s no way to do type constants, I guess?
15:25 aberrant: http://dpaste.com/18AEHBG
15:25 aberrant: like that
15:27 vtjnash: no. that declares a local variable inside the type (can be seen by the constructors and other methods in there)
15:27 vtjnash: instead, define `y(::Foo) = 6`
15:28 aberrant: is that mutable?
15:29 aberrant: hm, yeah, that’s not what I want though.
15:29 aberrant: but I guess I can use it.
15:30 vtjnash: not what you want, or not what other languages do?
15:31 vtjnash: multiple dispatch in julia allows you to collapse 4 or 5 or more different constructs needed by other OO languages into one abstraction
15:33 aberrant: oh, I see how it works.
15:33 aberrant: well, it’s a function and therefore more overhead
15:33 aberrant: basically, I want to set the “bitwidth” of an IPv4 address to be 32, and the “bitwidth” of an IPv6 address to be 128, and then be able to write a function that takes ::IPAddr and uses the appropriate bitwidth.
15:34 aberrant: I can do this with the function, but it seems like overhead to have a function return a constant.
15:35 aberrant: e.g., http://dpaste.com/3RXRCAG
15:36 vtjnash: don't assume that a function has more overhead
15:36 vtjnash: in this case, it would actually have less overhead
15:54 aberrant: wow, ok
15:54 aberrant: I don’t see how, but I’ll take your word for it :)
15:59 vtjnash: inlining
...
18:04 aberrant: there’s no way to associate a constant inside a type?
18:04 aberrant: it would make my life a whole lot easier.
18:04 mlubin: aberrant: t.constant or constant(t) is just a syntax difference
18:04 mlubin: julia uses the latter
18:05 aberrant: mlubin: the issue is that you have to instantiate the type first.
18:05 aberrant: mlubin: I need the equivalent of a class property.
18:05 mlubin: constant(::Type{Int32}) = 10
18:05 mlubin: constant(Int32)
18:05 aberrant: oh. wow.
18:06 Travisty: The only member of Type{T} is T, which is why this works
18:06 mlubin: mind=blown? ;)
18:06 aberrant: yeah
18:06 aberrant: that’s jacked up
18:07 aberrant: there’s NO WAY I would have ever thought of that on my own :(
18:07 mlubin: once you see it for the first time it becomes a lot more intuitive
18:07 aberrant: ipaddrwidth(::Type{IPv4}) = uint8(32)
18:08 aberrant: w00t
18:10 aberrant: can I do a const in front of that?
18:11 Travisty: I don’t think so, but usually when julia generates code, it should be inlined so that just the constant uint8(32) appears, instead of a function call

if you're looking for more examples: https://github.com/JuliaLang/julia/pull/7291

Net code reduction is nice from 6K C++ to about 1K Julia. Initial performance benchmarks show it's just under 2x the native C++.

how many languages do you know that can claim implementation of printf, from the operators (+-*/, Int, Float) to the numerical output formatting, in the language itself, and benchmark within a small margin of the libc versions? Maybe C++? C can't even claim this (it doesn't have operator overloading).

Python/MATLAB/etc. might be able to claim bits of this, but have you ever tried to use the string operators in MATLAB?

[tknopp]

This is a very interesting discussion and I like to add some points:

• What seems to be not really answered is how Julia compares to polymorphism in C++ speedwise. In C++ one can have static polymorphism through templates or dynamic polymorphism through OO. But the later requires a vtable and is always a little slower, especially when calling such things in a tight loop. For this reason one would never overload the array access operator since this would not be fast. Eigen for instance uses static polymorphism.
• In Julia there is no destination between a virtual and a non-virtual function and we still don't have performance penalties. My understanding is that the hugh advantage is that Julia can generate code "on the fly" and thus still generate fast code where in C++ runtime means one cannot modify the generated code anymore. Additionally in Julia we are able to do inlining at runtime. So no function call overhead. And in Julia we don't have to think about all this because it is done automatically
• One huge thing (again) is: generics without overhead and without all that pain C++ has (looking compile times) There are different methods how to implement generics and IMHO Julia uses the same model as C#, so this might not be entirely new (for Jeffs thesis). But still I think this could be something to discuss in the thesis.

[StefanKarpinski]

How hard it is to explain how effective multiple dispatch is, is itself quite interesting.

[timholy]

Glad we're on the same page now. Yes, our new SubArrays don't rely on having strides, nor do they require linear indexing (although they can use it, and do by default if it happens to be efficient). So they work efficiently for any AbstractArray, and support "non-strided" (Vector{Int} indexes) views as well.

To clarify, I presented our new SubArrays primarily as an example of staged functions, not multiple dispatch. That said, the current scheme would fall apart at the construction step without multiple dispatch: we call completely different constructor methods for slice(::AbstractArray, ::UnitRange{Int}, ::Int, ::Vector{Int}) and slice(::AbstractArray, ::Int, ::UnitRange{Int}, ::UnitRange{Int}), etc. Those constructors are generated by staged functions, but we need multiple dispatch for the system to work.

Without separate functions, a varargs constructor has to face the fact that looping over the entries in the indexes tuple is not type-stable: to process individual indexes in a loop, you have to assign them to a named scalar variable, and that variable is guaranteed not to be well-typed even though the input tuple can be. Of course, the staged function does oh-so-much more: almost all the "hard" decisions---are we dropping this as a sliced dimension, or keeping it, and what are the consequences for the internal representation of the object?---can be made at compile time, reducing the runtime constructor to a remarkably trivial operation. It just so happens that the particular trivial operation is different for each combination of input types, so to exploit this triviality and make construction fast you need to be calling different methods customized for each combination of input types. This is why you need both staged functions and multiple dispatch.

Personally, I doubt that it would be possible to do what we've done without the combination of staged functions and multiple dispatch, and that it's likely that julia now has the most flexible and efficient array views of any language. But I don't pretend to have undertaken any kind of actual study of this problem.

[timholy]

This makes me wonder if I should write up our SubArrays for publication---it does seem we have something new here.

[tknopp]

@StefanKarpinski: Lets be self-critical: Maybe we just have not written it up yet?

[tknopp]

@timholy I actually see two things here, the staged functions (which I have to admit not fully understood) and how Julia fits into C++ virtual functions vs templates universe (C# and Java are similar here).

For the C++ comparison it would be very interesting to do timings and experimentally proof that the we can reach compile time polymorphism (which I am sure we do)

Another thing I wanted to formulate for some time is how Julia compares to C++ with regards to C++ Concepts (lite). The subtyping we have is a large part of C++ Concepts.

[timholy]

Does one even need to write a bunch of C++ code and run timings to check this?

A = rand(3,5)
@code_llvm slice(A, 2, :)
@code_llvm slice(A, :, 2)

There's essentially nothing but load and store operations (i.e., memory access) in there.

[ViralBShah]

@timholy I do think it is worth writing this up. I believe (but have not verified) that many systems that implement Nd arrays have specialized implementations for the first few dimensions, and then generally fall back to something slow otherwise, and that shouldn't be the case for the new implementation.

Another related case is writing an efficient sparse N-d store (arrays being the common special case).

[tknopp]

Tim,

yes indeed looking at the resulting assembly code is a great tool. But still I think this has to be compared in a larger context, where (in the C++ implementation) virtual functions would have to be taken into account.

My guess is really that it is not multiple dispatch that is so fast but that we can inline at runtime and thus transform virtual function (which our generic functions are) into efficient instruction without function call indirection.

If this is true (please correct me if I am wrong @JeffBezanson @vtjnash @stevengj @timholy ) the multiple in multiple dispatch is not the reason of Julia beeing so fast but a neat side affect that certain code can be formulated nicer (where single dispatch is limiting)

[timholy]

I'm probably just not understanding, but I'm not sure of the distinction you're making. In Julia, "inline" and "runtime" don't really seem to go together; inlining is performed during type inference, i.e., at compile-time. Multiple dispatch is used to select the appropriate method for the inferred types.

[tknopp]

The confusion here is that "compile time" and "runtime" of C++ cannot be compared with that of Julia. Codegen can happen during "runtime" so yes I think when I do include("myscript.jl") inlining is performed. And even if "runtime" is not the right word for that from a C++ perspective it is "runtime".

And the dispatching on different types is like a vtable but more general, no?

[johnmyleswhite]

This makes me wonder if I should write up our SubArrays for publication---it does seem we have something new here.

It's a little far from standard topics, but you might consider submitting to JSS. We need more Julia articles.

[ViralBShah]

I would love it if Tim writes a blog post to describe this work to start with, as it will get a lot of people up to speed. JSS is a great idea, and perhaps some of the foundational work that had been done in data frames and distributions is worth writing up too? I certainly would enjoy reading it.

[timholy]

Well, what made me think of it is that much of it has been written up: http://docs.julialang.org/en/latest/devdocs/subarrays/. For a publication you'd want to go into a lot more detail, but this hits a fair number of the main big-picture points.

[ViralBShah]

To the question @stevengj raised about multiple dispatch, I would say that without multiple dispatch, it is pretty difficult to write our base library in julia. I am saying what is already known by everyone here, and wonder if this is not a compelling example for the reasons brought up here.

Details such as operations on numbers and the conversion/promotion are intricately tied to multiple dispatch. Since multiple dispatch is what essentially exposes type inference in the compiler to the way types are used in code, we are able to write a generic and fast numerical base library. To quote a statement you made - it helps separate the policy making out of the compiler and into libraries. For example, @JeffBezanson showed me once how the Scheme spec spends 1/3 of its space on numerical details.

Many interpreted systems often end up having general types and inspecting the types of their objects at runtime to make decisions on what code to execute. They often then have a separate implementation in C/C++/Fortran in the base library for each type, leading to a large and difficult to debug codebase. Often these are generated through an external macro system, but increasingly the use of C++ templates have avoided this specific issue. The issue of two languages and type inference still remains in these cases.

At some level, vectorization is how many scientific languages amortize the cost of determining types at runtime and doing the appropriate code selection. In Julia, with the combination of type inference, multiple dispatch and generic programming, our costs for both are significantly lower allowing us to write generic devectorized code - C without types.