SimonDanisch/FixedSizeArrays.jl

StaticArrays

Opened this issue ยท 39 comments

I think there should be only one package for fixed size arrays ;)
Of course I'm a bit biased towards FixedSizeArrays, since I'm more familiar with the code :D
@andyferris @c42f, could you explain to me what would speak for StaticArrays?
What were the initial motivations behind it?
Right now the biggest things seem to be inheritance from DenseArray and already up to dateness with 0.5.
These could be brought to FixedSizeArrays with a rewrite.

I'm trying to fix up FixedArrayBench with support for StaticArrays to further judge this.
Last time, something wasn't working with 0.5 and BenchmarkTools.jl.

Best,
Simon

HI Simon,

Thanks for inviting me to comment. Let me start by saying I've been using FixedSizeArrays for a while, and it has been quite suitable and fast and I really have appreciated all the effort you put into it!

So as a simple list, the complaints I had and things I was interested in overcoming are:

  • Non-compatibility with AbstractArray. This made any code that I want to use with FixedSizeArrays have to include signatures for Vec, FixedVector, etc. They were "special cased". If FixedSizeArrays inherits from AbstractArray, then people could write generic packages and users may switch back and forth between StaticArrays and FixedSizeArrays without a problem. Some new user who never used either but makes there own type which is effectively a FixedVectorNoTuple (I wonder how many implementations of Point there are in the wild?) could also use the package.
  • Non-compatibility with the "generic" algorithms in Base that work on AbstractArray. Mostly this is because there are no mutable fallbacks for similar.
  • Non-compatibility with the "external" algorithms in Base that work on DenseArray. These want a pointer to a dense structure to pass to BLAS or LAPACK.
  • An incomplete list of built-in types. StaticArrays offers arbitrary dimensional SArray and mutable MVector, MMatrix and MArray. There is very little code to implement each - the generic implementation looks after the rest. Soon I plan to implement pointer-based versions of these as well, which will be useful in certain circumstances (interacting with C/C++/etc and also for a hack I'm working on using llvmcall to get a mutable container on the stack).
  • I had trouble understanding the code (or codegen, I think). This is more my fault than yours, but like you say we tend to be motivated to our own code because we understand it. I think if you look at StaticArrays then each function should be more-or-less entirely self contained (provided you understand similar_type()), though matrix multiplication is super messy. The nested tuples in FixedSizeArrays didn't help me here.
  • General bugginess and missing methods. I don't blame you - there are many, many corner cases to cover (trust me, I know!) but the above made me reluctant to fix the problems I found and submit a PR. I think bugginess and completeness has improved somewhat this year, which is great.

To conclude, perhaps I should explain that when I started StaticArrays I wanted to understand how to make FixedSizeArrays better. The idea was to prototype a design and discuss porting it over to FixedSizeArrays with you. The funny "problem" however was that at some point I had achieved the first 4 things above, and it became more-or-less a complete implementation. I had achieved the "rewrite" you mentioned in your opening comment, but it seemed a bit too rude to say "Hey Simon, could you please change all your code to this, please?" (especially given some functionality (expm, etc) exists in FixedSizeArrays but not StaticArrays, StaticArrays has a different API for users, and only runs on Julia 0.5).

Given how things are now, I'm not sure what should be done? To me the major, long-term "organizational" goal should be to move these out of our personal accounts and have some JuliaX github page (like JuliaStats, JuliaBio, etc, but where X โ‰ˆ FastStaticallySizedArrayPackages). It would also be most efficient if we could start collaborating together.

Anyway, thanks again for starting this discussion and do let me know your thoughts :)

I'm trying to fix up FixedArrayBench with support for StaticArrays to further judge this.
Last time, something wasn't working with 0.5 and BenchmarkTools.jl.

Good benchmarks will be very useful no matter what happens in the future, so this is great :)

c42f commented

Hi Simon, thanks for starting the conversation, we definitely need to sort this out as fragmentation is happening - we now have three(!) different fixed size packages (counting ImmutableArrays which was about the best that could be implemented in julia-0.3, and which several packages still depend on.)

