How to construct a multi_array::index_gen at runtime

Question

In the following code, the ExtractSubArray function is totally generic, while the ExtractSubArrayCornerAndExtent requires knowledge of the dimensionality at the time the code is written (to construct the sequence of RangeType arguments). Is there any way to write a generic ExtractSubArrayCornerAndExtent (without SFINAE to use a different function for each Dimension (which would be annoying as well as would require a fixed set of possible dimensions).

#include <boost/multi_array.hpp>

template <unsigned int Dimension>
boost::multi_array<double, Dimension> ExtractSubArray(const boost::multi_array<double, Dimension>& array, const typename boost::detail::multi_array::index_gen<Dimension, Dimension>& indices)
{
    using ArrayType = boost::multi_array<double, Dimension>;
    using IndexType = boost::array<double, Dimension>;

    typename ArrayType::template const_array_view<Dimension>::type view = array[indices];

    IndexType subArraySize;
    for(size_t dim = 0 ; dim < Dimension; ++dim) {
        subArraySize[dim] = indices.ranges_[dim].finish_ - indices.ranges_[dim].start_;
    }

    ArrayType subArray = view;

    return subArray;
}

template <unsigned int Dimension>
boost::multi_array<double, Dimension> ExtractSubArrayCornerAndExtent(const boost::multi_array<double, Dimension>& array, const boost::array<double, Dimension>& corner,
                              const boost::array<double, Dimension>& subarraySize)
{
    using ArrayType = boost::multi_array<double, Dimension>;
    using RangeType = typename ArrayType::index_range;

    // Here I have assumed Dimension=3 to produce the second argument. How do you construct this second argument when you don't know Dimension ahead of time.
    return ExtractSubArray<Dimension>(array, boost::indices[RangeType(corner[0],subarraySize[0])][RangeType(corner[1],subarraySize[1])][RangeType(corner[2],subarraySize[2])]);
}

int main()
{
    using ArrayType = boost::multi_array<double, 3>;
    using IndexType = boost::array<double, 3>;

    ArrayType myArray(IndexType({{3,3,3}}));

    std::vector<double> data(9);
    for(size_t i = 0; i < data.size(); ++i) {
        data[i] = i;
    }

    boost::detail::multi_array::index_gen<3,3> indices = boost::indices[ArrayType::index_range(0,1)][ArrayType::index_range(0,1)][ArrayType::index_range(0,1)];

    ArrayType subArray = ExtractSubArray<3>(myArray, indices);

    IndexType corner = {0,0,0};
    IndexType subarraySize = {1,1,1};
    ArrayType subArray2 = ExtractSubArrayCornerAndExtent<3>(myArray, corner, subarraySize);

    return 0;
}

Answer 1

Compile-time iteration generally has to be implemented with recursion. You can build the index object using a recursive template. It can look something like this:

#include <boost/multi_array.hpp>

template <typename T, size_t Dimension>
boost::multi_array<T, Dimension> ExtractSubArray(
   const boost::multi_array<T, Dimension>& array,
   const typename boost::detail::multi_array::index_gen<Dimension, Dimension>& indices)
{
   using ArrayType = boost::multi_array<T, Dimension>;
   using IndexType = boost::array<size_t, Dimension>;

   typename ArrayType::template const_array_view<Dimension>::type view = array[indices];

   IndexType subArraySize;
   for(size_t dim = 0 ; dim < Dimension; ++dim) {
      subArraySize[dim] = indices.ranges_[dim].finish_ - indices.ranges_[dim].start_;
   }

   ArrayType subArray = view;

   return subArray;
}

// Helper functor to build indices.
template<typename RangeArrayType, size_t Dimension>
struct IndicesBuilder {
   // Recursively invoke the functor for the next lowest dimension and
   // add the next range.
   static auto build(const RangeArrayType& range)
      -> decltype(IndicesBuilder<RangeArrayType, Dimension - 1>::build(range)[range[Dimension - 1]]) {
      return IndicesBuilder<RangeArrayType, Dimension - 1>::build(range)[range[Dimension - 1]];
   }
};

// Helper functor specialization to terminate recursion.
template<typename RangeArrayType>
struct IndicesBuilder<RangeArrayType, 1> {
   static auto build(const RangeArrayType& range)
      -> decltype(boost::indices[range[0]]) {
      return boost::indices[range[0]];
   }
};

template <typename T, size_t Dimension>
boost::multi_array<T, Dimension> ExtractSubArrayCornerAndExtent(
   const boost::multi_array<T, Dimension>& array,
   const boost::array<size_t, Dimension>& corner,
   const boost::array<size_t, Dimension>& subarraySize)
{
    using ArrayType = boost::multi_array<T, Dimension>;
    using RangeType = typename ArrayType::index_range;

    // Build a random-access container with the ranges.
    std::vector<RangeType> range;
    for (size_t i = 0; i < Dimension; ++i)
       range.push_back(RangeType(corner[i], corner[i] + subarraySize[i]));

    // Use the helper functor to build the index object.
    const auto index = IndicesBuilder<decltype(range), Dimension>::build(range);

    return ExtractSubArray<T, Dimension>(array, index);
}

int main() {
   using ArrayType = boost::multi_array<double, 3>;
   using IndexType = boost::array<size_t, 3>;

   ArrayType myArray(IndexType({{3,3,3}}));

   IndexType corner = {0,0,0};
   IndexType subarraySize = {1,1,1};
   ArrayType subArray2 = ExtractSubArrayCornerAndExtent(myArray, corner, subarraySize);

   return 0;
};

The addition to your code is the IndicesBuilder template struct that has a static function to build one dimension of the index object. That static function calls itself (with different template arguments) to build all the dimensions. Recursion is terminated by a specialization of IndicesBuilder that returns the first dimension. Hooray for trailing return types!

I compiled and ran this code without fatal errors, but I did not do any further testing so there may be minor bugs. The general approach should work, though.