Rivet  1.8.3
VectorN.hh
1 #ifndef RIVET_MATH_VECTORN
2 #define RIVET_MATH_VECTORN
3 
4 #include "Rivet/Math/MathHeader.hh"
5 #include "Rivet/Math/MathUtils.hh"
6 
7 #include "Rivet/Math/eigen/vector.h"
8 
9 namespace Rivet {
10 
11 
12  template <size_t N>
13  class Vector;
14  template <size_t N>
15  class Matrix;
16 
17  template <size_t N>
18  Vector<N> multiply(const Matrix<N>& a, const Vector<N>& b);
19 
20 
22  template <size_t N>
23  class Vector {
24  template <size_t M>
25  friend Vector<M> multiply(const Matrix<M>& a, const Vector<M>& b);
26 
27  public:
28  Vector() { _vec.loadZero(); }
29 
30  Vector(const Vector<N>& other)
31  : _vec(other._vec) { }
32 
33  const double& get(const size_t index) const {
34  if (index >= N) {
35  throw std::runtime_error("Tried to access an invalid vector index.");
36  } else {
37  return _vec(index);
38  }
39  }
40 
42  const double& operator[](const size_t index) const {
43  return get(index);
44  }
45 
47  double& operator[](const size_t index) {
48  return get(index);
49  }
50 
52  Vector<N>& set(const size_t index, const double value) {
53  if (index >= N) {
54  throw std::runtime_error("Tried to access an invalid vector index.");
55  } else {
56  _vec[index] = value;
57  }
58  return *this;
59  }
60 
62  size_t size() const {
63  return N;
64  }
65 
67  bool isZero(double tolerance=1E-5) const {
68  for (size_t i=0; i < N; ++i) {
69  if (! Rivet::isZero(_vec[i], tolerance) ) return false;
70  }
71  return true;
72  }
73 
76  double mod2() const {
77  double mod2 = 0.0;
78  for (size_t i = 0; i < size(); ++i) {
79  const double element = get(i);
80  mod2 += element*element;
81  }
82  return mod2;
83  }
84 
87  double mod() const {
88  const double norm = mod2();
89  assert(norm >= 0);
90  return sqrt(norm);
91  }
92 
94  Vector<N> operator-() const {
95  Vector<N> rtn;
96  rtn._vec = -_vec;
97  return rtn;
98  }
99 
100  bool operator==(const Vector<N>& a) const {
101  return _vec == a._vec;
102  }
103 
104  bool operator!=(const Vector<N>& a) const {
105  return _vec != a._vec;
106  }
107 
108  bool operator<(const Vector<N>& a) const {
109  return _vec < a._vec;
110  }
111 
112  bool operator<=(const Vector<N>& a) const {
113  return _vec <= a._vec;
114  }
115 
116  bool operator>(const Vector<N>& a) const {
117  return _vec > a._vec;
118  }
119 
120  bool operator>=(const Vector<N>& a) const {
121  return _vec >= a._vec;
122  }
123 
124 
125  protected:
126  double& get(const size_t index) {
127  if (index >= N) {
128  throw std::runtime_error("Tried to access an invalid vector index.");
129  } else {
130  return _vec(index);
131  }
132  }
133 
135  Eigen::Vector<double,N> _vec;
136  };
137 
138 
139 
142  template <size_t N>
143  inline double mod2(const Vector<N>& v) {
144  return v.mod2();
145  }
146 
149  template <size_t N>
150  inline double mod(const Vector<N>& v) {
151  return v.mod();
152  }
153 
154 
156 
157 
159 
160 
162  template <size_t N>
163  inline const string toString(const Vector<N>& v) {
164  ostringstream out;
165  out << "(";
166  for (size_t i = 0; i < v.size(); ++i) {
167  out << (fabs(v[i]) < 1E-30 ? 0.0 : v[i]);
168  if (i < v.size()-1) out << ", ";
169  }
170  out << ")";
171  return out.str();
172  }
173 
175  template <size_t N>
176  inline std::ostream& operator<<(std::ostream& out, const Vector<N>& v) {
177  out << toString(v);
178  return out;
179  }
180 
182 
183 
185 
186 
188  template <size_t N>
189  inline bool fuzzyEquals(const Vector<N>& va, const Vector<N>& vb, double tolerance=1E-5) {
190  for (size_t i = 0; i < N; ++i) {
191  const double a = va.get(i);
192  const double b = vb.get(i);
193  if (!Rivet::fuzzyEquals(a, b, tolerance)) return false;
194  }
195  return true;
196  }
197 
198 
200  template <size_t N>
201  inline bool isZero(const Vector<N>& v, double tolerance=1E-5) {
202  return v.isZero(tolerance);
203  }
204 
205 
206 }
207 
208 #endif
Definition: MC_JetAnalysis.hh:9
bool isZero(double tolerance=1E-5) const
Check for nullness, allowing for numerical precision.
Definition: VectorN.hh:67
double mod2() const
Calculate the modulus-squared of a vector. .
Definition: VectorN.hh:76
double & operator[](const size_t index)
Direct access to vector elements by index.
Definition: VectorN.hh:47
General -dimensional mathematical matrix object.
Definition: MatrixN.hh:14
bool isZero(double val, double tolerance=1E-8)
Definition: MathUtils.hh:17
A minimal base class for -dimensional vectors.
Definition: VectorN.hh:13
size_t size() const
Vector dimensionality.
Definition: VectorN.hh:62
std::string toString(const AnalysisInfo &ai)
String representation.
Definition: AnalysisInfo.cc:234
const double & operator[](const size_t index) const
Direct access to vector elements by index.
Definition: VectorN.hh:42
Vector< N > operator-() const
Invert the vector.
Definition: VectorN.hh:94
bool fuzzyEquals(double a, double b, double tolerance=1E-5)
Compare two floating point numbers for equality with a degree of fuzziness.
Definition: MathUtils.hh:34
double mod() const
Calculate the modulus of a vector. .
Definition: VectorN.hh:87