A bit of back story - as of a few months ago, I was slowly fixing the various problems with FixedSizeArrays, eg trying to fix the type promotion problems with map, etc. In the process I got frustrated with several design problems but I was just planning to plug along fixing things as they came up since I take fragmentation and API compatibility fairly seriously. I guess I complained about the problems to Andy a bit, and he surprised me by writing an entire new package in a fit of productivity. With 0.5 imminent I've been persuaded to start using StaticArrays as a dependency in some places, simply because it fixes various API problems which I haven't had time to address in FixedSizeArrays, and also provides a cleaner implementation strategy for the future. So there was no fragmentation planned, I'm bummed that we seem to have fallen into it, and I'm not yet sure how to dig ourselves out! A common supertype would go a long way, but the packages as they stand are incompatible with doing that easily.

To summarize some frustrations I've had with FixedSizeArrays more specifically -

  • Obviously, the issue with subtyping AbstractArray. I started fixing this a while back, but I just ran out of steam :-( It seems doable though (Andy has done it!) so this could be pushed through if someone had time.
  • The API has various oddities which aren't consistent with AbstractArray in a duck-typing sense. A good example is the concatenating constructor function (::Type{FSA}){FSA <: FixedArray, T1 <: FixedArray}(a::T1, b...) - this makes sense in a glsl world, but is just different from Base.Array where you'd use vcat or something. Another example is the string parsing vector constructor; this is just weird to me. The nullary version of map() also seems strange; I preserved it through my map rewrite to avoid breaking people's code, but it's just an oddity. To fix these issues I think some deprecations are required to get a tight focus on compatibility with AbstractArray semantics.
  • The tuple-of-tuples representation is really a pain to work with internally. Unfortunately it seems completely necessary for good performance in 0.4 so it's nice to be able to make a clean break and ignore this in StaticArrays. Dunno what we could do here other than make a big break at some stage in the FixedArray implementation.
  • Perhaps a small point, but I really don't like FixedVectorNoTuple as a name. Andy has FieldVector which is better but perhaps there's something better still. Again, fixable of course at the cost of some API breakage.
  • Somewhat confusing code structure: the way that methods are separated into files, I still find it hard to find things. Totally fixable of course, it's just work for someone to go through and group functions more clearly.
  • Lack of clear documentation about where exactly the package is going and how to use the central features. (I'm super guilty of this in some of my own packages, but it is a thing.)

In general, the code has/had a fair bit of one-off patching to make things work with a small set of sizes or given types, without working in the general case. I think this is ok for graphics work because of the small number of shapes and types you typically use, but these things bite you all the time if you're doing general computation with FixedSizeArrays. Especially with funny types like Duals etc, and larger matrix sizes. I already fixed a lot of the type promotion problems, but I still find missing methods now and then.

Btw, some possibly useful contributor statistics on current number of lines:

$ git ls-tree -r -z --name-only HEAD -- *.jl | xargs -0 -n1 git blame --line-porcelain HEAD |grep  "^author " | sed -e 's/\(Sdanisch\|SimonDanisch\)/Simon Danisch/' | sort | uniq -c | sort -nr
    678 author Simon Danisch
    497 author Chris Foster
    302 author mschauer
     65 author Tim Holy
     45 author Jan Weidner
     34 author Andy Ferris
     21 author Steve Kelly
      4 author ScottPJones
      4 author Jens Schleede
      3 author Andrew Smith
      2 author David Sanders

@mschauer any thoughts? You're the next biggest contributor.

Thanks a lot, seems like reasonable concerns!
I see the way forward to actually eliminate most of the code, inherit from DenseArray and just have broadcast in Base do most of the work ;)
In that process, I wouldn't mind dropping 0.4 support at all.

Let me get familiar with StaticArrays codebase and procure the benchmarks! To get familiar, I could see how easy it'd be to integrate VecElement for better SIMD performance.

First thing striking me when looking over the code is, that there seem to be a few unnecessary staged functions ;)
Is this just an artifact? E.g. https://github.com/andyferris/StaticArrays.jl/blob/master/src/matrix_multiply.jl#L11

By the way, VecElement, broadcast, inheriting from DenseArray and proper integration with StructsOfArrays.jl would be the list of most important features for me!
@timholy, we could also warm up the discussion of making the color types a static array.
Btw. this is pretty slick:

