Julia: multi-threaded (@threads) dotcall/broadcast?

[stevengj]

It would be nice to be able to put @threads in front of a dot call, e.g. @threads X .= f.(Y), and have it call a multi-threaded version of broadcast (which would assume f is thread-safe).

[vchuravy]

I don't think this is something that should go into the @threads macro (since that is already too complex). @SimonDanisch has been working on doing this based on dispatch to a special array type in https://github.com/JuliaGPU/GPUArrays.jl/blob/master/src/backends/julia/julia.jl

The beauty of Simon's approach is that it is general and works for distributed arrays, GPU arrays, and native arrays. I consed that it might feel a bit unnatural for native arrays, and I wouldn't be opposed to add a macro that transforms X .= f.(Y) into the threaded dispatch version but I don't think that necessarily needs to be in Base.

[stevengj]

Regardless of the name of the macro, it would be nice to have something that just involved a decorator and didn't involve re-allocating or wrapping all of your arrays in some other type. This is the big distinction between threads and GPUs or distributed-memory — with threads, you don't need to decide in advance to put your data on a GPU or in another process.

[SimonDanisch]

I'd suggest having a macro similar to fastmath, that replaces the broadcast calls with e.g. threaded_broadcast. Would be nice to share the implementation between what I have and Base :)
I tried to create a minimal working version, but it seems like expand needs some extra attention to get the broadcast back.
Blocks of code expand into a Expr(:thunk, Toplevel LambdaInfo thunk), so I had to iterate through the sub expressions to make things easier.

[StefanKarpinski]

Would it be too crazy for broadcast to be implicitly parallel once threading is stable?

[martinholters]

Already, broadcast makes no guaranties concerning the order or exact number of times the given function is called (think sparse matrices), so passing non-pure functions is questionable anyway. That would be in favor of implicit parallelism. Con: When broadcasting over a collection with only few items, the complexity of the given function will decide whether it is worth any threading overhead, which cannot be easily decided by broadcast itself.

[ChrisRackauckas]

Could this get a 1.0 milestone? At least some kind of macro would be very useful, if not implicit parallelism (or a mixture: default implicit parallelism, which can be overridden with a macro). Since broadcasting is such a cool feature, this would really complete the story.

[Sacha0]

cc @lkuper

[vchuravy]

If anyone wants to tackle this (developing this outside of base at first is probably a good idea) my roadmap/ideas would be.

1. Thin wrapper type around Arrays
2. Macro that automates this wrapping
3. Co-existence with other types of parallelism. We currently have DistributedArrays as a seconder user of broadcast
   The nice thing about using a thin wrapper is that you can have DArray{Threaded{Array}} to combine distributed and threaded parallelism.
4. What happens in the case of mixed accelerated array types? e.g. A GPUArray meets a threaded array?
5. NUMA-aware parallelism. Communication across NUMA nodes is still quite expensive and memory should probably be pinned to the NUMA group that the threads belong to.

For other inspiration take a look at parallel collections in Scala http://docs.scala-lang.org/overviews/parallel-collections/overview.html (which were the inspiration for Java 8).

[SimonDanisch]

Just as an info: I'm already getting good speed ups out of @threads in my broadcast implementation in GPUArrays ;) It's still a bit crappy, but if I have some time I could polish it up and I can make a first PR for a tbroadcast or however we call it ;)

[RoyiAvital]

@ChrisRackauckas ,

I like the idea of letting the user choose.
Should be something like a Macro (Or the most efficient way to do it) with 3 states:

1. Auto (Default) - Adding heuristic to set On or OFF.
2. ON.
3. OFF.

Thank You.

[ChrisRackauckas]

It would be interesting if the heuristic could be applied to arbitrary loops via some macro as well.

[stevengj]

@ChrisRackauckas, parallelizing arbitrary loops is precisely what the @threads macro does, no? I'm simply proposing using it for broadcast (and dot calls) as well.

@vchuravy, I'm not sure I like the idea of a special array type, vs. just a @threads decorator on looping constructs. Not needing a special container type is one of the main attractions of shared-memory parallelism.

[ChrisRackauckas]

@ChrisRackauckas, parallelizing arbitrary loops is precisely what the @threads macro does, no? I'm simply proposing using it for broadcast (and dot calls) as well.

I was asking if there could be a way to apply whatever implicit parallelism heuristic to a loop. Essentially a macro for "multithread this if the size of the array is greater than x" or whatever is involved in the heuristic, and have the options be tweakable. Then broadcast would just be essentially applying that with the defaults. Your proposal just has a macro, but I'm wondering if implicit parallelism can be added as well.

Otherwise I could see applications wanting to have a bunch of conditionals to check if multithreading should be ran? That last part is dependent on the overhead of multithreading (which I found to be measurable in many small problems, but the benchmarks may be mixed up with https://github.com/JuliaLang/julia/issues/15276 issues).

[stevengj]

@ChrisRackauckas, it's not the size of the array, but the expense of the loop iterations that matters. There's also the issue of load balancing if the loop iterations have unequal cost. I agree that you want to automate this (both deciding how many threads to use and how to load balance) to the extent possible. My understanding is that Cilk (and subsequently OpenMP) mostly solved this issue.

Anyway, I see that as orthogonal to this issue.@threads f.(x) should use the same machinery as @threads for. Improvements in the latter should help the former.

[Sacha0]

Ref. https://github.com/JuliaLang/julia/pull/18278#issuecomment-279021388

[pabloferz]

Related #1802 (or is this a duplicate?)

[RoyiAvital]

@stevengj , There is a lot to consider whether or not to multi thread a loop.
Hence I just hope Julia will also allow the user to set it OFF or ON on his choice.

[jebej]

It would be nice to initially have a simple-case version of the macro as explained in the first post (for cases like cos.(x)) for all the people on the mailing list complaining that it is slower than MATLAB. This would make it really easy to answer. Right now people need to make a special function to see the difference:

```julia
function threadedcos(x::AbstractArray)
out = similar(x)
Threads.@threads for i in eachindex(x)
out[i]=cos(x[i])
end
return out
end

[vchuravy]

using GPUArrays you can automatically accelerate broadcast.

using GPUArrays
GPUArrays.init(:julia) # Otherwise an OpenCL or CUDA backend might be used
a = GPUArray(rand(Float32, 32, 32))
cos.(a)

[bramtayl]

I wonder if it wouldn't be possible to have some sort of option to switch out which function is used for . broadcasting. Something like ENV[:dotfunction] = threaded_broadcast

[Sacha0]

I wonder if it wouldn't be possible to have some sort of option to switch out which function is used for . broadcasting. Something like ENV[:dotfunction] = threaded_broadcast

Ref. the discussion around https://github.com/JuliaLang/julia/issues/16285#issuecomment-237320513. Best!

[yuyichao]

Do note that this should not be the default or a global setting, not before we require every functions to be threadsafe at least. No guarantee on execution order is a very weak requirement compare to thread safe.

[MasonProtter]

It should probably be noted here for those looking for this feature that such a macro has been implemented intoStrided.jl.