IMP-dev
Threads by month
- ----- 2024 -----
- November
- October
- September
- August
- July
- June
- May
- April
- March
- February
- January
- ----- 2023 -----
- December
- November
- October
- September
- August
- July
- June
- May
- April
- March
- February
- January
- ----- 2022 -----
- December
- November
- October
- September
- August
- July
- June
- May
- April
- March
- February
- January
- ----- 2021 -----
- December
- November
- October
- September
- August
- July
- June
- May
- April
- March
- February
- January
- ----- 2020 -----
- December
- November
- October
- September
- August
- July
- June
- May
- April
- March
- February
- January
- ----- 2019 -----
- December
- November
- October
- September
- August
- July
- June
- May
- April
- March
- February
- January
- ----- 2018 -----
- December
- November
- October
- September
- August
- July
- June
- May
- April
- March
- February
- January
- ----- 2017 -----
- December
- November
- October
- September
- August
- July
- June
- May
- April
- March
- February
- January
- ----- 2016 -----
- December
- November
- October
- September
- August
- July
- June
- May
- April
- March
- February
- January
- ----- 2015 -----
- December
- November
- October
- September
- August
- July
- June
- May
- April
- March
- February
- January
- ----- 2014 -----
- December
- November
- October
- September
- August
- July
- June
- May
- April
- March
- February
- January
- ----- 2013 -----
- December
- November
- October
- September
- August
- July
- June
- May
- April
- March
- February
- January
- ----- 2012 -----
- December
- November
- October
- September
- August
- July
- June
- May
- April
- March
- February
- January
- ----- 2011 -----
- December
- November
- October
- September
- August
- July
- June
- May
- April
- March
- February
- January
- ----- 2010 -----
- December
- November
- October
- September
- August
- July
- June
- May
- April
- March
- February
- January
- ----- 2009 -----
- December
- November
- October
- September
- August
- July
- June
- May
- April
- March
- February
- January
- ----- 2008 -----
- December
- November
- October
- September
- August
- July
- June
- May
- April
- March
- February
- January
- ----- 2007 -----
- December
- November
- October
August 2008
- 5 participants
- 23 discussions
OK, so given Daniel's and Keren's suggestions, here's what I'll do:
I will create three new modules: core, exp, and modeller. All developers
will have commit privileges to the core and exp modules. Modifications
to the kernel will still need to go through me. The other specific
modules will remain as they currently are - each with a given owner
(Keren for em and domino, me for modeller).
The core module will contain most of the code that most users and other
projects use. Thus, most of the concrete classes currently in the kernel
(headers in subdirectories, not the top-level directory - e.g.
ConnectivityRestraint, RestraintSet, Daniel's nonbonded lists) will move
into the core module. API changes, deletions from, or additions of new
restraints, etc. to core (but not simply new methods to existing
classes, or cleanups of internals, etc.) should be announced on imp-dev
*before* the code is committed. (A general outline of new features is
fine - I just don't want to wake up one day and find 20,000 lines of
code in some entirely new set of classes.)
Code in the core module can be at various levels of maturity. The only
requirement is that the test cases pass "most of the time". (The general
rule of thumb should be that the nightly builds do not break for more
than 3 consecutive days.) As Daniel suggested, immature code should be
marked in a Doxygen group so that users can be aware.
The exp or experimental module is for code that is (scientifically)
experimental in nature. I will put particle refiners, cover bonds,
tunnel singleton score, and wormlike chain in there. This code should
still have test cases, of course, and be marked if it is immature just
like core is, but won't be built or tested as part of the nightly
builds, so there is less requirement that the tests pass all the time.
If a bunch of classes are envisaged for a given project (e.g. NPC
dynamics, a new optimizer, etc.) then probably the best thing is a
distinct module for that, where you'd have more flexibility, or even
something entirely external (e.g. Assembler).
Finally, the modeller module will contain what's currently in
modeller.py, plus anything else that requires Modeller, so that the
kernel and core can be built and tested without needing Modeller around.
(It will also get some additional Modeller-specific restraints which
I'll need to entice users away from Modeller.)
I will decide what goes in each module. Movements between modules should
be exceedingly rare. It should certainly *not* be the case that somebody
uses exp simply to write dodgy code that doesn't work, and then moves it
to core when the bugs are fixed. If stuff really needs to be moved
between modules, *I* will do it.
I will have the nightly build system email imp-dev if the nightly builds
fail. The culprit will then hopefully fix the problem ;) (with my
assistance if necessary).
I will create a Bugzilla instance to keep track of IMP bugs and wishlist
items. Bug fixes should, of course, come with a testcase to ensure the
bug stays fixed, and to some extent document the nature of said bug.
I will also work with Keren to put Assembler into the automatic build
system.
I am reluctant to set too many rules, but there is an obvious
requirement that things rely reasonably stable, so that developers can
concentrate on the science and not on updating their code to work with
the most recent changes. It is my responsibility to keep things that
way, so I will revert patches (or users ;) if they cause extreme
breakage. I'll be relying on everyone to take care with their code and
its testcases. So while I cannot make everybody happy, I can do the next
best thing and make everybody equally unhappy (hopefully only moderately
so ;).
Obviously there is a bunch of stuff to do, and I am not going to spend
my holiday weekend working on it. Plus, I will be in Chicago for most of
next week. So probably things will be ready to move forward at the start
of the following week.
If anybody can point out an obvious problem with this solution, I may
have to think about (and delay) things some more. Otherwise, have a
relaxing and sunny Labor Day. ;)
Ben
--
ben(a)salilab.org http://salilab.org/~ben/
"It is a capital mistake to theorize before one has data."
- Sir Arthur Conan Doyle
2
4
Daniel Russel wrote:
> Just to check, you need to have the NBL tests and the attribute removal
> test, but no others?
Well, the two bipartite NBL tests, and the attribute removal test, yes.
(The crash actually seems to occur in set_is_optimized, not in
remove_attribute itself.) If I removed either one of the NBL tests, the
bug went away. You can have *more* tests, but obviously it takes longer
to run the test suite then.
> Now that I think about it, do we have any reasonable way of debugging
> with wine? I wouldn't expect gdb to work.
There is winedbg (I just installed it) but I don't know how well that
will work. Adding print statements is the time-honored fallback (not
sure I'd call that 'reasonable' though).
> Are you using CGAL with wine? If not, then bbox NBL is actually a grid
> NBL and I think there were a couple of bugs in that which I fixed but
> never made it into IMP svn. If you are using CGAL with wine, how do I
> do that?
No, CGAL isn't built for Wine on synth yet, so IMP is a non-CGAL build.
Ben
--
ben(a)salilab.org http://salilab.org/~ben/
"It is a capital mistake to theorize before one has data."
- Sir Arthur Conan Doyle
2
3
No api changes. The internal evaluation functions have been split in two
so they are more generally applicable. Also more tests for the
transformed distnace pair score
Index: kernel/src/pair_scores/DistancePairScore.cpp
===================================================================
--- kernel/src/pair_scores/DistancePairScore.cpp (revision 672)
+++ kernel/src/pair_scores/DistancePairScore.cpp (working copy)
@@ -16,10 +16,6 @@
DistancePairScore::DistancePairScore(UnaryFunction *f): f_(f){}
-struct Identity
-{
- Float operator()(Float t) const {return t;}
-};
Float DistancePairScore::evaluate(Particle *a, Particle *b,
DerivativeAccumulator *da) const
Index: kernel/src/pair_scores/TransformedDistancePairScore.cpp
===================================================================
--- kernel/src/pair_scores/TransformedDistancePairScore.cpp (revision 672)
+++ kernel/src/pair_scores/TransformedDistancePairScore.cpp (working copy)
@@ -52,46 +52,62 @@
}
+Vector3D TransformedDistancePairScore
+::get_transformed(const Vector3D &v) const {
+ return Vector3D(get_transformed(v,0),
+ get_transformed(v,1),
+ get_transformed(v,2));
+}
+
+
+Float TransformedDistancePairScore
+::get_transformed(const Vector3D &v, unsigned int i) const {
+ IMP_assert(i<3, "Bad coordinate");
+ return (v-c_)*r_[i] + tc_[i];
+}
+
+
+Vector3D TransformedDistancePairScore
+::get_back_rotated(const Vector3D &v) const {
+ return Vector3D(v*ri_[0],
+ v*ri_[1],
+ v*ri_[2]);
+}
+
+
/*
Compute R(x-c)+tc and the reverse
dt= R(d-c)+tc-> Rt(dt-tc)+c
*/
-struct TransformParticle
+struct TransformedDistancePairScore::TransformParticle
{
- const Vector3D *r_, *ri_;
- const Vector3D &tc_;
- const Vector3D &c_;
+ const TransformedDistancePairScore *trd_;
XYZDecorator d_;
- TransformParticle(const Vector3D *r,
- const Vector3D *ri,
- const Vector3D &tc,
- const Vector3D &c,
- Particle *p): r_(r), ri_(ri),
- tc_(tc), c_(c), d_(p){}
+ TransformParticle(const TransformedDistancePairScore *trd,
+ Particle *p): trd_(trd), d_(p){}
Float get_coordinate(unsigned int i) const {
- return (d_.get_vector()-c_)*r_[i] + tc_[i];
+ return trd_->get_transformed(d_.get_vector(), i);
}
void add_to_coordinates_derivative(const Vector3D& f,
DerivativeAccumulator &da) {
- Vector3D r(f*ri_[0],
- f*ri_[1],
- f*ri_[2]);
- d_.add_to_coordinates_derivative(r, da);
+
+ d_.add_to_coordinates_derivative(trd_->get_back_rotated(f), da);
}
};
Float TransformedDistancePairScore::evaluate(Particle *a, Particle *b,
DerivativeAccumulator *da) const
{
- TransformParticle tb(r_,ri_, tc_,c_,b);
- IMP_LOG(VERBOSE, "Transformed particle is "
+ TransformParticle tb(this,b);
+ /*IMP_LOG(VERBOSE, "Transformed particle is "
<< tb.get_coordinate(0) << " " << tb.get_coordinate(1)
- << " " << tb.get_coordinate(2) << std::endl);
+ << " " << tb.get_coordinate(2)
+ << " from " << XYZDecorator(b) << std::endl);*/
return internal::evaluate_distance_pair_score(XYZDecorator(a),
tb,
- da, f_.get(),
+ da, f_.get(),
boost::lambda::_1);
}
@@ -118,24 +134,55 @@
<< ri_[2] << std::endl);
}
-void TransformedDistancePairScore::set_translation(float t0, float t1, float t2)
+void TransformedDistancePairScore::set_translation(const Vector3D &t)
{
- tc_= Vector3D(t0, t1, t2) + c_;
+ tc_= t + c_;
}
-void TransformedDistancePairScore::set_center(float t0, float t1, float t2)
+void TransformedDistancePairScore::set_center(const Vector3D &t)
{
tc_= tc_-c_;
- c_= Vector3D(t0, t1, t2);
+ c_= t;
tc_= tc_+c_;
}
+void TransformedDistancePairScore::set_rotation(const Vector3D &v, Float theta)
+{
+ IMP_check(v.get_magnitude() <1.1 && v.get_magnitude() > .9,
+ "Vector must be normalized", ValueException);
+ IMP_LOG(TERSE, "Computing rotation from " << v << " and "
+ << theta << std::endl);
+ Float c= std::cos(theta);
+ Float s= std::sin(theta);
+ Float C= 1-c;
+ Vector3D vs = v*s;
+ Vector3D vC = v*C;
+
+ Float xyC = v[0]*vC[1];
+ Float yzC = v[1]*vC[2];
+ Float zxC = v[2]*vC[0];
+ set_rotation(v[0]*vC[0] + c, xyC - vs[2], zxC + vs[1],
+ xyC + vs[2], v[1]*vC[1]+c, yzC - vs[0],
+ zxC - vs[1], yzC + vs[0], v[2]*vC[2] + c);
+}
+
+
void TransformedDistancePairScore::show(std::ostream &out) const
{
out << "TransformedDistancePairScore using ";
f_->show(out);
+ out << "With rotation:\n";
+ out << r_[0] << std::endl;
+ out << r_[1] << std::endl;
+ out << r_[2] << std::endl;
+ out << "and inverse rotation:\n";
+ out << ri_[0] << std::endl;
+ out << ri_[1] << std::endl;
+ out << ri_[2] << std::endl;
+ out << "and center: " << c_ << std::endl;
+ out << "and translation: " << tc_-c_ << std::endl;
}
} // namespace IMP
Index: kernel/include/IMP/pair_scores/TransformedDistancePairScore.h
===================================================================
--- kernel/include/IMP/pair_scores/TransformedDistancePairScore.h (revision 672)
+++ kernel/include/IMP/pair_scores/TransformedDistancePairScore.h (working copy)
@@ -25,13 +25,28 @@
The second particle, x, is transformed as R*(x-center)+ translation+center
+ \note This score is asymetric, that is which point is passed first or
+ second makes a difference. If you don't know which direction the symetry
+ goes make sure you try both directions.
+
+ \note While it would be nice to have this apply a general pair
+ score to the transformed particle, this would require either creating
+ a temporary particle or rewriting and then restoring the coordinates
+ of the particle being transformed. Neither of which sound appealing.
+
\ingroup pairscore
*/
class IMPDLLEXPORT TransformedDistancePairScore : public PairScore
{
+ class TransformParticle;
+ friend class TransformParticle;
Pointer<UnaryFunction> f_;
Vector3D tc_, c_;
Vector3D r_[3], ri_[3];
+
+ Float get_transformed(const Vector3D &v, unsigned int i) const ;
+ // Transform a derivative back
+ Vector3D get_back_rotated(const Vector3D &v) const ;
public:
TransformedDistancePairScore(UnaryFunction *f);
virtual ~TransformedDistancePairScore(){}
@@ -39,11 +54,23 @@
DerivativeAccumulator *da) const;
virtual void show(std::ostream &out=std::cout) const;
+ //! Set the rotation from a rotation matrix
+ /** The matrix is multiplied by a vector from the right.
+ */
void set_rotation(float r00, float r01, float r02,
float r10, float r11, float r12,
float r20, float r21, float r22);
- void set_translation(float t0, float t1, float t2);
- void set_center(float t0, float t1, float t2);
+ //! Set the rotation from an axis and an angle (radians)
+ /** The axis should be close to a unit vector.
+ */
+ void set_rotation(const Vector3D &axis, Float angle);
+
+ void set_translation(const Vector3D &t);
+ //! Set the rotation center
+ void set_center(const Vector3D &t);
+
+ //! Apply the current transformation to the given vector
+ Vector3D get_transformed(const Vector3D &v) const ;
};
} // namespace IMP
Index: kernel/include/IMP/internal/evaluate_distance_pair_score.h
===================================================================
--- kernel/include/IMP/internal/evaluate_distance_pair_score.h (revision 672)
+++ kernel/include/IMP/internal/evaluate_distance_pair_score.h (working copy)
@@ -9,6 +9,8 @@
#define __IMP_EVALUATE_DISTANCE_PAIR_SCORE_H
#include "../Vector3D.h"
+#include <boost/tuple/tuple.hpp>
+#include <IMP/UnaryFunction.h>
namespace IMP
{
@@ -16,41 +18,51 @@
namespace internal
{
+template <class SD>
+Float compute_distance_pair_score(const Vector3D &delta,
+ const UnaryFunction *f,
+ Vector3D *d,
+ SD sd) {
+ static const Float MIN_DISTANCE = .00001;
+ Float distance= delta.get_magnitude();
+ Float shifted_distance = sd(distance);
+
+ // if needed, calculate the partial derivatives of the scores with respect
+ // to the particle attributes
+ // avoid division by zero if the distance is too small
+ Float score, deriv;
+ if (d && distance >= MIN_DISTANCE) {
+ boost::tie(score, deriv) = f->evaluate_with_derivative(shifted_distance);
+
+ *d= delta/distance *deriv;
+ } else {
+ // calculate the score based on the distance feature
+ score = f->evaluate(shifted_distance);
+ }
+ return score;
+}
+
+
template <class W0, class W1, class SD>
Float evaluate_distance_pair_score(W0 d0, W1 d1,
DerivativeAccumulator *da,
- UnaryFunction *f, SD sd)
+ const UnaryFunction *f, SD sd)
{
- static const Float MIN_DISTANCE = .00001;
IMP_CHECK_OBJECT(f);
- Float d2 = 0;
Vector3D delta;
- Float score;
for (int i = 0; i < 3; ++i) {
delta[i] = d0.get_coordinate(i) - d1.get_coordinate(i);
- d2 += square(delta[i]);
}
- Float distance = std::sqrt(d2);
+ Vector3D d;
+ Float score= compute_distance_pair_score(delta, f, (da? &d : NULL), sd);
- Float shifted_distance = sd(distance); //scale*(distance - offset);
- // if needed, calculate the partial derivatives of the scores with respect
- // to the particle attributes
- // avoid division by zero if the distance is too small
- if (da && distance >= MIN_DISTANCE) {
- Float deriv;
-
- score = f->evaluate_with_derivative(shifted_distance, deriv);
-
- Vector3D d= delta/distance *deriv;
+ if (da) {
d0.add_to_coordinates_derivative(d, *da);
d1.add_to_coordinates_derivative(-d, *da);
- } else {
- // calculate the score based on the distance feature
- score = f->evaluate(shifted_distance);
}
return score;
Index: kernel/test/pair_scores/test_transform.py
===================================================================
--- kernel/test/pair_scores/test_transform.py (revision 672)
+++ kernel/test/pair_scores/test_transform.py (working copy)
@@ -4,9 +4,10 @@
class DistanceTests(IMP.test.TestCase):
"""Test the symmetry restraint"""
- def test_symmetry(self):
+ def _test_symmetry(self):
"""Test the transform pair score basics"""
IMP.set_log_level(IMP.VERBOSE)
+ random.seed()
m= IMP.Model()
p0= IMP.Particle()
m.add_particle(p0)
@@ -14,7 +15,7 @@
p1= IMP.Particle()
m.add_particle(p1)
d1= IMP.XYZDecorator.create(p1)
- tps= IMP.TransformedDistancePairScore(IMP.Harmonic(0,1))
+ tps= IMP.dsr.TransformedDistancePairScore(IMP.Harmonic(0,1))
tps.set_translation(0,1,0)
tps.set_rotation(1, 0, 0,
0, 1, 0,
@@ -56,7 +57,7 @@
self.assert_(d1.get_coordinate_derivative(1) > 0)
self.assert_(d1.get_coordinate_derivative(0) == 0)
self.assert_(d1.get_coordinate_derivative(2) == 0)
- def test_symmetry2(self):
+ def _test_symmetry2(self):
"""Test the transform pair score optimization"""
IMP.set_log_level(IMP.VERBOSE)
m= IMP.Model()
@@ -70,27 +71,122 @@
d1.set_coordinates(IMP.Vector3D(20,20,40))
d0.set_coordinates_are_optimized(True)
d1.set_coordinates_are_optimized(True)
- tps= IMP.TransformedDistancePairScore(IMP.Harmonic(0,1))
- tps.set_translation(0,1,0)
+ tps= IMP.dsr.TransformedDistancePairScore(IMP.Harmonic(0,1))
+ cv= IMP.Vector3D(2,1,0)
+ tv= IMP.Vector3D(0,1,0)
+ tps.set_translation(tv[0], tv[1], tv[2])
+ tps.set_center(cv[0], cv[1], cv[2])
tps.set_rotation(1, 0, 0,
0, 0,-1,
0, 1, 0)
- pr= IMP.PairListRestraint(tps)
+ pr= IMP.dsr.PairListRestraint(tps)
pr.add_particle_pair(IMP.ParticlePair(p0, p1))
m.add_restraint(pr)
- cg= IMP.ConjugateGradients()
+ cg= IMP.SteepestDescent()
cg.set_model(m)
+ IMP.set_log_level(IMP.SILENT)
cg.optimize(100)
+ IMP.set_log_level(IMP.VERBOSE)
d0.show()
d1.show()
- vt= IMP.Vector3D(d1.get_vector()*IMP.Vector3D(1,0,0)+0,
- d1.get_vector()*IMP.Vector3D(0,0,-1)+1,
- d1.get_vector()*IMP.Vector3D(0,1,0)+0)
- print "trans"
- print str(vt[0]) + " " + str(vt[1])+" " + str(vt[2])
- self.assertEqual(vt[0], d0.get_coordinate(0))
- self.assertEqual(vt[1], d0.get_coordinate(1))
- self.assertEqual(vt[2], d0.get_coordinate(2))
+ ncv= d1.get_vector()- cv
+ vt= IMP.Vector3D(ncv*IMP.Vector3D(1,0,0),
+ ncv*IMP.Vector3D(0,0,-1),
+ ncv*IMP.Vector3D(0,1,0))
+ fvt= vt+ tv+cv
+ print "itrans: "
+ tps.get_transformed(d1.get_vector()).show()
+ print "trans: "
+ fvt.show()
+ self.assertInTolerance(fvt[0], d0.get_coordinate(0), .1)
+ self.assertInTolerance(fvt[1], d0.get_coordinate(1), .1)
+ self.assertInTolerance(fvt[2], d0.get_coordinate(2), .1)
+ def _setup_system(self):
+ IMP.set_log_level(IMP.VERBOSE)
+ m= IMP.Model()
+ random.seed(100)
+ p0= IMP.Particle()
+ m.add_particle(p0)
+ d0= IMP.XYZDecorator.create(p0)
+ p1= IMP.Particle()
+ m.add_particle(p1)
+ d1= IMP.XYZDecorator.create(p1)
+ d0.set_coordinates(IMP.random_vector_in_box(IMP.Vector3D(-10, -10, -10),
+ IMP.Vector3D(10,10,10)))
+ d1.set_coordinates(IMP.random_vector_in_box(IMP.Vector3D(-10, -10, -10),
+ IMP.Vector3D(10,10,10)))
+ d0.set_coordinates_are_optimized(True)
+ d1.set_coordinates_are_optimized(True)
+ tps= IMP.dsr.TransformedDistancePairScore(IMP.Harmonic(0,1))
+ cv= IMP.random_vector_in_box(IMP.Vector3D(-10, -10, -10),
+ IMP.Vector3D(10,10,10))
+ tv= IMP.random_vector_in_box(IMP.Vector3D(-10, -10, -10),
+ IMP.Vector3D(10,10,10))
+ tps.set_translation(tv)
+ tps.set_center(cv)
+ tps.set_rotation(IMP.random_vector_on_sphere(), 2.0*math.pi*random.random())
+ pr= IMP.dsr.PairListRestraint(tps)
+ pr.add_particle_pair(IMP.ParticlePair(p0, p1))
+ m.add_restraint(pr)
+ return (m, tps, d0, d1)
+
+
+ def _test_transoformish(self):
+ """Test optimizing with monte carlo to test the non-deriv part"""
+ IMP.set_log_level(IMP.VERBOSE)
+ (m, tps, d0, d1)= self._setup_system()
+
+ cg= IMP.dsr.MonteCarlo()
+ cg.set_model(m)
+ cg.set_temperature(0)
+ bm= IMP.dsr.BallMover(IMP.XYZDecorator.get_xyz_keys(),
+ 1.0, IMP.Particles(1, d0.get_particle()))
+ cg.add_mover(bm)
+ IMP.set_log_level(IMP.SILENT)
+ cg.optimize(2000)
+ IMP.set_log_level(IMP.VERBOSE)
+ d0.get_vector().show(); print
+ d1.get_vector().show(); print
+ self.assertInTolerance((d0.get_vector()
+ -tps.get_transformed(d1.get_vector())).get_magnitude(), 0,
+ 1)
+
+ def _test_transoformish2(self):
+ """Test optimizing with steepest descent for testing derivatives part"""
+ IMP.set_log_level(IMP.VERBOSE)
+ (m, tps, d0, d1)= self._setup_system()
+ tps.show()
+ cg= IMP.SteepestDescent()
+ cg.set_model(m)
+ IMP.set_log_level(IMP.SILENT)
+ cg.optimize(2000)
+ IMP.set_log_level(IMP.VERBOSE)
+ print "p0 is "
+ d0.get_vector().show(); print
+ print "p1 is "
+ d1.get_vector().show(); print
+ print "target is"
+ tps.get_transformed(d1.get_vector()).show(); print
+ self.assertInTolerance((d0.get_vector()
+ -tps.get_transformed(d1.get_vector())).get_magnitude(), 0,
+ 1)
+
+
+ def _test_set_rot(self):
+ """Testing setting the rotation in tps"""
+ (m, tps, d0, d1)= self._setup_system()
+ v= IMP.random_vector_on_sphere()
+ a= random.random()*2.0* math.pi
+ tps.set_rotation(v, a)
+ tv= tps.get_transformed(d0.get_vector())
+ tps.set_rotation(-v, -a)
+ tv2= tps.get_transformed(d0.get_vector())
+ self.assertInTolerance((tv-tv2).get_magnitude(), 0, .1);
+
+
+
+
+
if __name__ == '__main__':
unittest.main()
2
3
My guess is that r680 (Daniel's attribute table simplification) is to
blame here. A little pruning of the testcase tree (rm -rf of some
directories, disabling of some test methods) yields the following
minimal testset: particles/test_particles.py (test_remove_attributes)
and states/test_nonbonded_list.py (test_bipartite_quadratic and
test_bipartite_bbox):
synth% scons wine=true test
...
Test bipartite bbox NBL ... ok
Test bipartite quadratic NBL ... ok
Test that attributes can be removed ... wine: Unhandled page fault on
read access to 0x206fdbb0 at address 0x8d005a (thread 0009), starting
debugger...
Ben
-------- Original Message --------
Subject: Build system failure, 08/28/2008
Date: Thu, 28 Aug 2008 07:00:02 -0700
From: build-system(a)salilab.org
To: ben(a)salilab.org
Build system failed on 08/28/2008.
Please see https://salilab.org/internal/imp/tests.html for full details.
--
ben(a)salilab.org http://salilab.org/~ben/
"It is a capital mistake to theorize before one has data."
- Sir Arthur Conan Doyle
2
3
Note that this also exposed a minor bug in Modeller's own Python init
script, so to run the Modeller tests in IMP, update Modeller to r6455 or
later.
Ben
--
ben(a)salilab.org http://salilab.org/~ben/
"It is a capital mistake to theorize before one has data."
- Sir Arthur Conan Doyle
1
0
I don't feel that other people have really given serious thought to what
various proposals for reorganization of the repository would entail and
has included far to many wild assertions about what the outcomes of
various changes would be. I don't know who else people think is involved
in the project, but as far as I know, everyone involved is pretty
reasonable, able to follow directions, can communicate over email and
aware of their limitations.
There are three disparate issues to consider:
- the location of code in the svn repository: I propose that it be kept
as simple as possible. The only advantage gained by splitting things
into separate modules is a decrease in compilation time, which given the
size of the swig files, might be worthwhile. Many modules means more
things that have to be tested together, and makes it harder to reuse
code and harder to discover new code. My favoured solution would be an
IMP_kernel as discussed and an IMP_code which has all the current
restraints and things (including domino and em).
- access control to the repository: we are all pretty reasonable, able
to follow directions, can communicate over email and aware of our
limitations, so I don't see any reason to impose access restrictions on
anyone. As we saw in the last couple of days, adding layers of
complication to submission increases the chance of bugs being introduced
via the submission process. I haven't gone back and checked, but I am
pretty sure Ben has broken the repository by changing patches more times
than the number of bugs he has found by reviewing the code. In addition,
with access control, we run into problems like when Ben was on vacation
and Frido discovered a bug in the NonbondedLists and no one could commit
a patch.
- how bits of code are designated stable or unstable: I suggest that
this be through deoxygen groups. This makes it easy to move things from
unstable to stable without break people's code (which would break if it
involved changing modules). In addition, it ensures that, except in
unusual cases, there is only one version of a particular class, reducing
confusion.
All the setups which have been discussed would provide a set of clearly
marked stable functionality which doesn't change much. So raising issues
related to that is silly. In addition, conflating access to svn, the
location of code and their stability simply confuses the discussion.
Something should be done about these decisions sooner rather than later,
as there will only be more code that has to be updated later.
Anyway, I don't think I have anything more to say on the subject. I'm
off with my brother the next couple of days. If a solution which allows
me to reasonable commit changes to the code I contributed is not in
place by Tuesday, I think I will give up on IMP and stop submitting patches.
1
0
Most of the restraints and other useful functionality are moving out of
the kernel. Currently the idea is that they go into a module called dsr
(for lack of another name) which would be modifiable by everyone. Only
the angle restraints and a couple of other non-kernel objects would be
left in the kernel. I propose instead we call it "basic" and move all
existing restraints into it. That way we have less random split of code
between modules. After all, does anyone want to remember that
"DistanceRestraint" is in IMP but "DistancePairScore" or
"NonbondedListScoreState" is part of IMP.dsr?
In either case, the npc-specific code would stay elsewhere.
3
10
In order to avoid going back and forth a few times, here is a patch with
all the changed files
The contained patches are:
kernel/include/IMP/optimizers/movers/BallMover.h
kernel/src/IMP/optimizers/movers/BallMover.cpp
Somewhere between my copy of IMP and the general copy of IMP BallMover
became completely broken due to a missing &, showing why returning by
reference arguments is problematic :-)
----
kernel/pyext/IMP.i
For uniformity we should make all containers look the same. And Mover
has virtual methods so they
should get directors.
-----
kernel/include/IMP/optimizers/MonteCarlo.h
kernel/include/IMP/optimizers/Mover.h
Movers should probably be declared in Mover.h.
-----
kernel/include/IMP/utility.h
kernel/include/IMP/internal/AttributeTable.h
isnan function to make it clearer what we are testing for when we test
things
----
kernel/test/states/test_nonbonded_list.py
kernel/test/states/test_cover_bonds.py
kernel/pyext/IMP/test.py
kernel/src/Vector3D.cpp
kernel/src/decorators/XYZDecorator.cpp
kernel/src/SConscript
kernel/doc/examples/chain.py
kernel/include/IMP/decorators/macros.h
kernel/include/IMP/decorators/XYZDecorator.h
Move the random functions from XYZDecorator to Vector3D.
----
kernel/src/optimizers/MonteCarlo.cpp
don't try to get a null pointer
-----
kernel/test/pair_scores/test_transform.py
kernel/src/pair_scores/TransformedDistancePairScore.cpp
kernel/src/pair_scores/DistancePairScore.cpp
kernel/include/IMP/pair_scores/TransformedDistancePairScore.h
kernel/include/IMP/internal/evaluate_distance_pair_score.h
Clean up the evaluate distance pair score functions.
Index: kernel/test/states/test_nonbonded_list.py
===================================================================
--- kernel/test/states/test_nonbonded_list.py (revision 681)
+++ kernel/test/states/test_nonbonded_list.py (working copy)
@@ -86,7 +86,8 @@
for i in range(0,10):
score= m.evaluate(False)
for d in ds:
- d.randomize_in_sphere(d.get_vector(), 2.0)
+ d.set_coordinates(IMP.random_vector_in_sphere(d.get_vector(),
+ 2.0))
def do_test_bl(self, ss):
"""Test the bond decorator list"""
@@ -148,8 +149,8 @@
print "Index is " +str(p.get_index().get_index())
d=IMP.XYZDecorator.create(p)
d.set_coordinates_are_optimized(True)
- d.randomize_in_box(IMP.Vector3D(0,0,0),
- IMP.Vector3D(10,10,10));
+ d.set_coordinates(IMP.random_vector_in_box(IMP.Vector3D(0,0,0),
+ IMP.Vector3D(10,10,10)))
nps= IMP.Particles([p])
s.add_particles(nps)
score= m.evaluate(False)
@@ -193,10 +194,12 @@
score= m.evaluate(False)
for p in ps0:
d= IMP.XYZDecorator.cast(p)
- d.randomize_in_sphere(d.get_vector(), 1)
+ d.set_coordinates(IMP.random_vector_in_sphere(d.get_vector(),
+ 1))
for p in ps1:
d= IMP.XYZDecorator.cast(p)
- d.randomize_in_sphere(d.get_vector(), 1)
+ d.set_coordinates(IMP.random_vector_in_sphere(d.get_vector(),
+ 1))
def do_test_spheres(self, ss):
Index: kernel/test/states/test_cover_bonds.py
===================================================================
--- kernel/test/states/test_cover_bonds.py (revision 681)
+++ kernel/test/states/test_cover_bonds.py (working copy)
@@ -15,8 +15,8 @@
p= IMP.Particle()
m.add_particle(p)
d= IMP.XYZDecorator.create(p)
- d.randomize_in_box(IMP.Vector3D(0,0,0),
- IMP.Vector3D(10,10,10))
+ d.set_coordinates(IMP.random_vector_in_box(IMP.Vector3D(0,0,0),
+ IMP.Vector3D(10,10,10)))
ps.append(p)
bds= []
bb= IMP.BondedDecorator.create(ps[0])
Index: kernel/test/pair_scores/test_transform.py
===================================================================
--- kernel/test/pair_scores/test_transform.py (revision 681)
+++ kernel/test/pair_scores/test_transform.py (working copy)
@@ -1,6 +1,8 @@
import unittest
import IMP.utils
import IMP.test, IMP
+import math
+import random
class DistanceTests(IMP.test.TestCase):
"""Test the symmetry restraint"""
@@ -15,7 +17,7 @@
m.add_particle(p1)
d1= IMP.XYZDecorator.create(p1)
tps= IMP.TransformedDistancePairScore(IMP.Harmonic(0,1))
- tps.set_translation(0,1,0)
+ tps.set_translation(IMP.Vector3D(0,1,0))
tps.set_rotation(1, 0, 0,
0, 1, 0,
0, 0, 1)
@@ -23,10 +25,10 @@
d1.set_coordinates(IMP.Vector3D(2,2,4))
self.assertEqual(tps.evaluate(p0, p1, None), 0)
self.assert_(tps.evaluate(p1, p0, None) != 0)
- tps.set_center(10,13,4)
+ tps.set_center(IMP.Vector3D(10,13,4))
self.assertEqual(tps.evaluate(p0, p1, None), 0)
self.assert_(tps.evaluate(p1, p0, None) != 0)
- tps.set_center(0,0,0)
+ tps.set_center(IMP.Vector3D(0,0,0))
print "test rotation"
tps.set_rotation(0, 0,-1,
0, 1, 0,
@@ -34,7 +36,7 @@
d1.set_coordinates(IMP.Vector3D(4, 2, -2))
self.assertEqual(tps.evaluate(p0, p1, None), 0)
self.assert_(tps.evaluate(p1, p0, None) != 0)
- tps.set_translation(0,0,0)
+ tps.set_translation(IMP.Vector3D(0,0,0))
tps.set_rotation(0,-1, 0,
1, 0, 0,
0, 0, 1)
@@ -71,26 +73,121 @@
d0.set_coordinates_are_optimized(True)
d1.set_coordinates_are_optimized(True)
tps= IMP.TransformedDistancePairScore(IMP.Harmonic(0,1))
- tps.set_translation(0,1,0)
+ cv= IMP.Vector3D(2,1,0)
+ tv= IMP.Vector3D(0,1,0)
+ tps.set_translation(tv)
+ tps.set_center(cv)
tps.set_rotation(1, 0, 0,
0, 0,-1,
0, 1, 0)
pr= IMP.PairListRestraint(tps)
pr.add_particle_pair(IMP.ParticlePair(p0, p1))
m.add_restraint(pr)
- cg= IMP.ConjugateGradients()
+ cg= IMP.SteepestDescent()
cg.set_model(m)
+ IMP.set_log_level(IMP.SILENT)
cg.optimize(100)
+ IMP.set_log_level(IMP.VERBOSE)
d0.show()
d1.show()
- vt= IMP.Vector3D(d1.get_vector()*IMP.Vector3D(1,0,0)+0,
- d1.get_vector()*IMP.Vector3D(0,0,-1)+1,
- d1.get_vector()*IMP.Vector3D(0,1,0)+0)
- print "trans"
- print str(vt[0]) + " " + str(vt[1])+" " + str(vt[2])
- self.assertEqual(vt[0], d0.get_coordinate(0))
- self.assertEqual(vt[1], d0.get_coordinate(1))
- self.assertEqual(vt[2], d0.get_coordinate(2))
+ ncv= d1.get_vector()- cv
+ vt= IMP.Vector3D(ncv*IMP.Vector3D(1,0,0),
+ ncv*IMP.Vector3D(0,0,-1),
+ ncv*IMP.Vector3D(0,1,0))
+ fvt= vt+ tv+cv
+ print "itrans: "
+ tps.get_transformed(d1.get_vector()).show()
+ print "trans: "
+ fvt.show()
+ self.assertInTolerance(fvt[0], d0.get_coordinate(0), .1)
+ self.assertInTolerance(fvt[1], d0.get_coordinate(1), .1)
+ self.assertInTolerance(fvt[2], d0.get_coordinate(2), .1)
+ def _setup_system(self):
+ random.seed()
+ IMP.set_log_level(IMP.VERBOSE)
+ m= IMP.Model()
+ p0= IMP.Particle()
+ m.add_particle(p0)
+ d0= IMP.XYZDecorator.create(p0)
+ p1= IMP.Particle()
+ m.add_particle(p1)
+ d1= IMP.XYZDecorator.create(p1)
+ d0.set_coordinates(IMP.random_vector_in_box(IMP.Vector3D(-10, -10, -10),
+ IMP.Vector3D(10,10,10)))
+ d1.set_coordinates(IMP.random_vector_in_box(IMP.Vector3D(-10, -10, -10),
+ IMP.Vector3D(10,10,10)))
+ d0.set_coordinates_are_optimized(True)
+ d1.set_coordinates_are_optimized(True)
+ tps= IMP.TransformedDistancePairScore(IMP.Harmonic(0,1))
+ cv= IMP.random_vector_in_box(IMP.Vector3D(-10, -10, -10),
+ IMP.Vector3D(10,10,10))
+ tv= IMP.random_vector_in_box(IMP.Vector3D(-10, -10, -10),
+ IMP.Vector3D(10,10,10))
+ tps.set_translation(tv)
+ tps.set_center(cv)
+ tps.set_rotation(IMP.random_vector_on_sphere(), 2.0*math.pi*random.random())
+ pr= IMP.PairListRestraint(tps)
+ pr.add_particle_pair(IMP.ParticlePair(p0, p1))
+ m.add_restraint(pr)
+ return (m, tps, d0, d1)
+
+
+ def test_transoformish(self):
+ """Test optimizing with monte carlo to test the non-deriv part"""
+ IMP.set_log_level(IMP.VERBOSE)
+ (m, tps, d0, d1)= self._setup_system()
+
+ cg= IMP.MonteCarlo()
+ cg.set_model(m)
+ cg.set_temperature(0.0001)
+ bm= IMP.BallMover(IMP.XYZDecorator.get_xyz_keys(),
+ 1.0, IMP.Particles(1, d0.get_particle()))
+ cg.add_mover(bm)
+ IMP.set_log_level(IMP.SILENT)
+ cg.optimize(2000)
+ IMP.set_log_level(IMP.VERBOSE)
+ d0.get_vector().show(); print
+ d1.get_vector().show(); print
+ self.assertInTolerance((d0.get_vector()
+ -tps.get_transformed(d1.get_vector())).get_magnitude(), 0,
+ 1)
+
+ def test_transoformish2(self):
+ """Test optimizing with steepest descent for testing derivatives part"""
+ IMP.set_log_level(IMP.VERBOSE)
+ (m, tps, d0, d1)= self._setup_system()
+ tps.show()
+ cg= IMP.SteepestDescent()
+ cg.set_model(m)
+ IMP.set_log_level(IMP.SILENT)
+ cg.optimize(2000)
+ IMP.set_log_level(IMP.VERBOSE)
+ print "p0 is "
+ d0.get_vector().show(); print
+ print "p1 is "
+ d1.get_vector().show(); print
+ print "target is"
+ tps.get_transformed(d1.get_vector()).show(); print
+ self.assertInTolerance((d0.get_vector()
+ -tps.get_transformed(d1.get_vector())).get_magnitude(), 0,
+ 1)
+
+
+ def test_set_rot(self):
+ """Testing setting the rotation in tps"""
+ (m, tps, d0, d1)= self._setup_system()
+ v= IMP.random_vector_on_sphere()
+ a= random.random()*2.0* math.pi
+ tps.set_rotation(v, a)
+ tv= tps.get_transformed(d0.get_vector())
+ tps.set_rotation(-v, -a)
+ tv2= tps.get_transformed(d0.get_vector())
+ self.assertInTolerance((tv-tv2).get_magnitude(), 0, .1);
+
+
+
+
+
if __name__ == '__main__':
unittest.main()
Index: kernel/include/IMP/optimizers/MonteCarlo.h
===================================================================
--- kernel/include/IMP/optimizers/MonteCarlo.h (revision 681)
+++ kernel/include/IMP/optimizers/MonteCarlo.h (working copy)
@@ -15,8 +15,6 @@
namespace IMP
{
-typedef std::vector<Mover*> Movers;
-
//! A Monte Carlo optimizer.
/** The optimizer uses a set of Mover objects to propose steps. Currently
each Mover is called at each Monte Carlo iteration. This may change in
Index: kernel/include/IMP/optimizers/movers/BallMover.h
===================================================================
--- kernel/include/IMP/optimizers/movers/BallMover.h (revision 681)
+++ kernel/include/IMP/optimizers/movers/BallMover.h (working copy)
@@ -42,6 +42,8 @@
Float get_radius() const {
return radius_;
}
+
+ void show(std::ostream &out= std::cout) const;
protected:
/** \internal */
void generate_move(float a);
Index: kernel/include/IMP/optimizers/Mover.h
===================================================================
--- kernel/include/IMP/optimizers/Mover.h (revision 681)
+++ kernel/include/IMP/optimizers/Mover.h (working copy)
@@ -76,6 +76,8 @@
IMP_OUTPUT_OPERATOR(Mover);
+typedef std::vector<Mover*> Movers;
+
} // namespace IMP
#endif /* __IMP_MOVER_H */
Index: kernel/include/IMP/decorators/XYZDecorator.h
===================================================================
--- kernel/include/IMP/decorators/XYZDecorator.h (revision 681)
+++ kernel/include/IMP/decorators/XYZDecorator.h (working copy)
@@ -8,16 +8,13 @@
#ifndef __IMP_XYZ_DECORATOR_H
#define __IMP_XYZ_DECORATOR_H
-#include <vector>
-#include <deque>
-#include <limits>
-
-#include "../Particle.h"
-#include "../Model.h"
#include "../DecoratorBase.h"
#include "../Vector3D.h"
#include "utility.h"
+#include <vector>
+#include <limits>
+
namespace IMP
{
@@ -104,23 +101,23 @@
/** 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());
+ return Vector3D(get_x(), get_y(), get_z());
}
+ //! Get the vector of derivatives.
+ /** Somewhat suspect based on wanting a Point/Vector differentiation
+ but we don't have points */
+ Vector3D get_derivative_vector() const {
+ return Vector3D(get_coordinate_derivative(0),
+ get_coordinate_derivative(1),
+ get_coordinate_derivative(2));
+ }
+
//! Get a vector containing the keys for x,y,z
/** This is quite handy for initializing movers and things.
*/
- static const FloatKeys get_xyz_keys() {
- decorator_initialize_static_data();
- return key_;
- }
+ IMP_DECORATOR_GET_KEY(FloatKeys, xyz_keys, key_)
- //! Generate random coordinates in a sphere centered at the vector
- void randomize_in_sphere(const Vector3D ¢er, float radius);
-
- //! Generate random coordinates in a box defined by the vectors
- void randomize_in_box(const Vector3D &lower_corner,
- const Vector3D &upper_corner);
protected:
static FloatKey get_coordinate_key(unsigned int i) {
IMP_check(i <3, "Out of range coordinate",
Index: kernel/include/IMP/decorators/macros.h
===================================================================
--- kernel/include/IMP/decorators/macros.h (revision 681)
+++ kernel/include/IMP/decorators/macros.h (working copy)
@@ -44,7 +44,7 @@
} \
add_required; \
} \
- friend class DecoratorBase; \
+ friend class ::IMP::DecoratorBase; \
public: \
typedef Name This; \
/** \short The default constructor. This is used as a null value */ \
@@ -58,14 +58,14 @@
} \
/** Add the necessary attributes to p and return a decorator. */ \
static Name create(::IMP::Particle *p) { \
- return IMP::DecoratorBase::create<Name>(p); \
+ return ::IMP::DecoratorBase::create<Name>(p); \
} \
/** Check that p has the necessary attributes and return a decorator. \
\throws InvalidStateException if some required attributes are \
missing \
*/ \
static Name cast(::IMP::Particle *p) { \
- return IMP::DecoratorBase::cast<Name>(p); \
+ return ::IMP::DecoratorBase::cast<Name>(p); \
} \
static bool is_instance_of(::IMP::Particle *p) { \
decorator_initialize_static_data(); \
@@ -262,6 +262,17 @@
name##_data_.initialize();
+//! add a method to get a key
+/** One has to make sure to call the
+ decorator_initialize_static_data method first
+ */
+#define IMP_DECORATOR_GET_KEY(KeyType, key_name, variable_name)\
+ static KeyType get_##key_name() { \
+ decorator_initialize_static_data(); \
+ return variable_name; \
+ }
+
+
#define IMP_ATOM_TYPE_INDEX 8974343
#define IMP_RESIDUE_TYPE_INDEX 90784334
Index: kernel/include/IMP/pair_scores/TransformedDistancePairScore.h
===================================================================
--- kernel/include/IMP/pair_scores/TransformedDistancePairScore.h (revision 681)
+++ kernel/include/IMP/pair_scores/TransformedDistancePairScore.h (working copy)
@@ -25,13 +25,28 @@
The second particle, x, is transformed as R*(x-center)+ translation+center
+ \note This score is asymetric, that is which point is passed first or
+ second makes a difference. If you don't know which direction the symetry
+ goes make sure you try both directions.
+
+ \note While it would be nice to have this apply a general pair
+ score to the transformed particle, this would require either creating
+ a temporary particle or rewriting and then restoring the coordinates
+ of the particle being transformed. Neither of which sound appealing.
+
\ingroup pairscore
*/
class IMPDLLEXPORT TransformedDistancePairScore : public PairScore
{
+ class TransformParticle;
+ friend class TransformParticle;
Pointer<UnaryFunction> f_;
Vector3D tc_, c_;
Vector3D r_[3], ri_[3];
+
+ Float get_transformed(const Vector3D &v, unsigned int i) const ;
+ // Transform a derivative back
+ Vector3D get_back_rotated(const Vector3D &v) const ;
public:
TransformedDistancePairScore(UnaryFunction *f);
virtual ~TransformedDistancePairScore(){}
@@ -39,11 +54,23 @@
DerivativeAccumulator *da) const;
virtual void show(std::ostream &out=std::cout) const;
+ //! Set the rotation from a rotation matrix
+ /** The matrix is multiplied by a vector from the right.
+ */
void set_rotation(float r00, float r01, float r02,
float r10, float r11, float r12,
float r20, float r21, float r22);
- void set_translation(float t0, float t1, float t2);
- void set_center(float t0, float t1, float t2);
+ //! Set the rotation from an axis and an angle (radians)
+ /** The axis should be close to a unit vector.
+ */
+ void set_rotation(const Vector3D &axis, Float angle);
+
+ void set_translation(const Vector3D &t);
+ //! Set the rotation center
+ void set_center(const Vector3D &t);
+
+ //! Apply the current transformation to the given vector
+ Vector3D get_transformed(const Vector3D &v) const ;
};
} // namespace IMP
Index: kernel/include/IMP/internal/AttributeTable.h
===================================================================
--- kernel/include/IMP/internal/AttributeTable.h (revision 681)
+++ kernel/include/IMP/internal/AttributeTable.h (working copy)
@@ -47,7 +47,7 @@
}
static bool get_is_valid(Float f) {
if (std::numeric_limits<Float>::has_quiet_NaN) {
- return f==f;
+ return !is_nan(f);
} else {
return f != get_invalid();
}
Index: kernel/include/IMP/internal/evaluate_distance_pair_score.h
===================================================================
--- kernel/include/IMP/internal/evaluate_distance_pair_score.h (revision 681)
+++ kernel/include/IMP/internal/evaluate_distance_pair_score.h (working copy)
@@ -17,41 +17,51 @@
namespace internal
{
+template <class SD>
+Float compute_distance_pair_score(const Vector3D &delta,
+ const UnaryFunction *f,
+ Vector3D *d,
+ SD sd) {
+ static const Float MIN_DISTANCE = .00001;
+ Float distance= delta.get_magnitude();
+ Float shifted_distance = sd(distance);
+
+ // if needed, calculate the partial derivatives of the scores with respect
+ // to the particle attributes
+ // avoid division by zero if the distance is too small
+ Float score, deriv;
+ if (d && distance >= MIN_DISTANCE) {
+ boost::tie(score, deriv) = f->evaluate_with_derivative(shifted_distance);
+
+ *d= delta/distance *deriv;
+ } else {
+ if (d) *d= Vector3D(0,0,0);
+ // calculate the score based on the distance feature
+ score = f->evaluate(shifted_distance);
+ }
+ return score;
+}
+
+
template <class W0, class W1, class SD>
Float evaluate_distance_pair_score(W0 d0, W1 d1,
DerivativeAccumulator *da,
- UnaryFunction *f, SD sd)
+ const UnaryFunction *f, SD sd)
{
- static const Float MIN_DISTANCE = .00001;
IMP_CHECK_OBJECT(f);
- Float d2 = 0;
Vector3D delta;
- Float score;
for (int i = 0; i < 3; ++i) {
delta[i] = d0.get_coordinate(i) - d1.get_coordinate(i);
- d2 += square(delta[i]);
}
- Float distance = std::sqrt(d2);
+ Vector3D d;
+ Float score= compute_distance_pair_score(delta, f, (da? &d : NULL), sd);
- Float shifted_distance = sd(distance); //scale*(distance - offset);
-
- // if needed, calculate the partial derivatives of the scores with respect
- // to the particle attributes
- // avoid division by zero if the distance is too small
- if (da && distance >= MIN_DISTANCE) {
- Float deriv;
-
- boost::tie(score, deriv) = f->evaluate_with_derivative(shifted_distance);
-
- Vector3D d= delta/distance *deriv;
+ if (da) {
d0.add_to_coordinates_derivative(d, *da);
d1.add_to_coordinates_derivative(-d, *da);
- } else {
- // calculate the score based on the distance feature
- score = f->evaluate(shifted_distance);
}
return score;
Index: kernel/include/IMP/Vector3D.h
===================================================================
--- kernel/include/IMP/Vector3D.h (revision 681)
+++ kernel/include/IMP/Vector3D.h (working copy)
@@ -11,6 +11,7 @@
#include "IMP_config.h"
#include "base_types.h"
#include "macros.h"
+#include "exception.h"
#include <cmath>
@@ -22,7 +23,11 @@
*/
class IMPDLLEXPORT Vector3D
{
+ bool is_default() const {return false;}
public:
+ // public for swig
+ typedef Vector3D This;
+
//! Initialize the vector from separate x,y,z values.
Vector3D(Float x, Float y, Float z) {
vec_[0] = x;
@@ -33,6 +38,8 @@
//! Default constructor
Vector3D() {}
+ IMP_COMPARISONS_3(vec_[0], vec_[1], vec_[2]);
+
//! \return A single component of this vector (0-2).
Float operator[](unsigned int i) const {
IMP_assert(i < 3, "Invalid component of vector requested");
@@ -149,14 +156,6 @@
<< operator[](2) << ")";
}
- bool operator<(const Vector3D &o) const {
- for (unsigned int i=0; i< 3; ++i) {
- if (operator[](i) < o[i]) return true;
- else if (operator[](i) > o[i]) return false;
- }
- return false;
- }
-
private:
Float vec_[3];
};
@@ -171,6 +170,53 @@
o[2]*s);
}
+//! Generate a random vector in a box with uniform density
+/**
+ \ingroup uncommitted
+ */
+IMPDLLEXPORT Vector3D
+random_vector_in_box(const Vector3D &lb=Vector3D(0,0,0),
+ const Vector3D &ub=Vector3D(10,10,10));
+
+//! Generate a random vector in a sphere with uniform density
+/**
+ \ingroup uncommitted
+ */
+IMPDLLEXPORT Vector3D
+random_vector_in_sphere(const Vector3D ¢er=Vector3D(0,0,0),
+ Float radius=1);
+
+//! Generate a random vector on a sphere with uniform density
+/**
+ \ingroup uncommitted
+ */
+IMPDLLEXPORT Vector3D
+random_vector_on_sphere(const Vector3D ¢er=Vector3D(0,0,0),
+ Float radius=1);
+
+
+struct SpacesIO {
+ const Vector3D &v_;
+ SpacesIO(const Vector3D &v): v_(v){}
+};
+
+
+inline std::ostream &operator<<(std::ostream &out, const SpacesIO &s) {
+ out << s.v_[0] << " " << s.v_[1] << " " << s.v_[2];
+ return out;
+}
+
+struct CommasIO {
+ const Vector3D &v_;
+ CommasIO(const Vector3D &v): v_(v){}
+};
+
+
+inline std::ostream &operator<<(std::ostream &out, const CommasIO &s) {
+ out << s.v_[0] << ", " << s.v_[1] << ", " << s.v_[2];
+ return out;
+}
+
} // namespace IMP
#endif /* __IMP_VECTOR_3D_H */
Index: kernel/include/IMP/utility.h
===================================================================
--- kernel/include/IMP/utility.h (revision 681)
+++ kernel/include/IMP/utility.h (working copy)
@@ -8,8 +8,16 @@
#ifndef __IMP_UTILITY_H
#define __IMP_UTILITY_H
+#include "base_types.h"
#include "macros.h"
+#ifdef __GNUC__
+/** \todo Use an sconscript test for the isnan function directly
+ */
+#include <cmath>
+#endif
+
+
namespace IMP
{
@@ -20,6 +28,18 @@
return t*t;
}
+//! put in check for G++
+/** This should not be called isnan as C99 says isnan is a macro. It's
+ all a bit messed up, but it is nice to have something grepable.
+ */
+inline bool is_nan(Float a) {
+#ifdef __GNUC__
+ return std::isnan(a);
+#else
+ return a != a;
+#endif
+}
+
} // namespace IMP
#endif /* __IMP_UTILITY_H */
Index: kernel/doc/examples/chain.py
===================================================================
--- kernel/doc/examples/chain.py (revision 681)
+++ kernel/doc/examples/chain.py (working copy)
@@ -17,8 +17,8 @@
p= IMP.Particle()
pi= m.add_particle(p)
d= IMP.XYZDecorator.create(p)
- d.randomize_in_box(IMP.Vector3D(0,0,0),
- IMP.Vector3D(10,10,10))
+ d.set_coordinates(IMP.random_vector_in_box(IMP.Vector3D(0,0,0),
+ IMP.Vector3D(10,10,10)))
d.set_coordinates_are_optimized(True)
p.add_attribute(rk, radius, False)
chain.append(p)
Index: kernel/src/SConscript
===================================================================
--- kernel/src/SConscript (revision 681)
+++ kernel/src/SConscript (working copy)
@@ -22,7 +22,8 @@
files = ['base_types.cpp', 'Model.cpp',
'Particle.cpp', 'ScoreState.cpp', 'Object.cpp',
'OptimizerState.cpp', 'Log.cpp', 'Restraint.cpp', 'Optimizer.cpp',
- 'random.cpp', 'Key.cpp', 'exception.cpp', 'ParticleRefiner.cpp'
+ 'random.cpp', 'Key.cpp', 'exception.cpp', 'ParticleRefiner.cpp',
+ 'Vector3D.cpp'
] + decorators_files + restraints_files + optimizers_files \
+ unary_functions_files + pair_scores_files + singleton_scores_files \
+ triplet_scores_files + score_states_files + internal_files \
Index: kernel/src/optimizers/MonteCarlo.cpp
===================================================================
--- kernel/src/optimizers/MonteCarlo.cpp (revision 681)
+++ kernel/src/optimizers/MonteCarlo.cpp (working copy)
@@ -39,7 +39,7 @@
Float MonteCarlo::optimize(unsigned int max_steps)
{
IMP_CHECK_OBJECT(this);
- if (cg_.get() != NULL) {
+ if (cg_) {
IMP_CHECK_OBJECT(cg_.get());
IMP_check(cg_->get_model() == get_model(),
"The model used by the local optimizer does not match "\
@@ -58,7 +58,7 @@
(*it)->propose_move(probability_);
}
Float next_energy;
- if (cg_.get() != NULL && num_local_steps_!= 0) {
+ if (cg_ && num_local_steps_!= 0) {
IMP_LOG(VERBOSE,
"MC Performing local optimization "<< std::flush);
IMP_CHECK_OBJECT(cg_.get());
Index: kernel/src/optimizers/movers/BallMover.cpp
===================================================================
--- kernel/src/optimizers/movers/BallMover.cpp (revision 681)
+++ kernel/src/optimizers/movers/BallMover.cpp (working copy)
@@ -16,7 +16,7 @@
static void random_point_in_sphere(unsigned int D,
Float radius,
- std::vector<Float> v)
+ std::vector<Float>& v)
{
IMP_assert(radius > 0, "No volume there");
::boost::uniform_real<> rand(-radius, radius);
@@ -66,7 +66,21 @@
for (unsigned int j = 0; j < get_number_of_float_keys(); ++j) {
propose_value(i, j, npos[j]);
}
+ IMP_LOG(VERBOSE, "Proposing move to ");
+ for (unsigned int j=0; j< npos.size(); ++j) {
+ IMP_LOG(VERBOSE, npos[j] << " ");
+ }
+ IMP_LOG(VERBOSE, " for particle " << get_particle(i)->get_index()
+ << std::endl);
}
}
+
+
+void BallMover::show(std::ostream &out) const
+{
+ out << "Ball Mover with radius " << radius_ << std::endl;
+
+}
+
} // namespace IMP
Index: kernel/src/decorators/XYZDecorator.cpp
===================================================================
--- kernel/src/decorators/XYZDecorator.cpp (revision 681)
+++ kernel/src/decorators/XYZDecorator.cpp (working copy)
@@ -6,9 +6,7 @@
*/
#include "IMP/decorators/XYZDecorator.h"
-#include "IMP/random.h"
-#include <boost/random/uniform_real.hpp>
#include <cmath>
@@ -26,35 +24,6 @@
}
-void XYZDecorator::randomize_in_sphere(const Vector3D ¢er,
- float radius)
-{
- IMP_check(radius > 0, "Radius in randomize must be postive",
- ValueException);
- Vector3D min(center[0]-radius, center[1]-radius, center[2]-radius);
- Vector3D max(center[0]+radius, center[1]+radius, center[2]+radius);
- float norm;
- do {
- randomize_in_box(min, max);
- norm=0;
- for (int i=0; i< 3; ++i) {
- norm+= square(center[i]-get_coordinate(i));
- }
- norm = std::sqrt(norm);
- } while (norm > radius);
-}
-
-void XYZDecorator::randomize_in_box(const Vector3D &min,
- const Vector3D &max)
-{
- for (unsigned int i=0; i< 3; ++i) {
- IMP_check(min[i] < max[i], "Box for randomize must be non-empty",
- ValueException);
- ::boost::uniform_real<> rand(min[i], max[i]);
- set_coordinate(i, rand(random_number_generator));
- }
-}
-
IMP_DECORATOR_INITIALIZE(XYZDecorator, DecoratorBase,
{
key_[0] = FloatKey("x");
@@ -62,16 +31,9 @@
key_[2] = FloatKey("z");
})
-namespace {
- template <class T>
- T d(T a, T b){T d=a-b; return d*d;}
-}
-
Float distance(XYZDecorator a, XYZDecorator b)
{
- double d2= d(a.get_x(), b.get_x()) + d(a.get_y(), b.get_y())
- + d(a.get_z(), b.get_z());
- return std::sqrt(d2);
+ return (a.get_vector()-b.get_vector()).get_magnitude();
}
} // namespace IMP
Index: kernel/src/pair_scores/DistancePairScore.cpp
===================================================================
--- kernel/src/pair_scores/DistancePairScore.cpp (revision 681)
+++ kernel/src/pair_scores/DistancePairScore.cpp (working copy)
@@ -16,10 +16,6 @@
DistancePairScore::DistancePairScore(UnaryFunction *f): f_(f){}
-struct Identity
-{
- Float operator()(Float t) const {return t;}
-};
Float DistancePairScore::evaluate(Particle *a, Particle *b,
DerivativeAccumulator *da) const
Index: kernel/src/pair_scores/TransformedDistancePairScore.cpp
===================================================================
--- kernel/src/pair_scores/TransformedDistancePairScore.cpp (revision 681)
+++ kernel/src/pair_scores/TransformedDistancePairScore.cpp (working copy)
@@ -52,43 +52,59 @@
}
+Vector3D TransformedDistancePairScore
+::get_transformed(const Vector3D &v) const {
+ return Vector3D(get_transformed(v,0),
+ get_transformed(v,1),
+ get_transformed(v,2));
+}
+
+
+Float TransformedDistancePairScore
+::get_transformed(const Vector3D &v, unsigned int i) const {
+ IMP_assert(i<3, "Bad coordinate");
+ return (v-c_)*r_[i] + tc_[i];
+}
+
+
+Vector3D TransformedDistancePairScore
+::get_back_rotated(const Vector3D &v) const {
+ return Vector3D(v*ri_[0],
+ v*ri_[1],
+ v*ri_[2]);
+}
+
+
/*
Compute R(x-c)+tc and the reverse
dt= R(d-c)+tc-> Rt(dt-tc)+c
*/
-struct TransformParticle
+struct TransformedDistancePairScore::TransformParticle
{
- const Vector3D *r_, *ri_;
- const Vector3D &tc_;
- const Vector3D &c_;
+ const TransformedDistancePairScore *trd_;
XYZDecorator d_;
- TransformParticle(const Vector3D *r,
- const Vector3D *ri,
- const Vector3D &tc,
- const Vector3D &c,
- Particle *p): r_(r), ri_(ri),
- tc_(tc), c_(c), d_(p){}
+ TransformParticle(const TransformedDistancePairScore *trd,
+ Particle *p): trd_(trd), d_(p){}
Float get_coordinate(unsigned int i) const {
- return (d_.get_vector()-c_)*r_[i] + tc_[i];
+ return trd_->get_transformed(d_.get_vector(), i);
}
void add_to_coordinates_derivative(const Vector3D& f,
DerivativeAccumulator &da) {
- Vector3D r(f*ri_[0],
- f*ri_[1],
- f*ri_[2]);
- d_.add_to_coordinates_derivative(r, da);
+
+ d_.add_to_coordinates_derivative(trd_->get_back_rotated(f), da);
}
};
Float TransformedDistancePairScore::evaluate(Particle *a, Particle *b,
DerivativeAccumulator *da) const
{
- TransformParticle tb(r_,ri_, tc_,c_,b);
- IMP_LOG(VERBOSE, "Transformed particle is "
+ TransformParticle tb(this,b);
+ /*IMP_LOG(VERBOSE, "Transformed particle is "
<< tb.get_coordinate(0) << " " << tb.get_coordinate(1)
- << " " << tb.get_coordinate(2) << std::endl);
+ << " " << tb.get_coordinate(2)
+ << " from " << XYZDecorator(b) << std::endl);*/
return internal::evaluate_distance_pair_score(XYZDecorator(a),
tb,
da, f_.get(),
@@ -118,24 +134,55 @@
<< ri_[2] << std::endl);
}
-void TransformedDistancePairScore::set_translation(float t0, float t1, float t2)
+void TransformedDistancePairScore::set_translation(const Vector3D &t)
{
- tc_= Vector3D(t0, t1, t2) + c_;
+ tc_= t + c_;
}
-void TransformedDistancePairScore::set_center(float t0, float t1, float t2)
+void TransformedDistancePairScore::set_center(const Vector3D &t)
{
tc_= tc_-c_;
- c_= Vector3D(t0, t1, t2);
+ c_= t;
tc_= tc_+c_;
}
+void TransformedDistancePairScore::set_rotation(const Vector3D &v, Float theta)
+{
+ IMP_check(v.get_magnitude() <1.1 && v.get_magnitude() > .9,
+ "Vector must be normalized", ValueException);
+ IMP_LOG(TERSE, "Computing rotation from " << v << " and "
+ << theta << std::endl);
+ Float c= std::cos(theta);
+ Float s= std::sin(theta);
+ Float C= 1-c;
+ Vector3D vs = v*s;
+ Vector3D vC = v*C;
+
+ Float xyC = v[0]*vC[1];
+ Float yzC = v[1]*vC[2];
+ Float zxC = v[2]*vC[0];
+ set_rotation(v[0]*vC[0] + c, xyC - vs[2], zxC + vs[1],
+ xyC + vs[2], v[1]*vC[1]+c, yzC - vs[0],
+ zxC - vs[1], yzC + vs[0], v[2]*vC[2] + c);
+}
+
+
void TransformedDistancePairScore::show(std::ostream &out) const
{
out << "TransformedDistancePairScore using ";
f_->show(out);
+ out << "With rotation:\n";
+ out << r_[0] << std::endl;
+ out << r_[1] << std::endl;
+ out << r_[2] << std::endl;
+ out << "and inverse rotation:\n";
+ out << ri_[0] << std::endl;
+ out << ri_[1] << std::endl;
+ out << ri_[2] << std::endl;
+ out << "and center: " << c_ << std::endl;
+ out << "and translation: " << tc_-c_ << std::endl;
}
} // namespace IMP
Index: kernel/src/Vector3D.cpp
===================================================================
--- kernel/src/Vector3D.cpp (revision 0)
+++ kernel/src/Vector3D.cpp (revision 0)
@@ -0,0 +1,63 @@
+/**
+ * \file Vector3D.cpp \brief Classes to handle individual a vector.
+ *
+ * Copyright 2007-8 Sali Lab. All rights reserved.
+ *
+ */
+
+#include "IMP/Vector3D.h"
+#include "IMP/random.h"
+#include "IMP/internal/constants.h"
+#include "IMP/utility.h"
+
+#include <boost/random/uniform_real.hpp>
+
+namespace IMP
+{
+
+Vector3D random_vector_in_box(const Vector3D &min, const Vector3D &max)
+{
+ Vector3D ret;
+ for (unsigned int i=0; i< 3; ++i) {
+ IMP_check(min[i] < max[i], "Box for randomize must be non-empty",
+ ValueException);
+ ::boost::uniform_real<> rand(min[i], max[i]);
+ ret[i]=rand(random_number_generator);
+ }
+ return ret;
+}
+
+
+Vector3D random_vector_in_sphere(const Vector3D ¢er, Float radius) {
+ IMP_check(radius > 0, "Radius in randomize must be postive",
+ ValueException);
+ Vector3D min(center[0]-radius, center[1]-radius, center[2]-radius);
+ Vector3D max(center[0]+radius, center[1]+radius, center[2]+radius);
+ float norm;
+ Vector3D ret;
+ do {
+ ret=random_vector_in_box(min, max);
+ norm= (center- ret).get_magnitude();
+ } while (norm > radius);
+ return ret;
+}
+
+Vector3D random_vector_on_sphere(const Vector3D ¢er, Float radius) {
+ IMP_check(radius > 0, "Radius in randomize must be postive",
+ ValueException);
+ ::boost::uniform_real<> rand(-1,1);
+ Vector3D up;
+ up[2]= rand(random_number_generator);
+ ::boost::uniform_real<> trand(0, 2*internal::PI);
+ Float theta= trand(random_number_generator);
+ // radius of circle
+ Float r= std::sqrt(1-square(up[2]));
+ up[0]= std::sin(theta)*r;
+ up[1]= std::cos(theta)*r;
+ IMP_assert(std::abs(up.get_magnitude() -1) < .1,
+ "Error generating unit vector on sphere");
+ IMP_LOG(VERBOSE, "Random vector on sphere is " << up << std::endl);
+ return center+ up*radius;
+}
+
+} // namespace IMP
Index: kernel/pyext/IMP/test.py
===================================================================
--- kernel/pyext/IMP/test.py (revision 681)
+++ kernel/pyext/IMP/test.py (working copy)
@@ -84,7 +84,7 @@
p= IMP.Particle()
model.add_particle(p)
d= IMP.XYZDecorator.create(p)
- d.randomize_in_box(lbv, ubv)
+ d.set_coordinates(IMP.random_vector_in_box(lbv, ubv))
ps.append(p)
d.set_coordinates_are_optimized(True)
return ps
Index: kernel/pyext/IMP.i
===================================================================
--- kernel/pyext/IMP.i (revision 681)
+++ kernel/pyext/IMP.i (working copy)
@@ -96,14 +96,12 @@
return ret;
}
}
- IMP_CONTAINER_SWIG(Model, ScoreState, score_state)
- IMP_CONTAINER_SWIG(Model, Restraint, restraint)
- IMP_CONTAINER_SWIG(RestraintSet, Restraint, restraint)
- IMP_CONTAINER_SWIG(Optimizer, OptimizerState, optimizer_state)
- IMP_ADD_OBJECT(MonteCarlo, add_mover)
- IMP_ADD_OBJECTS(MonteCarlo, add_movers)
- IMP_ADD_OBJECT(NonbondedListScoreState, add_bonded_list)
- IMP_ADD_OBJECTS(NonbondedListScoreState, add_bonded_lists)
+ IMP_CONTAINER_SWIG(Model, ScoreState, score_state);
+ IMP_CONTAINER_SWIG(Model, Restraint, restraint);
+ IMP_CONTAINER_SWIG(RestraintSet, Restraint, restraint);
+ IMP_CONTAINER_SWIG(MonteCarlo, Mover, mover);
+ IMP_CONTAINER_SWIG(Optimizer, OptimizerState, optimizer_state);
+ IMP_CONTAINER_SWIG(NonbondedListScoreState, BondedListScoreState, bonded_list);
}
%feature("ref") Particle "$this->ref();"
@@ -124,6 +122,7 @@
%feature("director") IMP::TripletScore;
%feature("director") IMP::Optimizer;
%feature("director") IMP::ParticleRefiner;
+%feature("director") IMP::Mover;
%include "IMP/base_types.h"
%include "IMP/Object.h"
Index: SConstruct
===================================================================
--- SConstruct (revision 681)
+++ SConstruct (working copy)
@@ -52,7 +52,7 @@
#SConscript('domino/SConscript')
# bin script first requires kernel libraries to be built:
-env.Depends(bin, [src, pyext])
+env.Depends(bin, [pyext, src])
# Build the binaries by default:
env.Default(bin)
1
0
Here is the patch that was discussed a while ago which makes
UnaryFunction::evaluate_with_derivative return a std::pair instead of
using one pass by reference to return. It also removes the hack in IMP.i
to get the return value right.
Index: kernel/include/IMP/UnaryFunction.h
===================================================================
--- kernel/include/IMP/UnaryFunction.h (revision 669)
+++ kernel/include/IMP/UnaryFunction.h (working copy)
@@ -14,9 +14,13 @@
namespace IMP
{
+typedef std::pair<Float, Float> FloatPair;
+
//! Abstract single variable functor class for score functions.
/** These functors take a single feature value, and return a corresponding
score (and optionally also the first derivative).
+
+ \ingroup kernel
*/
class IMPDLLEXPORT UnaryFunction : public RefCountedObject
{
@@ -36,7 +40,7 @@
given feaure.
\return Score
*/
- virtual Float evaluate_with_derivative(Float feature, Float& deriv) const = 0;
+ virtual FloatPair evaluate_with_derivative(Float feature) const = 0;
virtual void show(std::ostream &out=std::cout) const = 0;
};
Index: kernel/include/IMP/optimizers/MonteCarlo.h
===================================================================
--- kernel/include/IMP/optimizers/MonteCarlo.h (revision 669)
+++ kernel/include/IMP/optimizers/MonteCarlo.h (working copy)
@@ -15,8 +15,6 @@
namespace IMP
{
-typedef std::vector<Mover*> Movers;
-
//! A Monte Carlo optimizer.
/** The optimizer uses a set of Mover objects to propose steps. Currently
each Mover is called at each Monte Carlo iteration. This may change in
Index: kernel/include/IMP/optimizers/Mover.h
===================================================================
--- kernel/include/IMP/optimizers/Mover.h (revision 669)
+++ kernel/include/IMP/optimizers/Mover.h (working copy)
@@ -76,6 +76,8 @@
IMP_OUTPUT_OPERATOR(Mover);
+typedef std::vector<Mover*> Movers;
+
} // namespace IMP
#endif /* __IMP_MOVER_H */
Index: kernel/include/IMP/internal/evaluate_distance_pair_score.h
===================================================================
--- kernel/include/IMP/internal/evaluate_distance_pair_score.h (revision 669)
+++ kernel/include/IMP/internal/evaluate_distance_pair_score.h (working copy)
@@ -9,6 +9,7 @@
#define __IMP_EVALUATE_DISTANCE_PAIR_SCORE_H
#include "../Vector3D.h"
+#include <boost/tuple/tuple.hpp>
namespace IMP
{
@@ -43,7 +44,7 @@
if (da && distance >= MIN_DISTANCE) {
Float deriv;
- score = f->evaluate_with_derivative(shifted_distance, deriv);
+ boost::tie(score, deriv) = f->evaluate_with_derivative(shifted_distance);
Vector3D d= delta/distance *deriv;
d0.add_to_coordinates_derivative(d, *da);
Index: kernel/include/IMP/unary_functions/ClosedCubicSpline.h
===================================================================
--- kernel/include/IMP/unary_functions/ClosedCubicSpline.h (revision 669)
+++ kernel/include/IMP/unary_functions/ClosedCubicSpline.h (working copy)
@@ -40,12 +40,10 @@
//! Calculate score and derivative with respect to the given feature.
/** \param[in] feature Value of feature being tested.
- \param[out] deriv Partial derivative of the score with respect to
- the feature value.
\exception ValueException Feature is out of defined range.
\return Score
*/
- virtual Float evaluate_with_derivative(Float feature, Float& deriv) const;
+ virtual FloatPair evaluate_with_derivative(Float feature) const;
void show(std::ostream &out=std::cout) const {
out << "Closed cubic spline of " << values_.size() << " values from "
Index: kernel/include/IMP/unary_functions/Linear.h
===================================================================
--- kernel/include/IMP/unary_functions/Linear.h (revision 669)
+++ kernel/include/IMP/unary_functions/Linear.h (working copy)
@@ -28,9 +28,8 @@
return (feature-offset_)*slope_;
}
- virtual Float evaluate_with_derivative(Float feature, Float& deriv) const {
- deriv= slope_;
- return evaluate(feature);
+ virtual FloatPair evaluate_with_derivative(Float feature) const {
+ return std::make_pair(evaluate(feature), slope_);
}
void set_slope(Float f) {
Index: kernel/include/IMP/unary_functions/WormLikeChain.h
===================================================================
--- kernel/include/IMP/unary_functions/WormLikeChain.h (revision 669)
+++ kernel/include/IMP/unary_functions/WormLikeChain.h (working copy)
@@ -62,10 +62,8 @@
//! Calculate the WormLikeChain energy given the length
/** \param[in] l Current length in angstroms
- \param[out] deriv force in kcal/angstrom mol
- \return Score
*/
- virtual Float evaluate_with_derivative(Float fl, Float& deriv) const {
+ virtual FloatPair evaluate_with_derivative(Float fl) const {
unit::Angstrom l(fl);
if (l < unit::Angstrom(0)) l=unit::Angstrom(0);
unit::Piconewton doubled;
@@ -81,9 +79,9 @@
// convert from picoNewton
unit::YoctoKilocaloriePerAngstrom du= unit::convert_J_to_Cal(doubled);
- deriv = (du*unit::ATOMS_PER_MOL).get_value();
+ Float deriv = (du*unit::ATOMS_PER_MOL).get_value();
//std::cout << "Which converts to " << d << std::endl;
- return evaluate(fl);
+ return std::make_pair(evaluate(fl), deriv);
}
void show(std::ostream &out=std::cout) const {
Index: kernel/include/IMP/unary_functions/Cosine.h
===================================================================
--- kernel/include/IMP/unary_functions/Cosine.h (revision 669)
+++ kernel/include/IMP/unary_functions/Cosine.h (working copy)
@@ -43,11 +43,9 @@
//! Calculate score and derivative with respect to the given feature.
/** \param[in] feature Value of feature being tested.
- \param[out] deriv Partial derivative of the score with respect to
- the feature value.
\return Score
*/
- virtual Float evaluate_with_derivative(Float feature, Float& deriv) const;
+ virtual FloatPair evaluate_with_derivative(Float feature) const;
void show(std::ostream &out=std::cout) const {
out << "Cosine function with force " << force_constant_
Index: kernel/include/IMP/unary_functions/Harmonic.h
===================================================================
--- kernel/include/IMP/unary_functions/Harmonic.h (revision 669)
+++ kernel/include/IMP/unary_functions/Harmonic.h (working copy)
@@ -56,20 +56,17 @@
\return Score
*/
virtual Float evaluate(Float feature) const {
- Float d;
- return evaluate_with_derivative(feature, d);
+ return evaluate_with_derivative(feature).first;
}
//! Calculate harmonic score and derivative with respect to the given feature.
/** \param[in] feature Value of feature being tested.
- \param[out] deriv Partial derivative of the score with respect to
- the feature value.
\return Score
*/
- virtual Float evaluate_with_derivative(Float feature, Float& deriv) const {
+ virtual FloatPair evaluate_with_derivative(Float feature) const {
Float e = (feature - mean_);
- deriv = k_ * e;
- return 0.5 * k_ * e * e;
+ Float deriv = k_ * e;
+ return std::make_pair(0.5 * k_ * e * e, deriv);
}
void show(std::ostream &out=std::cout) const {
Index: kernel/include/IMP/unary_functions/HarmonicLowerBound.h
===================================================================
--- kernel/include/IMP/unary_functions/HarmonicLowerBound.h (revision 669)
+++ kernel/include/IMP/unary_functions/HarmonicLowerBound.h (working copy)
@@ -37,16 +37,13 @@
//! Calculate lower-bound harmonic score and derivative for a feature.
/** If the feature is greater than or equal to the mean, the score is zero.
\param[in] feature Value of feature being tested.
- \param[out] deriv Partial derivative of the score with respect to
- the feature value.
\return Score
*/
- virtual Float evaluate_with_derivative(Float feature, Float& deriv) const {
+ virtual FloatPair evaluate_with_derivative(Float feature) const {
if (feature >= Harmonic::get_mean()) {
- deriv = 0.0;
- return 0.0;
+ return std::make_pair(0.0f, 0.0f);
} else {
- return Harmonic::evaluate_with_derivative(feature, deriv);
+ return Harmonic::evaluate_with_derivative(feature);
}
}
Index: kernel/include/IMP/unary_functions/OpenCubicSpline.h
===================================================================
--- kernel/include/IMP/unary_functions/OpenCubicSpline.h (revision 669)
+++ kernel/include/IMP/unary_functions/OpenCubicSpline.h (working copy)
@@ -25,7 +25,7 @@
\param[in] minrange Feature value at first spline point
\param[in] spacing Distance (in feature space) between points
*/
- OpenCubicSpline(const std::vector<Float> &values, Float minrange,
+ OpenCubicSpline(const Floats &values, Float minrange,
Float spacing);
virtual ~OpenCubicSpline() {}
@@ -39,12 +39,10 @@
//! Calculate score and derivative with respect to the given feature.
/** \param[in] feature Value of feature being tested.
- \param[out] deriv Partial derivative of the score with respect to
- the feature value.
\exception ValueException Feature is out of defined range.
\return Score
*/
- virtual Float evaluate_with_derivative(Float feature, Float& deriv) const;
+ virtual FloatPair evaluate_with_derivative(Float feature) const;
void show(std::ostream &out=std::cout) const {
out << "Open cubic spline of " << values_.size() << " values from "
Index: kernel/include/IMP/unary_functions/HarmonicUpperBound.h
===================================================================
--- kernel/include/IMP/unary_functions/HarmonicUpperBound.h (revision 669)
+++ kernel/include/IMP/unary_functions/HarmonicUpperBound.h (working copy)
@@ -37,16 +37,13 @@
//! Calculate upper-bound harmonic score and derivative for a feature.
/** If the feature is less than or equal to the mean, the score is zero.
\param[in] feature Value of feature being tested.
- \param[out] deriv Partial derivative of the score with respect to
- the feature value.
\return Score
*/
- virtual Float evaluate_with_derivative(Float feature, Float& deriv) const {
+ virtual FloatPair evaluate_with_derivative(Float feature) const {
if (feature <= Harmonic::get_mean()) {
- deriv = 0.0;
- return 0.0;
+ return std::make_pair(0.0f, 0.0f);
} else {
- return Harmonic::evaluate_with_derivative(feature, deriv);
+ return Harmonic::evaluate_with_derivative(feature);
}
}
Index: kernel/src/singleton_scores/AttributeSingletonScore.cpp
===================================================================
--- kernel/src/singleton_scores/AttributeSingletonScore.cpp (revision 669)
+++ kernel/src/singleton_scores/AttributeSingletonScore.cpp (working copy)
@@ -8,7 +8,7 @@
#include "IMP/singleton_scores/AttributeSingletonScore.h"
#include "IMP/UnaryFunction.h"
#include "IMP/Particle.h"
-
+#include <boost/tuple/tuple.hpp>
namespace IMP
{
@@ -20,8 +20,8 @@
DerivativeAccumulator *da) const
{
if (da) {
- Float d;
- float r= f_->evaluate_with_derivative(b->get_value(k_), d);
+ Float d, r;
+ boost::tie(d,r) = f_->evaluate_with_derivative(b->get_value(k_));
b->add_to_derivative(k_, d, *da);
return r;
} else {
Index: kernel/src/singleton_scores/TunnelSingletonScore.cpp
===================================================================
--- kernel/src/singleton_scores/TunnelSingletonScore.cpp (revision 669)
+++ kernel/src/singleton_scores/TunnelSingletonScore.cpp (working copy)
@@ -8,7 +8,7 @@
#include "IMP/singleton_scores/TunnelSingletonScore.h"
#include "IMP/decorators/XYZDecorator.h"
-
+#include <boost/tuple/tuple.hpp>
namespace IMP
{
@@ -58,7 +58,7 @@
// look below if changed
Float dist= -std::min(std::min(rd, hdu), hdd) - radius;
if (accum) {
- score= f_->evaluate_with_derivative(dist, deriv_scalar);
+ boost::tie(score, deriv_scalar)= f_->evaluate_with_derivative(dist);
} else {
score= f_->evaluate(dist);
}
Index: kernel/src/unary_functions/Cosine.cpp
===================================================================
--- kernel/src/unary_functions/Cosine.cpp (revision 669)
+++ kernel/src/unary_functions/Cosine.cpp (working copy)
@@ -18,11 +18,11 @@
- force_constant_ * std::cos(periodicity_ * feature + phase_);
}
-Float Cosine::evaluate_with_derivative(Float feature, Float& deriv) const
+FloatPair Cosine::evaluate_with_derivative(Float feature) const
{
- deriv = force_constant_ * periodicity_
- * std::sin(periodicity_ * feature + phase_);
- return evaluate(feature);
+ Float deriv = force_constant_ * periodicity_
+ * std::sin(periodicity_ * feature + phase_);
+ return std::make_pair(evaluate(feature), deriv);
}
} // namespace IMP
Index: kernel/src/unary_functions/OpenCubicSpline.cpp
===================================================================
--- kernel/src/unary_functions/OpenCubicSpline.cpp (revision 669)
+++ kernel/src/unary_functions/OpenCubicSpline.cpp (working copy)
@@ -66,8 +66,7 @@
* (spacing_ * spacing_) / 6.;
}
-Float OpenCubicSpline::evaluate_with_derivative(Float feature,
- Float& deriv) const
+FloatPair OpenCubicSpline::evaluate_with_derivative(Float feature) const
{
size_t lowbin = static_cast<size_t>((feature - minrange_) / spacing_);
// handle the case where feature ~= maxrange
@@ -79,11 +78,11 @@
Float a = 1. - b;
float sixthspacing = spacing_ / 6.;
- deriv = (values_[highbin] - values_[lowbin]) / spacing_
- - (3. * a * a - 1.) * sixthspacing * second_derivs_[lowbin]
- + (3. * b * b - 1.) * sixthspacing * second_derivs_[highbin];
+ Float deriv = (values_[highbin] - values_[lowbin]) / spacing_
+ - (3. * a * a - 1.) * sixthspacing * second_derivs_[lowbin]
+ + (3. * b * b - 1.) * sixthspacing * second_derivs_[highbin];
- return evaluate(feature);
+ return std::make_pair(evaluate(feature), deriv);
}
} // namespace IMP
Index: kernel/src/unary_functions/ClosedCubicSpline.cpp
===================================================================
--- kernel/src/unary_functions/ClosedCubicSpline.cpp (revision 669)
+++ kernel/src/unary_functions/ClosedCubicSpline.cpp (working copy)
@@ -82,8 +82,8 @@
* (spacing_ * spacing_) / 6.;
}
-Float ClosedCubicSpline::evaluate_with_derivative(Float feature,
- Float& deriv) const
+FloatPair
+ClosedCubicSpline::evaluate_with_derivative(Float feature) const
{
size_t lowbin = static_cast<size_t>((feature - minrange_) / spacing_);
size_t highbin = lowbin + 1;
@@ -99,11 +99,11 @@
Float a = 1. - b;
float sixthspacing = spacing_ / 6.;
- deriv = (values_[highbin] - values_[lowbin]) / spacing_
+ Float deriv = (values_[highbin] - values_[lowbin]) / spacing_
- (3. * a * a - 1.) * sixthspacing * second_derivs_[lowbin]
+ (3. * b * b - 1.) * sixthspacing * second_derivs_[highbin];
- return evaluate(feature);
+ return std::make_pair(evaluate(feature), deriv);
}
} // namespace IMP
Index: kernel/src/restraints/DihedralRestraint.cpp
===================================================================
--- kernel/src/restraints/DihedralRestraint.cpp (revision 669)
+++ kernel/src/restraints/DihedralRestraint.cpp (working copy)
@@ -15,6 +15,8 @@
#include "IMP/restraints/DihedralRestraint.h"
#include "IMP/decorators/XYZDecorator.h"
+#include <boost/tuple/tuple.hpp>
+
namespace IMP
{
@@ -80,7 +82,7 @@
if (accum) {
Float deriv;
- score = score_func_->evaluate_with_derivative(angle, deriv);
+ boost::tie(score, deriv) = score_func_->evaluate_with_derivative(angle);
// method for derivative calculation from van Schaik et al.
// J. Mol. Biol. 234, 751-762 (1993)
Index: kernel/src/triplet_scores/AngleTripletScore.cpp
===================================================================
--- kernel/src/triplet_scores/AngleTripletScore.cpp (revision 669)
+++ kernel/src/triplet_scores/AngleTripletScore.cpp (working copy)
@@ -8,7 +8,7 @@
#include "IMP/triplet_scores/AngleTripletScore.h"
#include "IMP/decorators/XYZDecorator.h"
#include "IMP/UnaryFunction.h"
-
+#include <boost/tuple/tuple.hpp>
namespace IMP
{
@@ -46,7 +46,7 @@
if (da) {
Float deriv;
- score = f_->evaluate_with_derivative(angle, deriv);
+ boost::tie(score, deriv) = f_->evaluate_with_derivative(angle);
Vector3D unit_rij = rij.get_unit_vector();
Vector3D unit_rkj = rkj.get_unit_vector();
Index: kernel/pyext/IMP.i
===================================================================
--- kernel/pyext/IMP.i (revision 669)
+++ kernel/pyext/IMP.i (working copy)
@@ -14,8 +14,17 @@
/* Return derivatives from unary functions */
%include "typemaps.i"
-%apply double &OUTPUT { IMP::Float& deriv };
+namespace IMP {
+ %typemap(out) std::pair<Float,Float> {
+ PyObject *tup= PyTuple_New(2);
+ PyTuple_SetItem(tup, 0, PyFloat_FromDouble($1.first));
+ PyTuple_SetItem(tup, 1, PyFloat_FromDouble($1.second));
+ $result= tup;
+ }
+}
+
+
%pythoncode %{
def check_particle(p, a):
if (not p.get_is_active()):
@@ -27,94 +36,30 @@
%}
namespace IMP {
- %pythonprepend Model::add_restraint %{
- args[1].thisown=0
- %}
- %pythonprepend Model::add_score_state %{
- args[1].thisown=0
- %}
- %pythonprepend Optimizer::add_optimizer_state %{
- args[1].thisown=0
- %}
- %pythonprepend RestraintSet::add_restraint %{
- args[1].thisown=0
- %}
- %pythonprepend NonbondedListScoreState::add_bonded_list %{
- args[1].thisown=0
- %}
- %pythonprepend DistanceRestraint::DistanceRestraint %{
- args[0].thisown=0
- %}
- %pythonprepend AngleRestraint::AngleRestraint %{
- args[0].thisown=0
- %}
- %pythonprepend DihedralRestraint::DihedralRestraint %{
- args[0].thisown=0
- %}
- %pythonprepend TorusRestraint::TorusRestraint %{
- args[3].thisown=0
- %}
- %pythonprepend NonbondedRestraint::NonbondedRestraint %{
- args[0].thisown=0
- %}
- %pythonprepend BondDecoratorRestraint::BondDecoratorRestraint %{
- args[0].thisown=0
- %}
- %pythonprepend SingletonListRestraint::SingletonListRestraint %{
- args[0].thisown=0
- %}
- %pythonprepend PairListRestraint::PairListRestraint %{
- args[0].thisown=0
- %}
- %pythonprepend TripletChainRestraint::TripletChainRestraint %{
- args[0].thisown=0
- %}
- %pythonprepend PairChainRestraint::PairChainRestraint %{
- args[0].thisown=0
- %}
- %pythonprepend ConnectivityRestraint::ConnectivityRestraint %{
- args[0].thisown=0
- %}
- %pythonprepend DistancePairScore::DistancePairScore %{
- args[0].thisown=0
- %}
- %pythonprepend TransformedDistancePairScore::TransformedDistancePairScore %{
- args[0].thisown=0
- %}
- %pythonprepend BondCoverPairScore::BondCoverPairScore %{
- args[0].thisown=0
- %}
- %pythonprepend SphereDistancePairScore::SphereDistancePairScore %{
- args[0].thisown=0
- %}
- %pythonprepend RefineOncePairScore::RefineOncePairScore %{
- args[0].thisown=0
- args[1].thisown=0
- %}
- %pythonprepend DistanceToSingletonScore::DistanceToSingletonScore %{
- args[0].thisown=0
- %}
- %pythonprepend AttributeSingletonScore::AttributeSingletonScore %{
- args[0].thisown=0
- %}
- %pythonprepend TunnelSingletonScore::TunnelSingletonScore %{
- args[0].thisown=0
- %}
- %pythonprepend AngleTripletScore::AngleTripletScore %{
- args[0].thisown=0
- %}
- %pythonprepend MonteCarlo::add_mover %{
- args[1].thisown=0
- %}
- %pythonprepend MonteCarlo::set_local_optimizer %{
- args[1].thisown=0
- %}
- %pythonprepend VRMLLogOptimizerState::add_particle_refiner %{
- args[1].thisown=0
- %}
- %pythonprepend TypedPairScore::set_pair_score %{
- args[1].thisown=0
- %}
+ // need to special case particle so can't add this to macro
+ IMP_OWN_FIRST_CONSTRUCTOR(DistanceRestraint)
+ IMP_OWN_FIRST_CONSTRUCTOR(AngleRestraint)
+ IMP_OWN_FIRST_CONSTRUCTOR(DihedralRestraint)
+ IMP_OWN_FIRST_CONSTRUCTOR(TorusRestraint)
+ IMP_OWN_FIRST_CONSTRUCTOR(NonbondedRestraint)
+ IMP_OWN_FIRST_CONSTRUCTOR(BondDecoratorRestraint)
+ IMP_OWN_FIRST_CONSTRUCTOR(SingletonListRestraint)
+ IMP_OWN_FIRST_CONSTRUCTOR(PairListRestraint)
+ IMP_OWN_FIRST_CONSTRUCTOR(TripletChainRestraint)
+ IMP_OWN_FIRST_CONSTRUCTOR(PairChainRestraint)
+ IMP_OWN_FIRST_CONSTRUCTOR(ConnectivityRestraint)
+ IMP_OWN_FIRST_CONSTRUCTOR(DistancePairScore)
+ IMP_OWN_FIRST_CONSTRUCTOR(TransformedDistancePairScore)
+ IMP_OWN_FIRST_CONSTRUCTOR(BondCoverPairScore)
+ IMP_OWN_FIRST_CONSTRUCTOR(SphereDistancePairScore)
+ IMP_OWN_FIRST_SECOND_CONSTRUCTOR(RefineOncePairScore)
+ IMP_OWN_FIRST_CONSTRUCTOR(DistanceToSingletonScore)
+ IMP_OWN_FIRST_CONSTRUCTOR(AttributeSingletonScore)
+ IMP_OWN_FIRST_CONSTRUCTOR(TunnelSingletonScore)
+ IMP_OWN_FIRST_CONSTRUCTOR(AngleTripletScore)
+ IMP_SET_OBJECT(MonteCarlo, set_local_optimizer)
+ IMP_SET_OBJECT(TypedPairScore, set_pair_score)
+
%pythonprepend Particle::get_value %{
check_particle(args[0], args[1])
%}
@@ -144,10 +89,19 @@
%}
- IMP_CONTAINER_SWIG(Model, Particle, particle);
+ // special case since particles are ref counted
+ %extend Model {
+ Particles get_particles() const {
+ IMP::Particles ret(self->particles_begin(), self->particles_end());
+ return ret;
+ }
+ }
IMP_CONTAINER_SWIG(Model, ScoreState, score_state);
IMP_CONTAINER_SWIG(Model, Restraint, restraint);
IMP_CONTAINER_SWIG(RestraintSet, Restraint, restraint);
+ IMP_CONTAINER_SWIG(MonteCarlo, Mover, mover);
+ IMP_CONTAINER_SWIG(Optimizer, OptimizerState, optimizer_state);
+ IMP_CONTAINER_SWIG(NonbondedListScoreState, BondedListScoreState, bonded_list);
}
%feature("ref") Particle "$this->ref();"
@@ -168,6 +122,7 @@
%feature("director") IMP::TripletScore;
%feature("director") IMP::Optimizer;
%feature("director") IMP::ParticleRefiner;
+%feature("director") IMP::Mover;
%include "IMP/Key.h"
%include "IMP/Object.h"
2
6
This patch moves the functions which create random vectors from
XYZDecorator to Vector3D.h so they can be used outside of the context of
an existing particle. XYZDecorator.randomize_in_box(lb, ub) goes to
XYZDecorator.set_coordinates(IMP.random_vector_in_box(lb,ub)).
Index: kernel/src/decorators/XYZDecorator.cpp
===================================================================
--- kernel/src/decorators/XYZDecorator.cpp (revision 672)
+++ kernel/src/decorators/XYZDecorator.cpp (working copy)
@@ -6,9 +6,7 @@
*/
#include "IMP/decorators/XYZDecorator.h"
-#include "IMP/random.h"
-#include <boost/random/uniform_real.hpp>
#include <cmath>
@@ -26,35 +24,6 @@
}
-void XYZDecorator::randomize_in_sphere(const Vector3D ¢er,
- float radius)
-{
- IMP_check(radius > 0, "Radius in randomize must be postive",
- ValueException);
- Vector3D min(center[0]-radius, center[1]-radius, center[2]-radius);
- Vector3D max(center[0]+radius, center[1]+radius, center[2]+radius);
- float norm;
- do {
- randomize_in_box(min, max);
- norm=0;
- for (int i=0; i< 3; ++i) {
- norm+= square(center[i]-get_coordinate(i));
- }
- norm = std::sqrt(norm);
- } while (norm > radius);
-}
-
-void XYZDecorator::randomize_in_box(const Vector3D &min,
- const Vector3D &max)
-{
- for (unsigned int i=0; i< 3; ++i) {
- IMP_check(min[i] < max[i], "Box for randomize must be non-empty",
- ValueException);
- ::boost::uniform_real<> rand(min[i], max[i]);
- set_coordinate(i, rand(random_number_generator));
- }
-}
-
IMP_DECORATOR_INITIALIZE(XYZDecorator, DecoratorBase,
{
key_[0] = FloatKey("x");
@@ -62,16 +31,9 @@
key_[2] = FloatKey("z");
})
-namespace {
- template <class T>
- T d(T a, T b){T d=a-b; return d*d;}
-}
-
Float distance(XYZDecorator a, XYZDecorator b)
{
- double d2= d(a.get_x(), b.get_x()) + d(a.get_y(), b.get_y())
- + d(a.get_z(), b.get_z());
- return std::sqrt(d2);
+ return (a.get_vector()-b.get_vector()).get_magnitude();
}
} // namespace IMP
Index: kernel/src/Vector3D.cpp
===================================================================
--- kernel/src/Vector3D.cpp (revision 0)
+++ kernel/src/Vector3D.cpp (revision 0)
@@ -0,0 +1,63 @@
+/**
+ * \file Vector3D.cpp \brief Classes to handle individual a vector.
+ *
+ * Copyright 2007-8 Sali Lab. All rights reserved.
+ *
+ */
+
+#include "IMP/Vector3D.h"
+#include "IMP/random.h"
+#include "IMP/internal/constants.h"
+#include "IMP/utility.h"
+
+#include <boost/random/uniform_real.hpp>
+
+namespace IMP
+{
+
+Vector3D random_vector_in_box(const Vector3D &min, const Vector3D &max)
+{
+ Vector3D ret;
+ for (unsigned int i=0; i< 3; ++i) {
+ IMP_check(min[i] < max[i], "Box for randomize must be non-empty",
+ ValueException);
+ ::boost::uniform_real<> rand(min[i], max[i]);
+ ret[i]=rand(random_number_generator);
+ }
+ return ret;
+}
+
+
+Vector3D random_vector_in_sphere(const Vector3D ¢er, Float radius) {
+ IMP_check(radius > 0, "Radius in randomize must be postive",
+ ValueException);
+ Vector3D min(center[0]-radius, center[1]-radius, center[2]-radius);
+ Vector3D max(center[0]+radius, center[1]+radius, center[2]+radius);
+ float norm;
+ Vector3D ret;
+ do {
+ ret=random_vector_in_box(min, max);
+ norm= (center- ret).get_magnitude();
+ } while (norm > radius);
+ return ret;
+}
+
+Vector3D random_vector_on_sphere(const Vector3D ¢er, Float radius) {
+ IMP_check(radius > 0, "Radius in randomize must be postive",
+ ValueException);
+ ::boost::uniform_real<> rand(-1,1);
+ Vector3D up;
+ up[2]= rand(random_number_generator);
+ ::boost::uniform_real<> trand(0, 2*internal::PI);
+ Float theta= trand(random_number_generator);
+ // radius of circle
+ Float r= std::sqrt(1-square(up[2]));
+ up[0]= std::sin(theta)*r;
+ up[1]= std::cos(theta)*r;
+ IMP_assert(std::abs(up.get_magnitude() -1) < .1,
+ "Error generating unit vector on sphere");
+ IMP_LOG(VERBOSE, "Random vector on sphere is " << up << std::endl);
+ return center+ up*radius;
+}
+
+} // namespace IMP
Index: kernel/include/IMP/decorators/XYZDecorator.h
===================================================================
--- kernel/include/IMP/decorators/XYZDecorator.h (revision 672)
+++ kernel/include/IMP/decorators/XYZDecorator.h (working copy)
@@ -8,16 +8,13 @@
#ifndef __IMP_XYZ_DECORATOR_H
#define __IMP_XYZ_DECORATOR_H
-#include <vector>
-#include <deque>
-#include <limits>
-
-#include "../Particle.h"
-#include "../Model.h"
#include "../DecoratorBase.h"
#include "../Vector3D.h"
#include "utility.h"
+#include <vector>
+#include <limits>
+
namespace IMP
{
@@ -67,19 +64,19 @@
return get_particle()->get_derivative(get_coordinate_key(i));
}
//! Add something to the derivative of the ith coordinate
- void add_to_coordinate_derivative(int i, Float v,
+ void add_to_coordinate_derivative(int i, Float v,
DerivativeAccumulator &d) {
get_particle()->add_to_derivative(get_coordinate_key(i), v, d);
}
//! Add something to the derivative of the coordinates
- void add_to_coordinates_derivative(const Vector3D& v,
+ void add_to_coordinates_derivative(const Vector3D& v,
DerivativeAccumulator &d) {
add_to_coordinate_derivative(0, v[0], d);
add_to_coordinate_derivative(1, v[1], d);
add_to_coordinate_derivative(2, v[2], d);
}
//! Get whether the coordinates are optimized
- /** \return true only if all of them are optimized.
+ /** \return true only if all of them are optimized.
*/
bool get_coordinates_are_optimized() const {
return get_particle()->get_is_optimized(get_coordinate_key(0))
@@ -104,23 +101,23 @@
/** 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());
+ return Vector3D(get_x(), get_y(), get_z());
}
+ //! Get the vector of derivatives.
+ /** Somewhat suspect based on wanting a Point/Vector differentiation
+ but we don't have points */
+ Vector3D get_derivative_vector() const {
+ return Vector3D(get_coordinate_derivative(0),
+ get_coordinate_derivative(1),
+ get_coordinate_derivative(2));
+ }
+
//! Get a vector containing the keys for x,y,z
/** This is quite handy for initializing movers and things.
*/
- static const FloatKeys get_xyz_keys() {
- decorator_initialize_static_data();
- return key_;
- }
+ IMP_DECORATOR_GET_KEY(FloatKeys, xyz_keys, key_)
- //! Generate random coordinates in a sphere centered at the vector
- void randomize_in_sphere(const Vector3D ¢er, float radius);
-
- //! Generate random coordinates in a box defined by the vectors
- void randomize_in_box(const Vector3D &lower_corner,
- const Vector3D &upper_corner);
protected:
static FloatKey get_coordinate_key(unsigned int i) {
IMP_check(i <3, "Out of range coordinate",
Index: kernel/include/IMP/Vector3D.h
===================================================================
--- kernel/include/IMP/Vector3D.h (revision 672)
+++ kernel/include/IMP/Vector3D.h (working copy)
@@ -11,6 +11,7 @@
#include "IMP_config.h"
#include "base_types.h"
#include "macros.h"
+#include "exception.h"
#include <cmath>
@@ -22,7 +23,11 @@
*/
class IMPDLLEXPORT Vector3D
{
+ bool is_default() const {return false;}
public:
+ // public for swig
+ typedef Vector3D This;
+
//! Initialize the vector from separate x,y,z values.
Vector3D(Float x, Float y, Float z) {
vec_[0] = x;
@@ -33,6 +38,8 @@
//! Default constructor
Vector3D() {}
+ IMP_COMPARISONS_3(vec_[0], vec_[1], vec_[2]);
+
//! \return A single component of this vector (0-2).
Float operator[](unsigned int i) const {
IMP_assert(i < 3, "Invalid component of vector requested");
@@ -60,23 +67,23 @@
//! product with scalar
Vector3D operator*(Float s) const {
- return Vector3D(operator[](0) * s,
+ return Vector3D(operator[](0) * s,
operator[](1) * s,
- operator[](2) * s);
+ operator[](2) * s);
}
//! divide by a scalar
Vector3D operator/(Float s) const {
- return Vector3D(operator[](0) / s,
+ return Vector3D(operator[](0) / s,
operator[](1) / s,
- operator[](2) / s);
+ operator[](2) / s);
}
//! negation
Vector3D operator-() const {
- return Vector3D(-operator[](0),
+ return Vector3D(-operator[](0),
-operator[](1),
- -operator[](2));
+ -operator[](2));
}
//! \return the vector product of two vectors.
@@ -149,14 +156,6 @@
<< operator[](2) << ")";
}
- bool operator<(const Vector3D &o) const {
- for (unsigned int i=0; i< 3; ++i) {
- if (operator[](i) < o[i]) return true;
- else if (operator[](i) > o[i]) return false;
- }
- return false;
- }
-
private:
Float vec_[3];
};
@@ -166,11 +165,58 @@
//! product with scalar
inline Vector3D operator*(Float s, const Vector3D &o) {
- return Vector3D(o[0]*s,
+ return Vector3D(o[0]*s,
o[1]*s,
- o[2]*s);
+ o[2]*s);
}
+//! Generate a random vector in a box with uniform density
+/**
+ \ingroup uncommitted
+ */
+IMPDLLEXPORT Vector3D
+random_vector_in_box(const Vector3D &lb=Vector3D(0,0,0),
+ const Vector3D &ub=Vector3D(10,10,10));
+
+//! Generate a random vector in a sphere with uniform density
+/**
+ \ingroup uncommitted
+ */
+IMPDLLEXPORT Vector3D
+random_vector_in_sphere(const Vector3D ¢er=Vector3D(0,0,0),
+ Float radius=1);
+
+//! Generate a random vector on a sphere with uniform density
+/**
+ \ingroup uncommitted
+ */
+IMPDLLEXPORT Vector3D
+random_vector_on_sphere(const Vector3D ¢er=Vector3D(0,0,0),
+ Float radius=1);
+
+
+struct SpacesIO {
+ const Vector3D &v_;
+ SpacesIO(const Vector3D &v): v_(v){}
+};
+
+
+inline std::ostream &operator<<(std::ostream &out, const SpacesIO &s) {
+ out << s.v_[0] << " " << s.v_[1] << " " << s.v_[2];
+ return out;
+}
+
+struct CommasIO {
+ const Vector3D &v_;
+ CommasIO(const Vector3D &v): v_(v){}
+};
+
+
+inline std::ostream &operator<<(std::ostream &out, const CommasIO &s) {
+ out << s.v_[0] << ", " << s.v_[1] << ", " << s.v_[2];
+ return out;
+}
+
} // namespace IMP
#endif /* __IMP_VECTOR_3D_H */
Index: kernel/test/states/test_nonbonded_list.py
===================================================================
--- kernel/test/states/test_nonbonded_list.py (revision 672)
+++ kernel/test/states/test_nonbonded_list.py (working copy)
@@ -86,7 +86,8 @@
for i in range(0,10):
score= m.evaluate(False)
for d in ds:
- d.randomize_in_sphere(d.get_vector(), 2.0)
+ d.set_coordinates(IMP.random_vector_in_sphere(d.get_vector(),
+ 2.0))
def do_test_bl(self, ss):
"""Test the bond decorator list"""
@@ -148,8 +149,8 @@
print "Index is " +str(p.get_index().get_index())
d=IMP.XYZDecorator.create(p)
d.set_coordinates_are_optimized(True)
- d.randomize_in_box(IMP.Vector3D(0,0,0),
- IMP.Vector3D(10,10,10));
+ d.set_coordinates(IMP.random_vector_in_box(IMP.Vector3D(0,0,0),
+ IMP.Vector3D(10,10,10)))
nps= IMP.Particles([p])
s.add_particles(nps)
score= m.evaluate(False)
@@ -193,10 +194,12 @@
score= m.evaluate(False)
for p in ps0:
d= IMP.XYZDecorator.cast(p)
- d.randomize_in_sphere(d.get_vector(), 1)
+ d.set_coordinates(IMP.random_vector_in_sphere(d.get_vector(),
+ 1))
for p in ps1:
d= IMP.XYZDecorator.cast(p)
- d.randomize_in_sphere(d.get_vector(), 1)
+ d.set_coordinates(IMP.random_vector_in_sphere(d.get_vector(),
+ 1))
def do_test_spheres(self, ss):
Index: kernel/test/states/test_cover_bonds.py
===================================================================
--- kernel/test/states/test_cover_bonds.py (revision 672)
+++ kernel/test/states/test_cover_bonds.py (working copy)
@@ -15,8 +15,8 @@
p= IMP.Particle()
m.add_particle(p)
d= IMP.XYZDecorator.create(p)
- d.randomize_in_box(IMP.Vector3D(0,0,0),
- IMP.Vector3D(10,10,10))
+ d.set_coordinates(IMP.random_vector_in_box(IMP.Vector3D(0,0,0),
+ IMP.Vector3D(10,10,10)))
ps.append(p)
bds= []
bb= IMP.BondedDecorator.create(ps[0])
Index: kernel/doc/examples/chain.py
===================================================================
--- kernel/doc/examples/chain.py (revision 672)
+++ kernel/doc/examples/chain.py (working copy)
@@ -17,8 +17,8 @@
p= IMP.Particle()
pi= m.add_particle(p)
d= IMP.XYZDecorator.create(p)
- d.randomize_in_box(IMP.Vector3D(0,0,0),
- IMP.Vector3D(10,10,10))
+ d.set_coordinates(IMP.random_vector_in_box(IMP.Vector3D(0,0,0),
+ IMP.Vector3D(10,10,10)))
d.set_coordinates_are_optimized(True)
p.add_attribute(rk, radius, False)
chain.append(p)
Index: kernel/pyext/IMP/test.py
===================================================================
--- kernel/pyext/IMP/test.py (revision 672)
+++ kernel/pyext/IMP/test.py (working copy)
@@ -84,7 +84,7 @@
p= IMP.Particle()
model.add_particle(p)
d= IMP.XYZDecorator.create(p)
- d.randomize_in_box(lbv, ubv)
+ d.set_coordinates(IMP.random_vector_in_box(lbv, ubv))
ps.append(p)
d.set_coordinates_are_optimized(True)
return ps
2
1