StructsOfArray(RGB, red, green, blue) # I think this actually already works...

Benchmark results (SVector vs FSA filtered for regressions):

7-element BenchmarkTools.BenchmarkGroup:
  ("vecvec",mean,3,Int64) => TrialJudgement(-32.69% => improvement)
  ("vecvec",mean,3,Float64) => TrialJudgement(-16.07% => improvement)
  ("vecvec",-,10,Int32) => TrialJudgement(-8.97% => improvement)
  ("vecvec",mean,3,Int32) => TrialJudgement(-23.08% => improvement)
  ("vecvec",-,3,Float64) => TrialJudgement(-9.44% => improvement)
  ("vecvec",sum,3,Int32) => TrialJudgement(+17.65% => regression)
  ("vecvec",sum,3,Float32) => TrialJudgement(+10.87% => regression)

More improvements than regressions... That's not too bad :)

c42f commented

the major, long-term "organizational" goal should be to move these out of our personal accounts and have some JuliaX github page

On that note, I think a future fixed size array type should be part of the stdlib, whenever that eventuates: it's a pretty core abstraction for geometric computation, among other things, and the correct choice of concrete fixed array container type should eventually be entirely obvious for users.

I think the inheritance from AbstractArray and mutable fallbacks are killer features. To be honest, in my new code I've switched entirely to StaticArrays, and I've not regretted it yet.

One could base RGB on SVector, but some platforms (OSX) I/O libraries return a 4-element color (RGB1 or RGB4) with one padding byte. I understand why they do that, but it also presents a challenge to generic julia code, where you want to pretend that element doesn't exist. However, that's pretty incompatible with a viewpoint based on static vectors. There's also the confusion about whether you mean machine ordering or "meaning" when it comes to a type like BGR. Bottom line is that I don't think Julia's types are expressive enough yet to migrate all of ColorTypes to some static array type. I think currently the closest thing we have to a viable strategy there is ImagesCore.jl (specifically the channelview and colorview containers) and the treatment of colors as if they are vectors in ColorVectorSpace.

I'm also mulling over adding SOVector and friends for static offset arrays (for which indexing does not necessarily start at 1).

I'd like to reply to everything but sorry for this mega-post!

@SimonDanisch said:

I see the way forward to actually eliminate most of the code, inherit from DenseArray and just have broadcast in Base do most of the work ;)
In that process, I wouldn't mind dropping 0.4 support at all.

Using broadcast is something I tried to design for in the package. However broadcast is Base is a lot slower than having an unrolled loop (and having the various dimensionalities figured out by an @generated function), but StaticArrays should use my specializatins of map and broadcast whenever appropriate, and when .+ etc are switched to broadcast we are ready.

@SimonDanisch Your description of the way forward sounds somewhat similar to dropping the FixedSizeArrays codebase and continuing from the StaticArrays code - is that right? The reason I would advocate for this is (as @c42f mentions) the shear number of API changes is overwhelming. It would take several rounds of deprecations in FixedSizeArrays, and IMHO it might be conceptually easier for users to transition to a new package. If the API of StaticArrays is close enough to AbstractArray (plus size on types) then there shouldn't be too much future API churn, so hopefully we wouldn't want a new package for v0.6 (even if we need a complete internal rewrite).

To get familiar, I could see how easy it'd be to integrate VecElement for better SIMD performance.

That would be really cool! My limited knowledge is this: LLVM already SIMDs the immutable containers (but does an even better job with -O3 switched on) but doesn't do so for the mutables (I guess this is stack vs. heap analysis). From what I read VecElement does nothing unless used in conjunction with llvmcall(), and LLVM only handles vectors of certain sizes and types, so we can't stick VecElement on things unless we know they are SIMD-friendly (the latter problem will improve in LLVM upstream over time). So short story: there might not be too much low-hanging fruit left, and e.g. 3D geometry won't get faster than just using SVector with -O3.

First thing striking me when looking over the code is, that there seem to be a few unnecessary staged functions ;)
Is this just an artifact? E.g. https://github.com/andyferris/StaticArrays.jl/blob/master/src/matrix_multiply.jl#L11

