A basic multi-dimensional array

Question

I have been lurking on the C++ side of stack-overflow only long enough to know there are a lot of beginners and intermediate programmers baffled by multi-dimensional arrays. I've seen a lot of monstrosities : three-star programmers, cache-unfriendly implementations, vectors of vectors of vectors ..., etc.

I thought to provide them with a basic multi-dimensional array:

• open for extensions,
• fit for production and experimentation,
• working without allocations,
• and cache-friendly.

Before I do so, I'd like some advice from the community.

1. What can I improve in the provided features of this class?
2. What new feature could one add, and why would it be useful?

#include <array>
#include <numeric>

namespace ysc
{
namespace _details
{
    template<class InputIt, class OutputIt>
    OutputIt partial_product(InputIt first, InputIt last, OutputIt output)
    { *output++ = 1; return partial_sum(first, last, output, std::multiplies<>{}); }

    // cache-friendly:
    // neighbor objects within the right-most coordinate are neighbors in memory
    template<class TDim, class TCoord>
    auto coordinates_to_index(TDim const& dimensions, TCoord const& coords)
    {
        std::array<std::size_t, dimensions.size()> dimension_product;
        using std::crbegin, std::crend, std::prev;
        partial_product(crbegin(dimensions), prev(crend(dimensions)), begin(dimension_product));
        return std::inner_product(cbegin(dimension_product), cend(dimension_product), crbegin(coords), 0);
    }
}

constexpr struct matrix_zero_t {} matrix_zero;

template<class T, std::size_t... Dimensions>
class matrix
{
template<class, std::size_t...> friend class matrix;

public:
    static constexpr std::size_t order      = sizeof...(Dimensions);
    static constexpr std::array  dimensions = { Dimensions... };

private:
    std::array<T, (Dimensions * ...)> _data;

public:
    friend void swap(matrix& lhs, matrix& rhs)
    {
        using std::swap;
        swap(lhs._data, rhs._data);
    }

public:
    matrix()                          = default;
    matrix(matrix&& other)            = default;
    matrix& operator=(matrix&& other) = default;

    matrix(matrix_zero_t) : _data({}) {}

    template<class U>
    matrix(matrix<U, Dimensions...> const& other) { std::copy(cbegin(other._data), cend(other._data), begin(_data)); }

    template<class U>
    matrix& operator=(matrix<U, Dimensions...> const& other)
    {
        matrix o{other};
        swap(*this, o);
        return *this;
    }

public:
    template<class... Args>
    T const& operator()(Args... coordinates) const
    { return _data[_details::coordinates_to_index(dimensions, std::array{coordinates...})]; }

    template<class... Args>
    T& operator()(Args... coordinates)
    { return _data[_details::coordinates_to_index(dimensions, std::array{coordinates...})]; }
};
} // namespace ysc

Usage demo: http://coliru.stacked-crooked.com/a/1652451c78275436

This is a C++17 implementation; this itself is not set in stone.

Answer 1

Here are some suggestions:

1. coordinates_to_index should return std::size_t instead of int because int may not be able to hold the required values.

   return std::inner_product(cbegin(dimension_product), cend(dimension_product),
                             crbegin(coords), std::size_t{0}); // instead of int
   

2. You don't need to compute the strides and store them. You can use Horner's rule directly.

   std::size_t result = 0;
   for (std::size_t i = 0; i < N; ++i)
       result = result * dimensions[i] + coords[i];
   return result;
   

3. Your calculation of the total size fails to consider the case where order == 0. Instead of

   (Dimensions * ...)                   // ill-formed if Dimensions is empty
   

   Use

   (Dimensions * ... * std::size_t{1})
   

   Moreover, I would suggest exposing the total size, along with the order and dimensions. Like this:

   static constexpr std::size_t size = (Dimensions * ... * std::size_t(1));
   

   And use it in place of ad-hoc computations further in code.

4. Your matrix is not Copyable. You declare move operations but not copy operations (the template doesn't make a difference), and the compiler synthesizes deleted copy operations which take precedence over the copying template. You should either default the copy operations explicitly (and maybe destructor as well), or omit the redundant default declarations directly.

5. Why use :_data({}) when :_data{} is more readable and less verbose?

6. Consider adding out-of-range detection for operator() in some form.

Update

11. Currently, the converting copy constructor first default-initializes the elements and then copy assign from the initializer matrix, which is kinda surprising. It is more reasonable to initialize directly (with some template hacks):

    namespace detail {
        template <class T, std::size_t N, class U, std::size_t... Is>
        std::array<T, N> construct(const std::array<U, N>& arr, std::index_sequence<Is...>)
        {
            return std::array<T, N>{static_cast<T>(arr[Is])...};
        }
    }
    
    template <class T, std::size_t... Dimensions>
    class matrix {
        // ...
    public:
        // ...
        template <class U>
        matrix(matrix<U, Dimensions...> const& other)
            : _data{detail::construct<T, N>(other._data, std::make_index_sequence<size>{})}
        {
        }
        // ...
    };
    

    The above implementation requires access to relies on access to _data of another matrix specialization, which can be granted with a friend declaration or an access function. It also relies on the C++17 guaranteed copy elision.

    Now, the converting copy assignment operator cannot delegate to the converting copy constructor:

    template <class U>
    matrix& operator=(const matrix<U, Dimensions...>& other)
    {
        std::copy_n(_data.begin(), size, other._data.begin());
        return *this;
    }
    

12. It makes sense to constrain the converting copy constructor and the converting copy assignment operator with SFINAE. Also, they should be made conditionally explicit and conditionally noexcept:

    template <class U, class = std::enable_if_t<
        std::is_constructible_v<T, const U&> && std::is_convertible_v<const U&, T>
    >> matrix(matrix<U, Dimensions...> const& other)
        noexcept(std::is_nothrow_constructible_v<T, const U&>);
    
    template <class U, class = std::enable_if_t<
        std::is_constructible_v<T, const U&> && !std::is_convertible_v<const U&, T>
    >> explicit matrix(matrix<U, Dimensions...> const& other)
        noexcept(std::is_nothrow_constructible_v<T, const U&>);
    

    Similarly for the converting copy assignment operator.

13. Consider support converting move construction and converting move assignment.

14. The class can be modified to support constexpr.