5
\$\begingroup\$

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.

Improve this question
\$\endgroup\$
10
  • \$\begingroup\$ Please do not forget to add c++ tag to make it visible on C++ page. \$\endgroup\$ – Incomputable Aug 9 '17 at 11:08
  • \$\begingroup\$ I wonder... do compilers vectorise this today on machines which are capabale of? \$\endgroup\$ – Maikel Aug 9 '17 at 11:39
  • \$\begingroup\$ I ask because simple loops are auto-vectorisable and this would clearly be preferable. \$\endgroup\$ – Maikel Aug 9 '17 at 11:45
  • \$\begingroup\$ Good point about the vectorization, I'll run some benchmarks once I've got some more stuff sorted out. My rationale for the fold expressions was actually that they are guaranteed to be compile time unrolled and never use temporary indices or values (afaik). \$\endgroup\$ – Stefan Aug 9 '17 at 14:59
  • 1
    \$\begingroup\$ @Maikel Getting back to you on the vectorization: The generated assembly for my test case is exactly the same, in both cases the main loop repeating the product function call is vectorized, nothing else (I tried bumping the vectors to 16D, no difference). The performance difference measured earlier must have been some fluke. \$\endgroup\$ – Stefan Aug 10 '17 at 18:28
2
\$\begingroup\$

  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).

Improve this answer
\$\endgroup\$
5
  • \$\begingroup\$ Another problem with multiply (specifically) is that it's not obvious whether it refers to dot-product or cross-product. In fact, an inspection of the code reveals that in this case it's neither - it seems to be element-wise products. I'm not sure what the use case is for that... \$\endgroup\$ – Toby Speight Aug 9 '17 at 11:55
  • \$\begingroup\$ @kraskevich: You have some good points. About 1: I was indeed after some separation. The operators supported by the iterators are simply the bare minimum for automatic forward iteration (that's enough for my use case). 2: As they are internal names I didn't give it too much thought yet, you're right that this could be clearer. 3: Thanks for pointing that out, it's a bug, plain and simple. \$\endgroup\$ – Stefan Aug 9 '17 at 15:03
  • \$\begingroup\$ @TobySpeight: I was actually not sure which to use as a default for operator*, do you have a suggestion or is there a standard? The element-wise multiplication I do use, but it doesn't need to be tied to the * operator. \$\endgroup\$ – Stefan Aug 9 '17 at 15:07
  • \$\begingroup\$ @Stefan - I don't have a suggestion (that's why I wrote a comment rather than my own review!) Whichever is more appropriate for the context in which you use it, I guess - my recommendation is "be clear", whichever interpretation you choose. Element-wise certainly has an advantage in that it has an obvious inverse with operator/ (subject to numeric precision and divide-by-zero). \$\endgroup\$ – Toby Speight Aug 9 '17 at 16:11
  • \$\begingroup\$ Math libraries usually differentiate between arrays and vectors for their multiplication behavior. As functions I suggest multiply and dot. \$\endgroup\$ – Maikel Aug 10 '17 at 5:27

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.