Initialization of an array of unknown dimension

Question

I was surprised that I could not find this question. I tried to generalize it (with some nice unchecked code) to the point that everyone can benefit.

Suppose I have a multidimensional Point :

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

Now I create a multidimensional array of them:

  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]; //... } 

As you can see, array is a multidimensional array of unknown dimension, presented in a 1D array.

---

My question is: how can I cleanly initialize this array for its multidimensional grid points?

Here is an example of the behavior I want in 1, 2 and 3 dimensions. Obviously, I do not want to write this for any possible dimension that I can use! The goal is to generalize this.

 //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]; } } } 

Again, my question is: summarize the above for any number of nested loops / measurements in pure form.

Note. I came up with several unpleasant ways and they are inefficient and slow. In particular, I want to avoid recursion, if possible, since this will often be called on multidimensional small data sets. Note. There are obvious parallels in C, so C or C ++ is fine. It was preferable to use C ++ 11.

Answer 1

Edit in response to comment and updated question

If you need performance and elegance, I would:

• abandon the multidimensional array approach and smooth it (i.e., one dimension of the array). No new , no pointer, use the modern C ++ approach with std::vector or std::array .
• provide an abstract container for your multidimensional array with handy methods like a nested loop generator
• Replace the variable arguments with a fixed-size array (since you know dims at compile time.

So, I found the following solution, which is quite consistent with your implementation and needs and try to stay simple.

I managed to rewrite the small MultiArray class into the "modern C ++ 11 way". I believe that the size of count may not be known at compile time, so now using std::vector . Of course, it is possible that a more general and compiled temporary code with std::array can be found in my original answer below.

 #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

As you can see in the results, data_ can be initialized as a whole for N measurements. If you need a higher-level abstract class, you can check below my initial answer, where you can perform some convenient operations (ie Access grid[{i,j,k}] to fill in the values).

Original answer

I needed a multidimensional mesh for my needs, and I asked for improvements to my code in the code review . A working example works here, where, of course, some functions may not be necessary for you ... My implementation relates to template and compile-time calculations. Note that the size size must be known at compile time.

In short, the class will look like this:

 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_; } 

And then you can build your multidimensional array and initialize it in the following ways:

 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. 

If you need to specify an algorithm for variable nested loops to fill in the values, you can look here and do it from the first answer:

 // 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)); 

If you need something that I missed in my implementation, feel free to ask. I would also be happy if someone can point out improvements.

Answer 2

Passing from X,Y,Z to a flattened array (F), we have the following equation

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

or

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

in an array of size 7 x 3 x 5, Z = 3, Y = 1, X = 2 will be 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 

in an array of size 7 x 3 x 5, Z = 4, Y = 2, X = 3 will be 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 

If we store cumulative products in an array, mult { 1, X, X*Y, X*Y*Z, ...} and all the points in the array, val

Points to a flat array:

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

Flat array to points:

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

Then we can iterate over F (flat array), reverse engineer from index to flat array, all points: X, Y, ... as indicated above, and perform initialization in a common loop:

Given mult as mult[0]=1; mult[d+1]=mult[d]*count[d]; 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

Here is my solution using C ++ 11 variation patterns and package extensions.

My "foobar" compiles as a constexpr function, so I would say pretty efficiently: p

As for elegance, you can be a judge, but I think that is pretty good.

Basically, the idea is to replace the for loops with a functional way of programming iteration, where we are simply trying to explicitly create the set that we want to iterate over first. This code can then be generalized to arbitrary dimension in a simpler way.

The code is completely self-contained, with the exception of the <array> header.

It compiles for me in the C ++ 11 standard with gcc 4.8.2 and clang 3.6.

Note. If you want to use the same technique for code where dimensions are a run-time parameter, basically you would override the Cartesian_Product metafile below as a run-time function using something like std::vector<std::vector<size_t>> , You can then create your own set of iterations by taking the dim number of Cartesian products and going through it, using a single loop to populate the result.

 #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); }