IMP logo
IMP Reference Guide  2.19.0
The Integrative Modeling Platform
VectorD.h
Go to the documentation of this file.
1 /**
2  * \file IMP/algebra/VectorD.h \brief Simple D vector class.
3  *
4  * Copyright 2007-2022 IMP Inventors. All rights reserved.
5  *
6  */
7 
8 #ifndef IMPALGEBRA_VECTOR_D_H
9 #define IMPALGEBRA_VECTOR_D_H
10 
11 #include <IMP/algebra/algebra_config.h>
12 #include "VectorBaseD.h"
13 #include <IMP/types.h>
14 #include <IMP/check_macros.h>
15 #include <IMP/exception.h>
16 #include <IMP/utility.h>
17 #include <IMP/InputAdaptor.h>
18 #include <IMP/random.h>
19 #include "algebra_macros.h"
20 #include <boost/random/variate_generator.hpp>
21 #include <boost/random/normal_distribution.hpp>
22 #include <boost/range.hpp>
23 #include "internal/vector.h"
24 
25 #include <limits>
26 #include <cmath>
27 #include <boost/random/normal_distribution.hpp>
28 #include <boost/static_assert.hpp>
29 
30 #if IMP_HAS_CHECKS >= IMP_INTERNAL
31 #define IMP_ALGEBRA_VECTOR_CHECK check_vector()
32 #else
33 #define IMP_ALGEBRA_VECTOR_CHECK
34 #endif
35 
36 #if IMP_HAS_CHECKS >= IMP_USAGE
37 #define IMP_ALGEBRA_VECTOR_CHECK_INDEX(i) check_index(i)
38 #define IMP_ALGEBRA_VECTOR_CHECK_COMPATIBLE(o) \
39  check_compatible_vector(o); \
40  o.check_vector()
41 #else
42 #define IMP_ALGEBRA_VECTOR_CHECK_INDEX(i)
43 #define IMP_ALGEBRA_VECTOR_CHECK_COMPATIBLE(o)
44 #endif
45 
46 IMPALGEBRA_BEGIN_NAMESPACE
47 //! A Cartesian vector in D-dimensions.
48 /** Store a vector of Cartesian coordinates. It supports all expected
49  mathematical operators, including using * for the dot product.
50  \see Vector3D
51  \see Vector2D
52 
53  \geometry
54  */
55 template <int D>
56 class VectorD : public VectorBaseD<D> {
57  /* implementing it via a specialization is in order to get swig
58  to only use the right constructors as well as C++. */
59  IMP_ALGEBRA_VECTOR_METHODS(D);
60 };
61 
62 template <>
63 class VectorD<-1> : public VectorBaseD<-1> {
64  std::vector<double> get_range(double x0, double x1, double x2, double x3,
65  double x4, double x5) {
66  IMP_USAGE_CHECK(x0 != std::numeric_limits<double>::max(),
67  "Bad init values");
68  std::vector<double> ret;
69  ret.push_back(x0);
70  if (x1 != std::numeric_limits<double>::max()) ret.push_back(x1);
71  if (x2 != std::numeric_limits<double>::max()) ret.push_back(x2);
72  if (x3 != std::numeric_limits<double>::max()) ret.push_back(x3);
73  if (x4 != std::numeric_limits<double>::max()) ret.push_back(x4);
74  if (x5 != std::numeric_limits<double>::max()) ret.push_back(x5);
75  return ret;
76  }
77 
78  public:
79  IMP_ALGEBRA_VECTOR_METHODS(-1);
80  explicit VectorD(double x0, double x1 = std::numeric_limits<double>::max(),
81  double x2 = std::numeric_limits<double>::max(),
82  double x3 = std::numeric_limits<double>::max(),
83  double x4 = std::numeric_limits<double>::max(),
84  double x5 = std::numeric_limits<double>::max())
85  : P(get_range(x0, x1, x2, x3, x4, x5)) {}
86 };
87 
88 template <>
89 class VectorD<1> : public VectorBaseD<1> {
90  public:
91  IMP_ALGEBRA_VECTOR_METHODS(1);
92 
93  //! Initialize the 1-vector from its value.
94  VectorD(double x) { P::operator[](0) = x; }
95 };
96 template <>
97 class VectorD<2> : public VectorBaseD<2> {
98  public:
99  IMP_ALGEBRA_VECTOR_METHODS(2);
100 
101  VectorD(double x, double y) {
102  P::operator[](0) = x;
103  P::operator[](1) = y;
104  }
105 };
106 template <>
107 class VectorD<3> : public VectorBaseD<3> {
108  public:
109  IMP_ALGEBRA_VECTOR_METHODS(3);
110 
111  VectorD(double x, double y, double z) {
112  P::operator[](0) = x;
113  P::operator[](1) = y;
114  P::operator[](2) = z;
115  }
116 };
117 template <>
118 class VectorD<4> : public VectorBaseD<4> {
119  public:
120  IMP_ALGEBRA_VECTOR_METHODS(4);
121 
122  VectorD(double x0, double x1, double x2, double x3) {
123  P::operator[](0) = x0;
124  P::operator[](1) = x1;
125  P::operator[](2) = x2;
126  P::operator[](3) = x3;
127  }
128 };
129 template <>
130 class VectorD<5> : public VectorBaseD<5> {
131  public:
132  IMP_ALGEBRA_VECTOR_METHODS(5);
133 
134  VectorD(double x0, double x1, double x2, double x3, double x4) {
135  P::operator[](0) = x0;
136  P::operator[](1) = x1;
137  P::operator[](2) = x2;
138  P::operator[](3) = x3;
139  P::operator[](4) = x4;
140  }
141 };
142 template <>
143 class VectorD<6> : public VectorBaseD<6> {
144  public:
145  IMP_ALGEBRA_VECTOR_METHODS(6);
146 
147  VectorD(double x0, double x1, double x2, double x3, double x4, double x5) {
148  P::operator[](0) = x0;
149  P::operator[](1) = x1;
150  P::operator[](2) = x2;
151  P::operator[](3) = x3;
152  P::operator[](4) = x4;
153  P::operator[](5) = x5;
154  }
155 };
156 
157 #ifndef IMP_DOXYGEN
158 
159 template <int D>
160 inline std::ostream &operator<<(std::ostream &out, const VectorD<D> &v) {
161  v.show(out);
162  return out;
163 }
164 
165 #ifndef SWIG
166 // SWIG 3.0.0 cannot parse operator>>
167 template <int D>
168 inline std::istream &operator>>(std::istream &in, VectorD<D> &v) {
169  for (unsigned int i = 0; i < D; ++i) {
170  in >> v[i];
171  }
172  return in;
173 }
174 #endif
175 
176 #endif
177 
178 //! lexicographic comparison of two vectors
179 /** Note that this is not very reliable and probably should not be used.
180  \see VectorD
181  */
182 template <int D>
183 inline int compare(const VectorD<D> &a, const VectorD<D> &b) {
184  IMP_USAGE_CHECK(a.get_dimension() == b.get_dimension(),
185  "Dimensions don't match.");
186  for (unsigned int i = 0; i < a.get_dimension(); ++i) {
187  if (a[i] < b[i])
188  return -1;
189  else if (a[i] > b[i])
190  return 1;
191  }
192  return 0;
193 }
194 
195 /** \see VectorD */
196 template <int D>
197 inline VectorD<D> operator*(double s, VectorD<D> o) {
198  return o *= s;
199 }
200 
201 //! Compute the squared distance between two vectors
202 /** \see VectorD
203  */
204 template <int D>
205 inline double get_squared_distance(const VectorD<D> &v1, const VectorD<D> &v2) {
206  return (v1 - v2).get_squared_magnitude();
207 }
208 
209 //! Compute the distance between two vectors
210 /** \see VectorD
211  */
212 template <int D>
213 inline double get_distance(const VectorD<D> &v1, const VectorD<D> &v2) {
214  return std::sqrt(get_squared_distance(v1, v2));
215 }
216 
217 //! Return the basis vector for the given coordinate
218 /** Return the unit vector pointing in the direction of the requested
219  coordinate. That is
220  \code
221  get_basis_vector_d<3>(2)== Vector3D(0,0,1);
222  \endcode
223  \see VectorD
224  */
225 template <int D>
226 inline VectorD<D> get_basis_vector_d(unsigned int coordinate) {
227  IMP_USAGE_CHECK(coordinate < D, "There are only " << D << " basis vectors");
228  double vs[D];
229  for (unsigned int i = 0; i < D; ++i) {
230  if (i == coordinate)
231  vs[i] = 1;
232  else
233  vs[i] = 0;
234  }
235  return VectorD<D>(vs, vs + D);
236 }
237 
238 //! Return the basis vector for the given coordinate
239 template <int D>
240 inline VectorD<D> get_basis_vector_kd(int Di, unsigned int coordinate) {
241  IMP_USAGE_CHECK(D == Di, "D must be equal");
242  IMP_UNUSED(Di);
243  return get_basis_vector_d<D>(coordinate);
244 }
245 
246 //! Return a dynamically sized basis vector
247 template <>
248 inline VectorD<-1> get_basis_vector_kd(int D, unsigned int coordinate) {
249  IMP_USAGE_CHECK(D > 0, "D must be positive");
250  IMP_USAGE_CHECK(coordinate < static_cast<unsigned int>(D),
251  "There are only " << D << " basis vectors");
252  boost::scoped_array<double> vs(new double[D]);
253  for (int i = 0; i < D; ++i) {
254  if (i == static_cast<int>(coordinate))
255  vs[i] = 1;
256  else
257  vs[i] = 0;
258  }
259  return VectorD<-1>(vs.get(), vs.get() + D);
260 }
261 
262 //! Return a dynamically sized basis vector
263 inline VectorD<-1> get_basis_vector_kd(int D, unsigned int coordinate) {
264  return get_basis_vector_kd<-1>(D, coordinate);
265 }
266 
267 //! Return a vector of zeros
268 template <int D>
270  IMP_USAGE_CHECK(D > 0, "D must be positive");
271  VectorD<D> ret;
272  for (int i = 0; i < D; ++i) {
273  ret[i] = 0;
274  }
275  return ret;
276 }
277 
278 //! Return a dynamically sized vector of zeros
279 template <int D>
281  IMP_USAGE_CHECK(D == Di, "D must be equal");
282  IMP_UNUSED(Di);
283  return get_zero_vector_d<D>();
284 }
285 
286 //! Return a dynamically sized vector of zeros
287 template<>
288 inline VectorD<-1> get_zero_vector_kd(int D) {
289  IMP_USAGE_CHECK(D > 0, "D must be positive");
290  Floats vs(D, 0);
291  return VectorD<-1>(vs.begin(), vs.end());
292 }
293 
294 //! Return a dynamically sized vector of zeros
295 inline VectorD<-1> get_zero_vector_kd(int D) {
296  return get_zero_vector_kd<-1>(D);
297 }
298 
299 //! Return a vector of ones (or another constant)
300 template <int D>
301 inline VectorD<D> get_ones_vector_d(double v = 1) {
302  IMP_USAGE_CHECK(D > 0, "D must be positive");
303  VectorD<D> ret;
304  for (unsigned int i = 0; i < D; ++i) {
305  ret[i] = v;
306  }
307  return ret;
308 }
309 
310 //! Return a vector of ones (or another constant)
311 /** Di must equal D. */
312 template <int D>
313 inline VectorD<D> get_ones_vector_kd(unsigned int Di, double v = 1) {
314  IMP_USAGE_CHECK(D == Di, "D must be equal");
315  IMP_UNUSED(Di);
316  return get_ones_vector_d<D>(v);
317 }
318 
319 //! Return a vector of ones (or another constant)
320 template <>
321 inline VectorD<-1> get_ones_vector_kd(unsigned int D, double v) {
322  IMP_USAGE_CHECK(D > 0, "D must be positive");
323  boost::scoped_array<double> vv(new double[D]);
324  for (unsigned int i = 0; i < D; ++i) {
325  vv[i] = v;
326  }
327  return VectorD<-1>(vv.get(), vv.get() + D);
328 }
329 
330 //! Return a dynamically sized vector of zeros
331 inline VectorD<-1> get_ones_vector_kd(unsigned int D, double v = 1) {
332  return get_ones_vector_kd<-1>(D, v);
333 }
334 
335 #ifndef SWIG
336 
337 /** \name Norms
338  We define a number of standard, \f$L^p\f$, norms on VectorD.
339  - \f$L^1\f$ is the Manhattan distance, the sum of the components
340  - \f$L^2\f$ is the standard Euclidean length
341  - \f$L^{\inf}\f$ is the maximum of the components
342  @{
343 */
344 
345 template <int D>
346 inline double get_l2_norm(const VectorD<D> &v) {
347  return v.get_magnitude();
348 }
349 
350 template <int D>
351 inline double get_l1_norm(const VectorD<D> &v) {
352  double n = std::abs(v[0]);
353  for (unsigned int i = 1; i < v.get_dimension(); ++i) {
354  n += std::abs(v[i]);
355  }
356  return n;
357 }
358 
359 template <int D>
360 inline double get_linf_norm(const VectorD<D> &v) {
361  double n = std::abs(v[0]);
362  for (unsigned int i = 1; i < v.get_dimension(); ++i) {
363  n = std::max(n, std::abs(v[i]));
364  }
365  return n;
366 }
367 
368 /** @} */
369 
370 #ifndef IMP_DOXYGEN
371 
372 template <int D>
373 struct SpacesIO {
374  const VectorD<D> &v_;
375  SpacesIO(const VectorD<D> &v) : v_(v) {}
376 };
377 
378 template <int D>
379 struct CommasIO {
380  const VectorD<D> &v_;
381  CommasIO(const VectorD<D> &v) : v_(v) {}
382 };
383 template <int D>
384 inline std::ostream &operator<<(std::ostream &out, const SpacesIO<D> &s) {
385  s.v_.show(out, " ", false);
386  return out;
387 }
388 template <int D>
389 inline std::ostream &operator<<(std::ostream &out, const CommasIO<D> &s) {
390  s.v_.show(out, ", ", false);
391  return out;
392 }
393 
394 //! Use this before outputting to delimited vector entries with a space
395 /** std::cout << spaces_io(v);
396  produces "1.0 2.0 3.0"
397  \see VectorD
398  */
399 template <int D>
400 inline SpacesIO<D> spaces_io(const VectorD<D> &v) {
401  return SpacesIO<D>(v);
402 }
403 
404 //! Use this before outputting to delimited vector entries with a comma
405 /** std::cout << commas_io(v);
406  produces "1.0, 2.0, 3.0"
407  \see VectorD
408  */
409 template <int D>
410 inline CommasIO<D> commas_io(const VectorD<D> &v) {
411  return CommasIO<D>(v);
412 }
413 #endif // doxygen
414 
415 #endif // swig
416 /** 1D vector typedef for swig */
418 /** 1D vectors typedef for swig */
420 /** 2D vector typedef for swig */
422 /** 2D vectors typedef for swig */
424 /** 3D vector typedef for swig */
426 /** 3D vectors typedef for swig */
428 /** 4D vector typedef for swig */
430 /** 4D vectors typedef for swig */
432 /** 5D vector typedef for swig */
434 /** 5D vectors typedef for swig */
436 /** 6D vector typedef for swig */
438 /** 6D vector typedef for swig */
440 /** KD vector typedef for swig */
441 typedef VectorD<-1> VectorKD;
442 /** KD vectors typedef for swig */
443 typedef Vector<VectorD<-1> > VectorKDs;
444 
445 #ifndef SWIG
446 /** \see VectorD \genericgeometry */
447 template <class C>
448 inline const VectorD<C::DIMENSION> &get_vector_geometry(const C &g) {
449  return g;
450 }
451 /** \see VectorD \genericgeometry */
452 template <class C, class E>
453 inline void set_vector_geometry(C &g, const E &v) {
454  g = v;
455 }
456 #endif
457 
458 //! Return the vector that is the elementwise product of the two.
459 /** \see VectorD */
460 template <int D>
462  const algebra::VectorD<D> &b) {
463  VectorD<D> ret(a);
464  for (unsigned int i = 0; i < ret.get_dimension(); ++i) {
465  ret[i] *= b[i];
466  }
467  return ret;
468 }
469 
470 //! Return the vector that is the elementwise product of the two.
471 /** \see VectorD */
472 template <int D>
474  const algebra::VectorD<D> &b) {
475  IMP_USAGE_CHECK(a.size() == b.get_dimension(), "Dimensions don't match,");
476  VectorD<D> ret(b);
477  for (unsigned int i = 0; i < ret.get_dimension(); ++i) {
478  ret[i] *= a[i];
479  }
480  return ret;
481 }
482 
483 IMPALGEBRA_END_NAMESPACE
484 
485 #endif /* IMPALGEBRA_VECTOR_D_H */
Basic types used by IMP.
VectorD< 6 > Vector6D
Definition: VectorD.h:437
VectorD< 2 > Vector2D
Definition: VectorD.h:421
VectorD< 5 > Vector5D
Definition: VectorD.h:433
Vector< VectorD< 1 > > Vector1Ds
Definition: VectorD.h:419
Vector< VectorD< 3 > > Vector3Ds
Definition: VectorD.h:427
double get_squared_distance(const VectorD< D > &v1, const VectorD< D > &v2)
Compute the squared distance between two vectors.
Definition: VectorD.h:205
Vector< VectorD< 5 > > Vector5Ds
Definition: VectorD.h:435
Exception definitions and assertions.
VectorD< 4 > Vector4D
Definition: VectorD.h:429
Vector< VectorD< 2 > > Vector2Ds
Definition: VectorD.h:423
const VectorD< C::DIMENSION > & get_vector_geometry(const C &g)
Definition: VectorD.h:448
Convenience class to accept multiple input types.
A Cartesian vector in D-dimensions.
Definition: VectorD.h:56
#define IMP_UNUSED(variable)
VectorD< D > get_elementwise_product(const Ints &a, const algebra::VectorD< D > &b)
Return the vector that is the elementwise product of the two.
Definition: VectorD.h:473
int compare(const VectorD< D > &a, const VectorD< D > &b)
lexicographic comparison of two vectors
Definition: VectorD.h:183
VectorD< D > get_zero_vector_d()
Return a vector of zeros.
Definition: VectorD.h:269
Vector< VectorD<-1 > > VectorKDs
Definition: VectorD.h:443
VectorD< D > get_basis_vector_d(unsigned int coordinate)
Return the basis vector for the given coordinate.
Definition: VectorD.h:226
VectorD< 1 > Vector1D
Definition: VectorD.h:417
void set_vector_geometry(C &g, const E &v)
Definition: VectorD.h:453
Vector< VectorD< 6 > > Vector6Ds
Definition: VectorD.h:439
For backwards compatibility.
Vector< VectorD< 4 > > Vector4Ds
Definition: VectorD.h:431
VectorD<-1 > get_ones_vector_kd(unsigned int D, double v)
Return a vector of ones (or another constant)
Definition: VectorD.h:321
Simple D vector class.
VectorD< D > get_ones_vector_d(double v=1)
Return a vector of ones (or another constant)
Definition: VectorD.h:301
Helper macros for throwing and handling exceptions.
VectorD< 3 > Vector3D
Definition: VectorD.h:425
VectorD<-1 > get_basis_vector_kd(int D, unsigned int coordinate)
Return a dynamically sized basis vector.
Definition: VectorD.h:248
#define IMP_USAGE_CHECK(expr, message)
A runtime test for incorrect usage of a class or method.
Definition: check_macros.h:168
VectorD< D > operator*(double s, VectorD< D > o)
Definition: VectorD.h:197
double get_distance(const VectorD< D > &v1, const VectorD< D > &v2)
Compute the distance between two vectors.
Definition: VectorD.h:213
Random number generators used by IMP.
A Cartesian vector in D-dimensions.
Definition: VectorBaseD.h:55
VectorD<-1 > VectorKD
Definition: VectorD.h:441
Various helper macros.
VectorD<-1 > get_zero_vector_kd(int D)
Return a dynamically sized vector of zeros.
Definition: VectorD.h:288