Concisely declare and initialize a multi-dimensional array in C++

Question

For example in 3 dimensions, I would normally do something like

vector<vector<vector<T>>> v(x, vector<vector<T>>(y, vector<T>(z, val)));

However this gets tedious for complex types and in large dimensions. Is it possible to define a type, say, tensor, whose usage would be like so:

tensor<T> t(x, y, z, val1);
t[i][j][k] = val2;

Any reason you can't just flatten it to 1D and use indexing like vec[i*rows+j]? — Tyler V
– Tyler V, Commented Dec 24, 2021 at 23:23
This is more of a theoretical question as to the cleanest way of doing this, and using that type of indexing is very error prone. Would it be possible to define a method that takes an index for each dimension and returns the element there (so you don't necessarily have to overload the [] operator)? — Simon Berens
– Simon Berens, Commented Dec 24, 2021 at 23:27
You can write a simple class that holds the 1D vector and overrides operator() with N arguments to handle the indexing easily enough. — Tyler V
– Tyler V, Commented Dec 24, 2021 at 23:30
A vector of a vector of a vector is very inefficient, as it gives triple indirection. Will this not just fit in an array? Anyhow these design questions are very old, and there are a number of libraries that have become very good at this, e.g. Boost UBLAS or Eigen. Also look at the new linear algebra proposal for C++. — JHBonarius
– JHBonarius, Commented Dec 25, 2021 at 7:18

Sandro J · Accepted Answer · 2021-12-24 23:36:43Z

3

It's possible with template metaprogramming.

Define a vector NVector

template<int D, typename T>
struct NVector : public vector<NVector<D - 1, T>> {
    template<typename... Args>
    NVector(int n = 0, Args... args) : vector<NVector<D - 1, T>>(n, NVector<D - 1, T>(args...)) {
    }
};

template<typename T>
struct NVector<1, T> : public vector<T> {
    NVector(int n = 0, const T &val = T()) : vector<T>(n, val) {
    }
};

You can use it like this

    const int n = 5, m = 5, k = 5;
    NVector<3, int> a(n, m, k, 0);
    cout << a[0][0][0] << '\n';

I think it's clear how it can be used. Let's still say NVector<# of dimensions, type> a(lengths of each dimension separated by coma (optional)..., default value (optional)).

edited Dec 24, 2021 at 23:36

answered Dec 24, 2021 at 23:29

Sandro J

5583 silver badges14 bronze badges

Sign up to request clarification or add additional context in comments.

Comments

Tyler V · Accepted Answer · 2021-12-25 19:55:46Z

The other answer shows a good way of making a vector of vectors with template metaprogramming. If you want a multidimensional array data structure with fewer allocations and contiguous storage underneath, here is an example of how to achieve that with an NDArray template class wrapping access to an underlying vector. This could be extended to define non-default operator=, copy operators, debug bounds checking per dimension, row-major or column-major storage, etc for extra convenience.

NDArray.h

#pragma once

#include <array>
#include <vector>

template<int N, typename ValueType>
class NDArray {
public:
    template<typename... Args>
    NDArray(Args... args)
    : dims({ args... }),
      offsets(compute_offsets(dims)),
      data(compute_size(dims), ValueType{})
    {
        static_assert(sizeof...(args) == N, 
            "Incorrect number of NDArray dimension arguments");
    }

    void fill(ValueType val) {
        std::fill(data.begin(), data.end(), val);
    }

    template<typename... Args>
    inline void resize(Args... args) {
        static_assert(sizeof...(args) == N,
            "Incorrect number of NDArray resize arguments");
        dims = { args... };
        offsets = compute_offsets(dims);
        data.resize(compute_size(dims));
        fill(ValueType{});
    }

    template<typename... Args>
    inline ValueType operator()(Args... args) const {
        static_assert(sizeof...(args) == N, 
            "Incorrect number of NDArray index arguments");
        return data[calc_index({ args... })];
    }

    template<typename... Args>
    inline ValueType& operator()(Args... args) {
        static_assert(sizeof...(args) == N, 
            "Incorrect number of NDArray index arguments");
        return data[calc_index({ args... })];
    }

    int length(int axis) const { return dims[axis]; }

    const int num_dims = N;

private:
    static std::array<int, N> compute_offsets(const std::array<int, N>& dims) {
        std::array<int, N> offsets{};
        offsets[0] = 1;
        for (int i = 1; i < N; ++i) {
            offsets[i] = offsets[i - 1] * dims[i - 1];
        }
        return offsets;
    }

    static int compute_size(const std::array<int, N>& dims) {
        int size = 1;
        for (auto&& d : dims) size *= d;
        return size;
    }

    inline int calc_index(const std::array<int, N>& indices) const {
        int idx = 0;
        for (int i = 0; i < N; ++i) idx += offsets[i] * indices[i];
        return idx;
    }

    std::array<int, N> dims;
    std::array<int, N> offsets;
    std::vector<ValueType> data;
};

This overrides the operator() with the correct number of arguments, and won't compile if the wrong number of arguments is given. Some example use

using Array2D = NDArray<2,double>;
using Array3D = NDArray<3,double>;

auto a = Array2D(3, 6);
a.fill(1.0);
a(2, 4) = 2.0;
//a(2,4,4) will not compile
std::cout << "a = " << std::endl << a << std::endl;

a.resize(2,2);
a(1,1) = 1.2;
std::cout << "a = " << std::endl << a << std::endl;

//auto b = Array3D(4, 4); // will not compile

auto b = Array3D(4, 3, 2);
b.fill(-1.0);
b(0, 0, 0) = 4.0;
b(1, 1, 1) = 2.0;
std::cout << "b = " << std::endl << b << std::endl;

(using helper output methods for 2D and 3D arrays)

std::ostream& operator<<(std::ostream& os, const Array2D& arr) {
    for (int i = 0; i < arr.length(0); ++i) {
        for (int j = 0; j < arr.length(1); ++j) {
            os << arr(i,j) << " ";
        }
        os << std::endl;
    }
    return os;
}

std::ostream& operator<<(std::ostream& os, const Array3D& arr) {
    for (int k = 0; k < arr.length(2); ++k) {
        os << "array(:,:,"<<k<<") = " << std::endl;
        for (int i = 0; i < arr.length(0); ++i) {
            os << "  ";
            for (int j = 0; j < arr.length(1); ++j) {
                os << arr(i, j, k) << " ";
            }
            os << std::endl;
        }
        os << std::endl;
    }
    return os;
}

Collectives™ on Stack Overflow

Concisely declare and initialize a multi-dimensional array in C++

2 Answers 2

Comments

Comments

Your Answer

Hot Network Questions

Collectives™ on Stack Overflow

2 Answers 2

Comments

Comments

Your Answer

Sign up or log in

Post as a guest

Related