IMP logo
IMP Reference Guide  2.20.0
The Integrative Modeling Platform
VectorBaseD.h
Go to the documentation of this file.
1 /**
2  * \file IMP/algebra/VectorBaseD.h \brief Simple D vector class.
3  *
4  * Copyright 2007-2022 IMP Inventors. All rights reserved.
5  *
6  */
7 
8 #ifndef IMPALGEBRA_VECTOR_BASE_D_H
9 #define IMPALGEBRA_VECTOR_BASE_D_H
10 
11 #include <IMP/algebra/algebra_config.h>
12 #include "GeometricPrimitiveD.h"
13 #include <IMP/check_macros.h>
14 #include <IMP/exception.h>
15 #include <IMP/random.h>
16 #include <IMP/utility.h>
17 #include <IMP/types.h>
18 #include <boost/random/variate_generator.hpp>
19 #include <boost/random/normal_distribution.hpp>
20 #include <boost/range.hpp>
21 #include <cereal/access.hpp>
22 #include "internal/vector.h"
23 
24 #include <limits>
25 #include <cmath>
26 #include <boost/random/normal_distribution.hpp>
27 #include <boost/static_assert.hpp>
28 
29 #if IMP_HAS_CHECKS >= IMP_INTERNAL
30 #define IMP_ALGEBRA_VECTOR_CHECK check_vector()
31 #else
32 #define IMP_ALGEBRA_VECTOR_CHECK
33 #endif
34 
35 #if IMP_HAS_CHECKS >= IMP_USAGE
36 #define IMP_ALGEBRA_VECTOR_CHECK_INDEX(i) check_index(i)
37 #define IMP_ALGEBRA_VECTOR_CHECK_COMPATIBLE(o) \
38  check_compatible_vector(o); \
39  o.check_vector()
40 #else
41 #define IMP_ALGEBRA_VECTOR_CHECK_INDEX(i)
42 #define IMP_ALGEBRA_VECTOR_CHECK_COMPATIBLE(o)
43 #endif
44 
45 IMPALGEBRA_BEGIN_NAMESPACE
46 //! A Cartesian vector in D-dimensions.
47 /** Store a vector of Cartesian coordinates. It supports all expected
48  mathematical operators, including using * for the dot product.
49  \see Vector3D
50  \see Vector2D
51 
52  \geometry
53  */
54 template <int D>
55 class VectorBaseD : public GeometricPrimitiveD<D> {
56  friend class cereal::access;
57 
58  template<class Archive>
59  void serialize(Archive &ar) {
60  ar(data_);
61  }
62 
63  void check_vector() const {
64  IMP_USAGE_CHECK(!data_.get_is_null(),
65  "Attempt to use uninitialized vector.");
66  }
67  template <int OD>
68  void check_compatible_vector(const VectorBaseD<OD> &o) const {
70  IMP_USAGE_CHECK(o.get_dimension() == get_dimension(),
71  "Dimensions don't match: " << get_dimension() << " vs "
72  << o.get_dimension());
73  }
74  void check_index(unsigned int i) const {
75 #if IMP_HAS_CHECKS < IMP_INTERNAL
76  IMP_UNUSED(i);
77 #endif
78  IMP_INTERNAL_CHECK(i < data_.get_dimension(),
79  "Invalid component of vector requested: "
80  << i << " of " << get_dimension());
81  }
82 
83  public:
84  /** Will accept a list of floats from Python. */
85  template <class Range>
86  explicit VectorBaseD(const Range &r) {
87  if (D != -1 && static_cast<int>(boost::distance(r)) != D) {
88  IMP_THROW("Expected " << D << " but got " << boost::distance(r),
90  }
92  for(double f : r) {
93  IMP_UNUSED(f);
94  IMP_USAGE_CHECK(!IMP::is_nan(f), "NaN passed to constructor");
95  }
96  }
97  data_.set_coordinates(boost::begin(r), boost::end(r));
98  }
99 
100 #ifndef SWIG
101  template <class R>
102  VectorBaseD<D> &operator=(const R &r) {
103  if (D != -1 && static_cast<int>(boost::distance(r)) != D) {
104  IMP_THROW("Expected " << D << " but got " << boost::distance(r),
106  }
108  for(double f : r) {
109  IMP_USAGE_CHECK(!IMP::is_nan(f), "NaN passed in equals");
110  }
111  }
112  data_.set_coordinates(boost::begin(r), boost::end(r));
113  }
114  //! Return the ith Cartesian coordinate.
115  /** In 3D use [0] to get the x coordinate etc. */
116  inline double operator[](unsigned int i) const {
117  IMP_ALGEBRA_VECTOR_CHECK_INDEX(i);
118  IMP_ALGEBRA_VECTOR_CHECK;
119  return data_.get_data()[i];
120  }
121  //! Return the ith Cartesian coordinate.
122  /** In 3D use [0] to get the x coordinate etc. */
123  inline double &operator[](unsigned int i) {
124  IMP_ALGEBRA_VECTOR_CHECK_INDEX(i);
125  return data_.get_data()[i];
126  }
127 
128 #endif
129  //! Default constructor
131 
132  double get_scalar_product(const VectorBaseD<D> &o) const {
133  IMP_ALGEBRA_VECTOR_CHECK_COMPATIBLE(o);
134  IMP_ALGEBRA_VECTOR_CHECK;
135  double ret = 0;
136  for (unsigned int i = 0; i < get_dimension(); ++i) {
137  ret += operator[](i) * o.operator[](i);
138  }
139  return ret;
140  }
141 
142  double get_squared_magnitude() const {
143  // Could be equivalently written as get_scalar_product(*this)
144  // This should be faster though since checks on 'o' can be skipped
145  IMP_ALGEBRA_VECTOR_CHECK;
146  double ret = 0;
147  const double *data = get_data();
148  for (unsigned int i = 0; i < get_dimension(); ++i) {
149  ret += data[i] * data[i];
150  }
151  return ret;
152  }
153 
154  double get_magnitude() const { return std::sqrt(get_squared_magnitude()); }
155 
156 #ifndef IMP_DOXYGEN
157  double operator*(const VectorBaseD<D> &o) const {
158  return get_scalar_product(o);
159  }
160 
161  VectorBaseD &operator+=(const VectorBaseD &o) {
162  IMP_ALGEBRA_VECTOR_CHECK_COMPATIBLE(o);
163  IMP_ALGEBRA_VECTOR_CHECK;
164  for (unsigned int i = 0; i < get_dimension(); ++i) {
165  operator[](i) += o[i];
166  }
167  return *this;
168  }
169 
170  VectorBaseD &operator-=(const VectorBaseD &o) {
171  IMP_ALGEBRA_VECTOR_CHECK_COMPATIBLE(o);
172  IMP_ALGEBRA_VECTOR_CHECK;
173  for (unsigned int i = 0; i < get_dimension(); ++i) {
174  operator[](i) -= o[i];
175  }
176  return *this;
177  }
178 
179  VectorBaseD &operator/=(double f) {
180  IMP_ALGEBRA_VECTOR_CHECK;
181  for (unsigned int i = 0; i < get_dimension(); ++i) {
182  operator[](i) /= f;
183  }
184  return *this;
185  }
186 
187  VectorBaseD &operator*=(double f) {
188  IMP_ALGEBRA_VECTOR_CHECK;
189  for (unsigned int i = 0; i < get_dimension(); ++i) {
190  operator[](i) *= f;
191  }
192  return *this;
193  }
194 
195  void show(std::ostream &out, std::string delim, bool parens = true) const {
196  IMP_ALGEBRA_VECTOR_CHECK;
197  if (parens) out << "(";
198  for (unsigned int i = 0; i < get_dimension(); ++i) {
199  out << operator[](i);
200  if (i != get_dimension() - 1) {
201  out << delim;
202  }
203  }
204  if (parens) out << ")";
205  }
206  IMP_SHOWABLE_INLINE(VectorBaseD, show(out, ", "););
207 #endif
208 
209 #ifndef SWIG
210  typedef double *iterator;
211  typedef const double *const_iterator;
212  iterator begin() { return data_.get_data(); }
213  iterator end() { return data_.get_data() + get_dimension(); }
214  const_iterator begin() const { return data_.get_data(); }
215  const_iterator end() const { return data_.get_data() + get_dimension(); }
216 
217  typedef double value_type;
218  typedef std::random_access_iterator_tag iterator_category;
219  typedef std::ptrdiff_t difference_type;
220  typedef double *pointer;
221  typedef double &reference;
222  typedef const double &const_reference;
223 
224  static const int DIMENSION = D;
225 #endif
226 
227 #ifndef SWIG
228  // For some reason, this method breaks IMP::atom::get_rmsd() in Python, so
229  // hide it from SWIG
230  Floats get_coordinates() const {
231  return Floats(begin(), end());
232  }
233 
234  //! Return a pointer to the data stored.
235  /** Useful for conversion to other types. */
236  const double *get_data() const { return data_.get_data(); }
237 #endif
238  unsigned int get_dimension() const { return data_.get_dimension(); }
239 
240  private:
241  internal::VectorData<double, D, false> data_;
242 };
243 
244 //! Returns a unit vector pointing at the same direction as this vector.
245 /**
246  @note If the magnitude of this vector is smaller than 1e-12
247  (an arbitrarily selected small number), returns a unit
248  vector pointing at a random direction.
249  */
250 template <class VT>
251 inline VT get_unit_vector(VT vt) {
252  static const double tiny_double =
253  256.0 * std::numeric_limits<double>::epsilon();
254  double mag = vt.get_magnitude();
255  if (mag > tiny_double) {
256  VT ret_value= vt/mag;
257  IMP_USAGE_CHECK(std::abs(ret_value.get_magnitude() - 1.0) < 256.0 * tiny_double,
258  "returned vector is not unit vector");
259  return ret_value;
260  } else {
261  // avoid division by zero - return random unit v
262  // NOTE: (1) avoids vector_generators / SphereD to prevent recursiveness
263  // (2) D might be -1, so use get_dimension()
264  static boost::variate_generator<RandomNumberGenerator,
265  boost::normal_distribution<> >
267  ::boost::normal_distribution<>(0, 1.0));
268  for (unsigned int i = 0; i < vt.get_dimension(); ++i) {
269  vt[i] = generator();
270  }
271  return get_unit_vector(vt);
272  }
273 }
274 
275 //! Returns the magnitude of vt and turns it to a unit vector in place.
276 /**
277  @note If the magnitude of this vector is smaller than 1e-12
278  (an arbitrarily selected small number), vt is turned into a unit
279  vector pointing at a random direction.
280  */
281 template <class VT>
283  const double tiny_double = 1e-12;
284  double mag = vt.get_magnitude();
285  if (mag > tiny_double) {
286  vt /= mag;
287  } else {
288  // avoid division by zero - return random unit v
289  // using get_unit_vector()
290  vt= get_unit_vector(vt);
291  }
292  return mag;
293 }
294 
295 IMPALGEBRA_END_NAMESPACE
296 
297 #endif /* IMPALGEBRA_VECTOR_BASE_D_H */
Base class for geometric types.
Basic types used by IMP.
#define IMP_USAGE_CHECK_VARIABLE(variable)
Definition: check_macros.h:190
#define IMP_SHOWABLE_INLINE(Name, how_to_show)
Declare the methods needed by an object that can be printed.
#define IMP_IF_CHECK(level)
Execute the code block if a certain level checks are on.
Definition: check_macros.h:104
VectorBaseD()
Default constructor.
Definition: VectorBaseD.h:130
IMP::Vector< Float > Floats
Standard way to pass a bunch of Float values.
Definition: types.h:46
VectorBaseD(const Range &r)
Definition: VectorBaseD.h:86
VT get_unit_vector(VT vt)
Returns a unit vector pointing at the same direction as this vector.
Definition: VectorBaseD.h:251
Exception definitions and assertions.
#define IMP_INTERNAL_CHECK(expr, message)
An assertion to check for internal errors in IMP. An IMP::ErrorException will be thrown.
Definition: check_macros.h:139
double get_magnitude_and_normalize_in_place(VT &vt)
Returns the magnitude of vt and turns it to a unit vector in place.
Definition: VectorBaseD.h:282
Base class for geometric types.
#define IMP_UNUSED(variable)
Various general useful functions for IMP.
std::ostream & show(Hierarchy h, std::ostream &out=std::cout)
Print the hierarchy using a given decorator to display each node.
#define IMP_THROW(message, exception_name)
Throw an exception with a message.
Definition: check_macros.h:50
Helper macros for throwing and handling exceptions.
#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
Random number generators used by IMP.
const double * get_data() const
Return a pointer to the data stored.
Definition: VectorBaseD.h:236
A Cartesian vector in D-dimensions.
Definition: VectorBaseD.h:55
An exception for an invalid value being passed to IMP.
Definition: exception.h:136
double operator[](unsigned int i) const
Return the ith Cartesian coordinate.
Definition: VectorBaseD.h:116
RandomNumberGenerator random_number_generator
A shared non-GPU random number generator.
double & operator[](unsigned int i)
Return the ith Cartesian coordinate.
Definition: VectorBaseD.h:123