C++ mathematical n-dimensional vector template using fold expressions

Question

I'm writing a template to unify a few mathematical vector types I was using. The goal is low overhead, constexpr where possible and generic utility functions built-in.

Normal usage would be like this:

typedef vector<float, 3> v3d_t;

constexpr v3d_t testv1 { 1.1f, 1.2f, 1.3f };
constexpr v3d_t testv2 { 1.4f, 1.5f, 1.6f };

constexpr auto testv = testv1 * testv2;

// This should compile to a precalculated value:
constexpr auto m = testv.magnitude();

The full source:

#pragma once

#include <cstddef>
#include <cmath>
#include <utility>

/**
 * Generic `S`-dimensional mathematical vector of type `T`
 * @param T Any type that can be parsed by `cmath`
 * @param S Vector dimensions (eg. 2 for a 2D point)
 */
template <typename T, size_t S>
class vector {
public:
    /** Create a new zero-initialized vector */
    constexpr vector() {}
    /** Create a vector based on provided values */
    template <typename...Ts>
    constexpr vector(Ts ...vals) : values{ vals... } {}
    ///TODO: Add clone / copy constructor:
    // constexpr vector(const vector<T, S>& other) : values{ ... } {}

    constexpr size_t size() const { return S; }

    constexpr T& operator[](const size_t index) { return values[index]; }
    constexpr T operator[](const size_t index) const { return values[index]; }

    //@{ Allow looping over values (for (auto&& dim : vector) {...})
    class iterator {
    public:
        constexpr explicit iterator(T *ptr): ptr(ptr) {}
        constexpr iterator operator++() { ++ptr; return *this; }
        constexpr bool     operator!=(const iterator& other) const { return ptr != other.ptr; }
        constexpr T&       operator* () { return *ptr; }
    private:
        T *ptr;
    };
    class citerator {
    public:
        constexpr explicit citerator(const T *ptr): ptr(ptr) {}
        constexpr citerator operator++() { ++ptr; return *this; }
        constexpr bool      operator!=(const citerator& other) const { return ptr != other.ptr; }
        constexpr T         operator* () const { return *ptr; }
    private:
        const T *ptr;
    };
    constexpr auto begin() { return iterator(values); }
    constexpr auto end()   { return iterator(values + S); }
    constexpr auto begin() const { return citerator(values); }
    constexpr auto end()   const { return citerator(values + S); }
    //@}

    template <typename Indices = std::make_index_sequence<S>>
    constexpr T magnitude() const {
        return magnitude_impl(Indices{});
    }

    template <typename Indices = std::make_index_sequence<S>>
    constexpr decltype(auto) direction() const {
        return division(magnitude(), Indices{});
    }

    // Equal sized vectors
    template <typename Indices = std::make_index_sequence<S>>
    constexpr decltype(auto) operator+(const vector<T, S> &other) const { return addition(other.values, Indices{}); }
    template <typename Indices = std::make_index_sequence<S>>
    constexpr decltype(auto) operator-(const vector<T, S> &other) const { return substraction(other.values, Indices{}); }
    template <typename Indices = std::make_index_sequence<S>>
    constexpr decltype(auto) operator*(const vector<T, S> &other) const { return multiplication(other.values, Indices{}); }
    template <typename Indices = std::make_index_sequence<S>>
    constexpr decltype(auto) operator/(const vector<T, S> &other) const { return division(other.values, Indices{}); }

    // Plain arrays with equal size
    template <typename Indices = std::make_index_sequence<S>>
    constexpr decltype(auto) operator+(const T (&other)[S]) const { return addition(other, Indices{}); }
    template <typename Indices = std::make_index_sequence<S>>
    constexpr decltype(auto) operator-(const T (&other)[S]) const { return substraction(other, Indices{}); }
    template <typename Indices = std::make_index_sequence<S>>
    constexpr decltype(auto) operator*(const T (&other)[S]) const { return multiplication(other, Indices{}); }
    template <typename Indices = std::make_index_sequence<S>>
    constexpr decltype(auto) operator/(const T (&other)[S]) const { return division(other, Indices{}); }

