IMP logo
IMP Reference Guide  2.11.0
The Integrative Modeling Platform
rigid_bodies.h
Go to the documentation of this file.
1 /**
2  * \file IMP/core/rigid_bodies.h
3  * \brief functionality for defining rigid bodies
4  *
5  * Copyright 2007-2019 IMP Inventors. All rights reserved.
6  */
7 
8 #ifndef IMPCORE_RIGID_BODIES_H
9 #define IMPCORE_RIGID_BODIES_H
10 
11 #include <IMP/core/core_config.h>
12 #include "internal/rigid_bodies.h"
13 
14 #include "XYZ.h"
15 #include "XYZR.h"
16 #include <IMP/SingletonContainer.h>
17 #include <IMP/SingletonModifier.h>
18 #include <IMP/Refiner.h>
19 #include <IMP/algebra/Vector3D.h>
20 #include <IMP/algebra/Rotation3D.h>
23 #include <Eigen/Dense>
24 
25 IMPCORE_BEGIN_NAMESPACE
26 
27 IMP_DECORATORS_DECL(RigidMember, RigidMembers);
28 IMP_DECORATORS_DECL(RigidBodyMember, RigidBodyMembers);
29 
30 //! A decorator for a rigid body
31 /** A rigid body particle describes a set of particles, known
32  as the members, which move rigidly together. The rigid body
33  is represented as an algebra::ReferenceFrame3D coupled
34  with local coordinates (RigidMember::get_internal_coordinates())
35  for the members expressed in that reference frame. The
36  global coordinates of the members are accessed, as with
37  other global coordinates, via the XYZ::get_coordinates().
38 
39  Since the
40  members are simply a set of particles which move together
41  they don't (necessarily) define a shape. For example,
42  the members of the rigid body made from a molecular hierarchy
43  would include particles corresponding to intermediate levels
44  of the hierarchy. As a result, methods
45  that use rigid bodies usually should simply take the list of
46  particles they are interested in and then check for rigid
47  bodies internally.
48 
49  The initial reference of the rigid body is computed from
50  the coordinates, masses and radii of the particles
51  passed to the constructor, based on diagonalizing the
52  inertial tensor (which is not stored, currently).
53 
54  The rigid body radius is the farthest point of any of its
55  members from the origin of its reference frame. For rigid
56  body members, this takes into account the radius of the
57  member.
58 
59  RigidBodies can be nested (that is, a RigidBody can have
60  another RigidBody as a member). This can be useful for
61  organizational reasons as well as for accelerating
62  computations since operations are affected by
63  the total number of children contained in the rigid body
64  being operated on. Examples of this include collision detection
65  where if you have multiple representations of geometry at
66  different resolutions it is faster to put each of them
67  in a separate rigid body and then create one rigid body
68  containing all of them.
69 
70  It is often desirable to randomize the orientation of a rigid
71  body:
72  \include randomize_rigid_body.py
73 
74  \usesconstraint
75 
76  \see RigidMember
77  \see NonRigidMember
78  \see RigidBodyMover
79  \see RigidClosePairsFinder
80  \see RigidBodyDistancePairScore
81  */
82 class IMPCOREEXPORT RigidBody : public XYZ {
83  private:
84  /* Computes the coordinates of p given its internal (local)
85  coordinates and the current position and orientation of the
86  rigid body.
87  */
89 
90  void add_member_internal(Particle *p,
91  const algebra::ReferenceFrame3D &rf);
92 
93  //! do updates to rigid body upon changes in its members
94  //! such as udpating the rigid body radius based on the
95  //! point/sphere distance of all of its point/sphere members
96  //! from its origin
97  void on_change();
98 
99  static void teardown_constraints(Particle *p);
100 
101  static ObjectKey get_constraint_key_0();
102 
103  static ObjectKey get_constraint_key_1();
104 
105  // setup rigid body attributes with particles in ps, using their
106  // center of mass, inertia tensor to initialize the reference frame
107  static void do_setup_particle(Model *m, ParticleIndex pi,
109 
110  // setup a rigid body with specified reference frame
111  static void do_setup_particle(Model *m, ParticleIndex pi,
112  const algebra::ReferenceFrame3D &rf);
113 
114  void setup_score_states();
115 
116  // add a member associated with xyz coords (if it has a ref frame,
117  // it is still being ignored)
118  void add_point_member(ParticleIndex pi);
119 
120  // add a member associated with a reference frame
121  void add_rigid_body_member(ParticleIndex pi);
122 
123  public:
124  RigidMembers get_rigid_members() const;
125 
126  //! Returns a list of all members that are not themselves decorated as
127  //! rigid bodies, in the form of particle indexes.
129  static ParticleIndexes empty;
130  if (get_model()->get_has_attribute(internal::rigid_body_data().members_,
131  get_particle_index())) {
132  return get_model()->get_attribute(internal::rigid_body_data().members_,
134  } else {
135  return empty;
136  }
137  }
138 
139  //! Get all members that are themselves decorated as rigid bodies,
140  //! as model particle indexes
142  static ParticleIndexes empty;
143  if (get_model()->get_has_attribute(
144  internal::rigid_body_data().body_members_, get_particle_index())) {
145  return get_model()->get_attribute(
146  internal::rigid_body_data().body_members_, get_particle_index());
147  } else {
148  return empty;
149  }
150  }
151 
152  //! Get the particle indexes of any member of this rigid body, regardless
153  //! of whether it is itself a rigid body or not
155  return get_member_particle_indexes() + get_body_member_particle_indexes();
156  }
157 
159 
160  /**
161  Create a rigid body for pi with the particle indexes ps as its members.
162  The coordinates of pi are set to the center of mass of ps and the rotation
163  of its reference frame is based on the diagonalized inertia tensor of ps.
164 
165  @note If size(ps)=1, then its reference frame is copied if it is a
166  rigid body, or its rotation is set to identity if it is not
167  a rigid body.
168  */
170 
171  /**
172  Create a rigid body with the passed reference frame as its initial
173  position.
174  */
176 
177  //! Make the rigid body no longer rigid.
178  static void teardown_particle(RigidBody rb);
179 
180  IMP_CXX11_DEFAULT_COPY_CONSTRUCTOR(RigidBody);
181  ~RigidBody();
182 
183  //! Return true if the particle is a rigid body
184  static bool get_is_setup(Model *m, ParticleIndex pi) {
185  return internal::get_has_required_attributes_for_body(m, pi);
186  }
187 
188  // swig doesn't support using, so the method is wrapped
189  //! Get the coordinates of the particle
190  //! (= translation from local to global rigid body coordinates)
192 
193  //! returns the rotation of the particle
194  //! (= rotation from local to global rigid body orientation)
196  return get_reference_frame().get_transformation_to().get_rotation();
197  }
198 
199  //! Get the reference frame of this rigid body, which enables
200  //! trnasformation between the local rigid body coordinates
201  //! global coordinates
204  get_model()->get_attribute(internal::rigid_body_data().quaternion_[0],
206  get_model()->get_attribute(internal::rigid_body_data().quaternion_[1],
208  get_model()->get_attribute(internal::rigid_body_data().quaternion_[2],
210  get_model()->get_attribute(internal::rigid_body_data().quaternion_[3],
211  get_particle_index()));
212  IMP_USAGE_CHECK_FLOAT_EQUAL(v.get_squared_magnitude(), 1,
213  "Rotation is not a unit vector: " << v);
214  /*if (v.get_squared_magnitude() > 0){
215  v = v.get_unit_vector();
216  } else {
217  v = algebra::VectorD<4>(1,0,0,0);
218  }*/
219  bool assume_normalized = true;
220  IMP::algebra::Rotation3D rot(v, assume_normalized);
223  }
224 
225  //! Set the current reference frame
226  /** All members of the rigid body will have their coordinates updated
227  immediately.
228  \see IMP::core::transform(RigidBody,const algebra::Transformation3D&)
229  \see set_reference_frame_lazy()
230  */
231  void set_reference_frame(const IMP::algebra::ReferenceFrame3D &tr);
232 
233  //! Change the reference, delay updating the members until evaluate
234  /** \see set_reference_frame()
235  */
236  inline void set_reference_frame_lazy
238  {
241  get_particle()->set_value(internal::rigid_body_data().quaternion_[0], v[0]);
242  get_particle()->set_value(internal::rigid_body_data().quaternion_[1], v[1]);
243  get_particle()->set_value(internal::rigid_body_data().quaternion_[2], v[2]);
244  get_particle()->set_value(internal::rigid_body_data().quaternion_[3], v[3]);
246  }
247 
248 #ifndef SWIG
249 #ifndef IMP_DOXYGEN
250  //! 'expert' method for setting the reference more quickly
251  //! use at own risk
252  inline void set_rotation_lazy_using_internal_tables
253  (const IMP::algebra::Rotation3D &rot,
254  double* quaternion_tables[])
255  {
257  rot.get_quaternion();
258  int pi=get_particle_index().get_index();
259  quaternion_tables[0][pi]=v[0];
260  quaternion_tables[1][pi]=v[1];
261  quaternion_tables[2][pi]=v[2];
262  quaternion_tables[3][pi]=v[3];
263  }
264 
265  //! 'expert' method for setting the reference more quickly
266  //! use at own risk
267  inline void apply_rotation_lazy_using_internal_tables
268  (const IMP::algebra::Rotation3D &rot,
269  double* quaternion_tables[])
270  {
271  int pi=get_particle_index().get_index();
273  ( quaternion_tables[0][pi],
274  quaternion_tables[1][pi],
275  quaternion_tables[2][pi],
276  quaternion_tables[3][pi] );
278  (cur_rot*rot).get_quaternion();;
279  quaternion_tables[0][pi]=v[0];
280  quaternion_tables[1][pi]=v[1];
281  quaternion_tables[2][pi]=v[2];
282  quaternion_tables[3][pi]=v[3];
283  }
284 
285 #endif // IMP_DOXYGEN
286 #endif // SWIG
287 
288 
289 
290 
291  /** Update the reference frame of the rigid body based on aligning
292  the current global coordinates of the passed rigid body members
293  onto their old local coordinates. Non-passed members are ignored.
294 
295  This method is useful for updating the rigid body after new
296  global coordinates were loaded for the members. The members are
297  passed explicitly since, typically, some are desired to just
298  move along with the newly loaded rigid body.
299 
300  \note This requires at least three members that are not colinear
301  to work.
302  */
303  void set_reference_frame_from_members(const ParticleIndexes &members);
304 
305  /** Update the translational and rotational derivatives
306  on the rigid body center of mass, using the Cartesian derivative
307  vector at a speicified location (the point where the force is
308  being applied).
309 
310  Updates both the quaternion derivatives and the torque.
311 
312  @param local_derivative The derivative vector in local rigid body coordinates
313  @param local_location The location where the derivative is taken in local
314  rigid body coordinates
315  @param da Accumulates the output derivative over the rigid body
316  center of mass (translation and rotation torque, quaternion)
317  */
318  inline void add_to_derivatives(const algebra::Vector3D &local_derivative,
319  const algebra::Vector3D &local_location,
320  DerivativeAccumulator &da);
321 
322  /** Faster version of the above, if all is cached.
323 
324  @param local_derivative The derivative vector in local rigid body coordinates
325  @param global_derivative The derivative vector in global coordinates
326  @param local_location The location where the derivative is taken in local
327  rigid body coordinates
328  @param rot_local_to_global Rotation matrix from local rigid body to
329  global coordinates
330  @param da Accumulates the output derivative over the rigid body
331  center of mass (translation and rotation torque, quaternion)
332  */
333  inline void add_to_derivatives(const algebra::Vector3D &local_derivative,
334  const algebra::Vector3D &global_derivative,
335  const algebra::Vector3D &local_location,
336  const algebra::Rotation3D &rot_local_to_global,
337  DerivativeAccumulator &da);
338 
339  /** Update the rotational derivatives from another body specified by the
340  rotation from the other body's local coordinates to this body's local
341  coordinates. The provided quaternion derivative on the other body are in
342  the reference frame of the other body.
343 
344  Updates only quaternion derivatives.
345 
346  @param other_qderiv The derivative on the quaternion taking the other body's
347  local coordinates to global.
348  @param rot_other_to_local Rotation taking the local coordinates of the other body
349  to this body's local coordinates.
350  @param rot_local_to_global Rotation taking this rigid body's local coordinates to
351  global coordinates.
352  @param da Accumulates the output derivatives.
353  */
354  inline void add_to_rotational_derivatives(const algebra::Vector4D &other_qderiv,
355  const algebra::Rotation3D &rot_other_to_local,
356  const algebra::Rotation3D &rot_local_to_global,
357  DerivativeAccumulator &da);
358 
359  /** Add to quaternion derivative of this rigid body
360  Note that this method does not update the torque.
361 
362  @param qderiv Derivative wrt to quaternion taking local coordinates to
363  global.
364  @param da Object for accumulating derivatives
365  */
366  inline void add_to_rotational_derivatives(const algebra::Vector4D &qderiv,
367  DerivativeAccumulator &da);
368 
369  /** Add torque to derivative table of this rigid body
370  Note that this method does not update the quaternion derivatives, so should
371  be used by optimizers that rely on torque only (e.g. BrownianDynamics)
372 
373  @param torque_local Torque vector in local reference frame,
374  in units of kCal/Mol/Radian
375  @param da Object for accumulating derivatives
376  */
377  inline void add_to_torque(const algebra::Vector3D &torque_local,
378  DerivativeAccumulator &da);
379 
380 
381  /** The units are kCal/Mol/Radian */
383  algebra::Vector3D ret;
384  for (unsigned int i = 0; i < 3; ++i) {
385  ret[i] = get_model()->get_derivative(
386  internal::rigid_body_data().torque_[i], get_particle_index());
387  }
388  return ret;
389  }
390 
391 #if !defined(SWIG) && !defined(IMP_DOXYGEN)
392  //! expert method for fast const-access to internal torque
393  //! of coordinate #i table
394  static double const* access_torque_i_data
395  (IMP::Model const* m, unsigned int i)
396  {
397  IMP_USAGE_CHECK(i<3,"torque is 3 dimensional");
398  FloatKey k=
399  internal::rigid_body_data().torque_[i];
400  double const* ret=m->access_derivative_data(k);
401  return ret;
402  }
403 
404  //! expert method for fast access to internal torque
405  //! of coordinate #i table
406  static double* access_torque_i_data
407  (IMP::Model* m, unsigned int i)
408  {
409  IMP_USAGE_CHECK(i<3,"torque is 3 dimensional");
410  FloatKey k=
411  internal::rigid_body_data().torque_[i];
412  double* ret=m->access_derivative_data(k);
413  return ret;
414  }
415 
416  //! expert method for fast const-access to internal quaternion coordinate #i table
417  static double const* access_quaternion_i_data
418  (IMP::Model const* m, unsigned int i)
419  {
420  IMP_USAGE_CHECK(i<4,"quaternion is 4 dimensional");
421  FloatKey k=
422  internal::rigid_body_data().quaternion_[i];
423  double const* ret=m->FloatAttributeTable::access_attribute_data(k);
424  return ret;
425  }
426 
427  //! expert method for fast access to internal quaternion coordinate #i table
428  static double* access_quaternion_i_data
429  (IMP::Model* m, unsigned int i)
430  {
431  IMP_USAGE_CHECK(i<4,"quaternion is 4 dimensional");
432  FloatKey k=
433  internal::rigid_body_data().quaternion_[i];
434  double* ret=m->FloatAttributeTable::access_attribute_data(k);
435  return ret;
436  }
437 
438 
439 #endif
440 
441  //! Returns true if the rigid body coordinates are flagged as
442  //! optimized for Optimizer objects
443  bool get_coordinates_are_optimized() const;
444 
445  //! Set whether the rigid body coordinates are flagged as optimized
446  //! for Optimizer objects
447  void set_coordinates_are_optimized(bool tf);
448 
449  //! Normalize the quaternion
450  void normalize_rotation();
451 
452  //! Update the global coordinates of the members based on
453  //! their local coordinates and this rigid body's reference frame
454  void update_members();
455 
456  //! Get the derivatives of the quaternion
457  algebra::VectorD<4> get_rotational_derivatives() const;
458 
459 #ifndef IMP_DOXYGEN
460  unsigned int get_number_of_members() const {
461  return get_body_member_particle_indexes().size() +
462  get_member_particle_indexes().size();
463  }
464 
465  RigidBodyMember get_member(unsigned int i) const;
466 #endif
467  //! Add a proper member that moves rigidly with this rigid body,
468  //! properly handling rigid bodies and XYZ particles.
469  /**
470  Add p to the list of members. If p is a valid RigidBody, it is added
471  as a rigid body member, otherwise it is added as a point member
472  (for which the rotation is not tracked). By default, p is considered
473  a strictly rigid member, in that its local coordinates are not expected
474  to change independently.
475 
476  The radius of the rigid body is updated to reflect the new member.
477 
478  \see add_non_rigid_member
479  */
480  void add_member(ParticleIndexAdaptor p);
481 
482  /** Add a non-rigid member, for which internal coordinates may change
483  independently.
484 
485  @note Currently RigidBody non-rigid members are not handled properly.
486  */
487  void add_non_rigid_member(ParticleIndexAdaptor p);
488 
489  /** Set whether a particular member is flagged as a rigid member
490  or as a non-rigid member. This affects the way the rigid body
491  updates the coordinates and / or reference frame of its members.
492 
493  The radius of the rigid body is updated to reflect this change.
494  */
495  void set_is_rigid_member(ParticleIndex pi, bool tf);
496 };
497 
498 #ifndef IMP_DOXYGEN
499 // inline implementation
500 void RigidBody::add_to_rotational_derivatives(const algebra::Vector4D &qderiv,
501  DerivativeAccumulator &da) {
502  for (unsigned int i = 0; i < 4; ++i) {
503  get_model()->add_to_derivative(internal::rigid_body_data().quaternion_[i],
504  get_particle_index(), qderiv[i], da);
505  }
506 }
507 
508 
509 // inline implementation
510 void RigidBody::add_to_torque(const algebra::Vector3D &torque_local,
511  DerivativeAccumulator &da) {
512  for (unsigned int i = 0; i < 3; ++i) {
513  get_model()->add_to_derivative(internal::rigid_body_data().torque_[i],
514  get_particle_index(), torque_local[i], da);
515  }
516 }
517 
518 
519 // inline implementation
520 void RigidBody::add_to_derivatives(const algebra::Vector3D &deriv_local,
521  const algebra::Vector3D &deriv_global,
522  const algebra::Vector3D &local,
523  const algebra::Rotation3D &rot_local_to_global,
524  DerivativeAccumulator &da) {
525  // IMP_LOG_TERSE( "Accumulating rigid body derivatives" << std::endl);
526  XYZ::add_to_derivatives(deriv_global, da);
527 
528  Eigen::RowVector4d q =
529  Eigen::RowVector3d(deriv_global.get_data()) *
530  rot_local_to_global.get_gradient(Eigen::Vector3d(local.get_data()));
531 
532  for (unsigned int i = 0; i < 4; ++i) {
533  get_model()->add_to_derivative(internal::rigid_body_data().quaternion_[i],
534  get_particle_index(), q[i], da);
535  }
536  algebra::Vector3D torque = algebra::get_vector_product(local, deriv_local);
537  add_to_torque(torque, da);
538 }
539 
540 // inline implementation
541 void RigidBody::add_to_derivatives(const algebra::Vector3D &deriv_local,
542  const algebra::Vector3D &local,
543  DerivativeAccumulator &da) {
544  algebra::Rotation3D rot_local_to_global =
545  get_reference_frame().get_transformation_to().get_rotation();
546  const algebra::Vector3D deriv_global = rot_local_to_global * deriv_local;
547  add_to_derivatives(deriv_local, deriv_global, local,
548  rot_local_to_global, da);
549 }
550 
551 // inline implementation
552 void RigidBody::add_to_rotational_derivatives(const algebra::Vector4D &other_qderiv,
553  const algebra::Rotation3D &rot_other_to_local,
554  const algebra::Rotation3D &rot_local_to_global,
555  DerivativeAccumulator &da) {
556  Eigen::MatrixXd derivs =
558  rot_local_to_global, rot_other_to_local);
559  Eigen::RowVector4d qderiv = Eigen::RowVector4d(other_qderiv.get_data()) * derivs;
560  for (unsigned int i = 0; i < 4; ++i) {
561  get_model()->add_to_derivative(internal::rigid_body_data().quaternion_[i],
562  get_particle_index(), qderiv[i], da);
563  }
564 }
565 #endif
566 
567 
568 /** It is often useful to store precalculated properties of the rigid body
569  for later use. These need to be cleared out when the rigid body changes.
570  To make sure this happens, register the key here.
571 */
572 void IMPCOREEXPORT add_rigid_body_cache_key(ObjectKey k);
573 
574 //! A member of a rigid body, it has internal (local) coordinates
575 class IMPCOREEXPORT RigidBodyMember : public XYZ {
577 
578  RigidBody get_rigid_body() const;
579 
580  //! Return the internal (local) coordinates of this member
581  /** These coordinates are relative to the reference frame of the
582  rigid body that owns it
583  */
585  return get_model()->get_internal_coordinates(get_particle_index());
586  }
587 
588  //! set the internal (local) coordinates for this member
590  get_model()->get_internal_coordinates(get_particle_index()) = v;
591  get_rigid_body().get_model()->clear_particle_caches(get_particle_index());
592  }
593 
594  //! Set the internal (local) coordinates of this member,
595  //! assuming it is a rigid body itself
596  /** Set the internal (local) coordinates of this rigid body
597  relative to the reference frame of the rigid body that owns it
598  */
601  get_model()->get_has_attribute(
602  internal::rigid_body_data().lquaternion_[0], get_particle_index()),
603  "Can only set the internal transformation if member is"
604  << " a rigid body itself.");
605  set_internal_coordinates(v.get_translation());
606 
607  get_model()->set_attribute(internal::rigid_body_data().lquaternion_[0],
609  v.get_rotation().get_quaternion()[0]);
610  get_model()->set_attribute(internal::rigid_body_data().lquaternion_[1],
612  v.get_rotation().get_quaternion()[1]);
613  get_model()->set_attribute(internal::rigid_body_data().lquaternion_[2],
615  v.get_rotation().get_quaternion()[2]);
616  get_model()->set_attribute(internal::rigid_body_data().lquaternion_[3],
618  v.get_rotation().get_quaternion()[3]);
619  get_rigid_body().get_model()->clear_particle_caches(get_particle_index());
620  }
621 
622  //! Return the internal (local) coordinates of this member,
623  //! assuming it is a rigid body itself
624  /** Return the internal (local) coordinates of this rigid body
625  relative to the reference frame of the rigid body that owns it
626  */
629  get_model()->get_has_attribute(
630  internal::rigid_body_data().lquaternion_[0], get_particle_index()),
631  "Can only set the internal transformation if member is a "
632  << "rigid body itself.");
633  algebra::Vector3D tr =
634  get_model()->get_internal_coordinates(get_particle_index());
636  get_model()->get_attribute(internal::rigid_body_data().lquaternion_[0],
638  get_model()->get_attribute(internal::rigid_body_data().lquaternion_[1],
640  get_model()->get_attribute(internal::rigid_body_data().lquaternion_[2],
642  get_model()->get_attribute(internal::rigid_body_data().lquaternion_[3],
643  get_particle_index()));
644  return algebra::Transformation3D(rot, tr);
645  }
646 
647  ~RigidBodyMember();
648  //! sets the global coordinates of this member using XYZ::set_coordinates()
649  // this is here since swig does like using statements
650  void set_coordinates(const algebra::Vector3D &center) {
651  XYZ::set_coordinates(center);
652  }
653 
654 #ifndef IMP_DOXYGEN
655  //! Set the global coordinates from the internal coordinates,
656  //! using tr as a reference frame
658  set_coordinates(tr.get_transformed(get_internal_coordinates()));
659  }
660 #endif
661  IMP_CXX11_DEFAULT_COPY_CONSTRUCTOR(RigidBodyMember);
662 
663  //! return true if it is a rigid member
665  return internal::get_has_required_attributes_for_member(m, p);
666  }
667 
668  static FloatKeys get_internal_coordinate_keys() {
669  return internal::rigid_body_data().child_keys_;
670  }
671 
672  static FloatKeys get_internal_rotation_keys() {
673  return internal::rigid_body_data().lquaternion_;
674  }
675 };
676 
677 //! A decorator for a particle that is part of a rigid body, and is
678 //! actually rigid
679 /**
680  RigidMember particles, as opposed to NonRigidMember particles, are
681  not expected to change their internal (local) coordinates or
682  reference frames, and their global coordinates are expected to
683  change only through setting the coordinates (or reference frame) of
684  the rigid body that owns them.
685 
686  \see RigidBody
687  */
688 class IMPCOREEXPORT RigidMember : public RigidBodyMember {
689  public:
691 
692  IMP_CXX11_DEFAULT_COPY_CONSTRUCTOR(RigidMember);
693  ~RigidMember();
694 
695  //! return true if it is a rigid member
697  return internal::get_has_required_attributes_for_rigid_member(m, p);
698  }
699 };
700 
701 //! A decorator for a particle that is part of a rigid body but not rigid
702 /** NonRigidMembers, like RigidMembers, have internal coordinates and move
703  along with the rigid body. However, it is expected that their internal
704  coordinates will change, and so they are not part of structures that
705  assume rigidity.
706 
707  \see RigidBody
708  */
709 class IMPCOREEXPORT NonRigidMember : public RigidBodyMember {
710  public:
712  IMP_CXX11_DEFAULT_COPY_CONSTRUCTOR(NonRigidMember);
713  ~NonRigidMember();
714 
715  //! return true if it is a rigid member
716  static bool get_is_setup(Model *m, ParticleIndex p) {
717  return internal::get_has_required_attributes_for_non_member(m, p);
718  }
719 
720  //! Add to derivatives of local coordinates.
721  /** @param deriv_parent The derivative vector in local coordinates of the
722  parent rigid body.
723  @param da Accumulates the derivative over the local translation.
724  */
726  DerivativeAccumulator &da) {
727  for (unsigned int i = 0; i < 3; ++i) {
728  get_model()->add_to_derivative(get_internal_coordinate_keys()[i],
729  get_particle_index(), deriv_parent[i], da);
730  }
731  }
732 
733  /** Update the rotational derivatives of the internal transformation.
734 
735  Updates only local quaternion derivatives.
736 
737  @param local_qderiv The derivative on the quaternion taking this non-rigid
738  body's local coordinates to global.
739  @param rot_local_to_parent Rotation taking the local coordinates of the non-rigid
740  body to its parent's.
741  @param rot_parent_to_global Rotation taking the parent rigid body's local coordinates
742  to global coordinates.
743  @param da Accumulates the output derivatives.
744  */
746  const algebra::Rotation3D &rot_local_to_parent,
747  const algebra::Rotation3D &rot_parent_to_global,
748  DerivativeAccumulator &da) {
749  Eigen::MatrixXd derivs =
751  rot_parent_to_global, rot_local_to_parent);
752  Eigen::RowVector4d qderiv = Eigen::RowVector4d(local_qderiv.get_data()) * derivs;
753  for (unsigned int i = 0; i < 4; ++i) {
754  get_model()->add_to_derivative(get_internal_rotation_keys()[i],
755  get_particle_index(), qderiv[i], da);
756  }
757  }
758 
759  /** Add to internal quaternion derivatives of this non-rigid body
760 
761  @param qderiv Derivative wrt to quaternion taking local coordinates to
762  parent.
763  @param da Object for accumulating derivatives
764  */
766  DerivativeAccumulator &da) {
768  get_model()->get_has_attribute(
769  get_internal_rotation_keys()[0], get_particle_index()),
770  "Can only set derivatives of internal rotation if member is a "
771  << "rigid body itself.");
772  for (unsigned int i = 0; i < 4; ++i) {
773  get_model()->add_to_derivative(get_internal_rotation_keys()[i],
774  get_particle_index(), qderiv[i], da);
775  }
776  }
777 
778 
779  //! Get derivatives wrt translation component of internal transformation.
781  algebra::Vector3D ret;
782  for (unsigned int i = 0; i < 3; ++i) {
783  ret[i] = get_model()->get_derivative(
784  get_internal_coordinate_keys()[i], get_particle_index());
785  }
786  return ret;
787  }
788 
789  //! Get derivatives wrt quaternion component of internal transformation.
791  algebra::Vector4D ret;
792  for (unsigned int i = 0; i < 4; ++i) {
793  ret[i] = get_model()->get_derivative(
794  get_internal_rotation_keys()[i], get_particle_index());
795  }
796  return ret;
797  }
798 };
799 
800 #ifndef IMP_DOXYGEN
801 
802 class IMPCOREEXPORT RigidMembersRefiner : public Refiner {
803  public:
804  RigidMembersRefiner(std::string name = "RigidMembersRefiner%d")
805  : Refiner(name) {}
806  virtual bool get_can_refine(Particle *) const IMP_OVERRIDE;
807 #ifndef SWIG
808  using Refiner::get_refined;
809 #endif
810  virtual const ParticlesTemp get_refined(Particle *) const
811  IMP_OVERRIDE;
812  virtual ModelObjectsTemp do_get_inputs(
813  Model *m, const ParticleIndexes &pis) const IMP_OVERRIDE;
814  IMP_OBJECT_METHODS(RigidMembersRefiner);
815 };
816 
817 namespace internal {
818 IMPCOREEXPORT RigidMembersRefiner *get_rigid_members_refiner();
819 }
820 #endif
821 
822 //! Transform a rigid body
823 /** The transformation is applied current conformation of the rigid
824  body, as opposed to replacing the current conformation, as in
825  RigidBody::set_reference_frame().
826 
827  \see RigidBody
828  \see algebra::Transformation3D
829 */
830 inline void transform(RigidBody a, const algebra::Transformation3D &tr) {
831  a.set_reference_frame(get_transformed(a.get_reference_frame(), tr));
832 }
833 
834 /** Compute the rigid body reference frame given a set of input particles.
835  */
836 IMPCOREEXPORT algebra::ReferenceFrame3D get_initial_reference_frame(
837  Model *m, const ParticleIndexes &pis);
838 
839 inline algebra::ReferenceFrame3D get_initial_reference_frame(
840  const ParticlesTemp &ps) {
841  if (ps.empty()) {
842  return algebra::ReferenceFrame3D();
843  }
844  return get_initial_reference_frame(ps[0]->get_model(),
845  get_indexes(ps));
846 }
847 
848 /** Create a set of rigid bodies that are bound together for efficiency.
849  These rigid bodies cannot nest or have other dependencies among them.
850 
851  All rigid bodies have the default reference frame.
852 
853  \note Do not use this with DOMINO as all the rigid bodies use the same
854  ScoreState and so will be considered inter-dependent.
855 */
856 IMPCOREEXPORT ParticlesTemp create_rigid_bodies(Model *m,
857  unsigned int n,
858  bool no_members =
859  false);
860 
861 IMP_DECORATORS_DEF(RigidMember, RigidMembers);
862 IMP_DECORATORS(RigidBody, RigidBodies, XYZs);
863 
864 /** Show the rigid body hierarchy rooted at passed body. */
865 IMPCOREEXPORT void show_rigid_body_hierarchy(RigidBody rb,
866  TextOutput out =
867  TextOutput(std::cout));
868 
869 //! Return the index of the outer-most rigid body containing the member.
870 /** Use this to, for example, group particles into rigid bodies. */
871 IMPCOREEXPORT ParticleIndex get_root_rigid_body(RigidMember m);
872 
873 IMPCORE_END_NAMESPACE
874 
875 #endif /* IMPCORE_RIGID_BODIES_H */
void add_to_derivatives(const algebra::Vector3D &v, DerivativeAccumulator &d)
Add the vector v to the derivative vector of the x,y,z coordinates.
Definition: XYZ.h:82
void set_internal_coordinates(const algebra::Vector3D &v) const
set the internal (local) coordinates for this member
Definition: rigid_bodies.h:589
A Modifier on ParticlesTemp.
Simple 3D transformation class.
ParticleIndex get_particle_index() const
Returns the particle index decorated by this decorator.
Definition: Decorator.h:188
algebra::Vector3D get_internal_derivatives() const
Get derivatives wrt translation component of internal transformation.
Definition: rigid_bodies.h:780
algebra::Vector4D get_internal_rotational_derivatives() const
Get derivatives wrt quaternion component of internal transformation.
Definition: rigid_bodies.h:790
Key< 0 > FloatKey
The type used to identify float attributes in the Particles.
Definition: base_types.h:32
A member of a rigid body, it has internal (local) coordinates.
Definition: rigid_bodies.h:575
A container for Singletons.
ParticleIndex get_root_rigid_body(RigidMember m)
Return the index of the outer-most rigid body containing the member.
#define IMP_USAGE_CHECK_FLOAT_EQUAL(expra, exprb, message)
Definition: check_macros.h:178
algebra::Vector3D get_coordinates() const
Definition: rigid_bodies.h:191
IMP::algebra::Rotation3D get_rotation() const
Definition: rigid_bodies.h:195
#define IMP_DECORATOR_SETUP_1(Name, FirstArgumentType, first_argument_name)
#define IMP_OBJECT_METHODS(Name)
Define the basic things needed by any Object.
Definition: object_macros.h:25
void add_to_internal_rotational_derivatives(const algebra::Vector4D &local_qderiv, const algebra::Rotation3D &rot_local_to_parent, const algebra::Rotation3D &rot_parent_to_global, DerivativeAccumulator &da)
Definition: rigid_bodies.h:745
Model * get_model() const
Returns the Model containing the particle.
Definition: Decorator.h:191
void add_to_internal_rotational_derivatives(const algebra::Vector4D &qderiv, DerivativeAccumulator &da)
Definition: rigid_bodies.h:765
ParticlesTemp create_rigid_bodies(Model *m, unsigned int n, bool no_members=false)
Index< ParticleIndexTag > ParticleIndex
Definition: base_types.h:154
Eigen::MatrixXd get_gradient_of_composed_with_respect_to_second(const Rotation3D &q, const Rotation3D &p, bool projected=true)
Get gradient of quaternion product with respect to second quaternion.
bool get_coordinates_are_optimized() const
Get whether the coordinates are optimized.
Definition: XYZ.h:89
void add_rigid_body_cache_key(ObjectKey k)
A more IMP-like version of the std::vector.
Definition: Vector.h:39
const Vector4D & get_quaternion() const
Return the quaternion so that it can be stored.
Definition: Rotation3D.h:207
Simple XYZ decorator.
VectorD< 4 > Vector4D
Definition: VectorD.h:399
A reference frame in 3D.
void show_rigid_body_hierarchy(RigidBody rb, TextOutput out=TextOutput(std::cout))
void clear_particle_caches(ParticleIndex pi)
Clear all the cache attributes of a given particle.
IMP::Vector< IMP::WeakPointer< ModelObject > > ModelObjectsTemp
Definition: base_types.h:82
Class for storing model, its restraints, constraints, and particles.
Definition: Model.h:72
const Rotation3D & get_rotation() const
Return the rotation associated with this transformation.
algebra::Vector3D get_torque() const
Definition: rigid_bodies.h:382
static bool get_is_setup(Model *m, ParticleIndexAdaptor p)
return true if it is a rigid member
Definition: rigid_bodies.h:696
virtual bool get_can_refine(Particle *) const
Return true if this refiner can refine that particle.
Definition: Refiner.h:51
Vector3D get_vector_product(const Vector3D &p1, const Vector3D &p2)
Return the vector product (cross product) of two vectors.
Definition: Vector3D.h:31
void set_reference_frame(const IMP::algebra::ReferenceFrame3D &tr)
Set the current reference frame.
A Cartesian vector in D-dimensions.
Definition: VectorD.h:52
Refine a particle into a list of particles.
void set_coordinates(const algebra::Vector3D &v)
set all coordinates from a vector
Definition: XYZ.h:62
const Transformation3D & get_transformation_to() const
Return transformation from local to global coordinates.
static bool get_is_setup(Model *m, ParticleIndex pi)
Return true if the particle is a rigid body.
Definition: rigid_bodies.h:184
const ParticleIndexes & get_member_particle_indexes() const
Definition: rigid_bodies.h:128
static bool get_is_setup(Model *m, ParticleIndexAdaptor p)
return true if it is a rigid member
Definition: rigid_bodies.h:664
Represent an XYZR particle with a sphere.
void set_attribute(TypeKey attribute_key, ParticleIndex particle, Type value)
set the value of particle attribute with the specified key
void add_to_internal_derivatives(const algebra::Vector3D &deriv_parent, DerivativeAccumulator &da)
Add to derivatives of local coordinates.
Definition: rigid_bodies.h:725
A decorator for a particle with x,y,z coordinates.
Definition: XYZ.h:30
const algebra::Vector3D & get_coordinates() const
Convert it to a vector.
Definition: XYZ.h:109
virtual const ParticlesTemp get_refined(Particle *a) const =0
Refine the passed particle into a set of particles.
Vector3D get_transformed(const Vector3D &o) const
Transform.
A decorator for a particle that is part of a rigid body but not rigid.
Definition: rigid_bodies.h:709
Simple 3D rotation class.
void set_internal_transformation(const algebra::Transformation3D &v)
Definition: rigid_bodies.h:599
Particle * get_particle() const
Returns the particle decorated by this decorator.
Definition: Decorator.h:171
Key< 4 > ObjectKey
The type used to identify an Object attribute.
Definition: base_types.h:48
3D rotation class.
Definition: Rotation3D.h:46
void set_coordinates_are_optimized(bool tf) const
Set whether the coordinates are optimized.
Definition: XYZ.h:95
void set_coordinates(const algebra::Vector3D &center)
sets the global coordinates of this member using XYZ::set_coordinates()
Definition: rigid_bodies.h:650
#define IMP_DECORATOR_METHODS(Name, Parent)
VectorD< 3 > Vector3D
Definition: VectorD.h:395
Abstract class to implement hierarchical methods.
Definition: Refiner.h:34
Simple 3D vector class.
Class to handle individual particles of a Model object.
Definition: Particle.h:41
#define IMP_USAGE_CHECK(expr, message)
A runtime test for incorrect usage of a class or method.
Definition: check_macros.h:168
#define IMP_DECORATORS(Name, PluralName, Parent)
Define the types for storing sets of decorators.
void transform(RigidBody a, const algebra::Transformation3D &tr)
Transform a rigid body.
Definition: rigid_bodies.h:830
Eigen::MatrixXd get_gradient_of_composed_with_respect_to_first(const Rotation3D &q, const Rotation3D &p, bool projected=true)
Get gradient of quaternion product with respect to first quaternion.
A reference frame in 3D.
ParticleIndexes get_member_indexes() const
Definition: rigid_bodies.h:154
const ParticleIndexes & get_body_member_particle_indexes() const
Definition: rigid_bodies.h:141
A decorator for a rigid body.
Definition: rigid_bodies.h:82
const Vector3D & get_translation() const
Return the translation vector associated with this transformation.
algebra::Transformation3D get_internal_transformation() const
Definition: rigid_bodies.h:627
Type get_attribute(TypeKey attribute_key, ParticleIndex particle)
get the value of the particle attribute with the specified key
Decorator for a sphere-like particle.
ParticleIndexes get_indexes(const ParticlesTemp &ps)
static bool get_is_setup(Model *m, ParticleIndex p)
return true if it is a rigid member
Definition: rigid_bodies.h:716
#define IMP_OVERRIDE
Cause a compile error if this method does not override a parent method.
Class for adding derivatives from restraints to the model.
const algebra::Vector3D & get_internal_coordinates() const
Return the internal (local) coordinates of this member.
Definition: rigid_bodies.h:584
IMP::algebra::ReferenceFrame3D get_reference_frame() const
Definition: rigid_bodies.h:202