初始化一个未知维度的数组

Question

我很惊讶我找不到这个问题。我试图将它（使用一些不错的未经测试的代码）概括为每个人都可以从中受益的东西。

假设我有一个多维Point：

template <int dims> class Point { public: double data[dims]; };

现在我为它们创建一个多维数组：

 template <int dims> void foobar(int count0, ...) {
    //Using variadic function.  Could also use variadic templates in C++ (arguably better)
    int counts[dims], total_count=count0; counts[0]=count0;
    va_list args; va_start(args,count0);
    for (int i=1;i<dims;++i) {
        int count = va_arg(args,int);
        counts[i] = count;
        total_count *= count;
    }
    va_end(args);

    Point<dims>* array = new Point<dims>[total_count];

    //...
}

如您所见，array 是一个维度未知的多维数组，以一维数组表示。

---

我的问题：我如何干净地将此数组初始化为其多维网格点？

这是我想要的 1、2 和 3 维的示例行为。显然，我不想为我可能想要使用的每个可能的维度都写这个！目标是概括这个。

//Example: dim==1
for (int x=0; x<counts[0]; ++x) {
    Point<1>& point = array[x];
    point.data[0] = (x+0.5) / (double)counts[0];
}

//Example: dim==2
for (int y=0; y<counts[1]; ++y) {
    for (int x=0; x<counts[0]; ++x) {
        Point<2>& point = array[y*counts[0]+x];
        point.data[0] = (x+0.5) / (double)counts[0];
        point.data[1] = (y+0.5) / (double)counts[1];
    }
}

//Example: dim==3
for (int z=0; z<counts[2]; ++z) {
    for (int y=0; y<counts[1]; ++y) {
        for (int x=0; x<counts[0]; ++x) {
            Point<3>& point = array[(z*counts[1]+y)*counts[0]+x];
            point.data[0] = (x+0.5) / (double)counts[0];
            point.data[1] = (y+0.5) / (double)counts[1];
            point.data[2] = (z+0.5) / (double)counts[2];
        }
    }
}

再次，我的问题：以简洁的方式将上述内容概括为任意数量的嵌套循环/维度。

注意：我想出了一些讨厌的方法，它们既不优雅又很慢。特别是，如果可能的话，我想避免递归，因为这将在高维小型数据集上非常频繁地调用。 注意：C 中有明显的相似之处，因此 C 或 C++ 都可以。首选 C++11。

Answer 1

根据评论和更新问题进行编辑

如果您需要性能和“优雅”，我会：

• 放弃您的多维数组方法并将其展平（即一个数组维度）。没有new，没有指针，使用带有std::vector 或std::array 的C++ 现代方法。
• 为您的多维度数组提供一个抽象容器，使用方便的方法，例如通用嵌套循环“生成器”
• 用固定大小的数组替换可变参数（因为您在编译时就知道dims。

因此，我找到了以下解决方案，该解决方案与您的实施和需求非常一致，并尽量保持简单。

我已经设法以“现代 C++11 方式”重写了一个小的 MultiArray 类。我在这里认为count 维度在编译时可能不知道，因此现在使用std::vector。当然可以使用std::array 获得更通用的编译时间代码，请参阅下面的原始答案。

#include <iostream>
#include <array>
#include <vector>
#include <numeric>

template<size_t DIMS>
class MultiArray {
public:
    // Point here is just an array
    using Point = std::array<double,DIMS>;

    // fill data_ with an init array
    // not that count is just a fix sized array here no variadic arguments needed
    MultiArray(const std::array<size_t, DIMS>& count) 
        : data_{init_array(count)} {}

private:   
    // the init functions are used for the constructor
    void init_point(Point& point, const std::array<size_t,DIMS>& coord, const std::array<size_t, DIMS>& count) {
        std::cout << " -> { ";
        for (size_t i = 0; i < DIMS; i ++) {
            point[i] = (coord[i] + 0.5) / count[i];
            std::cout << point[i] << ";";
        }
        std::cout << " }\n";
    }

    std::vector<Point> init_array(const std::array<size_t, DIMS>& count) {
        std::vector<Point> data(std::accumulate(count.begin(), count.end(), 1, std::multiplies<int>())); // accumulate computes the prod of DIMS total_count
        std::array<size_t, DIMS> current{};
        size_t i=0;
        do {
             for (size_t i = 0; i < DIMS; i ++)
                 std::cout << current[i] << ";";
            init_point(data[i++],current,count);
        } while (increment(current, count));
        return data;
    }