Yes, it is an artifact - an artifact of Julia! Inference is not allowed to "approximate" (or "widen") the input types to a function since the generator could perform arbitrary logic depending on these types. This can also happen for return types, hence everything is also inlined. Standard functions do not always handle static arrays with large numbers of elements smoothly (I honestly can't predict when it works and when it doesn't!) since it can result in NTuples with unknown N. On this note, I've made a (surprisingly small!) list of compiler features that would completely eliminate all need for @generated functions, including

  • inference support of arbitrary length NTuple (on Jeff's radar)
  • an @unroll macro that inserts an :unroll meta or similar so that loops can be unrolled (from constants like N - not literals, which a macro could do)
  • inlining during inference (until convergence, also on Jeff's radar, though I think this might be awaiting a final decision)
  • some improved support for constant propagation and associated branch elimination (I think improving these is already a goal).

If these were done, then I would think we're in a good place to bring static arrays to the standard library, as @c42f says. Incidentally, this would also eliminate all the generated functions in TypedTables (assuming@unroll works over tuples of mixed type, e.g. occurs during inference). @timholy what's your thoughts on how feasible this "wish list" is? And how would you see a statically-sized array fitting into Base's array framework in the future?

@timholy said:

... and mutable fallbacks are killer features

The mutable static arrays are more than fallbacks - I've found them also extremely useful on there own. They are faster than Base.Array in many cases but are slower at some things (which I don't understand yet - maybe Array gets some special love in allocation or SIMD or something. It will be easier to tell once Buffer is implemented). If I would add a 5th item to my wishlist it would be stack-allocated mutables.

but some platforms (OSX) I/O libraries return a 4-element color (RGB1 or RGB4) with one padding byte. I understand why they do that, but it also presents a challenge to generic julia code, where you want to pretend that element doesn't exist.

@timholy You could kind-of solve this problem with FieldVector, an inner constructor that only takes 3 values (the 4th value could never be set), and also overloading size(::Type{RGB4}}) = 3 so you have 4 fields yet 3 components. Conversions to other 3-vectors will ignore the 4th field (and any code in StaticArrays will ignore the 4th field, or if you overload getindex, you could also ignore the 1st field or whatever).

I'm also mulling over adding SOVector and friends for static offset arrays (for which indexing does not necessarily start at 1).

LOL. Yes I realized I was being more specific than strictly necessary, but there seemed to be enough complication to begin with :)

Using broadcast is something I tried to design for in the package. However broadcast is Base is a lot slower than having an unrolled loop (and having the various dimensionalities figured out by an @generated function), but StaticArrays should use my specializatins of map and broadcast whenever appropriate, and when .+ etc are switched to broadcast we are ready.

Yeah I meant the last part of that ;)

Yes, it is an artifact - an artifact of Julia!

Well but in the cited line, the generated function is not needed at all, is it?

Well but in the cited line, the generated function is not needed at all, is it?

If anyone has a proof of that I'm all for simplifying it! (please!)

c42f commented

Well but in the cited line, the generated function is not needed at all, is it?

To clarify, I think the apparent need for this was what Andy was referring to in his big post above, in the section containing

Inference is not allowed to "approximate" (or "widen") the input types to a function since the generator could perform arbitrary logic depending on these types

That is, having @generated there is (might be?) forcing the compiler to ignore some of its heuristics. Or possibly, it's just working around a bug in inference, it's hard to say without a test case. @andyferris do you have a minimal test case showing a concrete example of the problem you were working around there?

do you have a minimal test case showing a concrete example of the problem you were working around there?

Nope, only complicated cases that scared me sufficiently that I became paranoid! :)

I did of course try to make a minimal test case at a few points, but couldn't pin down the problems (simple stuff mostly works). This is something that needs to be researched further (or we ask Jeff or someone to explain the specifics). Any help here is greatly appreciated. Simplifying the code is definitely a goal, but not at the cost of abandoning performant, large static arrays.

There's even a chain of @pure functions that was crashing Julia earlier, but I couldn't come up with a simple example of that either. One of the similar_type functions is unnecessarily a @generated function instead of a @pure function because of this.

Oh I see, so my ignorance then ;)
So that's a sad corner case of Julia then...

So that's a sad corner case of Julia then...

