Faster initialisation on a gpu

[richfitz]

If we had jump coefficients for 2^128, 2^129, ... (up to some reasonable number 128 + m - 1) then we can create up to 2^m starting points very quickly by jumping (or not) through each of these coefficients:

• let i be the thread index on 0..k-1 for k threads
• encode i in binary with ceiling(log2(k)), e.g. 0110 for i=6, k=8
• for each bit in our series advance using a decreasing series of our jump coefficients (2^131, 2^130, 2^129, 2^128 with k =8 so that the thread 6 jumps with 2^130 and 2^131)

This gives in log(k) steps (walltime) a set of evenly spaced starting points exactly as if we'd jumped in series.

To pull this off we need to precompute some number of coefficients. Getting 20 of them would be nice as that gets us up to a million particles.

Resources:

• https://peteroupc.github.io/jump.html
• Seeking comments and review on several articles posted here peteroupc/peteroupc.github.io#6
• https://doi.org/10.1016/j.cam.2016.11.006
• https://groups.google.com/g/prng/c/1azP-3QVqeE
• https://groups.google.com/g/prng/c/nHTEt73CFMQ
• https://web.archive.org/web/20180909014435/http://xoroshiro.di.unimi.it/xorshift-1.2.tgz

Other:

• https://libntl.org/