Best practice on loop on cartesian set A1 x … x An

[tobiasBora]

I need to loop over a cartesian set of an "unknown" (i.e. specified in a variable) number of copies (for now) of a set. Like [0, …, 9] x …(repeat N times)… x [0, ..., 9].

For now I do it this way, i.e. loop over integers from 0 to 10^N-1, and compute the "digits" of each number like:

def ith_digit basis (digit_pos:i64) (x:i64) = (x / basis**digit_pos) % basis
def list_from_number (basis:i64) (size_list:i64) (x:i64) = (iota size_list) |> map (\i -> ith_digit basis i x)
def multi_cartesian_product_int_version (basis:i64) (size_list:i64) = (iota (basis**size_list)) |> map (list_from_number basis size_list)

this seems to work fairly well, and should be quite efficient, but I've a few questions as I'm new to GPU computations and futhark:

• is there a more idiomatic way to deal with this? Like what if I have different sets (specified in a list) ? (I guess I can still decompose it similarly, but the math is a bit more ugly.)
• then I need to do a sequential operation for each element $x$ of this cartesian set that works on an array of size like max(x). My understanding is that futhark works better if all arrays in each branch have the same size, so does it mean that the space of my computation will be quite large since each branch will have the largest possible space size?
• related question: for some branches (possibly a lot), I know that they are not worth exploring (I will then take the "maximum" number of each branch and math helps me to prune some branches). But my understanding is that filter is quite expensive, so it may not be a good idea to use filter to remove these branches, yet adapting my number-to-digit operation to account for these holes might be quite complex… Yet if these branches still use a large array, my memory will quickly explode… So should I still use a filter and hope that this will be negligible compared to the sequential operation that runs after? Or can I somehow let futhark know that this branch is not worth exploring and let it fill this space with a more interesting thread?
• what if the GPU (CUDA) cannot store all branches in memory? Will futhark automatically send them in chunks that fit in memory, or do I need to take care of this with an external program?

[athas]

is there a more idiomatic way to deal with this? Like what if I have different sets (specified in a list) ? (I guess I can still decompose it similarly, but the math is a bit more ugly.)

The neat thing about niche languages is that you can be part of deciding what is idiomatic. The code you have written looks fine to my eyes. You can optimise the computation of all those ith_digits by using a scan to avoid potentially expensive ** operations (remember that these are essentially loops). E.g. instead of

tabulate 20 (\i -> 1.1**f64.i64 (i+1))

you can write

scan (*) 1 (replicate 20 1.1)

then I need to do a sequential operation for each element x of this cartesian set that works on an array of size like max(x). My understanding is that futhark works better if all arrays in each branch have the same size, so does it mean that the space of my computation will be quite large since each branch will have the largest possible space size?

This is unbalanced work, which is intrinsically difficult. And yes, Futhark will compute the maximum memory requirement and multiply that with the number of threads. (At one point we summed the memory requirements instead, but that results in far more complicated access patterns, and was in practice pretty slow.) There is a trick that you maybe can use but it requires somewhat careful (and occasionally impractical) programming.

So should I still use a filter and hope that this will be negligible compared to the sequential operation that runs after?

This is also a difficult tradeoff. The Futhark compiler will not help you do anything clever (it's actually not a very smart compiler, it is just very thorough and consistent). In practice I sometimes see people do filter-based pruning after doing multiple steps of their algorithm, with the number of steps between pruning being some heuristically determined parameter. It depends a lot on your algorithm.

what if the GPU (CUDA) cannot store all branches in memory? Will futhark automatically send them in chunks that fit in memory, or do I need to take care of this with an external program?

You have to do this manually, although CUDA will also sometimes do transparent paging between GPU memory and CPU memory.

[tobiasBora]

Thank you very much for your answer.

The code you have written looks fine to my eyes.

Ok good then. I was hoping to see a bit of a simpler construct as cartesian sets are fairly common (I guess), but at least it is very clear what the code is doing.

You can optimise the computation of all those ith_digits by using a scan to avoid potentially expensive ** operations (remember that these are essentially loops)

Really? Do this (at least) use stuff like square and multiply (in O(log(n))) to implement this instead of the naive O(n) operation? (if so, wouldn't it equal the complexity of scan that, I guess, will also involve a log(n) depth… I can't find any info in here) The number I deal with are fairly small so I would not expect this to take much time compared to the rest of the computation, and I'm afraid that the additional cost of spawning a thread might overcome the benefits… This brings another related question: is the cost of spawning new threads noticeable in these scenarios? Said differently, can I always safely use a parallel implementation instead of a sequential one, or can thread spawn create additional costs?

This is unbalanced work, which is intrinsically difficult.

I see, reading about memory interleaving helps to understand difficulties here.

There is a trick that you maybe can use

Wait, I'm not sure to understand which trick you refer to. I see one trick which is to define a bigger array above the loop, basically containing all sub-arrays (this is an interesting trick that may work for simpler settings (I guess you won't have interleaving here, so you must deal with this manually if this is important for you), but I should then be able to find a mapping from the thread ID to the position in the big array, which may also not be trivial to obtain/compute…), and I see another trick mentioning scan but without further details… What is this trick exactly?

In practice I sometimes see people do filter-based pruning after doing multiple steps of their algorithm, with the number of steps between pruning being some heuristically determined parameter. It depends a lot on your algorithm.

Ok thanks, I'll experiment!

You have to do this manually, although CUDA will also sometimes do transparent paging between GPU memory and CPU memory.

Oh good to know, thanks!