    // Singleton values
    template <typename Indices = std::make_index_sequence<S>>
    constexpr decltype(auto) operator+(const T &other) const { return addition(other, Indices{}); }
    template <typename Indices = std::make_index_sequence<S>>
    constexpr decltype(auto) operator-(const T &other) const { return substraction(other, Indices{}); }
    template <typename Indices = std::make_index_sequence<S>>
    constexpr decltype(auto) operator*(const T &other) const { return multiplication(other, Indices{}); }
    template <typename Indices = std::make_index_sequence<S>>
    constexpr decltype(auto) operator/(const T &other) const { return division(other, Indices{}); }

    ///TODO: add +=, -=, etc operators

    // https://en.wikipedia.org/wiki/Euclidean_vector#Equality
    template <typename Indices = std::make_index_sequence<S>>
    constexpr bool operator==(const vector &other) const {
        return equality(other.values, Indices{});
    }
    template <typename Indices = std::make_index_sequence<S>>
    constexpr bool operator!=(const vector &other) const {
        return !equality(other.values, Indices{});
    }

    //@{ These compare the magnitude of two vectors, might be non-standard behaviour, beware!
    constexpr bool operator< (const vector &other) const { return magnitude() <  other.magnitude(); }
    constexpr bool operator<=(const vector &other) const { return magnitude() <= other.magnitude(); }
    constexpr bool operator> (const vector &other) const { return magnitude() >  other.magnitude(); }
    constexpr bool operator>=(const vector &other) const { return magnitude() >= other.magnitude(); }
    // template <typename U> constexpr bool operator==(const U other) const { return magnitude() == other; }
    // template <typename U> constexpr bool operator!=(const U other) const { return magnitude() != other; }
    template <typename U> constexpr bool operator< (const U other) const { return magnitude() <  other; }
    template <typename U> constexpr bool operator<=(const U other) const { return magnitude() <= other; }
    template <typename U> constexpr bool operator> (const U other) const { return magnitude() >  other; }
    template <typename U> constexpr bool operator>=(const U other) const { return magnitude() >= other; }
    //@}
private:
    T values[S] {0};

    template <size_t... I>
    constexpr vector<T, S> addition(const T (&other)[S], std::index_sequence<I...>)  const {
        return { (values[I] + other[I])... };
    }
    template <size_t... I>
    constexpr vector<T, S> addition(const T &other, std::index_sequence<I...>)  const {
        return { (values[I] + other)... };
    }
    template <size_t... I>
    constexpr vector<T, S> substraction(const T (&other)[S], std::index_sequence<I...>)  const {
        return { (values[I] - other[I])... };
    }
    template <size_t... I>
    constexpr vector<T, S> substraction(const T &other, std::index_sequence<I...>)  const {
        return { (values[I] - other)... };
    }
    template <size_t... I>
    constexpr vector<T, S> multiplication(const T (&other)[S], std::index_sequence<I...>)  const {
        return { (values[I] * other[I])... };
    }
    template <size_t... I>
    constexpr vector<T, S> multiplication(const T &other, std::index_sequence<I...>)  const {
        return { (values[I] * other)... };
    }
    template <size_t... I>
    constexpr vector<T, S> division(const T (&other)[S], std::index_sequence<I...>)  const {
        return { (values[I] / other[I])... };
    }
    template <size_t... I>
    constexpr vector<T, S> division(const T &other, std::index_sequence<I...>)  const {
        return { (values[I] / other)... };
    }

    template <size_t... I>
    constexpr T magnitude_impl(std::index_sequence<I...>) const {
        return sqrt(((values[I] * values[I]) + ...));
    }

    template <size_t... I>
    constexpr bool equality(const T (&other)[S], std::index_sequence<I...>) const {
        //TODO: Verify if fold expressions can actually short circuit.
        return ((values[I] == other[I]) && ...);
    }
};

It's a work in progress (still some TODOs here and there). I'd appreciate some early feedback.

Answer 1

1. You can return a raw pointer instead of creating custom iterators. Your current version doesn't support some operations, like postfix increment. A raw pointer is already a random-access iterator. It'll make the code much simpler. The only advantage of your approach is that it makes vector iterators a separate entity so that one cannot accidentally assign a raw pointer to them (or vice versa).

2. I don't think the names like addition, multiplication and so on are good. Addition is a noun representing an operation. But it actually returns a sum. Multiplication returns a product. I'd call them sum, product, difference and so on. add, multiply ... would also do fine.

3. A citerator operator * should return a const reference, not a copy. It would be more efficient that way (it doesn't really matter for types like int as they're easy to copy, though).