Unfortunately it's necessary that inference makes approximations. Inference will work it's way through e.g. recursive function calls, and if you make a mistake and your code is non-terminating, compilation would be toast. Currently you run the code and quickly get a stack-overflow error and a simplified traceback.

What I've suggested to Jeff Bezanson & co. is that we might want a :try_harder meta or something that tells the compiler not to be a wimp. :)

Getting back to the original issue, where are we here?

I'm pretty set on switching to StaticArrays. Just need some time to play around with it and do the switch!

OK, great, thanks Simon. It will be great to have you on board.

I think we should now try decide where to put StaticArrays more centrally.

Maybe it's time for JuliaArrays?
Tim Holy alone has around 10000 packages around Arrays ;)
[1, 2, 3, 4, 5, ..., 10000]
I'm sure we could fill that org quickly!

...with another 10000 planned ๐Ÿ˜„. Sure, I'd be happy to deposit a lot of stuff in some org.

done! ;)

Can I be an owner, please? I can't create a repo and transfer it...

All organizations need an icon... Any preferences/ideas?

julia_arrays_3x3

julia_arrays_4x4

julia_arrays_5x5

c42f commented

How about something a little more squareish? A banded array of julia colors maybe?

juliaarrays

juliaarrayssmall2

juliaarrayssmall

juliaarrays2

I'd be for the squares :) Sorry, didn't realize I invited you only as members! Should be fixed.

I took the liberty to chose the icon... We can of course change it! ;)

c42f commented

Thanks Simon! Ah, too bad that github decides to resize it, I'll try a 100px version of the same. Gotta get the important things right - someone else may have decided the color of the bikeshed, but we get to decide how many doors.

By the way, github makes all organization members private by default, and as far as I know it's up to members to make themselves public. So you may want to do that :-)

I moved StaticArrays over now. You should each be invited for admin access.

Squares are definitely better :)

Great! Thanks everyone :)

OK, now StaticArrays looks less lonely over at JuliaArrays.

c42f commented

Yikes, that's a whole bunch of array packages :-)

Can we move ImmutableArrays to JuliaArrays? There are a few packages that use ImmutableArrays that can use some updates as well.

Well, I didn't move FixedSizeArrays to make the statement that it's not the official fixed size array package anymore! The same should apply to ImmutableArrays. Can we just make them not depend on it anymore? Should be easy with adding a few typealiases?

c42f commented

I think it might make make sense to move ImmutableArrays and FixedSizeArrays, but just put a big deprecation notice in their documentation.

If we also release a new 0.5-only version of both these packages, and depwarn a prominent part of their interfaces with a message which recommends StaticArrays, that would really encourage migration: people using 0.5 would see the depwarn, but people on older julia versions would be able to continue using the older packages without fuss. A few typealias might make sense for migration, but I wouldn't want to encourage too much syntactical divergence.

Speaking of which, I rather liked writing Vec and Mat - can/should we recover these for use in StaticArrays? SVector and SMatrix are admittedly more systematic.

I kinda agree that the packages could be amalgamated at JuliaArrays with prominent plans for the future advertised.

Speaking of which, I rather liked writing Vec and Mat - can/should we recover these for use in StaticArrays? SVector and SMatrix are admittedly more systematic.

Make an issue there, and we can try hash it out. I guess my main issue is that I knew I would have arbitrary dimension arrays (Arr?) and multiple types of each (currently immutable and mutable, but I'm planning a pointer version). But I'm open to sensible reform :)

SVector and SMatrix are admittedly more systematic.

I forgot to mention: ArrayFire.jl has AFArray and DistributedArrays.jl has DArray so it's becoming a little idiomatic. We need to be careful not to have a package race/war to see who can claim Vec and Mat first (would we really be the package that should claim these? Well, our arrays are smaller, so the name is smaller... hmm....)

c42f commented

As far as I know, FixedSizeArrays was the first to claim Vec quite some time ago, though SIMD.jl has also overloaded the name now.

+1 for moving ImmutableArrays and FixedSizeArrays into JuliaArrays and adding prominent deprecations warnings in the README

I'll put this here - I made a relavant announcement about StaticArrays on julia-users

https://groups.google.com/forum/#!topic/julia-users/wxtAnqk3FWk