Cartesian product for multiple sets at compile time

Question

I am struggling with an implementation of the Cartesian product for multiple indices with a given range 0,...,n-1.

The basic idea is to have a function:

cartesian_product<std::size_t range, std::size_t sets>()

with an output array that contains tuples that hold the different products

[(0,..,0), (0,...,1), (0,...,n-1),...., (n-1, ..., n-1)]

An simple example would be the following:

auto result = cartesian_product<3, 2>();

with the output type std::array<std::tuple<int, int>, (3^2)>:

[(0,0), (0,1), (0,2), (1,0), (1,1), (1,2), (2,0), (2,1), (2,2)]

My main problem is that my version of the Cartesian product is slow and creates a stack overflow if you choose to have more than 5 sets. I believe that my code has too many recursions and temporary variables.

My implementation (C++17) can be found here: cartesian_product

#include <stdio.h>
#include <iostream>

#include <tuple>


template<typename T, std::size_t ...is>
constexpr auto flatten_tuple_i(T tuple, std::index_sequence<is...>) {

    return std::tuple_cat(std::get<is>(tuple)...);

}

template<typename T>
constexpr auto flatten_tuple(T tuple) {
    return flatten_tuple_i(tuple, std::make_index_sequence<std::tuple_size<T>::value>{});
}

template<std::size_t depth, typename T>
constexpr auto recursive_flatten_tuple(T tuple){

    if constexpr(depth <= 1){
        return tuple;
    }else{
        return recursive_flatten_tuple<depth-1>(flatten_tuple(tuple));
    }
}

template<std::size_t depth, typename T, std::size_t ...is>
constexpr auto wdh(T&& tuple, std::index_sequence<is...>){

    if constexpr (depth == 0) {
        return tuple;
    }else{
        //return (wdh<depth-1>(std::tuple_cat(tuple, std::make_tuple(is)),std::make_index_sequence<sizeof...(is)>{})...);
        return std::make_tuple(wdh<depth-1>(std::tuple_cat(tuple, std::make_tuple(is)), std::make_index_sequence<sizeof...(is)>{})...);
    }
}

template<std::size_t sets, typename T, std::size_t ...is>
constexpr auto to_array(T tuple, std::index_sequence<is...>){

    if constexpr (sets == 0){
        auto t = (std::make_tuple(std::get<is>(tuple)),...);
        std::array<decltype(t), sizeof...(is)> arr = {std::make_tuple(std::get<is>(tuple))...};
        //decltype(arr)::foo = 1;
        return arr;
    }else{
        auto t = ((std::get<is>(tuple)),...);
        std::array<decltype(t), sizeof...(is)> arr = {std::get<is>(tuple)...};
        return arr;
    }
}

template<std::size_t sets, std::size_t ...is>
constexpr auto ct_i(std::index_sequence<is...>){

    if constexpr (sets == 0){

        auto u = std::tuple_cat(wdh<sets>(std::make_tuple(is), std::make_index_sequence<sizeof...(is)>{})...);
        auto arr = to_array<sets>(u, std::make_index_sequence<std::tuple_size<decltype(u)>::value>{});

        return arr;

    }else {

        auto u = std::tuple_cat(wdh<sets>(std::make_tuple(is), std::make_index_sequence<sizeof...(is)>{})...);

        auto r = recursive_flatten_tuple<sets>(u);

        auto d = to_array<sets>(r, std::make_index_sequence<std::tuple_size<decltype(r)>::value>{});


        return d;
    }

}

template<std::size_t range, std::size_t sets>
constexpr auto cartesian_product(){

    static_assert( (range > 0), "lowest input must be cartesian<1,1>" );
    static_assert( (sets > 0), "lowest input must be cartesian<1,1>" );
    return ct_i<sets-1>(std::make_index_sequence<range>{});
}


int main()
{
    constexpr auto crt = cartesian_product<3, 2>();

    for(auto&& ele : crt){

        std::cout << std::get<0>(ele) << " " << std::get<1>(ele) << std::endl;

    }

    return 0;
}

Answer 1

With Boost.Mp11, this is... alright, it's not a one-liner, but it's still not so bad:

template <typename... Lists>
using list_product = mp_product<mp_list, Lists...>;

template <typename... Ts>
constexpr auto list_to_tuple(mp_list<Ts...>) {
    return std::make_tuple(int(Ts::value)...);
}

template <typename... Ls>
constexpr auto list_to_array(mp_list<Ls...>) {
    return std::array{list_to_tuple(Ls{})...};
}

template <size_t R, size_t N>
constexpr auto cartesian_product()
{
    using L = mp_repeat_c<mp_list<mp_iota_c<R>>, N>; 
    return list_to_array(mp_apply<list_product, L>{});
}

With C++20, you can declare the two helper function templates as lambdas inside of cartesian_product, which makes this read nicer (top to bottom instead of bottom to top).

Explanation of what's going on, based on the OP example of cartesian_product<3, 2>:

• mp_iota_c<R> gives us the list [0, 1, 2] (but as integral constant types)
• mp_repeat_c<mp_list<mp_iota_c<R>>, N> gives us [[0, 1, 2], [0, 1, 2]]. We just repeat the list, but we want a list of lists (hence the extra mp_list in the middle).
• mp_apply<list_product, L> does mp_product, which is a cartesian product of all the lists you pass in... sticking the result in an mp_list. This gives you [[0, 0], [0, 1], [0, 2], ..., [2, 2]], but as an mp_list of mp_list of integral constants.
• At this point the hard part is over, we just have to convert the result back to an array of tuples. list_to_tuple takes an mp_list of integral constants and turns that into a tuple<int...> with the right values. And list_to_array takes an mp_list of mp_lists of integral constants and turns that into an std::array of tuples.

---

A slightly different approach using just the single helper function:

template <template <typename...> class L,
    typename... Ts, typename F>
constexpr auto unpack(L<Ts...>, F f) {
    return f(Ts{}...);
}

template <size_t R, size_t N>
constexpr auto cartesian_product()
{
    using P = mp_apply_q<
        mp_bind_front<mp_product_q, mp_quote<mp_list>>,
        mp_repeat_c<mp_list<mp_iota_c<R>>, N>>;

    return unpack(P{},
        [](auto... lists){
            return std::array{
                unpack(lists, [](auto... values){
                    return std::make_tuple(int(values)...);
                })...
            };
        });
}

This approach is harder to read though, but it's the same algorithm.