First, I changed vector so that it is get_magnitude rather than just magnitude to be consistent with the rest of IMP. I also added operator* for the dot product and operator- (and difference for python) since that is pretty standard and added a method to XYZDecorator to get a vector for the current position.
Second, I added a check for NaN derivative values in Particle.
Third, I added a comment about the internal checks in nonbonded list test.
Index: kernel/src/restraints/DihedralRestraint.cpp =================================================================== --- kernel/src/restraints/DihedralRestraint.cpp (revision 457) +++ kernel/src/restraints/DihedralRestraint.cpp (working copy) @@ -58,7 +58,7 @@ Vector3D v1 = rij.vector_product(rkj); Vector3D v2 = rkj.vector_product(rkl); Float scalar_product = v1.scalar_product(v2); - Float mag_product = v1.magnitude() * v2.magnitude(); + Float mag_product = v1.get_magnitude() * v2.get_magnitude();
// avoid division by zero Float cosangle = std::abs(mag_product) > 1e-12 @@ -86,9 +86,9 @@ // J. Mol. Biol. 234, 751-762 (1993) Vector3D vijkj = rij.vector_product(rkj); Vector3D vkjkl = rkj.vector_product(rkl); - Float sijkj2 = vijkj.squared_magnitude(); - Float skjkl2 = vkjkl.squared_magnitude(); - Float skj = rkj.magnitude(); + Float sijkj2 = vijkj.get_squared_magnitude(); + Float skjkl2 = vkjkl.get_squared_magnitude(); + Float skj = rkj.get_magnitude(); Float rijkj = rij.scalar_product(rkj); Float rkjkl = rkj.scalar_product(rkl);
Index: kernel/include/IMP/decorators/XYZDecorator.h =================================================================== --- kernel/include/IMP/decorators/XYZDecorator.h (revision 457) +++ kernel/include/IMP/decorators/XYZDecorator.h (working copy) @@ -82,6 +82,13 @@ b.get_coordinate(2) - get_coordinate(2)); }
+ //! Convert it to a vector. + /** Somewhat suspect based on wanting a Point/Vector differentiation + but we don't have points*/ + Vector3D get_vector() const { + return Vector3D(get_x(), get_y(), get_z()); + } + //! Get a vector containing the keys for x,y,z /** This is quite handy for initializing movers and things. */ Index: kernel/include/IMP/Vector3D.h =================================================================== --- kernel/include/IMP/Vector3D.h (revision 457) +++ kernel/include/IMP/Vector3D.h (working copy) @@ -52,6 +52,11 @@ return vec_[0] * vec2.vec_[0] + vec_[1] * vec2.vec_[1] + vec_[2] * vec2.vec_[2]; } + + //! scalar product + Float operator*(const Vector3D &o) const { + return scalar_product(o); + }
//! \return the vector product of two vectors. /** \param[in] vec2 The other vector to use in the product. @@ -63,23 +68,33 @@ }
//! \return The square of the magnitude of this vector. - Float squared_magnitude() const { + Float get_squared_magnitude() const { return vec_[0] * vec_[0] + vec_[1] * vec_[1] + vec_[2] * vec_[2]; }
//! \return The magnitude of this vector. - Float magnitude() const { - return std::sqrt(squared_magnitude()); + Float get_magnitude() const { + return std::sqrt(get_squared_magnitude()); }
//! \return This vector normalized to unit length. Vector3D get_unit_vector() const { - Float mag = magnitude(); + Float mag = get_magnitude(); // avoid division by zero mag = std::max(mag, static_cast<Float>(1e-12)); return Vector3D(vec_[0] / mag, vec_[1] / mag, vec_[2] / mag); }
+ //! for C++ + Vector3D operator-(const Vector3D &o) const { + return Vector3D(operator[](0)-o[0], + operator[](1)-o[1], + operator[](2)-o[2]); + } + //! for python + Vector3D difference(const Vector3D &o) const { + return operator-(o); + } private: Float vec_[3]; };
Index: kernel/include/IMP/Particle.h =================================================================== --- kernel/include/IMP/Particle.h (revision 457) +++ kernel/include/IMP/Particle.h (working copy) @@ -336,6 +336,7 @@ { IMP_check(get_is_active(), "Do not touch inactive particles", InactiveParticleException()); + IMP_assert(value==value, "Can't add NaN to derivative"); float_indexes_.get_value(name).derivative+= da(value); }
Index: kernel/test/states/test_nonbonded_list.py =================================================================== --- kernel/test/states/test_nonbonded_list.py (revision 457) +++ kernel/test/states/test_nonbonded_list.py (working copy) @@ -119,6 +119,8 @@ self.assertEqual(score, 25, "Wrong score") def test_spheres2(self): """Test the nonbonded list of spheres (num pairs)""" + # This uses the internal checks of the nonbonded list to + # verify correctness m= IMP.Model() rk= IMP.FloatKey("radius") for i in range(0,100):
Daniel Russel wrote: > First, I changed vector so that it is get_magnitude rather than just > magnitude to be consistent with the rest of IMP. I also added operator* > for the dot product and operator- (and difference for python) since that > is pretty standard and added a method to XYZDecorator to get a vector > for the current position.
Doesn't compile. But I fixed it for you. ;)
There is also no need for a difference method. Python can override operators too (except for =, of course, since everything in Python is a reference). And we shouldn't be putting "only for Python" junk into the C++ API anyway - you can use %extend in the SWIG interface for that sort of thing. I'll fix this in a later patch.
I also intentionally did not put in operator* for dot product before (I do know what an overloaded operator is!), so that people wouldn't get confused and think it was cross product. But what do others think?
> Third, I added a comment about the internal checks in nonbonded list test.
Was that comment intended for test_frido_spheres() ? That's the test case which seems to make no sense to me - not test_spheres2().
Ben
> Doesn't compile. But I fixed it for you. ;) > Thanks. I must have left out some file. > There is also no need for a difference method. Python can override > operators too (except for =, of course, since everything in Python is a > reference). And we shouldn't be putting "only for Python" junk into the > C++ API anyway - you can use %extend in the SWIG interface for that sort > of thing. I'll fix this in a later patch. > Fair enough. I don't know swig well enough to mess with such things. > I also intentionally did not put in operator* for dot product before (I > do know what an overloaded operator is!), so that people wouldn't get > confused and think it was cross product. But what do others think? > Most other libraries I have used support it. There is no real ambiguity in C++ as the return types have to be different (so it is hard to write an expressions which compiles but where you thought it was the cross product). You would get errors too in python, but not until later on in the code. > Was that comment intended for test_frido_spheres() ? That's the test > case which seems to make no sense to me - not test_spheres2(). > I must have scrolled too far :-) Yes, it should be the frido test.
participants (2)
-
Ben Webb -
Daniel Russel