1 files changed, 113 insertions, 0 deletions
diff --git a/lib/wolff_models/vector.hpp b/lib/wolff_models/vector.hpp
new file mode 100644
index 0000000..4106e62
--- /dev/null
+++ b/lib/wolff_models/vector.hpp
@@ -0,0 +1,113 @@
+
+#ifndef WOLFF_MODELS_VECTOR_H
+#define WOLFF_MODELS_VECTOR_H
+
+#include <cmath>
+#include <array>
+#include <iostream>
+
+namespace wolff {
+
+  template <unsigned q, class T>
+  class vector_t : public std::array<T, q> {
+    public: 
+
+      vector_t() {
+        this->fill((T)0);
+        (*this)[0] = (T)1;
+      }
+
+      vector_t(const T *x) {
+        for (unsigned i = 0; i < q; i++) {
+          (*this)[i] = x[i];
+        }
+      }
+
+      typedef vector_t <q, T> M_t;
+      typedef vector_t <q, double> F_t;
+
+      template <class U>
+      inline vector_t<q, T>& operator+=(const vector_t<q, U> &v) {
+        for (unsigned i = 0; i < q; i++) {
+          (*this)[i] += (T)v[i];
+        }
+        return *this;
+      }
+
+      template <class U>
+      inline vector_t<q, T>& operator-=(const vector_t<q, U> &v) {
+        for (unsigned i = 0; i < q; i++) {
+          (*this)[i] -= (T)v[i];
+        }
+        return *this;
+      }
+
+      inline vector_t<q, T> operator*(unsigned x) const {
+        vector_t<q, T> result;
+        for (unsigned i = 0; i < q; i++) {
+          result[i] = x * (*this)[i];
+        }
+
+        return result;
+      }
+
+      inline vector_t<q, double> operator*(double x) const {
+        vector_t<q, double> result;
+        for (unsigned i = 0; i < q; i++) {
+          result[i] = x * (*this)[i];
+        }
+
+        return result;
+      }
+
+      inline vector_t<q, T> operator-(const vector_t<q, T>& v) const {
+        vector_t<q, T> diff = *this;
+        diff -= v;
+        return diff;
+      }
+
+      inline T operator*(const vector_t<q, T>& v) const {
+        double prod = 0;
+
+        for (unsigned i = 0; i < q; i++) {
+          prod += v[i] * (*this)[i];
+        }
+
+        return prod;
+      }
+
+      template <class U>
+      inline vector_t<q, T> operator/(U a) const {
+        vector_t<q, T> result;
+        for (unsigned i = 0; i < q; i++) {
+          result[i] = (*this)[i] / a;
+        }
+
+        return result;
+      }
+  };
+
+  template<unsigned q, class T>
+  inline vector_t<q, T> operator*(unsigned a, const vector_t<q, T>&v) {
+    return v * a;
+  }
+
+  template<unsigned q, class T>
+  inline vector_t<q, double> operator*(double a, const vector_t<q, T>&v) {
+    return v * a;
+  }
+
+  template<unsigned q, class T>
+  std::ostream& operator<<(std::ostream& os, const vector_t<q, T>&v) {
+    os << "( ";
+    for (T vi : v) {
+      os << vi << " ";
+    }
+    os << ")";
+    return os;
+  }
+
+}
+
+#endif
+