    // the following function allows to imitate the nested loop by incrementing multidim coordinates
    bool increment( std::array<size_t, DIMS>& v, const std::array<size_t, DIMS>& upper) {
        for (auto i = v.size(); i-- != 0; ) {
            ++v[i];
            if (v[i] != upper[i]) {
                return true;
            }
            v[i] = 0;
        }
        return false;
    }   
private:
    std::vector<Point> data_; // A flatten multi dim vector of points
};


int main() {
   std::array<size_t, 3> count{{4,5,3}};
   MultiArray<3> test{count};
}

Live on Coliru

正如您在结果中看到的那样，data_ 可以针对N 维度进行一般初始化。如果您需要更高级别的抽象类，可以查看我的原始答案，您可以在其中执行一些方便的操作（即访问grid[{i,j,k}] 以填充值）。

原答案

我需要一个多维网格来满足我的需要，并且碰巧在code review 上询问我的代码改进。这里是一个有效的example，当然你可能不需要某些功能......我的实现与模板和编译时间计算有关。请注意，尺寸大小必须在编译时已知。

简单地说，这个类应该是这样的：

template<typename T, size_t... DIMS> // variadic template here for the dimensions size
class MultiGrid {
   // Access from regular idx such as grid[64]
   T& operator[] (size_type idx)       { return values_[idx]; };
   // Access from multi dimensional coordinates such as grid[{6,3,4}]
   T& operator[] (const std::array<size_t,sizeof...(DIMS)>& coord)       { // can give code for runtime here };
private:
  std::array<T,sizeof...(DIMS)> data_;
}

然后你可以构造你的多维数组并用这些方式初始化它：

MultiGrid<float,DIM1,DIM2,DIM3> data; // 3D
// MultiGrid<float,DIM1,DIM2,DIM3,DIM4> data; // 4D
// etc...

// initialize it like this with nested arrays
for (size_t z=0; z < DIM3; z ++)
  for (size_t y=0; y < DIM2; y ++)
    for (size_t x=0; x < DIM1; x ++)
      data[{x,y,z}] = [...] // whatever

// or like this in C++11/14 way
for (auto &x : data) x = [...] // this is convenient to provide a container like approach since no nested arrays are needed here.

如果您需要为可变参数嵌套循环指定一个算法来填写值，您可以查看here，并在第一个答案中这样做：

// here lower_bound is 0-filled vector
std::vector<int> current = lower_bound;
do {
   data[current] = [...] // fill in where current is a coordinate
} while (increment(current, lower_bound, upper_bound));

如果您需要我在实施过程中遗漏的东西，请随时提出。如果有人能指出改进，我也会很高兴。

Answer 2

从X,Y,Z 到扁平数组 (F) 我们有以下等式

F=(Z*DimY+y)*DimX+X

或

F=Z*DimY*DimX+Y*DimX+X

X = F % DimX
Y = F % DimX*DimY/DimX
Z = F % DimX*DimY*DimZ/DimX*DimY

在 7 x 3 x 5 数组中，Z=3, Y=1, X=2 将位于 3*3*5 + 1*5 + 2= 45+5+2= 52

X = `52` % 5 = 2
Y = `52` % (5 * 3) / 5 = 7 / 5 = 1
Z = `52` % (7 * 5 * 3)/15 = 52/15 = 3

在 7 x 3 x 5 数组中，Z=4, Y=2, X=3 将位于 4*3*5 + 2*5 + 3= 60+10+3= 73

X = `73` % 5 = 3
Y = `73` % (5 * 3) / 5 = 13 / 5 = 2
Z = `73` % (7 * 5 * 3)/15 = 73/15 = 4

如果我们将累积乘积保存在数组中，mult{ 1, X, X*Y, X*Y*Z, ...} 和数组中的所有点，val

指向平面数组：

F=sum(mult[i]*val[i]);

平面数组到点：

i[0]=F%mult[1]/mult[0];
i[1]=F%mult[2]/mult[1];
...

然后我们可以迭代 F（平面数组），从索引逆向工程到平面数组中的所有点：X,Y,... 如上所述，并在通用循环中进行您想要的初始化：

给定mult 为mult[0]=1; mult[d+1]=mult[d]*count[d];

for (int i = 0; i < total_count; ++i) {
    for (int d=0; d < dims; ++d) {
      int dim=(i%mult[d+1])/mult[d];
      point.data[d] = (dim+0.5) / (double)counts[d];
  }
}

Answer 3

这是我的解决方案，使用 C++11 可变参数模板和包扩展。

我的“foobar”编译为 constexpr 函数，所以我认为这非常有效：p

至于优雅，你可以评判，但我觉得还不错。

基本上，这个想法是用函数式编程方式替换for 循环进行迭代，我们只是寻求显式地构建我们想要首先迭代的集合。然后可以以更直接的方式将该代码推广到任意维度。

除了<array> 标头之外，代码完全独立。

它使用 gcc 4.8.2 和 clang 3.6 按照 C++11 标准为我编译。

注意：如果您想对维度是运行时参数的代码使用相同的技术，基本上您要做的就是使用类似 std::vector<std::vector<size_t>> 的东西将下面的 Cartesian_Product 元函数重新实现为运行时函数.然后，您可以通过获取dim 笛卡尔积的数量来构建您的迭代集，并使用单个 for 循环对其进行迭代以填充结果。

#include <array>
#include <iostream>

/////////////////////////////////////////////
// The point class from problem statement
/////////////////////////////////////////////

template <size_t dims> class Point { public: double data[dims]; };

// Some basic type list and meta function objects to work with

// List of indices
template<size_t... sl>
struct StList {
    static constexpr size_t length = sizeof...(sl);
};

// Metafunction to compute a volume
template <typename T>
struct Compute_Volume;

template<size_t s>
struct Compute_Volume<StList<s>> {
    static constexpr size_t value = s;
};

template <size_t s1, size_t s2, size_t... sl>
struct Compute_Volume<StList<s1, s2, sl...>> {
    static constexpr size_t value = s1 * Compute_Volume<StList<s2, sl...>>::value;
};

// Concatenate
template<typename TL, typename TR>
struct Concat;

template<size_t... SL, size_t... SR>
struct Concat<StList<SL...>, StList<SR...>> {
    typedef StList<SL..., SR...> type;
};

template<typename TL, typename TR>
using Concat_t = typename Concat<TL, TR>::type;

// Meta function to check if a typelist is component-wise less than another typelist
// For testing purposes
template <typename TL, typename TR>
struct all_less_than;

template <size_t l1, size_t... SL, size_t r1, size_t... SR>
struct all_less_than<StList<l1, SL...>, StList<r1, SR...>> {
    static constexpr bool value = (l1 < r1) && all_less_than<StList<SL...>, StList<SR...>>::value;
};

template<>
struct all_less_than<StList<>, StList<>> {
    static constexpr bool value = true;
};

/////////////////////////////////////////////////////////////////////////////
// constexpr template function for the point initializations you described
/////////////////////////////////////////////////////////////////////////////
template <typename index, typename dims>
struct point_maker;

template <size_t... idx, size_t... dim>
struct point_maker<StList<idx...>, StList<dim...>> {
    static_assert(sizeof...(idx) == sizeof...(dim), "misuse of 'point_maker' template, idx and dim must have same number of coordinates");
    static_assert(all_less_than<StList<idx...>, StList<dim...>>::value, "misuse of 'point_maker' template, idx is out of bounds");
    static constexpr Point<sizeof...(idx)> make_point() {
        return {{ ((idx + 0.5) / static_cast<double>(dim))... }};
    }
};

//////////////////////////////////////////////////////////////////////////////////////////
// Now we need to abstract the for loops. We need a little more infrastructure for this.
//////////////////////////////////////////////////////////////////////////////////////////

// A basic typelist object
template <typename... tl>
struct TypeList {
    static constexpr size_t length = sizeof...(tl);
};

// Specialization for the Concat metafunction
template<typename... TL, typename... TR>
struct Concat<TypeList<TL...>, TypeList<TR...>> {
    typedef TypeList<TL..., TR...> type;
};


// Metafunction to compute the cartesian product of two lists of lists, and evaluate the `Concat` metafunction for each pair.
template <typename S, typename T>
struct Cartesian_Product;

template <typename s1, typename s2, typename... sl, typename T>
struct Cartesian_Product<TypeList<s1, s2, sl...>, T> {
    typedef Concat_t<typename Cartesian_Product<TypeList<s1>, T>::type, typename Cartesian_Product<TypeList<s2, sl...>, T>::type> type;
};

template <typename s, typename t1, typename t2, typename... tl>
struct Cartesian_Product<TypeList<s>, TypeList<t1, t2, tl...>> {
    typedef Concat_t<typename Cartesian_Product<TypeList<s>, TypeList<t1>>::type, typename Cartesian_Product<TypeList<s>, TypeList<t2, tl...>>::type> type;
};

template <typename s, typename t>
struct Cartesian_Product<TypeList<s>, TypeList<t>> {
    typedef TypeList<Concat_t<s, t>> type;
};

template <typename S, typename T>
using Cartesian_Product_t = typename Cartesian_Product<S, T>::type;

// Some unit tests for the above :)

static_assert( std::is_same<Cartesian_Product_t<TypeList<StList<1>, StList<2>>, TypeList<StList<3>, StList<4>>>, TypeList<StList<1,3>, StList<1,4>, StList<2,3>, StList<2,4>>>::value , "cartesian product not working");

static_assert( std::is_same<Cartesian_Product_t<TypeList<StList<1,5>, StList<2>>, TypeList<StList<3>, StList<4>>>, TypeList<StList<1,5,3>, StList<1,5,4>, StList<2,3>, StList<2,4>>>::value , "cartesian product not working");

static_assert( std::is_same<Cartesian_Product_t<TypeList<StList<1,5>, StList<2>>, TypeList<StList<>>>, TypeList<StList<1,5>, StList<2>>>::value , "cartesian product not working");

// Metafunction to count from 0 to n

// Count from zero to n-1, and make a sequence of singleton sets containing the numbers
template<size_t s>
struct Count {
    typedef Concat_t<typename Count<s-1>::type, TypeList<StList<s-1>>> type;
};

template<>
struct Count<0> {
    typedef TypeList<> type;
};

template<size_t s>
using Count_t = typename Count<s>::type;


// Metafunction to abstract a series of for loops, generating the list of all the index tuples the collection of loops would generate

template <typename T>
struct Volume_Maker;

template <>
struct Volume_Maker<StList<>> {
    typedef TypeList<StList<>> type;
};

template <size_t s, size_t... sl>
struct Volume_Maker<StList<s, sl...>> {
    typedef Cartesian_Product_t<Count_t<s>, typename Volume_Maker<StList<sl...>>::type> type;
};

template <typename T>
using Volume_t = typename Volume_Maker<T>::type;

// Some more quick unit tests

template <typename T>
struct Volume_Test {
    static_assert( Volume_t<T>::length == Compute_Volume<T>::value, "volume computation mismatch");
};

Volume_Test<StList<1,2,3>> test1{};
Volume_Test<StList<1,1,1>> test2{};
Volume_Test<StList<1>> test3{};
Volume_Test<StList<7,6,8>> test4{};

/////////////////
// Grand finale
/////////////////

template <typename dim_list, typename idx_list = Volume_t<dim_list>>
struct foobar_helper;

template <size_t... dim, typename... idx>
struct foobar_helper<StList<dim...>, TypeList<idx...>> {
    typedef StList<dim...> dim_list;
    typedef std::array<Point<sizeof...(dim)>, sizeof...(idx)> array_type;

    static constexpr array_type make_array() {
        return {{ point_maker<idx, dim_list>::make_point()... }};
    }
};

template <size_t... dim>
constexpr std::array<Point<sizeof...(dim)>, Compute_Volume<StList<dim...>>::value> foobar() {
    return foobar_helper<StList<dim...>>::make_array();
}

/////////////////
// Example usage
/////////////////

template <size_t ndim>
std::ostream & operator << (std::ostream & s, const Point<ndim> & p) {
    s << "[ ";
    for (size_t i = 0; i < ndim; ++i) {
        if (i) { s << ", "; }
        s << p.data[i];
    }
    s << " ]";
    return s;
}

template<size_t ndim, size_t N>
void print_array(std::ostream & s, const std::array<Point<ndim>, N> & a) {
    for (size_t i = 0; i < N; ++i) {
        s << a[i] << std::endl;
    }
}

int main() {
    constexpr auto array1 = foobar<2,2,2>();
    print_array(std::cout, array1);

    constexpr auto array2 = foobar<3,3,3,3>();
    print_array(std::cout, array2);
}