Multithreading

[robertdj]

@tknopp: I have finally gotten around to parallelizing parts of NFFT. With the new Threads.@threads the convolve! function uses multithreading -- check my branch.

Compare a profiling of the old code (The purple part on the same level as the highlighted convolve! is fft.)

With the new code (with export JULIA_NUM_THREADS=2):

The (regular) FFT call now consumes a significant proportion of the time.

Recalling stevengj's comments in #17 FFTW 3.3.5 is needed to multithread fft and this version is used in Julia v0.6. But I can't find info on how to use multithreading for fft.

[tknopp]

nice, so funny seeing this stuff landing in this package. I worked on the original multithreading code in Julia JuliaLang/julia#6741

Regarding FFTW: You just have to call FFTW.set_num_threads.

[robertdj]

I didn't realize that multithreading in Julia was so old!

Would it make sense to add

FFTW.set_num_threads(Threads.nthread())

in NFFT?

[tknopp]

Yes this sounds good. Although this probably needs a more general solution since its always a little bit bas if we change global state in a user library. Maybe

if Threads.nthread() > 1
  FFTW.set_num_threads(Threads.nthread())
end

is a little safer since we would not set it to 1 if a user runs Julia in serial mode but has set FFTW to be parallel

[robertdj]

Good point with the global state.

There's not much documentation for FFTW.set_num_threads, but do you know if it can be set only for the NFFT module?

[tknopp]

No thats not possible. Its a global value in libfftw. There seems to be currently even no getter so that we cannot know to what is had been set.

But I would not care to much. More interesting is how good your speed ups are. Parallel nfft is quite nice and that often seen in the wild... ;-)

[robertdj]

I'm currently on an old computer with only 2 threads and here it's about 30% faster for a 2D example I'm testing.

[robertdj]

More timing suggests that the threading constructs take quite some time. E.g. for 2D case in nfft_performance:

• 4 threads: 2.5 times faster
• 8 threads: 3.6 times faster

[tknopp]

Hey, that not that bad. In multithreading its not that simple to achieve full speedups. So this is not necessary an issue with the Julia implementation of threads.

The issue in multithreading is usually memory accesses to similar regions from different threads. This will lead fighting CPU caches.

Do we already sort the k-space nodes? Might be a good idea...

[robertdj]

Are the k-nodes part of the NFFTPlan?

[tknopp]

I meant x

https://github.com/tknopp/NFFT.jl/blob/master/src/NFFT.jl#L56

I usually use NFFT for MRI reconstruction where x is k (Fourier space is named k-space)

[robertdj]

We don't currently sort x. I haven't thought about this before, but indexing the input and window arrays would be more cache efficient with a suitable sorting of x.
Let me check how sorting affects the timing.

[robertdj]

One problem with sorting in the plan creation is that fHat must be sorted similarly...

[tknopp]

Yes indeed. One might for with an indirect index. Might make sense to look how the solved this in the NFFT C library:

https://github.com/NFFT/nfft/blob/9fbc5ab42fe6c690d2b82b971b3566be4e6df165/kernel/nfft/nfft.c#L111

[tknopp]

I have merged your multithreaded convolve into master