Index: kernel/test/particles/test_refcount.py
===================================================================
--- kernel/test/particles/test_refcount.py (revision 649)
+++ kernel/test/particles/test_refcount.py (working copy)
@@ -12,6 +12,8 @@
IMP.set_log_level(IMP.VERBOSE)
self.basenum= IMP.RefCountedObject.get_number_of_live_objects()
print "The base number of objects is " + str(self.basenum)
+ # set this to true if the model is ref counted
+ self.count_model=0
def _check_number(self, expected):
print "Expected "+str(expected)\
@@ -19,7 +21,7 @@
- self.basenum)
self.assertEqual(IMP.RefCountedObject.get_number_of_live_objects() - self.basenum,
expected,
- "wrong number of particles")
+ "wrong number of objects")
def __test_simple(self):
"""Check that ref counting of particles works within python"""
@@ -47,11 +49,12 @@
m= IMP.Model()
p= IMP.Particle()
pi= m.add_particle(p)
- self._check_number(1)
+ self._check_number(1+self.count_model)
m.remove_particle(pi)
- self._check_number(1)
self.assert_(not p.get_is_active(), "Removed particle is still active")
p=1
+ self._check_number(0+ self.count_model)
+ m=1
self._check_number(0)
# This test does not work since swig refcounting is broken
@@ -123,17 +126,17 @@
print "Add particle to mc"
mc.add_particle(p)
# also have the score state now
- self._check_number(2)
+ self._check_number(2+ self.count_model)
print "Remove from model"
m.remove_particle(pi)
- self._check_number(2)
+ self._check_number(2+ self.count_model)
p=1
- self._check_number(2)
+ self._check_number(2+ self.count_model)
print "Remove from mc"
mc.clear_particles()
- self._check_number(1)
+ self._check_number(1+ self.count_model)
mc=0
- self._check_number(0)
+ self._check_number(0+ self.count_model)
def test_skip(self):
"""Check that removed particles are skipped"""
@@ -148,6 +151,7 @@
self.assertEqual(len(ps), 0, "Should no particles particle")
def test_restraints(self):
+ """Check reference counting of restraints"""
m= IMP.Model()
r= IMP.ConstantRestraint(1)
s= IMP.RestraintSet()
@@ -160,7 +164,27 @@
m=0
self._check_number(0)
+ def __test_model(self):
+ """Check reference counting of models"""
+ # disabled due to the usual issues with swig
+ self._check_number(0)
+ m= IMP.Model()
+ o= IMP.ConjugateGradients()
+ self._check_number(1)
+ o.set_model(m)
+ self._check_number(1)
+ o=[]
+ self._check_number(1)
+ m=[]
+ self._check_number(0)
+ m= IMP.Model()
+ o= IMP.ConjugateGradients()
+ self._check_number(1)
+ o.set_model(m)
+ self._check_number(1)
+ m=[]
+ self._check_number(1)
+ o=[]
-
if __name__ == '__main__':
unittest.main()
Index: kernel/test/distance/test_distance.py
===================================================================
--- kernel/test/distance/test_distance.py (revision 649)
+++ kernel/test/distance/test_distance.py (working copy)
@@ -37,6 +37,7 @@
IMP.HarmonicLowerBound(mean, 0.1),
IMP.Harmonic(mean, 0.1)):
r = IMP.DistanceRestraint(sf, p1, p2)
+ r.set_was_owned(True)
self.rsrs.append(r)
def _make_restraints(self):
@@ -56,6 +57,7 @@
IMP.HarmonicLowerBound(3.0, 0.1), IMP.Harmonic(3.0, 0.1)):
r = IMP.DistanceRestraint(fs, self.particles[1],
self.particles[0])
+ r.set_was_owned(True)
self.rsrs.append(r)
# exceed lower bound
@@ -63,6 +65,7 @@
IMP.HarmonicLowerBound(5.0, 0.1), IMP.Harmonic(5.0, 0.1)):
r = IMP.DistanceRestraint(fs, self.particles[1],
self.particles[2])
+ r.set_was_owned(True)
self.rsrs.append(r)
# exceed upper bound
@@ -70,6 +73,7 @@
IMP.HarmonicLowerBound(4.0, 0.1), IMP.Harmonic(4.0, 0.1)):
r = IMP.DistanceRestraint(fs, self.particles[0],
self.particles[2])
+ r.set_was_owned(True)
self.rsrs.append(r)
def test_show(self):
@@ -77,6 +81,7 @@
r = IMP.DistanceRestraint( IMP.Harmonic(0.0, 0.1),
self.particles[1],
self.particles[0])
+ r.set_was_owned(True)
r.show()
def test_distance(self):
Index: kernel/test/restraints/test_angles.py
===================================================================
--- kernel/test/restraints/test_angles.py (revision 649)
+++ kernel/test/restraints/test_angles.py (working copy)
@@ -31,7 +31,7 @@
m.add_particle(p)
d= IMP.XYZDecorator.create(p)
l0.append(p)
- self.randomize_particles(l0, 20)
+ IMP.test.randomize_particles(l0, 20)
return l0
def create_angle_r(self, s, ps):
@@ -53,21 +53,15 @@
IMP.set_log_level(IMP.VERBOSE)
m=IMP.Model()
l0= self.create_particles(m, 3)
- l1= self.create_particles(m, 10)
- l2= self.create_particles(m, 2)
l= One()
t= IMP.AngleTripletScore(l)
r= IMP.TripletChainRestraint(t)
- r.add_chain(l0)
- r.add_chain(l1)
- r.add_chain(l2)
+ r.set_particles(l0)
m.add_restraint(r)
s= IMP.RestraintSet("angle restraints")
m.add_restraint(s)
print "creating angle restraints"
self.create_angle_r(s, l0)
- self.create_angle_r(s, l1)
- self.create_angle_r(s, l2)
ss= s.evaluate(None)
rs= r.evaluate(None)
diff = ss-rs
Index: kernel/test/restraints/test_lowestn.py
===================================================================
--- kernel/test/restraints/test_lowestn.py (revision 0)
+++ kernel/test/restraints/test_lowestn.py (revision 0)
@@ -0,0 +1,63 @@
+import unittest
+import IMP.utils
+import IMP.test, IMP
+import math
+import random
+
+class DistanceTests(IMP.test.TestCase):
+ """Test the symmetry restraint"""
+
+
+ def test_symmetry3(self):
+ """Test lowest n restraint"""
+ IMP.set_log_level(IMP.VERBOSE)
+ m= IMP.Model()
+ ps= IMP.Particles()
+ ds= IMP.XYZDecorators()
+ coords = [
+ [33.2415, 42.2329, 36.6919],
+ [44.6913, 32.324, 39.0969],
+ [58.5871, 33.3032, 36.0686],
+ [69.121, 45.868, 39.4236],
+ [67.4134, 59.7567, 40.5254],
+ [56.1621, 67.1539, 37.9005],
+ [41.6894, 68.0878, 40.6074],
+ [31.9436, 55.7674, 39.9397],
+ ]
+
+ for i in range(0,8):
+ p= IMP.Particle()
+ m.add_particle(p)
+ ps.append(p)
+ d= IMP.XYZDecorator.create(p)
+
+ d.set_coordinates(IMP.Vector3D(coords[i][0], coords[i][1], coords[i][2]))
+ d.set_coordinates_are_optimized(True)
+ ds.append(d)
+
+ ds[0].set_coordinates_are_optimized(False)
+ #ds[4].set_coordinates_are_optimized(False)
+ pps= IMP.ParticlePairs()
+ # want all pairs since score is asymetric
+ for i in range(0,8):
+ for j in range(0,8):
+ pps.append(IMP.ParticlePair(ps[i], ps[j]))
+ tps= IMP.TransformedDistancePairScore(IMP.Harmonic(0,1))
+ cv= IMP.Vector3D(50, 50, 50)
+ tps.set_center(cv)
+ tps.set_rotation(IMP.Vector3D(0,0,1), 2*math.pi/8.0)
+ pr= IMP.LowestNPairListRestraint(tps, 8, pps)
+ m.add_restraint(pr)
+ cg= IMP.SteepestDescent()
+ cg.set_model(m)
+ m.evaluate(False)
+ IMP.set_log_level(IMP.SILENT)
+ cg.optimize(1000)
+ IMP.set_log_level(IMP.VERBOSE)
+ for d in ds:
+ print str(d.get_x()) + " " + str(d.get_y()) + " " + str(d.get_z())
+
+
+
+if __name__ == '__main__':
+ unittest.main()
Index: kernel/test/restraints/test_pairchain.py
===================================================================
--- kernel/test/restraints/test_pairchain.py (revision 649)
+++ kernel/test/restraints/test_pairchain.py (working copy)
@@ -26,20 +26,14 @@
p= IMP.Particle()
m.add_particle(p)
ps0.append(p)
- ps1= IMP.Particles()
- for i in range(0,10):
- p= IMP.Particle()
- m.add_particle(p)
- ps1.append(p)
os= IndexDiff()
s= IMP.PairChainRestraint(os)
- s.add_chain(ps0)
- s.add_chain(ps1)
+ s.set_particles(ps0)
m.add_restraint(s)
score= m.evaluate(False)
print str(score)
- self.assertEqual(score, 18, "Wrong score")
- s.clear_chains()
+ self.assertEqual(score, 9, "Wrong score")
+ s.clear_particles()
self.assertEqual(m.evaluate(False), 0, "Should be no terms")
Index: kernel/test/particle_refiners/test_refine_once_ps.py
===================================================================
--- kernel/test/particle_refiners/test_refine_once_ps.py (revision 649)
+++ kernel/test/particle_refiners/test_refine_once_ps.py (working copy)
@@ -9,7 +9,7 @@
def test_rops(self):
- """Make sure that bond cover coordinates are correct"""
+ """Test refine one pair score"""
IMP.set_log_level(IMP.VERBOSE)
m= IMP.Model()
pp= IMP.Particle()
Index: kernel/test/particle_refiners/test_bond_pr.py
===================================================================
--- kernel/test/particle_refiners/test_bond_pr.py (revision 649)
+++ kernel/test/particle_refiners/test_bond_pr.py (working copy)
@@ -10,8 +10,8 @@
def _set_up_stuff(self, n):
# return [model, particles, bonds]
m= IMP.Model()
- ps= self.create_particles_in_box(m)
- bds= self.create_chain(ps, 10)
+ ps= IMP.test.create_particles_in_box(m)
+ bds= IMP.test.create_chain(ps, 10)
bl= IMP.BondDecoratorListScoreState(ps)
ss= IMP.CoverBondsScoreState(bl, rk)
m.add_score_state(bl)
@@ -56,7 +56,7 @@
print vol*(n+1.0)/n
print n
# crap check
- self.assert_(.5*vol* (n-1.0)/n < 1 and 2*vol *(n+1.0)/n > 1)
+ self.assert_(.5*vol* (n-1.0)/n < 2 and 2*vol *(n+1.0)/n > .5)
Index: kernel/test/modeller/test_pdb_read.py
===================================================================
--- kernel/test/modeller/test_pdb_read.py (revision 649)
+++ kernel/test/modeller/test_pdb_read.py (working copy)
@@ -33,7 +33,7 @@
self.assertEqual(i_num_atom_type, f_num_atom_type, "too many atom types")
self.assertEqual(1377, hc.get_count(),
"Wrong number of particles created")
- rd= IMP.molecular_hierarchy_get_residue(mp, 29)
+ rd= IMP.get_residue(mp, 29)
self.assertEqual(rd.get_index(), 29);
def test_bonds(self):
@@ -41,8 +41,8 @@
m = IMP.Model()
mp= IMP.pdb.read_pdb('modeller/single_protein.pdb', m)
#mp= IMP.MolecularHierarchyDecorator.cast(p)
- all_atoms= IMP.molecular_hierarchy_get_by_type(mp,
- IMP.MolecularHierarchyDecorator.ATOM);
+ all_atoms= IMP.get_by_type(mp,
+ IMP.MolecularHierarchyDecorator.ATOM);
self.assertEqual(1221, all_atoms.size(),
"Wrong number of atoms found in protein")
Index: kernel/test/states/test_nonbonded_list.py
===================================================================
--- kernel/test/states/test_nonbonded_list.py (revision 649)
+++ kernel/test/states/test_nonbonded_list.py (working copy)
@@ -11,16 +11,16 @@
def test_quadratic(self):
"""Test quadratic NBL"""
- ss='IMP.AllNonbondedListScoreState(md, self.rk, IMP.AllNonbondedListScoreState.QUADRATIC)'
+ ss='IMP.AllNonbondedListScoreState(md)'
#eval(ss)
- self.do_test_all(ss)
+ self.do_test_all(ss, IMP.AllNonbondedListScoreState.QUADRATIC)
def test_grid(self):
"""Test grid NBL"""
#IMP.set_log_level(IMP.VERBOSE)
- ss='IMP.AllNonbondedListScoreState(md, self.rk, IMP.AllNonbondedListScoreState.GRID)'
- self.do_test_all(ss)
+ ss='IMP.AllNonbondedListScoreState(md)'
+ self.do_test_all(ss,IMP.AllNonbondedListScoreState.GRID )
def test_bbox(self):
"""Test bbox NBL"""
@@ -29,18 +29,18 @@
except:
print "No CGAL support"
else:
- ss= 'IMP.AllNonbondedListScoreState(md, self.rk, IMP.AllNonbondedListScoreState.BBOX)'
- self.do_test_all(ss)
+ ss= 'IMP.AllNonbondedListScoreState(md)'
+ self.do_test_all(ss, IMP.AllNonbondedListScoreState.BBOX)
def test_default(self):
"""Test default NBL"""
- ss= 'IMP.AllNonbondedListScoreState(md, self.rk)'
- self.do_test_all(ss)
+ ss= 'IMP.AllNonbondedListScoreState(md)'
+ self.do_test_all(ss, IMP.AllNonbondedListScoreState.DEFAULT)
def test_bipartite_quadratic(self):
"""Test bipartite quadratic NBL"""
- ss='IMP.BipartiteNonbondedListScoreState(md, self.rk, IMP.BipartiteNonbondedListScoreState.QUADRATIC)'
- self.do_test_all_bi(ss)
+ ss='IMP.BipartiteNonbondedListScoreState(md)'
+ self.do_test_all_bi(ss, IMP.BipartiteNonbondedListScoreState.QUADRATIC)
def test_bipartite_bbox(self):
"""Test bipartite bbox NBL"""
@@ -49,27 +49,27 @@
except:
print "No CGAL support"
else:
- ss= 'IMP.BipartiteNonbondedListScoreState(md, self.rk, IMP.BipartiteNonbondedListScoreState.BBOX)'
- self.do_test_all_bi(ss)
+ ss= 'IMP.BipartiteNonbondedListScoreState(md)'
+ self.do_test_all_bi(ss, IMP.BipartiteNonbondedListScoreState.BBOX)
- def do_test_all(self, ss):
- self.do_test_bl(ss)
- self.do_test_distfilt(ss)
- self.do_test_spheres(ss)
+ def do_test_all(self, ss, algo):
+ self.do_test_bl(ss, algo)
+ self.do_test_distfilt(ss, algo)
+ self.do_test_spheres(ss, algo)
- def do_test_all_bi(self, ss):
- self.do_test_bi(ss)
- self.do_test_bi_update(ss)
+ def do_test_all_bi(self, ss, algo):
+ self.do_test_bi(ss, algo)
+ self.do_test_bi_update(ss, algo)
def make_spheres(self, m, num, lbv, ubv, minr, maxr):
- ps=self.create_particles_in_box(m, num, lbv, ubv)
+ ps=IMP.test.create_particles_in_box(m, num, lbv, ubv)
for p in ps:
p.add_attribute(self.rk, random.uniform(minr, maxr), False)
return ps
- def do_test_update(self, ss):
+ def do_test_update(self, ss, algo):
m= IMP.Model()
self.make_spheres(m, 20, [0,0,0], [10,10,10], .1, 1)
ds=[]
@@ -77,6 +77,7 @@
ds.append(IMP.XYZDecorator.cast(p))
md=5
s=eval(ss)
+ s.set_algorithm(algo)
s.set_particles(m.get_particles())
m.add_score_state(s)
o= IMP.test.ConstPairScore(1)
@@ -88,7 +89,7 @@
for d in ds:
d.randomize_in_sphere(d.get_vector(), 2.0)
- def do_test_bl(self, ss):
+ def do_test_bl(self, ss, algo):
"""Test the bond decorator list"""
m= IMP.Model()
bds=[]
@@ -99,6 +100,7 @@
IMP.custom_bond(bds[i-1], bds[i], 1, .1)
md= 20
s=eval(ss)
+ s.set_algorithm(algo)
s.set_particles(pts)
b= IMP.BondDecoratorListScoreState(pts)
s.add_bonded_list(b)
@@ -114,17 +116,18 @@
print "Score with bonds is " + str(score)
self.assertEqual(score, 190+1900-19, "Wrong score")
- def do_test_distfilt(self, ss):
+ def do_test_distfilt(self, ss, algo):
"""Test filtering based on distance in nonbonded list"""
#IMP.set_log_level(IMP.TERSE)
m= IMP.Model()
ps=IMP.Particles()
- ps= self.create_particles_in_box(m, 20, [0,0,0], [10,10,10])
- pts= self.create_particles_in_box(m, 20, [160,160,160], [170,170,170])
+ ps= IMP.test.create_particles_in_box(m, 20, [0,0,0], [10,10,10])
+ pts= IMP.test.create_particles_in_box(m, 20, [160,160,160], [170,170,170])
for p in pts:
ps.append(p)
md=15
s=eval(ss)
+ s.set_algorithm(algo)
s.set_particles(ps)
m.add_score_state(s)
o= IMP.test.ConstPairScore(1)
@@ -148,8 +151,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)
@@ -157,13 +160,14 @@
self.assertEqual(score, 2*190, "Wrong score after insert")
- def do_test_bi(self, ss):
+ def do_test_bi(self, ss, algo):
"""Test the bipartite nonbonded list and restraint which uses it"""
m= IMP.Model()
ps0=self.make_spheres(m, 5, [0,0,0], [10,10,10], .1, 1)
ps1=self.make_spheres(m, 5, [0,0,0], [10,10,10], .1, 1)
md=15
s=eval(ss)
+ s.set_algorithm(algo)
#IMP.QuadraticBipartiteNonbondedListScoreState(IMP.FloatKey(),
# ps0, ps1)
s.set_particles(ps0, ps1)
@@ -175,13 +179,14 @@
print score
self.assertEqual(score, 25, "Wrong score")
- def do_test_bi_update(self, ss):
+ def do_test_bi_update(self, ss, algo):
"""Test the bipartite nonbonded list and restraint which uses it"""
m= IMP.Model()
ps0=self.make_spheres(m, 5, [0,0,0], [10,10,10], .1, 1)
ps1=self.make_spheres(m, 5, [0,0,0], [10,10,10], .1, 1)
md=5
s=eval(ss)
+ s.set_algorithm(algo)
#IMP.QuadraticBipartiteNonbondedListScoreState(IMP.FloatKey(),
# ps0, ps1)
s.set_particles(ps0, ps1)
@@ -193,19 +198,20 @@
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):
+ def do_test_spheres(self, ss, algo):
"""Test the nonbonded list of spheres (collision detection)"""
m= IMP.Model()
ps= self.make_spheres(m, 20, IMP.Vector3D(0,0,0),
IMP.Vector3D(10,10,10), 0, 100)
md= 1
s=eval(ss)
+ s.set_algorithm(algo)
s.set_particles(m.get_particles())
#IMP.AllNonbondedListScoreState(rk, m.get_particles())
m.add_score_state(s)
Index: kernel/test/states/test_cover_bonds.py
===================================================================
--- kernel/test/states/test_cover_bonds.py (revision 649)
+++ 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/decorators/test_yaml.py
===================================================================
--- kernel/test/decorators/test_yaml.py (revision 0)
+++ kernel/test/decorators/test_yaml.py (revision 0)
@@ -0,0 +1,47 @@
+import unittest
+import IMP
+import IMP.test
+
+class YamlTests(IMP.test.TestCase):
+ def _create_model(self):
+ IMP.set_log_level(IMP.VERBOSE)
+ m = IMP.Model()
+ p= IMP.Particle()
+ m.add_particle(p)
+ IMP.XYZDecorator.create(p)
+ p= IMP.Particle()
+ m.add_particle(p)
+ IMP.XYZDecorator.create(p)
+ return m
+ def test_yaml(self):
+ """Check writing to yaml """
+ m= self._create_model()
+ IMP.write_yaml(m)
+ def test_read(self):
+ """Check reading from yaml"""
+ m= self._create_model()
+ st="""particle: 0
+ float-attributes:
+ x: 0
+ y: 1
+ z: 2
+ int-attributes:
+ string-attributes:
+ particle-attributes:
+particle: 1
+ float-attributes:
+ x: 3
+ y: 4
+ z: 5
+ int-attributes:
+ string-attributes:
+ particle-attributes:
+"""
+ IMP.read_yaml(st, m)
+ IMP.write_yaml(m)
+
+
+
+
+if __name__ == '__main__':
+ unittest.main()
Index: kernel/test/optimizers/test_sd_optimizer.py
===================================================================
--- kernel/test/optimizers/test_sd_optimizer.py (revision 649)
+++ kernel/test/optimizers/test_sd_optimizer.py (working copy)
@@ -70,7 +70,7 @@
# Start off with all particles in close proximity (but not actually
# colocated, as the derivative of zero distance is zero):
- self.randomize_particles(self.particles, .01)
+ IMP.test.randomize_particles(self.particles, .01)
self.steepest_descent.optimize(50)
Index: kernel/test/pair_scores/test_typed_pair_score.py
===================================================================
--- kernel/test/pair_scores/test_typed_pair_score.py (revision 649)
+++ kernel/test/pair_scores/test_typed_pair_score.py (working copy)
@@ -28,6 +28,11 @@
# The ordering of the particles should not matter:
self.assertEqual(ps.evaluate(pa, pb, da), 5.0)
self.assertEqual(ps.evaluate(pb, pa, da), 5.0)
+ print "ps"
+ ps=[]
+ print "cps"
+ cps= []
+ print "done"
def test_invalid_type(self):
"""Check TypedPairScore behavior with invalid particle types"""
Index: kernel/test/pair_scores/test_transform.py
===================================================================
--- kernel/test/pair_scores/test_transform.py (revision 649)
+++ kernel/test/pair_scores/test_transform.py (working copy)
@@ -1,12 +1,18 @@
import unittest
import IMP.utils
import IMP.test, IMP
+import math
+import random
+def transform(pt, r, c, t):
+ return (pt-c)*r
+
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)
@@ -56,7 +62,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()
@@ -71,26 +77,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[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.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.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)
+ 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/test/pair_scores/test_refine_high.py
===================================================================
--- kernel/test/pair_scores/test_refine_high.py (revision 0)
+++ kernel/test/pair_scores/test_refine_high.py (revision 0)
@@ -0,0 +1,45 @@
+import unittest
+import IMP.utils
+import IMP.test, IMP
+
+class DistanceTests(IMP.test.TestCase):
+ """Test the refine close pair score"""
+ def _build_tree(self, m, levels):
+ p= IMP.Particle()
+ m.add_particle(p)
+ d=IMP.XYZRDecorator.create(p)
+ h= IMP.HierarchyDecorator.create(p)
+ if levels > 0:
+ cs= IMP.Particles()
+ for i in range(0,6):
+ hc= self._build_tree(m, levels-1)
+ h.add_child(hc)
+ cs.append(hc.get_particle())
+ IMP.set_enclosing_sphere(cs, d)
+ else:
+ d.set_coordinates(IMP.random_vector_in_box(IMP.Vector3D(0,0,0),
+ IMP.Vector3D(10,10,10)))
+ d.set_radius(0)
+ return h
+ def test_close(self):
+ """Test the refine close pair score"""
+ m= IMP.Model()
+ h0= self._build_tree(m, 4)
+ h1= self._build_tree(m, 4)
+ pr= IMP.ChildrenParticleRefiner()
+ hlb= IMP.HarmonicLowerBound(0,1)
+ sd= IMP.SphereDistancePairScore(hlb)
+ rps= IMP.RefineHighPairScore(sd, pr, 0)
+ score= rps.evaluate(h0.get_particle(), h1.get_particle(), None)
+ l0= IMP.get_leaves(h0)
+ l1= IMP.get_leaves(h1)
+ ms= 0
+ for la in l0:
+ for lb in l1:
+ s= sd.evaluate(la, lb, None)
+ if s > 0:
+ ms= ms+s
+ self.assertInTolerance(ms, score, .1)
+
+if __name__ == '__main__':
+ unittest.main()
Index: kernel/test/connectivity/test_connectivity.py
===================================================================
--- kernel/test/connectivity/test_connectivity.py (revision 649)
+++ kernel/test/connectivity/test_connectivity.py (working copy)
@@ -26,8 +26,6 @@
IMP.set_log_level(IMP.VERBOSE)
m = IMP.Model()
- rk= IMP.FloatKey("radius")
-
p0=IMP.Particles()
p1=IMP.Particles()
p2=IMP.Particles()
@@ -35,28 +33,28 @@
for i in range(13):
p= IMP.Particle()
m.add_particle(p)
- d= IMP.XYZDecorator.create(p)
+ d= IMP.XYZRDecorator.create(p)
d.set_coordinates_are_optimized(True)
if i % 3 == 0:
p0.append(p)
- p.add_attribute(rk, 1);
+ d.set_radius(1);
elif i % 3 == 1:
p1.append(p)
- p.add_attribute(rk, 2);
+ d.set_radius(2);
else:
p2.append(p)
- p.add_attribute(rk, 3);
+ d.set_radius(3);
p= IMP.Particle()
m.add_particle(p)
- d= IMP.XYZDecorator.create(p)
+ d= IMP.XYZRDecorator.create(p)
d.set_coordinates_are_optimized(True)
p3.append(p)
- p.add_attribute(rk, 1);
+ d.set_radius(1);
o = IMP.ConjugateGradients()
o.set_threshold(1e-4)
o.set_model(m)
- self.randomize_particles(m.get_particles(), 50.0)
+ IMP.test.randomize_particles(m.get_particles(), 50.0)
# add connectivity restraints
Index: kernel/include/IMP/singleton_scores/SConscript
===================================================================
--- kernel/include/IMP/singleton_scores/SConscript (revision 649)
+++ kernel/include/IMP/singleton_scores/SConscript (working copy)
@@ -1,10 +1,10 @@
+Import('env')
import os.path
-Import('env')
-
-files = ['DistanceToSingletonScore.h', 'AttributeSingletonScore.h',
- 'TunnelSingletonScore.h']
-
-# Install the include files:
-includedir = os.path.join(env['includedir'], 'IMP', 'singleton_scores')
+files=[
+ 'AttributeSingletonScore.h',
+ 'DistanceToSingletonScore.h',
+ 'TunnelSingletonScore.h',
+ ]
+includedir = os.path.join(env['includedir'], 'IMP', 'singleton_scores' )
inst = env.Install(includedir, files)
env.Alias('install', inst)
Index: kernel/include/IMP/singleton_scores/TunnelSingletonScore.h
===================================================================
--- kernel/include/IMP/singleton_scores/TunnelSingletonScore.h (revision 649)
+++ kernel/include/IMP/singleton_scores/TunnelSingletonScore.h (working copy)
@@ -14,6 +14,7 @@
#include "../Vector3D.h"
#include "../UnaryFunction.h"
#include "../Pointer.h"
+#include "../decorators/XYZRDecorator.h"
#include "../internal/kernel_version_info.h"
namespace IMP
@@ -42,7 +43,8 @@
Pointer<UnaryFunction> f_;
FloatKey rk_;
public:
- TunnelSingletonScore(UnaryFunction* f, FloatKey r);
+ TunnelSingletonScore(UnaryFunction* f,
+ FloatKey r= XYZRDecorator::get_radius_key());
void set_center(Vector3D c) {
center_=c;
Index: kernel/include/IMP/base_types.h
===================================================================
--- kernel/include/IMP/base_types.h (revision 649)
+++ kernel/include/IMP/base_types.h (working copy)
@@ -14,7 +14,6 @@
#include <string>
#include <vector>
-#include <map>
namespace IMP
{
Index: kernel/include/IMP/Model.h
===================================================================
--- kernel/include/IMP/Model.h (revision 649)
+++ kernel/include/IMP/Model.h (working copy)
@@ -28,6 +28,7 @@
/** All attribute data for particles is stored through indexing in the
model_data_ structure.
Currently no suport for constraints (e.g. rigid bodies).
+
\ingroup kernel
*/
class IMPDLLEXPORT Model: public Object
Index: kernel/include/IMP/UnaryFunction.h
===================================================================
--- kernel/include/IMP/UnaryFunction.h (revision 649)
+++ 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/score_states/ManualBondDecoratorListScoreState.h
===================================================================
--- kernel/include/IMP/score_states/ManualBondDecoratorListScoreState.h (revision 0)
+++ kernel/include/IMP/score_states/ManualBondDecoratorListScoreState.h (revision 0)
@@ -0,0 +1,50 @@
+/**
+ * \file BondDecoratorListScoreState.h
+ * \brief Allow iteration through pairs of a set of atoms.
+ *
+ * Copyright 2007-8 Sali Lab. All rights reserved.
+ */
+
+#ifndef __IMP_MANUAL_BOND_DECORATOR_LIST_SCORE_STATE_H
+#define __IMP_MANUAL_BOND_DECORATOR_LIST_SCORE_STATE_H
+
+#include "BondedListScoreState.h"
+#include "../decorators/bond_decorators.h"
+
+#include <vector>
+
+namespace IMP
+{
+
+class ManualBondDecoratorListScoreState;
+typedef Index<ManualBondDecoratorListScoreState> BondDecoratorListIndex;
+
+//! Keep track of particles that are connected by BondDecorator bonds.
+/** This class just keeps track of the bonds it has been told about.
+ Use the non-manual version if you want auto discovery of bonds.
+ \ingroup bond
+ */
+class IMPDLLEXPORT ManualBondDecoratorListScoreState:
+ public BondedListScoreState
+{
+ public:
+ //! Find bonds amongst the following points.
+ /** \param [in] ps The set of particles to use.
+ */
+ ManualBondDecoratorListScoreState(const BondDecorators &ps= BondDecorators());
+ virtual ~ManualBondDecoratorListScoreState(){}
+
+ virtual void set_particles(const Particles &ps) {}
+ virtual void add_particles(const Particles &ps) {}
+
+ IMP_LIST(public, BondDecorator, bond_decorator, BondDecorator);
+
+ virtual bool are_bonded(Particle *a, Particle *b) const;
+
+protected:
+ virtual void do_before_evaluate() {}
+};
+
+} // namespace IMP
+
+#endif /* __IMP_BOND_DECORATOR_LIST_SCORE_STATE_H */
Index: kernel/include/IMP/score_states/BondDecoratorListScoreState.h
===================================================================
--- kernel/include/IMP/score_states/BondDecoratorListScoreState.h (revision 649)
+++ kernel/include/IMP/score_states/BondDecoratorListScoreState.h (working copy)
@@ -8,25 +8,21 @@
#ifndef __IMP_BOND_DECORATOR_LIST_SCORE_STATE_H
#define __IMP_BOND_DECORATOR_LIST_SCORE_STATE_H
-#include "BondedListScoreState.h"
+#include "ManualBondDecoratorListScoreState.h"
#include "../decorators/bond_decorators.h"
#include <vector>
namespace IMP
{
-
-class BondDecoratorListScoreState;
-typedef Index<BondDecoratorListScoreState> BondDecoratorListIndex;
-
//! Keep track of particles that are connected by BondDecorator bonds.
/** We also may want to add lazy rescanning of bonds rather than doing
it every update call and a faster lookup of bonds.
\ingroup bond
*/
-class IMPDLLEXPORT BondDecoratorListScoreState: public BondedListScoreState
+class IMPDLLEXPORT BondDecoratorListScoreState:
+ public ManualBondDecoratorListScoreState
{
- std::vector<BondDecorator> bonds_;
Particles ps_;
public:
//! Find bonds amongst the following points.
@@ -36,24 +32,8 @@
virtual ~BondDecoratorListScoreState(){}
virtual void set_particles(const Particles &ps);
+ virtual void add_particles(const Particles &ps);
- virtual bool are_bonded(Particle *a, Particle *b) const;
-
- //! This iterates through the pairs of bonded particles
- /** \note update() must be called first for this to be valid.
- */
- typedef std::vector<BondDecorator>::const_iterator BondIterator;
- BondIterator bonds_begin() const {
- return bonds_.begin();
- }
- BondIterator bonds_end() const {
- return bonds_.end();
- }
-
- unsigned int get_number_of_bonds() const {
- return bonds_.size();
- }
-
protected:
virtual void do_before_evaluate();
};
Index: kernel/include/IMP/score_states/NonbondedListScoreState.h
===================================================================
--- kernel/include/IMP/score_states/NonbondedListScoreState.h (revision 649)
+++ kernel/include/IMP/score_states/NonbondedListScoreState.h (working copy)
@@ -9,6 +9,7 @@
#define __IMP_NONBONDED_LIST_SCORE_STATE_H
#include "../ScoreState.h"
+#include "MaxChangeScoreState.h"
#include "../internal/kernel_version_info.h"
#include "../decorators/XYZDecorator.h"
#include "BondedListScoreState.h"
@@ -16,7 +17,6 @@
#include <boost/iterator/filter_iterator.hpp>
#include <vector>
-#include <limits>
namespace IMP
{
@@ -51,6 +51,7 @@
*/
unsigned int max_nbl_size_;
ParticlePairs nbl_;
+ Pointer<MaxChangeScoreState> mcr_;
struct NBLTooLargeException{};
@@ -184,12 +185,22 @@
BoxesOverlap boxes_overlap_object(float cut) const {
return BoxesOverlap(get_radius_object(), cut);
}
+
+ Float get_max_radius_change() const;
+
+ void reset_max_radius_change();
+
+ void update_max_radius_change(unsigned int i);
+
+ void add_particles_for_max_radius_change(const Particles &ps);
+ void clear_particles_for_max_radius_change();
+
public:
- NonbondedListScoreState(Float cutoff, FloatKey rk);
+ NonbondedListScoreState(Float cutoff);
~NonbondedListScoreState();
FloatKey get_radius_key() const {return rk_;}
- void set_radius_key(FloatKey rk) {rk_=rk;}
+ void set_radius_key(FloatKey rk);
//! Set the maximum allowable size for the NBL
/** The NBL will keep reducing the slack and trying to
Index: kernel/include/IMP/score_states/SConscript
===================================================================
--- kernel/include/IMP/score_states/SConscript (revision 649)
+++ kernel/include/IMP/score_states/SConscript (working copy)
@@ -1,12 +1,16 @@
+Import('env')
import os.path
-Import('env')
-
-files = ['BondedListScoreState.h', 'NonbondedListScoreState.h',
- 'BipartiteNonbondedListScoreState.h', 'MaxChangeScoreState.h',
- 'BondDecoratorListScoreState.h', 'AllNonbondedListScoreState.h',
- 'GravityCenterScoreState.h', 'CoverBondsScoreState.h']
-
-# Install the include files:
-includedir = os.path.join(env['includedir'], 'IMP', 'score_states')
+files=[
+ 'AllNonbondedListScoreState.h',
+ 'BipartiteNonbondedListScoreState.h',
+ 'BondDecoratorListScoreState.h',
+ 'BondedListScoreState.h',
+ 'CoverBondsScoreState.h',
+ 'GravityCenterScoreState.h',
+ 'ManualBondDecoratorListScoreState.h',
+ 'MaxChangeScoreState.h',
+ 'NonbondedListScoreState.h',
+ ]
+includedir = os.path.join(env['includedir'], 'IMP', 'score_states' )
inst = env.Install(includedir, files)
env.Alias('install', inst)
Index: kernel/include/IMP/score_states/CoverBondsScoreState.h
===================================================================
--- kernel/include/IMP/score_states/CoverBondsScoreState.h (revision 649)
+++ kernel/include/IMP/score_states/CoverBondsScoreState.h (working copy)
@@ -9,7 +9,7 @@
#define __IMP_COVER_BONDS_SCORE_STATE_H
#include "../ScoreState.h"
-#include "BondDecoratorListScoreState.h"
+#include "ManualBondDecoratorListScoreState.h"
#include "../internal/kernel_version_info.h"
#include "../Pointer.h"
@@ -30,13 +30,13 @@
*/
class IMPDLLEXPORT CoverBondsScoreState: public ScoreState
{
- Pointer<BondDecoratorListScoreState> bl_;
+ Pointer<ManualBondDecoratorListScoreState> bl_;
FloatKey rk_;
public:
/** Get the list of bonds from the BondDecoratorListScoreState. This list is
not owned and update is not called on this list automatically.
*/
- CoverBondsScoreState(BondDecoratorListScoreState *bl,
+ CoverBondsScoreState(ManualBondDecoratorListScoreState *bl,
FloatKey rk=FloatKey("radius"));
~CoverBondsScoreState();
Index: kernel/include/IMP/score_states/BipartiteNonbondedListScoreState.h
===================================================================
--- kernel/include/IMP/score_states/BipartiteNonbondedListScoreState.h (revision 649)
+++ kernel/include/IMP/score_states/BipartiteNonbondedListScoreState.h (working copy)
@@ -20,7 +20,9 @@
/** To iterate through the list of pairs use the NonbondedListScoreState::begin,
NonbondedListScoreState::end functions.
- The radius key can be the default key.
+ By default the XYZRDecorator radius is used as the radius. The radius key
+ can be set to FloatKey to force all particles to be point particles.
+ Particles without radius attribute are treated as points.
Changes in coordinates and radii are handled properly.
@@ -43,7 +45,7 @@
};
protected:
Algorithm a_;
- Pointer<MaxChangeScoreState> mc0_, mc1_, mcr_;
+ Pointer<MaxChangeScoreState> mc0_, mc1_;
void process_sets(const Particles &p0,
const Particles &p1);
@@ -54,22 +56,15 @@
public:
/** \param[in] cutoff The distance cutoff to use.
- \param[in] radius The key to use to get the radius
- \param[in] a Which algorithm to use. The default is the best.
*/
- BipartiteNonbondedListScoreState(Float cutoff, FloatKey radius,
- Algorithm a=DEFAULT);
+ BipartiteNonbondedListScoreState(Float cutoff);
/** \param[in] cutoff The distance cutoff to use.
- \param[in] radius The key to use to get the radius
\param[in] ps0 The first set.
\param[in] ps1 The second set.
- \param[in] a Which algorithm to use. The default is the best.
*/
BipartiteNonbondedListScoreState(Float cutoff,
- FloatKey radius,
const Particles &ps0,
- const Particles &ps1,
- Algorithm a=DEFAULT);
+ const Particles &ps1);
~BipartiteNonbondedListScoreState();
IMP_SCORE_STATE(internal::kernel_version_info);
@@ -78,6 +73,12 @@
void add_particles_0(const Particles &ps);
//! Add the particles to the second set
void add_particles_1(const Particles &ps);
+
+ //! Add the particles to the first set
+ void add_particle_0(Particle* ps);
+ //! Add the particles to the second set
+ void add_particle_1(Particle* ps);
+
//! Remove all particles
void clear_particles();
//! Replace the set of particles
Index: kernel/include/IMP/score_states/AllNonbondedListScoreState.h
===================================================================
--- kernel/include/IMP/score_states/AllNonbondedListScoreState.h (revision 649)
+++ kernel/include/IMP/score_states/AllNonbondedListScoreState.h (working copy)
@@ -26,7 +26,9 @@
/** To iterate through the list of pairs use the NonbondedListScoreState::begin,
NonbondedListScoreState::end functions.
- The radius key can be the default key.
+ By default the XYZRDecorator radius is used as the radius. The radius key
+ can be set to FloatKey to force all particles to be point particles.
+ Particles without radius attribute are treated as points.
Changes in coordinates and radii are handled properly unless grid is used,
then changes in radii may not be handled properly.
@@ -50,7 +52,7 @@
};
protected:
Algorithm a_;
- Pointer<MaxChangeScoreState> mc_, mcr_;
+ Pointer<MaxChangeScoreState> mc_;
void check_nbl() const;
@@ -71,21 +73,11 @@
/** \param[in] cutoff The distance cutoff to use.
- \param[in] radius The key to use to get the radius
- \param[in] a Which algorithm to use. The default is the best.
- */
- AllNonbondedListScoreState(Float cutoff,
- FloatKey radius,
- Algorithm a= DEFAULT);
- /** \param[in] cutoff The distance cutoff to use.
- \param[in] radius The key to use to get the radius
\param[in] ps A list of particles to use.
- \param[in] a Which algorithm to use. The default is the best.
*/
AllNonbondedListScoreState(Float cutoff,
- FloatKey radius,
- const Particles &ps,
- Algorithm a= DEFAULT);
+ const Particles &ps= Particles());
+
~AllNonbondedListScoreState();
IMP_SCORE_STATE(internal::kernel_version_info);
@@ -97,6 +89,11 @@
//! Replace the set of particles
void set_particles(const Particles &ps);
+ //! hack
+ void add_particle(Particle *p) {
+ add_particles(Particles(1,p));
+ }
+
//! If there is CGAL support, a more efficient algorithm BBOX can be used
void set_algorithm(Algorithm a);
};
Index: kernel/include/IMP/score_states/MaxChangeScoreState.h
===================================================================
--- kernel/include/IMP/score_states/MaxChangeScoreState.h (revision 649)
+++ kernel/include/IMP/score_states/MaxChangeScoreState.h (working copy)
@@ -20,13 +20,14 @@
//! Keeps track of the maximum change of a set of attributes.
/** The score state maintains a list of particle and a list of
float attribute keys and keeps track of the maximum amount
- any of these have changed since the last time reset was called.
-
+ any of these have changed since the last time reset was called.
*/
class IMPDLLEXPORT MaxChangeScoreState: public ScoreState
{
FloatKeys keys_;
- FloatKeys origkeys_;
+ typedef internal::AttributeTable<internal::FloatAttributeTableTraits>
+ Table;
+ std::vector<Table> old_values_;
float max_change_;
public:
//! Track the changes with the specified keys.
Index: kernel/include/IMP/restraints/PairChainRestraint.h
===================================================================
--- kernel/include/IMP/restraints/PairChainRestraint.h (revision 649)
+++ kernel/include/IMP/restraints/PairChainRestraint.h (working copy)
@@ -20,31 +20,28 @@
namespace IMP
{
-//! Restrain each pair of consecutive particles in each chain.
+//! Restrain each pair of consecutive particles in chain.
/** \ingroup restraint
*/
class IMPDLLEXPORT PairChainRestraint : public Restraint
{
public:
- //! Create the pair restraint.
+ //! Create the chain restraint.
/** \param[in] pair_score Pair score to apply.
*/
PairChainRestraint(PairScore* pair_score);
virtual ~PairChainRestraint(){}
- IMP_RESTRAINT(internal::kernel_version_info)
+ IMP_RESTRAINT(internal::kernel_version_info);
- //! Add a chain of particles
- /** Each two successive particles are restrained.
- */
- void add_chain(const Particles &ps);
+ using Restraint::add_particles;
+ using Restraint::set_particles;
+ using Restraint::clear_particles;
- //! Clear all the stored chains
- void clear_chains();
+ virtual ParticlesList get_interacting_particles() const;
protected:
Pointer<PairScore> ts_;
- std::vector<unsigned int> chain_splits_;
};
} // namespace IMP
Index: kernel/include/IMP/restraints/SingletonListRestraint.h
===================================================================
--- kernel/include/IMP/restraints/SingletonListRestraint.h (revision 649)
+++ kernel/include/IMP/restraints/SingletonListRestraint.h (working copy)
@@ -40,6 +40,7 @@
using Restraint::clear_particles;
using Restraint::set_particles;
+ virtual ParticlesList get_interacting_particles() const;
protected:
Pointer<SingletonScore> ss_;
};
Index: kernel/include/IMP/restraints/LowestNPairListRestraint.h
===================================================================
--- kernel/include/IMP/restraints/LowestNPairListRestraint.h (revision 0)
+++ kernel/include/IMP/restraints/LowestNPairListRestraint.h (revision 0)
@@ -0,0 +1,46 @@
+/**
+ * \file LowestNPairListRestraint.h
+ * \brief Apply a PairScore to each particle pair in a list.
+ *
+ * Copyright 2007-8 Sali Lab. All rights reserved.
+ *
+ */
+
+#ifndef __IMP_BESTNPAIR_RESTRAINT_H
+#define __IMP_BESTNPAIR_RESTRAINT_H
+
+#include "PairListRestraint.h"
+
+namespace IMP
+{
+
+//! Computes the best N of the pairs and uses the score for them.
+/** The restraint finds the n edges with the lowest score and uses
+ them to evaluate the derivative and score.
+
+ \ingroup restraint
+ */
+class IMPDLLEXPORT LowestNPairListRestraint : public PairListRestraint
+{
+ unsigned int n_;
+public:
+ //! Create the restraint.
+ /** \param[in] ss The function to apply to each particle.
+ \param[in] n The number of pairs to use.
+ \param[in] ps The list of particle pairs to use in the restraints
+ */
+ LowestNPairListRestraint(PairScore *ss, unsigned int n,
+ const ParticlePairs &ps=ParticlePairs());
+ virtual ~LowestNPairListRestraint();
+
+ IMP_RESTRAINT(internal::kernel_version_info);
+
+ //! Set the number of "edges" to use
+ void set_n(unsigned int n) {
+ n_=n;
+ }
+};
+
+} // namespace IMP
+
+#endif /* __IMP_PAIR_LIST_RESTRAINT_H */
Index: kernel/include/IMP/restraints/SConscript
===================================================================
--- kernel/include/IMP/restraints/SConscript (revision 649)
+++ kernel/include/IMP/restraints/SConscript (working copy)
@@ -1,15 +1,21 @@
+Import('env')
import os.path
-Import('env')
-
-files = ['ConnectivityRestraint.h',
- 'DistanceRestraint.h', 'AngleRestraint.h',
- 'DihedralRestraint.h', 'RestraintSet.h',
- 'NonbondedRestraint.h', 'BondDecoratorRestraint.h',
- 'SingletonListRestraint.h', 'PairListRestraint.h',
- 'TripletChainRestraint.h', 'PairChainRestraint.h',
- 'ConstantRestraint.h']
-
-# Install the include files:
-includedir = os.path.join(env['includedir'], 'IMP', 'restraints')
+files=[
+ 'AngleRestraint.h',
+ 'BondDecoratorRestraint.h',
+ 'ConnectivityRestraint.h',
+ 'ConstantRestraint.h',
+ 'DihedralRestraint.h',
+ 'DistanceRestraint.h',
+ 'LowestNPairListRestraint.h',
+ 'NonbondedRestraint.h',
+ 'PairChainRestraint.h',
+ 'PairListRestraint.h',
+ 'PairRestraint.h',
+ 'RestraintSet.h',
+ 'SingletonListRestraint.h',
+ 'TripletChainRestraint.h',
+ ]
+includedir = os.path.join(env['includedir'], 'IMP', 'restraints' )
inst = env.Install(includedir, files)
env.Alias('install', inst)
Index: kernel/include/IMP/restraints/NonbondedRestraint.h
===================================================================
--- kernel/include/IMP/restraints/NonbondedRestraint.h (revision 649)
+++ kernel/include/IMP/restraints/NonbondedRestraint.h (working copy)
@@ -16,12 +16,11 @@
#include "../internal/kernel_version_info.h"
#include "../Pointer.h"
#include "../score_states/NonbondedListScoreState.h"
+#include "../PairScore.h"
namespace IMP
{
-class PairScore;
-
//! Apply a PairScore to all nonbonded pairs of particles
/**
\ingroup restraint
@@ -38,6 +37,8 @@
NonbondedRestraint(PairScore *ps, NonbondedListScoreState *nbl);
virtual ~NonbondedRestraint(){}
+ ParticlesList get_interacting_particles() const;
+
IMP_RESTRAINT(internal::kernel_version_info)
protected:
Index: kernel/include/IMP/restraints/BondDecoratorRestraint.h
===================================================================
--- kernel/include/IMP/restraints/BondDecoratorRestraint.h (revision 649)
+++ kernel/include/IMP/restraints/BondDecoratorRestraint.h (working copy)
@@ -16,12 +16,11 @@
#include "../internal/kernel_version_info.h"
#include "../Pointer.h"
#include "../UnaryFunction.h"
+#include "../score_states/ManualBondDecoratorListScoreState.h"
namespace IMP
{
-class BondDecoratorListScoreState;
-
//! Restrain all pairs of non-bonded particles
/** This restraint currently only works for bonds which have their
length set explicitly. Eventually we should add a table for standard
@@ -40,15 +39,18 @@
\param[in] bl The BondDecoratorListScoreState to use to get the list
of bonds.
*/
- BondDecoratorRestraint(UnaryFunction *f, BondDecoratorListScoreState *bl);
+ BondDecoratorRestraint(UnaryFunction *f,
+ ManualBondDecoratorListScoreState *bl);
virtual ~BondDecoratorRestraint(){}
IMP_RESTRAINT(internal::kernel_version_info)
void set_function(UnaryFunction *f) {f_=f;}
+ ParticlesList get_interacting_particles() const;
+
protected:
- BondDecoratorListScoreState *bl_;
+ Pointer<ManualBondDecoratorListScoreState> bl_;
Pointer<UnaryFunction> f_;
};
Index: kernel/include/IMP/restraints/TripletChainRestraint.h
===================================================================
--- kernel/include/IMP/restraints/TripletChainRestraint.h (revision 649)
+++ kernel/include/IMP/restraints/TripletChainRestraint.h (working copy)
@@ -20,7 +20,7 @@
namespace IMP
{
-//! Restrain each triplet of consecutive particles in each chain.
+//! Restrain each triplet of consecutive particles in the chain.
/** \ingroup restraint
*/
class IMPDLLEXPORT TripletChainRestraint : public Restraint
@@ -28,23 +28,22 @@
public:
//! Create the triplet restraint.
/** \param[in] trip_score Triplet score to apply.
+ \param[in] chain The chain to restrain.
*/
- TripletChainRestraint(TripletScore* trip_score);
+ TripletChainRestraint(TripletScore* trip_score,
+ const Particles &chain=Particles());
virtual ~TripletChainRestraint(){}
+ virtual ParticlesList get_interacting_particles() const;
+ using Restraint::add_particle;
+ using Restraint::set_particles;
+ using Restraint::clear_particles;
+ using Restraint::add_particles;
+
IMP_RESTRAINT(internal::kernel_version_info)
- //! Add a chain of particles
- /** Each three successive particles are restrained.
- */
- void add_chain(const Particles &ps);
-
- //! Clear all the stored chains
- void clear_chains();
-
protected:
Pointer<TripletScore> ts_;
- std::vector<unsigned int> chain_splits_;
};
} // namespace IMP
Index: kernel/include/IMP/restraints/PairListRestraint.h
===================================================================
--- kernel/include/IMP/restraints/PairListRestraint.h (revision 649)
+++ kernel/include/IMP/restraints/PairListRestraint.h (working copy)
@@ -39,7 +39,12 @@
void add_particle_pair(ParticlePair p);
void clear_particle_pairs();
void add_particle_pairs(const ParticlePairs &ps);
+ void set_particle_pairs(const ParticlePairs &ps) {
+ clear_particle_pairs();
+ add_particle_pairs(ps);
+ }
+ virtual ParticlesList get_interacting_particles() const;
protected:
Pointer<PairScore> ss_;
};
Index: kernel/include/IMP/restraints/DistanceRestraint.h
===================================================================
--- kernel/include/IMP/restraints/DistanceRestraint.h (revision 649)
+++ kernel/include/IMP/restraints/DistanceRestraint.h (working copy)
@@ -10,7 +10,7 @@
#include "../IMP_config.h"
#include "../pair_scores/DistancePairScore.h"
-#include "../Restraint.h"
+#include "PairRestraint.h"
#include "../internal/kernel_version_info.h"
#include <iostream>
@@ -23,7 +23,7 @@
\note If the particles are closer than a certain distance, then
the contributions to the derivatives are set to 0.
*/
-class IMPDLLEXPORT DistanceRestraint : public Restraint
+class IMPDLLEXPORT DistanceRestraint : public PairRestraint
{
public:
//! Create the distance restraint.
@@ -35,11 +35,7 @@
Particle* p1, Particle* p2);
virtual ~DistanceRestraint() {}
- IMP_RESTRAINT(internal::kernel_version_info)
-
-protected:
- //! scoring function for this restraint
- DistancePairScore dp_;
+ void show(std::ostream &out=std::cout) const;
};
} // namespace IMP
Index: kernel/include/IMP/restraints/PairRestraint.h
===================================================================
--- kernel/include/IMP/restraints/PairRestraint.h (revision 0)
+++ kernel/include/IMP/restraints/PairRestraint.h (revision 0)
@@ -0,0 +1,47 @@
+/**
+ * \file PairRestraint.h \brief Pair restraint between two particles.
+ *
+ * Copyright 2007-8 Sali Lab. All rights reserved.
+ *
+ */
+
+#ifndef __IMP_PAIR_RESTRAINT_H
+#define __IMP_PAIR_RESTRAINT_H
+
+#include "../IMP_config.h"
+#include "../PairScore.h"
+#include "../Restraint.h"
+#include "../Pointer.h"
+#include "../internal/kernel_version_info.h"
+
+#include <iostream>
+
+namespace IMP
+{
+
+//! Pair restraint between two particles
+/**
+ Apply a pair score to a pair of particles.
+ */
+class IMPDLLEXPORT PairRestraint : public Restraint
+{
+public:
+ //! Create the pair restraint.
+ /** \param[in] score_func Scoring function for the restraint.
+ \param[in] p1 Pointer to first particle in distance restraint.
+ \param[in] p2 Pointer to second particle in distance restraint.
+ */
+ PairRestraint(PairScore* score_func,
+ Particle* p1, Particle* p2);
+ virtual ~PairRestraint() {}
+
+ IMP_RESTRAINT(internal::kernel_version_info)
+
+protected:
+ //! scoring function for this restraint
+ Pointer<PairScore> dp_;
+};
+
+} // namespace IMP
+
+#endif /* __IMP_PAIR_RESTRAINT_H */
Index: kernel/include/IMP/restraints/RestraintSet.h
===================================================================
--- kernel/include/IMP/restraints/RestraintSet.h (revision 649)
+++ kernel/include/IMP/restraints/RestraintSet.h (working copy)
@@ -42,6 +42,8 @@
//! Get weight for all restraints contained by this set.
Float get_weight() const { return weight_; }
+ ParticlesList get_interacting_particles() const;
+
protected:
//! Weight for all restraints.
Index: kernel/include/IMP/Restraint.h
===================================================================
--- kernel/include/IMP/Restraint.h (revision 649)
+++ kernel/include/IMP/Restraint.h (working copy)
@@ -27,6 +27,8 @@
namespace IMP
{
+typedef std::vector<Particles> ParticlesList;
+
class Model;
/** \defgroup restraint General purpose restraints
Classes to define and help in defining restraints. The restraints
@@ -51,6 +53,13 @@
\note Physical restraints should use the units of kcal/mol for restraint
values and kcal/mol/A for derivatives.
+
+ \note Restraints will print a warning message if they are destroyed
+ without ever having been added to a model as this is an easy mistake
+ to make. To disable this warning for a particular restraint, call
+ set_was_owned(true).
+
+ \ingroup kernel
*/
class IMPDLLEXPORT Restraint : public RefCountedObject
{
@@ -94,18 +103,43 @@
Model *get_model() const {
IMP_assert(model_,
"get_model() called before set_model()");
- return model_.get();
+ return model_;
}
+ /** A warning is printed if a restraint is destroyed
+ without ever having belonged to a restraint set or a model.
+ */
+ void set_was_owned(bool tf) {
+ was_owned_=tf;
+ }
+
+ //! Return a list of sets of particles that are restrained by this restraint
+ /** This function returns a list of sets of particles that are
+ interacting within this restraint. Particles can appear in more
+ than one set. However, if two particles never appear in the same
+ set, then changing one of the particles should not change the
+ derivatives of the other particle.
+
+ The default implementation returns all particles stored in this restraint.
+ */
+ virtual ParticlesList get_interacting_particles() const;
+
+
IMP_LIST(protected, Particle, particle, Particle*)
private:
- Pointer<Model> model_;
+ // not ref counted as the model has a pointer to this
+ Model *model_;
/* True if restraint has not been deactivated.
If it is not active, evaluate should not be called
*/
bool is_active_;
+
+ /* keep track of whether the restraint ever was in a model.
+ Give warnings on destruction if it was not.
+ */
+ bool was_owned_;
};
IMP_OUTPUT_OPERATOR(Restraint);
Index: kernel/include/IMP/ParticleRefiner.h
===================================================================
--- kernel/include/IMP/ParticleRefiner.h (revision 649)
+++ kernel/include/IMP/ParticleRefiner.h (working copy)
@@ -11,6 +11,8 @@
#include "base_types.h"
#include "VersionInfo.h"
#include "RefCountedObject.h"
+#include "Pointer.h"
+#include "Particle.h"
namespace IMP
{
@@ -56,6 +58,54 @@
IMP_OUTPUT_OPERATOR(ParticleRefiner);
+//! A class to ensure that cleanup is called with a particle refiner
+/** The guard makes sure that all calls of its refine method are paired
+ with a cleanup.
+
+ It is refcounted so it can be used in a vector. Move semantics will make
+ that unnecessary eventually.
+ */
+class ParticleRefinerGuard
+{
+ Pointer<ParticleRefiner> pr_;
+ Pointer<Particle> p_;
+ Particles ps_;
+ void cleanup() {
+ if (pr_ && !ps_.empty() && get_was_refined()) {
+ pr_->cleanup_refined(p_, ps_);
+ }
+ pr_= NULL;
+ p_= NULL;
+ ps_.clear();
+ }
+public:
+ ~ParticleRefinerGuard() {
+ cleanup();
+ }
+ const Particles& refine(ParticleRefiner *pr, Particle *p) {
+ cleanup();
+ pr_=pr;
+ p_=p;
+ if (pr->get_can_refine(p)) {
+ ps_=pr->get_refined(p);
+ } else {
+ ps_.push_back(p);
+ }
+ return ps_;
+ }
+ const Particles &get_particles() const {
+ return ps_;
+ }
+ Particle* get_particle(unsigned int i) const {
+ return ps_[i];
+ }
+ unsigned int get_number_of_particles() const {
+ return ps_.size();
+ }
+ bool get_was_refined() const {
+ return ps_[0] != p_;
+ }
+};
} // namespace IMP
#endif /* __IMP_PARTICLE_REFINER_H */
Index: kernel/include/IMP/optimizers/states/SConscript
===================================================================
--- kernel/include/IMP/optimizers/states/SConscript (revision 649)
+++ kernel/include/IMP/optimizers/states/SConscript (working copy)
@@ -1,10 +1,10 @@
+Import('env')
import os.path
-Import('env')
-
-files = ['VRMLLogOptimizerState.h', 'CMMLogOptimizerState.h',
- 'VelocityScalingOptimizerState.h']
-
-# Install the include files:
-includedir = os.path.join(env['includedir'], 'IMP', 'optimizers', 'states')
+files=[
+ 'CMMLogOptimizerState.h',
+ 'VelocityScalingOptimizerState.h',
+ 'VRMLLogOptimizerState.h',
+ ]
+includedir = os.path.join(env['includedir'], 'IMP', 'optimizers', 'states' )
inst = env.Install(includedir, files)
env.Alias('install', inst)
Index: kernel/include/IMP/optimizers/states/VRMLLogOptimizerState.h
===================================================================
--- kernel/include/IMP/optimizers/states/VRMLLogOptimizerState.h (revision 649)
+++ kernel/include/IMP/optimizers/states/VRMLLogOptimizerState.h (working copy)
@@ -31,6 +31,9 @@
x,y,z coordinates (as defined by XYZDecorator) and,
optionally have a radius.
+ In addition, the derivatives of particles added to the
+ derivatives list will be drawn.
+
Documentation for the VRML file format can be found at
http://www.cgl.ucsf.edu/chimera/docs/UsersGuide/bild.html
@@ -39,7 +42,9 @@
class IMPDLLEXPORT VRMLLogOptimizerState : public OptimizerState
{
public:
- VRMLLogOptimizerState(std::string filename,
+ typedef Particles Derivatives;
+
+ VRMLLogOptimizerState(std::string filename,
const Particles &pis=Particles());
virtual ~VRMLLogOptimizerState(){}
@@ -72,16 +77,19 @@
void set_color(int c, Vector3D v);
IMP_LIST(public, Particle, particle, Particle*);
- IMP_CONTAINER(ParticleRefiner, particle_refiner, ParticleRefinerIndex);
+ IMP_LIST(public, Derivative, derivative, Particle*);
+ void set_particle_refiner(ParticleRefiner *pr) {
+ pr_=pr;
+ }
+
//! Force it to write the next file
void write_next_file();
void write(std::string name) const;
-protected:
- //! A helper function to just write a list of particles to a file
- void write(std::ostream &out, const Particles &ps) const;
+ protected:
+
std::string filename_;
int file_number_;
int call_number_;
@@ -89,6 +97,7 @@
FloatKey radius_;
IntKey color_;
std::map<int, Vector3D > colors_;
+ Pointer<ParticleRefiner> pr_;
};
} // namespace IMP
Index: kernel/include/IMP/optimizers/ConjugateGradients.h
===================================================================
--- kernel/include/IMP/optimizers/ConjugateGradients.h (revision 649)
+++ kernel/include/IMP/optimizers/ConjugateGradients.h (working copy)
@@ -35,11 +35,11 @@
void set_max_change(Float t) { max_change_ = t; }
private:
- Float get_score(std::vector<FloatIndex> float_indices,
+ Float get_score(const std::vector<FloatIndex> &float_indices,
std::vector<Float> &x, std::vector<Float> &dscore);
bool line_search(std::vector<Float> &x, std::vector<Float> &dx,
- float &alpha, const std::vector<FloatIndex> &float_indices,
- int &ifun, float &f, float &dg, float &dg1,
+ Float &alpha, const std::vector<FloatIndex> &float_indices,
+ int &ifun, Float &f, Float &dg, Float &dg1,
int max_steps, const std::vector<Float> &search,
const std::vector<Float> &estimate);
Float threshold_;
Index: kernel/include/IMP/optimizers/SConscript
===================================================================
--- kernel/include/IMP/optimizers/SConscript (revision 649)
+++ kernel/include/IMP/optimizers/SConscript (working copy)
@@ -1,12 +1,16 @@
+Import('env')
import os.path
-Import('env')
-
-files = ['ConjugateGradients.h', 'SteepestDescent.h', 'MolecularDynamics.h',
- 'MoverBase.h', 'MonteCarlo.h', 'Mover.h', 'BrownianDynamics.h']
-
-# Install the include files:
-includedir = os.path.join(env['includedir'], 'IMP', 'optimizers')
+files=[
+ 'BrownianDynamics.h',
+ 'ConjugateGradients.h',
+ 'MolecularDynamics.h',
+ 'MonteCarlo.h',
+ 'MoverBase.h',
+ 'Mover.h',
+ 'SteepestDescent.h',
+ ]
+includedir = os.path.join(env['includedir'], 'IMP', 'optimizers' )
inst = env.Install(includedir, files)
env.Alias('install', inst)
-SConscript('states/SConscript')
-SConscript('movers/SConscript')
+SConscript('movers/SConscript' )
+SConscript('states/SConscript' )
Index: kernel/include/IMP/optimizers/MoverBase.h
===================================================================
--- kernel/include/IMP/optimizers/MoverBase.h (revision 649)
+++ kernel/include/IMP/optimizers/MoverBase.h (working copy)
@@ -75,6 +75,11 @@
void propose_value(unsigned int i, unsigned int j, Float t) {
if (get_particle(i)->get_is_optimized(get_float_key(j))) {
get_particle(i)->set_value(get_float_key(j), t);
+ } else {
+ IMP_LOG(SILENT, "Trying to change unoptimized attribute "
+ << get_float_key(j) << " in particle "
+ << get_particle(i)->get_value(get_float_key(j))
+ << std::endl);
}
}
//! Propose a value
Index: kernel/include/IMP/optimizers/movers/SConscript
===================================================================
--- kernel/include/IMP/optimizers/movers/SConscript (revision 649)
+++ kernel/include/IMP/optimizers/movers/SConscript (working copy)
@@ -1,8 +1,9 @@
+Import('env')
import os.path
-Import('env')
-
-files = ['BallMover.h', 'NormalMover.h']
-
+files=[
+ 'BallMover.h',
+ 'NormalMover.h',
+ ]
includedir = os.path.join(env['includedir'], 'IMP', 'optimizers', 'movers' )
inst = env.Install(includedir, files)
env.Alias('install', inst)
Index: kernel/include/IMP/optimizers/movers/BallMover.h
===================================================================
--- kernel/include/IMP/optimizers/movers/BallMover.h (revision 649)
+++ kernel/include/IMP/optimizers/movers/BallMover.h (working copy)
@@ -42,6 +42,9 @@
Float get_radius() const {
return radius_;
}
+
+ virtual void show(std::ostream &out=std::cout) const;
+
protected:
/** \internal */
void generate_move(float a);
Index: kernel/include/IMP/optimizers/BrownianDynamics.h
===================================================================
--- kernel/include/IMP/optimizers/BrownianDynamics.h (revision 649)
+++ kernel/include/IMP/optimizers/BrownianDynamics.h (working copy)
@@ -36,7 +36,7 @@
class IMPDLLEXPORT BrownianDynamics : public Optimizer
{
public:
- //! \param[in] radkey The key for the Stokes radius
+ //! \param[in] dkey The key for diffusion constant
BrownianDynamics(FloatKey dkey= FloatKey("D"));
virtual ~BrownianDynamics();
Index: kernel/include/IMP/decorators/XYZDecorator.h
===================================================================
--- kernel/include/IMP/decorators/XYZDecorator.h (revision 649)
+++ 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/XYZRDecorator.h
===================================================================
--- kernel/include/IMP/decorators/XYZRDecorator.h (revision 0)
+++ kernel/include/IMP/decorators/XYZRDecorator.h (revision 0)
@@ -0,0 +1,63 @@
+/**
+ * \file XYZRDecorator.h \brief Simple xyz decorator.
+ *
+ * Copyright 2007-8 Sali Lab. All rights reserved.
+ *
+ */
+
+#ifndef __IMP_XYZR_DECORATOR_H
+#define __IMP_XYZR_DECORATOR_H
+
+#include "XYZDecorator.h"
+
+#include <vector>
+#include <deque>
+#include <limits>
+
+namespace IMP
+{
+
+//! A a decorator for a particle with x,y,z coordinates.
+/** \ingroup helper
+ */
+class IMPDLLEXPORT XYZRDecorator: public XYZDecorator
+{
+ IMP_DECORATOR(XYZRDecorator, XYZDecorator,
+ {
+ return p->has_attribute(radius_);
+ },
+ { p->add_attribute(radius_,0);
+ });
+
+protected:
+ static FloatKey radius_;
+
+public:
+ IMP_DECORATOR_GET_SET(radius, radius_, Float, Float);
+ IMP_DECORATOR_GET_KEY(FloatKey, radius_key, radius_);
+};
+
+IMP_OUTPUT_OPERATOR(XYZRDecorator);
+
+//! Compute the distance between a pair of particles
+/** \ingroup helper
+ */
+IMPDLLEXPORT Float distance(XYZRDecorator a, XYZRDecorator b);
+
+//! Set the coordinates and radius of the last to enclose the list
+/** \ingroup helper
+ \param[in] v The vector of XYZR objects to enclose
+ \param[out] b The one whose values should be set
+
+ \note this takes a Particles object rather than a vector of
+ something else since you can't easily cast vectors of
+ different things to one another. Well, you can cast vectors
+ of decorators in C++, but you can't in python without adding
+ explicit support.
+ */
+IMPDLLEXPORT void set_enclosing_sphere(const Particles &v,
+ XYZRDecorator b);
+
+} // namespace IMP
+
+#endif /* __IMP_XYZR_DECORATOR_H */
Index: kernel/include/IMP/decorators/HierarchyDecorator.h
===================================================================
--- kernel/include/IMP/decorators/HierarchyDecorator.h (revision 649)
+++ kernel/include/IMP/decorators/HierarchyDecorator.h (working copy)
@@ -19,6 +19,7 @@
#include <limits>
#include <vector>
#include <deque>
+#include <memory>
namespace IMP
{
@@ -135,14 +136,14 @@
\ingroup hierarchy
*/
IMPDLLEXPORT
-void breadth_first_traversal(HierarchyDecorator d, HierarchyVisitor &v);
+void breadth_first_traversal(HierarchyDecorator d, HierarchyVisitor *v);
//! Depth first traversal of the hierarchy
/** See breadth_first_traversal and HierarchyVisitor for more information
\ingroup hierarchy
*/
IMPDLLEXPORT
-void depth_first_traversal(HierarchyDecorator d, HierarchyVisitor &v);
+void depth_first_traversal(HierarchyDecorator d, HierarchyVisitor *v);
//! Apply functor F to each particle, traversing the hierarchy breadth first.
@@ -313,23 +314,35 @@
};
+/** This returns an auto_ptr so that the objects are cleaned up properly.
+ I have to pass a non-const ref (or a pointer) to the algorithms which
+ are getting the gather object. Returning an autoptr ensures that the
+ object is alive until after the function call has completed (if gather
+ is called in the parens for the function).
+
+ \internal
+*/
+template <class F, class Out>
+std::auto_ptr<Gather<F, Out> > gather(F f, Out out) {
+ return std::auto_ptr<Gather<F, Out> >(new Gather<F, Out>(f, out));
+}
+
} // namespace internal
//! Gather all the Particle* in the hierarchy which meet some criteria
/** \ingroup hierarchy
*/
-template <class Out, class F>
-Out hierarchy_gather(HierarchyDecorator h, F f, Out out)
+template <class F>
+Particles hierarchy_get(HierarchyDecorator h, F f)
{
- internal::Gather<F,Out> gather(f,out);
- depth_first_traversal(h, gather);
- return gather.get_out();
+ Particles ret;
+ depth_first_traversal(h, internal::gather(f,std::back_inserter(ret)).get());
+ return ret;
}
namespace internal
{
-
template <class K, class V>
struct MatchAttribute
{
@@ -341,20 +354,19 @@
else return o->get_value(k_) == v_;
}
};
+}
-} // namespace internal
-
//! Gather all the Particle* in the hierarchy which match on an attribute
/** \ingroup hierarchy
*/
-template <class Out, class K, class V>
-Out hierarchy_gather_by_attribute(HierarchyDecorator h, K k, V v, Out out)
+template <class K, class V>
+Particles hierarchy_get_by_attribute(HierarchyDecorator h, K k,
+ V v)
{
- internal::Gather<internal::MatchAttribute<K, V>,Out>
- gather(internal::MatchAttribute<K,V>(k,v),
- out);
- depth_first_traversal(h, gather);
- return gather.get_out();
+ Particles ret;
+ depth_first_traversal(h, gather(internal::MatchAttribute<K,V>(k,v),
+ std::back_inserter(ret)).get());
+ return ret;
}
@@ -380,20 +392,27 @@
}
};
+template <class K0, class V0, class K1, class V1>
+MatchAttributes<K0, V0, K1, V1>
+match_attributes(K0 k0, V0 v0,
+ K1 k1, V1 v1) {
+ return MatchAttributes<K0, V0, K1, V1>(k0, v0, k1, v1);
+}
+
} // namespace internal
//! Gather all the Particle* in the hierarchy which match on two attributes
/** \ingroup hierarchy
*/
-template <class Out, class K0, class V0, class K1, class V1>
-Out hierarchy_gather_by_attributes(HierarchyDecorator h, K0 k0,
- V0 v0, K1 k1, V1 v1, Out out)
+template <class K0, class V0, class K1, class V1>
+Particles hierarchy_get_by_attributes(HierarchyDecorator h, K0 k0,
+ V0 v0, K1 k1, V1 v1)
{
- internal::Gather<internal::MatchAttributes<K0, V0, K1, V1>,Out>
- gather(internal::MatchAttributes<K0,V0, K1, V1>(k0,v0, k1, v1),
- out);
- depth_first_traversal(h, gather);
- return gather.get_out();
+ Particles ret;
+ depth_first_traversal(h,
+ gather(internal::match_attributes(k0,v0, k1, v1),
+ std::back_inserter(ret)).get());
+ return ret;
}
@@ -427,15 +446,15 @@
//! Get all the leaves of the bit of hierarchy
IMPDLLEXPORT Particles
-hierarchy_get_leaves(HierarchyDecorator mhd);
+get_leaves(HierarchyDecorator mhd);
//! Get the bonds internal to this tree
IMPDLLEXPORT BondDecorators
-hierarchy_get_internal_bonds(HierarchyDecorator mhd);
+get_internal_bonds(HierarchyDecorator mhd);
//! Get all the particles in the subtree
IMPDLLEXPORT Particles
-hierarchy_get_all_descendants(HierarchyDecorator mhd);
+get_all_descendants(HierarchyDecorator mhd);
namespace internal
Index: kernel/include/IMP/decorators/SConscript
===================================================================
--- kernel/include/IMP/decorators/SConscript (revision 649)
+++ kernel/include/IMP/decorators/SConscript (working copy)
@@ -1,11 +1,18 @@
+Import('env')
import os.path
-Import('env')
-
-files = ['HierarchyDecorator.h', 'NameDecorator.h', 'utility.h',
- 'XYZDecorator.h', 'ResidueDecorator.h', 'AtomDecorator.h',
- 'MolecularHierarchyDecorator.h', 'bond_decorators.h', 'macros.h']
-
-# Install the include files:
-includedir = os.path.join(env['includedir'], 'IMP', 'decorators')
+files=[
+ 'AtomDecorator.h',
+ 'bond_decorators.h',
+ 'HierarchyDecorator.h',
+ 'macros.h',
+ 'MolecularHierarchyDecorator.h',
+ 'NameDecorator.h',
+ 'ResidueDecorator.h',
+ 'utility.h',
+ 'XYZDecorator.h',
+ 'XYZRDecorator.h',
+ 'yaml.h',
+ ]
+includedir = os.path.join(env['includedir'], 'IMP', 'decorators' )
inst = env.Install(includedir, files)
env.Alias('install', inst)
Index: kernel/include/IMP/decorators/yaml.h
===================================================================
--- kernel/include/IMP/decorators/yaml.h (revision 0)
+++ kernel/include/IMP/decorators/yaml.h (revision 0)
@@ -0,0 +1,49 @@
+#ifndef IMP_DECORATORS_YAML_H
+#define IMP_DECORATORS_YAML_H
+#include "../Particle.h"
+#include "../Model.h"
+
+namespace IMP
+{
+
+ //! Write the particle to a stream as a YAML stream
+ /** \param[in] p The particle to write
+ \param[in] out The stream to write to
+ \param[in] indent The base level of indentation
+ */
+ IMPDLLEXPORT void write_yaml(Particle *p,
+ std::ostream &out= std::cout,
+ std::string indent="");
+
+//! Write the model to a stream as a YAML stream
+ /** \param[in] m The model to write
+ \param[in] out The stream to write to
+ \param[in] indent The base level of indentation
+ */
+ IMPDLLEXPORT void write_yaml(Model *m,
+ std::ostream &out= std::cout,
+ std::string indent="");
+
+ //! Read the particle from a YAML stream
+ /** The model must already have particles matching all read particles.
+ Currently the particles must already have the same attributes
+ as are being read, but this probably should change due to lists
+ stored in attributes.
+
+ The intentended usage model is that the model (with restraints) is
+ initialized. Optimization is performed and then the model is written
+ out to a file. Some time later, when you want to reload the model,
+ you can reuse the initialization code to set up the restraints,
+ and then read in the values for the attributes.
+
+ \note The base indent is determined from the first line. */
+ IMPDLLEXPORT void read_yaml(std::istream &in,
+ Model *m);
+
+ //! primarily for debuggin*/
+ IMPDLLEXPORT void read_yaml(std::string contents,
+ Model *m);
+
+}
+
+#endif
Index: kernel/include/IMP/decorators/AtomDecorator.h
===================================================================
--- kernel/include/IMP/decorators/AtomDecorator.h (revision 649)
+++ kernel/include/IMP/decorators/AtomDecorator.h (working copy)
@@ -14,9 +14,6 @@
#include "utility.h"
#include "XYZDecorator.h"
-#include <vector>
-#include <deque>
-
namespace IMP
{
Index: kernel/include/IMP/decorators/macros.h
===================================================================
--- kernel/include/IMP/decorators/macros.h (revision 649)
+++ kernel/include/IMP/decorators/macros.h (working copy)
@@ -261,5 +261,15 @@
#define IMP_DECORATOR_ARRAY_INIT(DecoratorType, name) \
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; \
+ }
+
#endif /* __IMP_DECORATOR_MACROS_H */
Index: kernel/include/IMP/decorators/MolecularHierarchyDecorator.h
===================================================================
--- kernel/include/IMP/decorators/MolecularHierarchyDecorator.h (revision 649)
+++ kernel/include/IMP/decorators/MolecularHierarchyDecorator.h (working copy)
@@ -156,8 +156,8 @@
\ingroup hierarchy
*/
IMPDLLEXPORT Particles
-molecular_hierarchy_get_by_type(MolecularHierarchyDecorator mhd,
- MolecularHierarchyDecorator::Type t);
+get_by_type(MolecularHierarchyDecorator mhd,
+ MolecularHierarchyDecorator::Type t);
class ResidueDecorator;
@@ -174,8 +174,8 @@
\ingroup hierarchy
*/
IMPDLLEXPORT ResidueDecorator
-molecular_hierarchy_get_residue(MolecularHierarchyDecorator mhd,
- unsigned int index);
+get_residue(MolecularHierarchyDecorator mhd,
+ unsigned int index);
//! Create a fragment containing the specified nodes
Index: kernel/include/IMP/decorators/bond_decorators.h
===================================================================
--- kernel/include/IMP/decorators/bond_decorators.h (revision 649)
+++ kernel/include/IMP/decorators/bond_decorators.h (working copy)
@@ -73,6 +73,9 @@
Float, -1);
IMP_DECORATOR_GET_SET_OPT(stiffness, internal::bond_stiffness_key_, Float,
Float, -1);
+
+ //! Get the key used to store the ith bonded particle
+ static ParticleKey get_bonded_key(unsigned int i);
};
IMP_OUTPUT_OPERATOR(BondDecorator);
Index: kernel/include/IMP/OptimizerState.h
===================================================================
--- kernel/include/IMP/OptimizerState.h (revision 649)
+++ kernel/include/IMP/OptimizerState.h (working copy)
@@ -60,11 +60,11 @@
Optimizer *get_optimizer() const {
IMP_assert(optimizer_,
"Must call set_optimizer before get_optimizer on state");
- return optimizer_.get();
+ return optimizer_;
}
protected:
- //! Stored optimizer
- Pointer<Optimizer> optimizer_;
+ //! Stored optimizer. Non-ref counted
+ Optimizer *optimizer_;
};
IMP_OUTPUT_OPERATOR(OptimizerState);
Index: kernel/include/IMP/pair_scores/SphereDistancePairScore.h
===================================================================
--- kernel/include/IMP/pair_scores/SphereDistancePairScore.h (revision 649)
+++ kernel/include/IMP/pair_scores/SphereDistancePairScore.h (working copy)
@@ -11,6 +11,7 @@
#include "../PairScore.h"
#include "../Pointer.h"
#include "../UnaryFunction.h"
+#include "../decorators/XYZRDecorator.h"
namespace IMP
{
@@ -24,7 +25,7 @@
FloatKey radius_;
public:
SphereDistancePairScore(UnaryFunction *f,
- FloatKey radius=FloatKey("radius"));
+ FloatKey radius= XYZRDecorator::get_radius_key());
virtual ~SphereDistancePairScore(){}
virtual Float evaluate(Particle *a, Particle *b,
DerivativeAccumulator *da) const;
Index: kernel/include/IMP/pair_scores/RefineHighPairScore.h
===================================================================
--- kernel/include/IMP/pair_scores/RefineHighPairScore.h (revision 0)
+++ kernel/include/IMP/pair_scores/RefineHighPairScore.h (revision 0)
@@ -0,0 +1,54 @@
+/**
+ * \file RefineHighPairScore.h
+ * \brief Refine particles at most high with a ParticleRefiner.
+ *
+ * Copyright 2007-8 Sali Lab. All rights reserved.
+ */
+
+#ifndef __IMP_REFINE_HIGH_PAIR_SCORE_H
+#define __IMP_REFINE_HIGH_PAIR_SCORE_H
+
+#include "../PairScore.h"
+#include "../UnaryFunction.h"
+#include "../ParticleRefiner.h"
+#include "../Pointer.h"
+
+namespace IMP
+{
+
+//! Refine the input particles based on a score.
+/** This pair score can be used to implement adaptive hierarchical collision
+ detection. The score for the pair score is the sum of the scores of all
+ pairs of leaves whose score is higher than a threshold. The pair score
+ assumes that if the score of a pair is below the threshold, the score
+ of all child pairs is also below the threshold.
+
+ \note It might make sense to add a second PairScore to control the
+ refinement separately from the final computed score (it could
+ default to being the same).
+
+ \ingroup pairscore
+ */
+class IMPDLLEXPORT RefineHighPairScore : public PairScore
+{
+ Pointer<PairScore> f_;
+ Pointer<ParticleRefiner> pr_;
+ Float threshold_;
+public:
+ /** \param[in] f The pair score to apply to the generated pairs
+ \param[in] pr The ParticleRefiner to use
+ \param[in] threshold The distance below which to attempt to refine.
+ */
+ RefineHighPairScore(PairScore *f,
+ ParticleRefiner *pr,
+ Float threshold);
+ virtual ~RefineHighPairScore(){}
+ virtual Float evaluate(Particle *a, Particle *b,
+ DerivativeAccumulator *da) const;
+ virtual void show(std::ostream &out=std::cout) const;
+
+};
+
+} // namespace IMP
+
+#endif /* __IMP_REFINE_HIGH_PAIR_SCORE_H */
Index: kernel/include/IMP/pair_scores/SConscript
===================================================================
--- kernel/include/IMP/pair_scores/SConscript (revision 649)
+++ kernel/include/IMP/pair_scores/SConscript (working copy)
@@ -1,11 +1,13 @@
+Import('env')
import os.path
-Import('env')
-
-files = ['DistancePairScore.h', 'SphereDistancePairScore.h',
- 'RefineOncePairScore.h', 'TypedPairScore.h',
- 'TransformedDistancePairScore.h']
-
-# Install the include files:
-includedir = os.path.join(env['includedir'], 'IMP', 'pair_scores')
+files=[
+ 'DistancePairScore.h',
+ 'RefineHighPairScore.h',
+ 'RefineOncePairScore.h',
+ 'SphereDistancePairScore.h',
+ 'TransformedDistancePairScore.h',
+ 'TypedPairScore.h',
+ ]
+includedir = os.path.join(env['includedir'], 'IMP', 'pair_scores' )
inst = env.Install(includedir, files)
env.Alias('install', inst)
Index: kernel/include/IMP/pair_scores/DistancePairScore.h
===================================================================
--- kernel/include/IMP/pair_scores/DistancePairScore.h (revision 649)
+++ kernel/include/IMP/pair_scores/DistancePairScore.h (working copy)
@@ -27,6 +27,9 @@
virtual Float evaluate(Particle *a, Particle *b,
DerivativeAccumulator *da) const;
virtual void show(std::ostream &out=std::cout) const;
+ UnaryFunction *get_unary_function() const {
+ return f_;
+ }
};
} // namespace IMP
Index: kernel/include/IMP/pair_scores/TransformedDistancePairScore.h
===================================================================
--- kernel/include/IMP/pair_scores/TransformedDistancePairScore.h (revision 649)
+++ kernel/include/IMP/pair_scores/TransformedDistancePairScore.h (working copy)
@@ -25,13 +25,35 @@
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 {
+ IMP_assert(i<3, "Bad coordinate");
+ return (v-c_)*r_[i] + tc_[i];
+ }
+ // Transform a derivative back
+ Vector3D get_back_rotated(const Vector3D &v) const {
+ return Vector3D(v*ri_[0],
+ v*ri_[1],
+ v*ri_[2]);
+ }
public:
TransformedDistancePairScore(UnaryFunction *f);
virtual ~TransformedDistancePairScore(){}
@@ -39,11 +61,27 @@
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 {
+ return Vector3D(get_transformed(v,0),
+ get_transformed(v,1),
+ get_transformed(v,2));
+ }
};
} // namespace IMP
Index: kernel/include/IMP/internal/graph_base.h
===================================================================
--- kernel/include/IMP/internal/graph_base.h (revision 649)
+++ kernel/include/IMP/internal/graph_base.h (working copy)
@@ -32,6 +32,9 @@
node_keys_[1]=ParticleKey((P::get_prefix()+" node 1").c_str());
P::initialize();
}
+ ParticleKey get_node_key(unsigned int i) const {
+ return node_keys_[i];
+ }
ParticleKey node_keys_[2];
};
Index: kernel/include/IMP/internal/SConscript
===================================================================
--- kernel/include/IMP/internal/SConscript (revision 649)
+++ kernel/include/IMP/internal/SConscript (working copy)
@@ -1,14 +1,24 @@
+Import('env')
import os.path
-Import('env')
-
-files = ['AttributeTable.h', 'graph_base.h', 'Grid3D.h', 'Vector.h',
- 'ref_counting.h', 'ObjectContainer.h', 'ParticleGrid.h',
- 'kernel_version_info.h', 'constants.h', 'units.h',
- 'utility.h', 'bbox_nbl_helpers.h',
- 'ArrayOnAttributesHelper.h', 'Unit.h', 'ExponentialNumber.h',
- 'evaluate_distance_pair_score.h']
-
-# Install the include files:
-includedir = os.path.join(env['includedir'], 'IMP', 'internal')
+files=[
+ 'ArrayOnAttributesHelper.h',
+ 'AttributeTable.h',
+ 'bbox_nbl_helpers.h',
+ 'constants.h',
+ 'evaluate_distance_pair_score.h',
+ 'ExponentialNumber.h',
+ 'graph_base.h',
+ 'Grid3D.h',
+ 'kernel_version_info.h',
+ 'MinimalSet.h',
+ 'ObjectContainer.h',
+ 'ParticleGrid.h',
+ 'ref_counting.h',
+ 'Unit.h',
+ 'units.h',
+ 'utility.h',
+ 'Vector.h',
+ ]
+includedir = os.path.join(env['includedir'], 'IMP', 'internal' )
inst = env.Install(includedir, files)
env.Alias('install', inst)
Index: kernel/include/IMP/internal/MinimalSet.h
===================================================================
--- kernel/include/IMP/internal/MinimalSet.h (revision 0)
+++ kernel/include/IMP/internal/MinimalSet.h (revision 0)
@@ -0,0 +1,74 @@
+/**
+ * \file MinimalSet.h \brief Various useful constants.
+ *
+ * Copyright 2007-8 Sali Lab. All rights reserved.
+ *
+ */
+
+#ifndef __IMP_MINIMAL_SET_H
+#define __IMP_MINIMAL_SET_H
+
+#include <utility>
+#include <algorithm>
+#include <limits>
+
+namespace IMP
+{
+
+namespace internal
+{
+
+
+/** Store the lowest n items seen so far.
+
+ */
+template <class Score, class Data>
+class MinimalSet {
+ unsigned int n_;
+ typedef std::pair<Score, Data> MP;
+ typedef std::vector<MP> Vector;
+ Vector data_;
+
+ struct CompareFirst {
+ template <class T>
+ bool operator()(const T &a, const T &b) const {
+ return a.first < b.first;
+ }
+ };
+
+public:
+ MinimalSet(unsigned int n): n_(n){}
+
+ bool can_insert(Score s) const {
+ if (data_.size() < n_) return true;
+ else return data_.back().first > s;
+ }
+
+ void insert(Score s, const Data &d) {
+ IMP_assert(can_insert(s), "Invalid insert");
+ std::pair<Score, Data> pair(s, d);
+ typename Vector::iterator it=
+ std::upper_bound(data_.begin(), data_.end(), pair,
+ CompareFirst());
+ data_.insert(it, pair);
+ if (size() > n_) {
+ data_.pop_back();
+ }
+ }
+
+ unsigned int size() const {
+ return data_.size();
+ }
+
+ const MP &operator[](unsigned int i) const {
+ return data_[i];
+ }
+
+};
+
+
+} // namespace internal
+
+} // namespace IMP
+
+#endif /* __IMP_CONSTANTS_H */
Index: kernel/include/IMP/internal/units.h
===================================================================
--- kernel/include/IMP/internal/units.h (revision 649)
+++ kernel/include/IMP/internal/units.h (working copy)
@@ -119,6 +119,12 @@
typedef Multiply<Multiply<Centimeter, Centimeter>::type,
Centimeter>::type CubicCentimeter;
typedef Divide<Gram, CubicCentimeter>::type GramPerCubicCentimeter;
+typedef Shift<Second, -9>::type Nanosecond;
+typedef Shift<Meter, -9>::type Nanometer;
+typedef Multiply<Multiply<Nanometer, Nanometer>::type,
+ Nanometer>::type CubicNanometer;
+typedef Shift<Second, -15>::type Femtosecond;
+typedef Divide<Piconewton, Nanometer>::type PiconewtonPerNanometer;
Index: kernel/include/IMP/internal/ObjectContainer.h
===================================================================
--- kernel/include/IMP/internal/ObjectContainer.h (revision 649)
+++ kernel/include/IMP/internal/ObjectContainer.h (working copy)
@@ -39,7 +39,7 @@
std::vector<int> free_;
struct OK {
bool operator()(const O*a) const {
- return a != NULL;
+ return a;
}
};
Index: kernel/include/IMP/internal/AttributeTable.h
===================================================================
--- kernel/include/IMP/internal/AttributeTable.h (revision 649)
+++ kernel/include/IMP/internal/AttributeTable.h (working copy)
@@ -42,9 +42,9 @@
}
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();
+ return f != get_invalid();
}
}
};
@@ -183,6 +183,11 @@
map_[k.get_index()] = v;
}
+ void set_values(Value v) {
+ for (unsigned int i=0; i< size_; ++i) {
+ map_[i]=v;
+ }
+ }
void insert(Key k, Value v) {
IMP_assert(!contains(k),
@@ -212,7 +217,7 @@
}
- void show(std::ostream &out, const char *prefix="") const;
+ void show(std::ostream &out, const std::string prefix="") const;
std::vector<Key> get_keys() const;
@@ -273,7 +278,7 @@
template <class Traits>
inline void AttributeTable<Traits>::show(std::ostream &out,
- const char *prefix) const
+ const std::string prefix) const
{
for (unsigned int i=0; i< size_; ++i) {
if (Traits::get_is_valid(map_[i])) {
Index: kernel/include/IMP/internal/Grid3D.h
===================================================================
--- kernel/include/IMP/internal/Grid3D.h (revision 649)
+++ kernel/include/IMP/internal/Grid3D.h (working copy)
@@ -139,6 +139,16 @@
}
};
+
+
+
+
+
+
+
+
+
+
//! Represent a real cell in a grid
/** These indexes represent an actual cell in the grid with no mapping.
They can only be constructed by the grid.
@@ -160,6 +170,16 @@
}
};
+
+
+
+
+
+
+
+
+
+
/** \brief A voxel grid in 3D space.
VT can be any class.
\internal
@@ -340,12 +360,12 @@
}
//! Get the data in a particular cell
- VoxelData& get_voxel(Index gi) {
+ typename std::vector<VT>::reference get_voxel(Index gi) {
return data_[index(gi)];
}
//! Get the data in a particular cell
- const VoxelData& get_voxel(Index gi) const {
+ typename std::vector<VT>::const_reference get_voxel(Index gi) const {
return data_[index(gi)];
}
Index: kernel/include/IMP/internal/ref_counting.h
===================================================================
--- kernel/include/IMP/internal/ref_counting.h (revision 649)
+++ kernel/include/IMP/internal/ref_counting.h (working copy)
@@ -27,7 +27,6 @@
template <class O>
static void eval(O* o) {
BOOST_STATIC_ASSERT((!boost::is_base_of<RefCountedObject, O >::value));
- IMP_LOG(VERBOSE, "Not refing particle " << o << std::endl);
}
};
@@ -36,8 +35,6 @@
{
template <class O>
static void eval(O* o) {
- IMP_LOG(VERBOSE, "Refing particle " << o->get_index()
- << o->get_ref_count() << std::endl);
o->assert_is_valid();
o->ref();
}
@@ -49,7 +46,6 @@
template <class O>
static void eval(O* o) {
BOOST_STATIC_ASSERT((!boost::is_base_of<RefCountedObject, O >::value));
- IMP_LOG(VERBOSE, "Not Unrefing object " << o << std::endl);
}
};
@@ -58,8 +54,6 @@
{
template <class O>
static void eval(O *o) {
- IMP_LOG(VERBOSE, "Unrefing particle " << o->get_index()
- << " " << o->get_ref_count() << std::endl);
o->assert_is_valid();
o->unref();
if (!o->get_has_ref()) {
@@ -104,9 +98,6 @@
template <class O>
void disown(O* o)
{
- /*IMP_LOG(VERBOSE, "Disown called with "
- << (boost::is_base_of<RefCountedObject, O >::value)
- << " for " << o << " " << o->get_ref_count() << std::endl);*/
o->unref();
if (!o->get_has_ref()) {
delete o;
@@ -118,10 +109,6 @@
template <class O>
void own(O* o)
{
- /*IMP_LOG(VERBOSE, "Own called with "
- << (boost::is_base_of<RefCountedObject, O >::value)
- << " for " << o
- << " " << o->get_ref_count() << std::endl);*/
if (boost::is_base_of<RefCountedObject, O >::value) {
// no checks
} else {
@@ -133,6 +120,9 @@
}
+
+
+
} // namespace internal
} // namespace IMP
Index: kernel/include/IMP/internal/evaluate_distance_pair_score.h
===================================================================
--- kernel/include/IMP/internal/evaluate_distance_pair_score.h (revision 649)
+++ 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
{
@@ -16,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 {
+ // 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/include/IMP/Key.h
===================================================================
--- kernel/include/IMP/Key.h (revision 649)
+++ kernel/include/IMP/Key.h (working copy)
@@ -70,7 +70,7 @@
typedef std::vector<Name> Name##s
-/** This must occur in exactly one .o in the internal namespace. Should
+/** This must occur in exactly one .o in the application. Should
be used in the IMP namespace.*/
#define IMP_DEFINE_KEY_TYPE(Name, Tag) \
namespace internal { \
@@ -249,8 +249,8 @@
{
IMP_assert(static_cast<unsigned int>(i)
< get_rmap().size(),
- "Corrupted " << " Key " << i
- << " vs " << get_rmap().size());
+ "Corrupted " << " Key Table asking for key " << i
+ << " with a table of size " << get_rmap().size());
return get_rmap()[i];
}
Index: kernel/include/IMP/Vector3D.h
===================================================================
--- kernel/include/IMP/Vector3D.h (revision 649)
+++ kernel/include/IMP/Vector3D.h (working copy)
@@ -171,6 +171,21 @@
o[2]*s);
}
+//! Generate a random vector in a box with uniform density
+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
+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
+IMPDLLEXPORT Vector3D random_vector_on_sphere(const Vector3D ¢er
+ =Vector3D(0,0,0),
+ Float radius=1);
+
} // namespace IMP
#endif /* __IMP_VECTOR_3D_H */
Index: kernel/include/IMP/exception.h
===================================================================
--- kernel/include/IMP/exception.h (revision 649)
+++ kernel/include/IMP/exception.h (working copy)
@@ -20,6 +20,7 @@
namespace IMP
{
+
//! The general base class for IMP exceptions
/** This way we can catch IMP exceptions without getting memory allocation
errors and everything. And it enforces having a description.
@@ -65,6 +66,7 @@
}
};
+
//! A general exception for an error in IMP.
/** \ingroup assert
*/
@@ -147,6 +149,10 @@
*/
IMPDLLEXPORT void check_fail(const char *msg);
+
+//! Control this through functions in Log.h
+extern bool print_exceptions;
+
} // namespace internal
} // namespace IMP
Index: kernel/include/IMP/unary_functions/ClosedCubicSpline.h
===================================================================
--- kernel/include/IMP/unary_functions/ClosedCubicSpline.h (revision 649)
+++ 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 649)
+++ 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 649)
+++ 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 649)
+++ 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 649)
+++ 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 649)
+++ 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/SConscript
===================================================================
--- kernel/include/IMP/unary_functions/SConscript (revision 649)
+++ kernel/include/IMP/unary_functions/SConscript (working copy)
@@ -1,11 +1,15 @@
+Import('env')
import os.path
-Import('env')
-
-files = ['Harmonic.h', 'HarmonicLowerBound.h', 'HarmonicUpperBound.h',
- 'OpenCubicSpline.h', 'ClosedCubicSpline.h', 'Cosine.h', 'Linear.h',
- 'WormLikeChain.h']
-
-# Install the include files:
-includedir = os.path.join(env['includedir'], 'IMP', 'unary_functions')
+files=[
+ 'ClosedCubicSpline.h',
+ 'Cosine.h',
+ 'Harmonic.h',
+ 'HarmonicLowerBound.h',
+ 'HarmonicUpperBound.h',
+ 'Linear.h',
+ 'OpenCubicSpline.h',
+ 'WormLikeChain.h',
+ ]
+includedir = os.path.join(env['includedir'], 'IMP', 'unary_functions' )
inst = env.Install(includedir, files)
env.Alias('install', inst)
Index: kernel/include/IMP/unary_functions/OpenCubicSpline.h
===================================================================
--- kernel/include/IMP/unary_functions/OpenCubicSpline.h (revision 649)
+++ 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 649)
+++ 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/include/IMP/SConscript
===================================================================
--- kernel/include/IMP/SConscript (revision 649)
+++ kernel/include/IMP/SConscript (working copy)
@@ -1,28 +1,44 @@
+Import('env')
import os.path
-Import('env')
-
-files = ['base_types.h', 'random.h', 'Index.h', 'Model.h',
- 'Particle.h', 'ScoreState.h', 'OptimizerState.h', 'IMP_config.h',
- 'log.h', 'DerivativeAccumulator.h',
- 'Key.h', 'utility.h', 'Restraint.h', 'Optimizer.h',
- 'DecoratorBase.h', 'Vector3D.h',
- 'UnaryFunction.h', 'PairScore.h', 'SingletonScore.h', 'macros.h',
- 'TripletScore.h', 'exception.h', 'VersionInfo.h',
- 'Object.h', 'Pointer.h', 'RefCountedObject.h']
-
-# Install the include files:
-includedir = os.path.join(env['includedir'], 'IMP')
+files=[
+ 'base_types.h',
+ 'DecoratorBase.h',
+ 'DerivativeAccumulator.h',
+ 'exception.h',
+ 'IMP_config.h',
+ 'Index.h',
+ 'Key.h',
+ 'log.h',
+ 'macros.h',
+ 'Model.h',
+ 'Object.h',
+ 'Optimizer.h',
+ 'OptimizerState.h',
+ 'PairScore.h',
+ 'Particle.h',
+ 'ParticleRefiner.h',
+ 'Pointer.h',
+ 'random.h',
+ 'RefCountedObject.h',
+ 'Restraint.h',
+ 'ScoreState.h',
+ 'SingletonScore.h',
+ 'TripletScore.h',
+ 'UnaryFunction.h',
+ 'utility.h',
+ 'Vector3D.h',
+ 'VersionInfo.h',
+ ]
+includedir = os.path.join(env['includedir'], 'IMP' )
inst = env.Install(includedir, files)
env.Alias('install', inst)
-
-# Subdirectories
-SConscript('restraints/SConscript')
-SConscript('optimizers/SConscript')
-SConscript('decorators/SConscript')
-SConscript('internal/SConscript')
-SConscript('unary_functions/SConscript')
-SConscript('pair_scores/SConscript')
-SConscript('singleton_scores/SConscript')
-SConscript('triplet_scores/SConscript')
-SConscript('score_states/SConscript')
-SConscript('particle_refiners/SConscript')
+SConscript('decorators/SConscript' )
+SConscript('internal/SConscript' )
+SConscript('optimizers/SConscript' )
+SConscript('pair_scores/SConscript' )
+SConscript('particle_refiners/SConscript' )
+SConscript('restraints/SConscript' )
+SConscript('score_states/SConscript' )
+SConscript('singleton_scores/SConscript' )
+SConscript('triplet_scores/SConscript' )
+SConscript('unary_functions/SConscript' )
Index: kernel/include/IMP/Particle.h
===================================================================
--- kernel/include/IMP/Particle.h (revision 649)
+++ kernel/include/IMP/Particle.h (working copy)
@@ -98,7 +98,7 @@
/** \return all particle data in the model.
*/
Model* get_model() const {
- return model_.get();
+ return model_;
}
//! Add a Float attribute to this particle.
@@ -369,7 +369,8 @@
// Set pointer to model particle data.
void set_model(Model *md, ParticleIndex pi);
- Pointer<Model> model_;
+ // not ref counted as the model has a pointer to this
+ Model *model_;
// true if particle is active
bool is_active_;
@@ -421,6 +422,7 @@
{
IMP_check(get_is_active(), "Do not touch inactive particles",
InactiveParticleException);
+ IMP_assert(!is_nan(value), "Can't set attribute value to NaN " << *this);
floats_.set_value(name, value);
}
@@ -453,7 +455,8 @@
{
IMP_check(get_is_active(), "Do not touch inactive particles",
InactiveParticleException);
- IMP_assert(value==value, "Can't add NaN to derivative in particle " << *this);
+ IMP_assert(!is_nan(value),
+ "Can't add NaN to derivative in particle " << *this);
derivatives_.set_value(name, derivatives_.get_value(name)+ da(value));
}
Index: kernel/include/IMP/utility.h
===================================================================
--- kernel/include/IMP/utility.h (revision 649)
+++ 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,15 @@
return t*t;
}
+//! put in check for G++
+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/include/IMP/particle_refiners/BondCoverParticleRefiner.h
===================================================================
--- kernel/include/IMP/particle_refiners/BondCoverParticleRefiner.h (revision 649)
+++ kernel/include/IMP/particle_refiners/BondCoverParticleRefiner.h (working copy)
@@ -8,12 +8,14 @@
#ifndef __IMP_BOND_COVER_PARTICLE_REFINER_H
#define __IMP_BOND_COVER_PARTICLE_REFINER_H
+#include "../ParticleRefiner.h"
#include "../internal/kernel_version_info.h"
-#include "../ParticleRefiner.h"
namespace IMP
{
+class Particle;
+
//! Cover a bond with a constant volume set of spheres.
/** Perhaps I want to add various custom bond types so that
this will only expand some custom bonds. Currently any
@@ -23,9 +25,17 @@
{
FloatKey rk_;
FloatKey vk_;
+ IntKey tk_;
public:
+ /** \param[in] rk The key to use to store the radius
+ \param[in] vk The key to use to extract the volume from the bond
+ \param[in] tk A key whose value is propagated from the bond to the
+ created particles (if non-default) in order to set their type
+
+ */
BondCoverParticleRefiner(FloatKey rk,
- FloatKey vk);
+ FloatKey vk,
+ IntKey tk=IntKey());
virtual ~BondCoverParticleRefiner() {}
Index: kernel/include/IMP/particle_refiners/SConscript
===================================================================
--- kernel/include/IMP/particle_refiners/SConscript (revision 649)
+++ kernel/include/IMP/particle_refiners/SConscript (working copy)
@@ -1,9 +1,9 @@
Import('env')
import os.path
-
-files=['BondCoverParticleRefiner.h', 'ChildrenParticleRefiner.h']
-
-# Install the include files
-includedir = os.path.join(env['includedir'], 'IMP', 'particle_refiners')
+files=[
+ 'BondCoverParticleRefiner.h',
+ 'ChildrenParticleRefiner.h',
+ ]
+includedir = os.path.join(env['includedir'], 'IMP', 'particle_refiners' )
inst = env.Install(includedir, files)
env.Alias('install', inst)
Index: kernel/include/IMP/triplet_scores/SConscript
===================================================================
--- kernel/include/IMP/triplet_scores/SConscript (revision 649)
+++ kernel/include/IMP/triplet_scores/SConscript (working copy)
@@ -1,9 +1,8 @@
-import os.path
Import('env')
-
-files = ['AngleTripletScore.h']
-
-# Install the include files:
-includedir = os.path.join(env['includedir'], 'IMP', 'triplet_scores')
+import os.path
+files=[
+ 'AngleTripletScore.h',
+ ]
+includedir = os.path.join(env['includedir'], 'IMP', 'triplet_scores' )
inst = env.Install(includedir, files)
env.Alias('install', inst)
Index: kernel/include/IMP/ScoreState.h
===================================================================
--- kernel/include/IMP/ScoreState.h (revision 649)
+++ kernel/include/IMP/ScoreState.h (working copy)
@@ -36,6 +36,7 @@
logging is TERSE the restraint should print out only a constant number
of lines per update call.
+ \ingroup kernel
*/
class IMPDLLEXPORT ScoreState : public RefCountedObject
{
@@ -82,7 +83,7 @@
Model *get_model() const {
IMP_assert(model_,
"Must call set_model before get_model on state");
- return model_.get();
+ return model_;
}
protected:
@@ -120,8 +121,8 @@
unsigned int update_iteration_;
unsigned int after_iteration_;
- // all of the particle data
- Pointer<Model> model_;
+ // not ref counted or destructed
+ Model *model_;
std::string name_;
};
Index: kernel/include/IMP/log.h
===================================================================
--- kernel/include/IMP/log.h (revision 649)
+++ kernel/include/IMP/log.h (working copy)
@@ -183,6 +183,20 @@
+//! Print error messages
+/** If this is true, then failures of the IMP_check
+ macro will print their messages to std::cerr. This
+ should be true if the code is used from C++ and false
+ if used from python.
+ */
+IMPDLLEXPORT bool get_print_exception_messages();
+
+
+//! Turn on and off printing of error messages
+IMPDLLEXPORT void set_print_exception_messages(bool tf);
+
+
+
} // namespace IMP
#ifndef IMP_DISABLE_LOGGING
Index: kernel/include/SConscript
===================================================================
--- kernel/include/SConscript (revision 649)
+++ kernel/include/SConscript (working copy)
@@ -1,11 +1,9 @@
Import('env')
-
-files = ['IMP.h']
-
-# Install the include files:
+import os.path
+files=[
+ 'IMP.h',
+ ]
includedir = env['includedir']
inst = env.Install(includedir, files)
env.Alias('install', inst)
-
-# Subdirectories
-SConscript('IMP/SConscript')
+SConscript('IMP/SConscript' )
Index: kernel/include/IMP.h
===================================================================
--- kernel/include/IMP.h (revision 649)
+++ kernel/include/IMP.h (working copy)
@@ -1,94 +1,119 @@
/**
- * \file IMP.h \brief IMP, an Integrative Modeling Platform.
- *
- * Copyright 2007-8 Sali Lab. All rights reserved.
- *
- */
-
+* ile IMP.h rief IMP, an Integrative Modeling Platform.
+*
+* Copyright 2007-8 Sali Lab. All rights reserved.
+*
+*/
#ifndef __IMP_H
#define __IMP_H
-
-#include "IMP/IMP_config.h"
-#include "IMP/log.h"
-#include "IMP/random.h"
-#include "IMP/base_types.h"
-#include "IMP/Particle.h"
-#include "IMP/Optimizer.h"
-#include "IMP/Restraint.h"
-#include "IMP/exception.h"
-#include "IMP/Object.h"
-#include "IMP/RefCountedObject.h"
-#include "IMP/Pointer.h"
-#include "IMP/UnaryFunction.h"
-#include "IMP/unary_functions/Harmonic.h"
-#include "IMP/unary_functions/HarmonicLowerBound.h"
-#include "IMP/unary_functions/HarmonicUpperBound.h"
-#include "IMP/unary_functions/OpenCubicSpline.h"
-#include "IMP/unary_functions/ClosedCubicSpline.h"
-#include "IMP/unary_functions/Cosine.h"
-#include "IMP/unary_functions/Linear.h"
-#include "IMP/unary_functions/WormLikeChain.h"
-#include "IMP/Model.h"
-#include "IMP/PairScore.h"
-#include "IMP/SingletonScore.h"
-#include "IMP/TripletScore.h"
-#include "IMP/Vector3D.h"
-#include "IMP/VersionInfo.h"
-#include "IMP/ParticleRefiner.h"
-#include "IMP/particle_refiners/BondCoverParticleRefiner.h"
-#include "IMP/particle_refiners/ChildrenParticleRefiner.h"
-#include "IMP/decorators/HierarchyDecorator.h"
-#include "IMP/decorators/MolecularHierarchyDecorator.h"
-#include "IMP/decorators/NameDecorator.h"
-#include "IMP/decorators/AtomDecorator.h"
-#include "IMP/decorators/ResidueDecorator.h"
-#include "IMP/decorators/XYZDecorator.h"
-#include "IMP/decorators/bond_decorators.h"
-#include "IMP/optimizers/SteepestDescent.h"
-#include "IMP/optimizers/ConjugateGradients.h"
-#include "IMP/optimizers/MolecularDynamics.h"
-#include "IMP/optimizers/BrownianDynamics.h"
-#include "IMP/optimizers/MonteCarlo.h"
-#include "IMP/optimizers/Mover.h"
-#include "IMP/optimizers/MoverBase.h"
-#include "IMP/optimizers/movers/BallMover.h"
-#include "IMP/optimizers/movers/NormalMover.h"
-#include "IMP/optimizers/states/VRMLLogOptimizerState.h"
-#include "IMP/optimizers/states/CMMLogOptimizerState.h"
-#include "IMP/optimizers/states/VelocityScalingOptimizerState.h"
-#include "IMP/pair_scores/DistancePairScore.h"
-#include "IMP/pair_scores/SphereDistancePairScore.h"
-#include "IMP/pair_scores/RefineOncePairScore.h"
-#include "IMP/pair_scores/TypedPairScore.h"
-#include "IMP/pair_scores/TransformedDistancePairScore.h"
-#include "IMP/singleton_scores/DistanceToSingletonScore.h"
-#include "IMP/singleton_scores/AttributeSingletonScore.h"
-#include "IMP/singleton_scores/TunnelSingletonScore.h"
-#include "IMP/triplet_scores/AngleTripletScore.h"
-#include "IMP/restraints/RestraintSet.h"
-#include "IMP/restraints/ConstantRestraint.h"
-#include "IMP/restraints/DistanceRestraint.h"
-#include "IMP/restraints/AngleRestraint.h"
-#include "IMP/restraints/DihedralRestraint.h"
-#include "IMP/restraints/ConnectivityRestraint.h"
-#include "IMP/restraints/NonbondedRestraint.h"
-#include "IMP/restraints/BondDecoratorRestraint.h"
-#include "IMP/restraints/SingletonListRestraint.h"
-#include "IMP/restraints/PairListRestraint.h"
-#include "IMP/restraints/TripletChainRestraint.h"
-#include "IMP/restraints/PairChainRestraint.h"
-#include "IMP/score_states/BipartiteNonbondedListScoreState.h"
-#include "IMP/score_states/MaxChangeScoreState.h"
-#include "IMP/score_states/NonbondedListScoreState.h"
-#include "IMP/score_states/BondedListScoreState.h"
-#include "IMP/score_states/BondDecoratorListScoreState.h"
-#include "IMP/score_states/AllNonbondedListScoreState.h"
-#include "IMP/score_states/GravityCenterScoreState.h"
-#include "IMP/score_states/CoverBondsScoreState.h"
-
-
-/**
- \namespace IMP The IMP namespace.
- */
-
+#include <IMP/base_types.h>
+#include <IMP/DecoratorBase.h>
+#include <IMP/DerivativeAccumulator.h>
+#include <IMP/exception.h>
+#include <IMP/IMP_config.h>
+#include <IMP/Index.h>
+#include <IMP/Key.h>
+#include <IMP/log.h>
+#include <IMP/macros.h>
+#include <IMP/Model.h>
+#include <IMP/Object.h>
+#include <IMP/Optimizer.h>
+#include <IMP/OptimizerState.h>
+#include <IMP/PairScore.h>
+#include <IMP/Particle.h>
+#include <IMP/ParticleRefiner.h>
+#include <IMP/Pointer.h>
+#include <IMP/random.h>
+#include <IMP/RefCountedObject.h>
+#include <IMP/Restraint.h>
+#include <IMP/ScoreState.h>
+#include <IMP/SingletonScore.h>
+#include <IMP/TripletScore.h>
+#include <IMP/UnaryFunction.h>
+#include <IMP/utility.h>
+#include <IMP/Vector3D.h>
+#include <IMP/VersionInfo.h>
+#include <IMP/decorators/AtomDecorator.h>
+#include <IMP/decorators/bond_decorators.h>
+#include <IMP/decorators/HierarchyDecorator.h>
+#include <IMP/decorators/macros.h>
+#include <IMP/decorators/MolecularHierarchyDecorator.h>
+#include <IMP/decorators/NameDecorator.h>
+#include <IMP/decorators/ResidueDecorator.h>
+#include <IMP/decorators/utility.h>
+#include <IMP/decorators/XYZDecorator.h>
+#include <IMP/decorators/XYZRDecorator.h>
+#include <IMP/decorators/yaml.h>
+#include <IMP/internal/ArrayOnAttributesHelper.h>
+#include <IMP/internal/AttributeTable.h>
+#include <IMP/internal/bbox_nbl_helpers.h>
+#include <IMP/internal/constants.h>
+#include <IMP/internal/evaluate_distance_pair_score.h>
+#include <IMP/internal/ExponentialNumber.h>
+#include <IMP/internal/graph_base.h>
+#include <IMP/internal/Grid3D.h>
+#include <IMP/internal/kernel_version_info.h>
+#include <IMP/internal/MinimalSet.h>
+#include <IMP/internal/ObjectContainer.h>
+#include <IMP/internal/ParticleGrid.h>
+#include <IMP/internal/ref_counting.h>
+#include <IMP/internal/Unit.h>
+#include <IMP/internal/units.h>
+#include <IMP/internal/utility.h>
+#include <IMP/internal/Vector.h>
+#include <IMP/optimizers/BrownianDynamics.h>
+#include <IMP/optimizers/ConjugateGradients.h>
+#include <IMP/optimizers/MolecularDynamics.h>
+#include <IMP/optimizers/MonteCarlo.h>
+#include <IMP/optimizers/MoverBase.h>
+#include <IMP/optimizers/Mover.h>
+#include <IMP/optimizers/SteepestDescent.h>
+#include <IMP/pair_scores/DistancePairScore.h>
+#include <IMP/pair_scores/RefineHighPairScore.h>
+#include <IMP/pair_scores/RefineOncePairScore.h>
+#include <IMP/pair_scores/SphereDistancePairScore.h>
+#include <IMP/pair_scores/TransformedDistancePairScore.h>
+#include <IMP/pair_scores/TypedPairScore.h>
+#include <IMP/particle_refiners/BondCoverParticleRefiner.h>
+#include <IMP/particle_refiners/ChildrenParticleRefiner.h>
+#include <IMP/restraints/AngleRestraint.h>
+#include <IMP/restraints/BondDecoratorRestraint.h>
+#include <IMP/restraints/ConnectivityRestraint.h>
+#include <IMP/restraints/ConstantRestraint.h>
+#include <IMP/restraints/DihedralRestraint.h>
+#include <IMP/restraints/DistanceRestraint.h>
+#include <IMP/restraints/LowestNPairListRestraint.h>
+#include <IMP/restraints/NonbondedRestraint.h>
+#include <IMP/restraints/PairChainRestraint.h>
+#include <IMP/restraints/PairListRestraint.h>
+#include <IMP/restraints/PairRestraint.h>
+#include <IMP/restraints/RestraintSet.h>
+#include <IMP/restraints/SingletonListRestraint.h>
+#include <IMP/restraints/TripletChainRestraint.h>
+#include <IMP/score_states/AllNonbondedListScoreState.h>
+#include <IMP/score_states/BipartiteNonbondedListScoreState.h>
+#include <IMP/score_states/BondDecoratorListScoreState.h>
+#include <IMP/score_states/BondedListScoreState.h>
+#include <IMP/score_states/CoverBondsScoreState.h>
+#include <IMP/score_states/GravityCenterScoreState.h>
+#include <IMP/score_states/ManualBondDecoratorListScoreState.h>
+#include <IMP/score_states/MaxChangeScoreState.h>
+#include <IMP/score_states/NonbondedListScoreState.h>
+#include <IMP/singleton_scores/AttributeSingletonScore.h>
+#include <IMP/singleton_scores/DistanceToSingletonScore.h>
+#include <IMP/singleton_scores/TunnelSingletonScore.h>
+#include <IMP/triplet_scores/AngleTripletScore.h>
+#include <IMP/unary_functions/ClosedCubicSpline.h>
+#include <IMP/unary_functions/Cosine.h>
+#include <IMP/unary_functions/Harmonic.h>
+#include <IMP/unary_functions/HarmonicLowerBound.h>
+#include <IMP/unary_functions/HarmonicUpperBound.h>
+#include <IMP/unary_functions/Linear.h>
+#include <IMP/unary_functions/OpenCubicSpline.h>
+#include <IMP/unary_functions/WormLikeChain.h>
+#include <IMP/optimizers/movers/BallMover.h>
+#include <IMP/optimizers/movers/NormalMover.h>
+#include <IMP/optimizers/states/CMMLogOptimizerState.h>
+#include <IMP/optimizers/states/VelocityScalingOptimizerState.h>
+#include <IMP/optimizers/states/VRMLLogOptimizerState.h>
#endif /* __IMP_H */
Index: kernel/doc/examples/chain.py
===================================================================
--- kernel/doc/examples/chain.py (revision 649)
+++ kernel/doc/examples/chain.py (working copy)
@@ -9,18 +9,17 @@
#IMP.set_log_level(IMP.VERBOSE)
np=20
radius =1.0
-rk= IMP.FloatKey("radius")
m= IMP.Model()
# The particles in the chain
chain= IMP.Particles()
for i in range(0,np):
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= IMP.XYZRDecorator.create(p)
+ 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)
+ d.set_radius(radius)
chain.append(p)
# create a bond between successive particles
@@ -36,15 +35,14 @@
p.show()
# Set up the nonbonded list
-nbl= IMP.AllNonbondedListScoreState(1, rk, chain)
+nbl= IMP.AllNonbondedListScoreState(1, chain)
nbli= m.add_score_state(nbl)
# This ScoreState uses the bonds constructed above to restrain
bl= IMP.BondDecoratorListScoreState(chain)
bli= nbl.add_bonded_list(bl)
# Set up excluded volume
-ps= IMP.SphereDistancePairScore(IMP.HarmonicLowerBound(0,1),
- rk)
+ps= IMP.SphereDistancePairScore(IMP.HarmonicLowerBound(0,1))
evr= IMP.NonbondedRestraint(ps, nbl)
evri= m.add_restraint(evr)
@@ -55,7 +53,7 @@
# Just for fun to make the chain straight (angles in radians)
ats= IMP.AngleTripletScore(IMP.Harmonic(3.1415, .1))
ar= IMP.TripletChainRestraint(ats)
-ar.add_chain(chain)
+ar.set_particles(chain)
ari= m.add_restraint(ar)
# Tie the ends of the chain
@@ -63,7 +61,7 @@
p= IMP.ParticlePair(chain[0], chain[-1])
pps= IMP.ParticlePairs()
pps.append(p)
-cr= IMP.PairListRestraint(IMP.SphereDistancePairScore(IMP.Harmonic(3,1), rk),
+cr= IMP.PairListRestraint(IMP.SphereDistancePairScore(IMP.Harmonic(3,1)),
pps)
cri=m.add_restraint(cr)
@@ -74,7 +72,6 @@
# Write the progression of states as the system is optimized to
# the files state.000.vrml, state.001.vrml etc.
vrml= IMP.VRMLLogOptimizerState("state.%03d.vrml", chain)
-vrml.set_radius_key(rk)
vrml.update()
vrml.set_skip_steps(100)
o.add_optimizer_state(vrml)
Index: kernel/src/exception.cpp
===================================================================
--- kernel/src/exception.cpp (revision 649)
+++ kernel/src/exception.cpp (working copy)
@@ -10,6 +10,12 @@
namespace IMP
{
+namespace internal {
+ // The error message is already in the exception
+ bool print_exceptions=true;
+}
+
+
static CheckLevel check_mode =
#ifdef NDEBUG
CHEAP;
@@ -30,12 +36,9 @@
namespace internal
{
-// The error message is already in the exception
-bool print_exceptions=false;
-
void assert_fail(const char *msg)
{
- if (print_exceptions) {
+ if (internal::print_exceptions) {
IMP_ERROR(msg);
}
throw ErrorException(msg);
@@ -43,11 +46,14 @@
void check_fail(const char *msg)
{
- if (print_exceptions) {
+ if (internal::print_exceptions) {
IMP_ERROR(msg);
}
}
} // namespace internal
+
+
+
} // namespace IMP
Index: kernel/src/singleton_scores/AttributeSingletonScore.cpp
===================================================================
--- kernel/src/singleton_scores/AttributeSingletonScore.cpp (revision 649)
+++ kernel/src/singleton_scores/AttributeSingletonScore.cpp (working copy)
@@ -9,6 +9,8 @@
#include "IMP/UnaryFunction.h"
#include "IMP/Particle.h"
+#include <boost/tuple/tuple.hpp>
+
namespace IMP
{
@@ -20,8 +22,8 @@
DerivativeAccumulator *da) const
{
if (da) {
- Float d;
- float r= f_->evaluate_with_derivative(b->get_value(k_), d);
+ Float r,d;
+ boost::tie( r,d)= 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 649)
+++ kernel/src/singleton_scores/TunnelSingletonScore.cpp (working copy)
@@ -9,6 +9,8 @@
#include "IMP/singleton_scores/TunnelSingletonScore.h"
#include "IMP/decorators/XYZDecorator.h"
+#include <boost/tuple/tuple.hpp>
+
namespace IMP
{
@@ -58,7 +60,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/singleton_scores/SConscript
===================================================================
--- kernel/src/singleton_scores/SConscript (revision 649)
+++ kernel/src/singleton_scores/SConscript (working copy)
@@ -1,7 +1,7 @@
Import('env')
-
-files = ['DistanceToSingletonScore.cpp', 'AttributeSingletonScore.cpp',
- 'TunnelSingletonScore.cpp']
-
-files = [File(x) for x in files]
+files=[
+ File( 'AttributeSingletonScore.cpp' ),
+ File( 'DistanceToSingletonScore.cpp' ),
+ File( 'TunnelSingletonScore.cpp' ),
+ ]
Return('files')
Index: kernel/src/Model.cpp
===================================================================
--- kernel/src/Model.cpp (revision 649)
+++ kernel/src/Model.cpp (working copy)
@@ -79,8 +79,8 @@
for (RestraintIterator it = restraints_begin();
it != restraints_end(); ++it) {
IMP_CHECK_OBJECT(*it);
- IMP_LOG(TERSE, "Evaluate restraint "
- << std::endl << **it);
+ IMP_LOG(TERSE, "Evaluate restraint: "
+ << std::distance(restraints_begin(), it) << std::endl);
Float tscore=0;
if ((*it)->get_is_active()) {
tscore = (*it)->evaluate(accpt);
Index: kernel/src/Particle.cpp
===================================================================
--- kernel/src/Particle.cpp (revision 649)
+++ kernel/src/Particle.cpp (working copy)
@@ -15,6 +15,7 @@
Particle::Particle()
{
+ model_=NULL;
is_active_ = true;
}
@@ -43,17 +44,15 @@
void Particle::zero_derivatives()
{
- for (FloatKeyIterator it= float_keys_begin(); it != float_keys_end(); ++it) {
- derivatives_.set_value(*it, 0);
- }
+ derivatives_.set_values(0);
}
void Particle::show(std::ostream& out) const
{
- const char* inset = " ";
+ const std::string inset(" ");
out << std::endl;
- out << "--" << get_index() << "--" << std::endl;
+ out << "Particle: " << get_index() << std::endl;
if (is_active_) {
out << inset << inset << "active";
} else {
@@ -62,18 +61,24 @@
out << std::endl;
if (get_model() != NULL) {
- out << inset << inset << "float attributes:" << std::endl;
- floats_.show(out, " ");
+ out << inset << "float attributes:" << std::endl;
+ floats_.show(out, inset+inset);
- out << inset << inset << "int attributes:" << std::endl;
- ints_.show(out, " ");
+ out << inset << "float derivatives:" << std::endl;
+ derivatives_.show(out, inset+inset);
- out << inset << inset << "string attributes:" << std::endl;
- strings_.show(out, " ");
+ out << inset << "optimizeds:" << std::endl;
+ optimizeds_.show(out, inset+inset);
- out << inset << inset << "particle attributes:" << std::endl;
- particles_.show(out, " ");
+ out << inset << "int attributes:" << std::endl;
+ ints_.show(out, inset+inset);
+ out << inset << "string attributes:" << std::endl;
+ strings_.show(out, inset+inset);
+
+ out << inset << "particle attributes:" << std::endl;
+ particles_.show(out, inset+inset);
+
}
}
Index: kernel/src/unary_functions/Cosine.cpp
===================================================================
--- kernel/src/unary_functions/Cosine.cpp (revision 649)
+++ 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 649)
+++ 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/SConscript
===================================================================
--- kernel/src/unary_functions/SConscript (revision 649)
+++ kernel/src/unary_functions/SConscript (working copy)
@@ -1,6 +1,7 @@
Import('env')
-
-files = ['OpenCubicSpline.cpp', 'ClosedCubicSpline.cpp', 'Cosine.cpp']
-
-files = [File(x) for x in files]
+files=[
+ File( 'ClosedCubicSpline.cpp' ),
+ File( 'Cosine.cpp' ),
+ File( 'OpenCubicSpline.cpp' ),
+ ]
Return('files')
Index: kernel/src/unary_functions/ClosedCubicSpline.cpp
===================================================================
--- kernel/src/unary_functions/ClosedCubicSpline.cpp (revision 649)
+++ 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/SConscript
===================================================================
--- kernel/src/SConscript (revision 649)
+++ 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/Restraint.cpp
===================================================================
--- kernel/src/Restraint.cpp (revision 649)
+++ kernel/src/Restraint.cpp (working copy)
@@ -17,6 +17,7 @@
Restraint::Restraint()
{
+ model_=NULL;
is_active_ = true; // active by default
}
@@ -24,6 +25,12 @@
//! Destructor
Restraint::~Restraint()
{
+ if (!was_owned_) {
+ // can't use virtual functions in the destructor
+ std::cerr << "Restraint " << this << " is being destroyed "
+ << "without ever having been added to a model."
+ << std::endl;
+ }
}
@@ -43,9 +50,10 @@
{
IMP_assert(model==NULL || get_number_of_particles()==0
|| model == get_particle(0)->get_model()
- || (model_ && model_.get() == model),
+ || (model_ && model_ == model),
"Model* different from Particle Model*");
model_=model;
+ was_owned_=true;
}
void Restraint::show(std::ostream& out) const
@@ -59,6 +67,10 @@
get_version_info().show(out);
}
+ParticlesList Restraint::get_interacting_particles() const {
+ return ParticlesList(1, Particles(particles_begin(),
+ particles_end()));
+}
// The index line is to disable a warning
IMP_LIST_IMPL(Restraint, Particle, particle,Particle*, {
Index: kernel/src/ScoreState.cpp
===================================================================
--- kernel/src/ScoreState.cpp (revision 649)
+++ kernel/src/ScoreState.cpp (working copy)
@@ -17,6 +17,7 @@
//! Constructor
ScoreState::ScoreState(std::string name) : name_(name)
{
+ model_=NULL;
update_iteration_= std::numeric_limits<unsigned int>::max();
after_iteration_= std::numeric_limits<unsigned int>::max();
IMP_LOG(VERBOSE, "ScoreState constructed " << name << std::endl);
Index: kernel/src/score_states/BondDecoratorListScoreState.cpp
===================================================================
--- kernel/src/score_states/BondDecoratorListScoreState.cpp (revision 649)
+++ kernel/src/score_states/BondDecoratorListScoreState.cpp (working copy)
@@ -21,7 +21,7 @@
{
IMP_LOG(TERSE, "Updating BondDecoratorList for "
<< ps_.size() << " particles" << std::endl);
- bonds_.clear();
+ std::vector<BondDecorator> bonds;
for (unsigned int i=0; i< ps_.size(); ++i) {
if (!ps_[i]->get_is_active()) continue;
BondedDecorator di(ps_[i]);
@@ -36,31 +36,26 @@
}
if (di < dj) {
IMP_LOG(VERBOSE, "Found bond " << di.get_bond(j) << std::endl);
- bonds_.push_back(di.get_bond(j));
+ bonds.push_back(di.get_bond(j));
}
}
}
- IMP_LOG(TERSE, "Found " << bonds_.size() << " bonds"<< std::endl);
+ ManualBondDecoratorListScoreState::set_bond_decorators(bonds);
+ IMP_LOG(TERSE, "Found " << bonds.size() << " bonds"<< std::endl);
}
void BondDecoratorListScoreState::set_particles(const Particles &ps)
{
ps_=ps;
std::sort(ps_.begin(), ps_.end());
- bonds_.clear();
+ clear_bond_decorators();
}
-
-bool BondDecoratorListScoreState::are_bonded(Particle *a, Particle *b) const
+void BondDecoratorListScoreState::add_particles(const Particles &ps)
{
- try {
- BondedDecorator da= BondedDecorator::cast(a);
- BondedDecorator db= BondedDecorator::cast(b);
- return get_bond(da, db) != BondDecorator();
- } catch (...) {
- IMP_LOG(VERBOSE, "Exception thrown in are_bonded"<< std::endl);
- }
- return false;
+ ps_.insert(ps_.end(), ps.begin(), ps.end());
+ std::sort(ps_.begin(), ps_.end());
+ clear_bond_decorators();
}
} // namespace IMP
Index: kernel/src/score_states/NonbondedListScoreState.cpp
===================================================================
--- kernel/src/score_states/NonbondedListScoreState.cpp (revision 649)
+++ kernel/src/score_states/NonbondedListScoreState.cpp (working copy)
@@ -6,7 +6,7 @@
*/
#include "IMP/score_states/NonbondedListScoreState.h"
-#include "IMP/decorators/XYZDecorator.h"
+#include "IMP/decorators/XYZRDecorator.h"
#include "IMP/internal/Grid3D.h"
#include "IMP/score_states/MaxChangeScoreState.h"
@@ -32,16 +32,15 @@
NonbondedListScoreState
-::NonbondedListScoreState(Float cut,
- FloatKey rk): rk_(rk),
- cutoff_(cut),
- nbl_is_valid_(false)
+::NonbondedListScoreState(Float cut): cutoff_(cut),
+ nbl_is_valid_(false)
{
slack_=cutoff_;
number_of_updates_=1;
number_of_rebuilds_=0;
number_of_overflows_=0;
max_nbl_size_= std::numeric_limits<unsigned int>::max();
+ set_radius_key(XYZRDecorator::get_radius_key());
}
@@ -51,6 +50,48 @@
{
}
+void NonbondedListScoreState::set_radius_key(FloatKey rk) {
+ rk_=rk;
+ Particles old_particles;
+ if (mcr_) {
+ old_particles= mcr_->get_particles();
+ }
+ if (rk != FloatKey()) {
+ FloatKeys ks;
+ ks.push_back(rk);
+ mcr_= new MaxChangeScoreState(ks);
+ } else {
+ mcr_= new MaxChangeScoreState(FloatKeys());
+ }
+ mcr_->add_particles(old_particles);
+}
+
+
+Float NonbondedListScoreState::get_max_radius_change() const {
+ if (rk_ != FloatKey()) {
+ return mcr_->get_max();
+ } else {
+ return 0;
+ }
+}
+
+void NonbondedListScoreState::reset_max_radius_change() {
+ if (mcr_) mcr_->reset();
+}
+
+void NonbondedListScoreState::update_max_radius_change(unsigned int i) {
+ if (mcr_) mcr_->before_evaluate(i);
+}
+
+void NonbondedListScoreState
+::add_particles_for_max_radius_change(const Particles &ps) {
+ mcr_->add_particles(particles_with_radius(ps));
+}
+
+void NonbondedListScoreState::clear_particles_for_max_radius_change() {
+ mcr_->clear_particles();
+}
+
void NonbondedListScoreState::show_statistics(std::ostream &out) const
{
out << "Nonbonded list averaged "
Index: kernel/src/score_states/CoverBondsScoreState.cpp
===================================================================
--- kernel/src/score_states/CoverBondsScoreState.cpp (revision 649)
+++ kernel/src/score_states/CoverBondsScoreState.cpp (working copy)
@@ -14,8 +14,9 @@
{
-CoverBondsScoreState::CoverBondsScoreState(BondDecoratorListScoreState *bl,
- FloatKey rk): bl_(bl), rk_(rk)
+CoverBondsScoreState
+::CoverBondsScoreState(ManualBondDecoratorListScoreState *bl,
+ FloatKey rk): bl_(bl), rk_(rk)
{
}
@@ -25,8 +26,9 @@
void CoverBondsScoreState::do_before_evaluate()
{
- for (BondDecoratorListScoreState::BondIterator it= bl_->bonds_begin();
- it != bl_->bonds_end(); ++it) {
+ for (BondDecoratorListScoreState::BondDecoratorIterator it
+ = bl_->bond_decorators_begin();
+ it != bl_->bond_decorators_end(); ++it) {
BondDecorator bd= *it;
BondedDecorator pa= bd.get_bonded(0);
BondedDecorator pb= bd.get_bonded(1);
@@ -62,8 +64,9 @@
void CoverBondsScoreState::after_evaluate(DerivativeAccumulator *dva)
{
if (dva) {
- for (BondDecoratorListScoreState::BondIterator it= bl_->bonds_begin();
- it != bl_->bonds_end(); ++it) {
+ for (BondDecoratorListScoreState::BondDecoratorIterator it
+ = bl_->bond_decorators_begin();
+ it != bl_->bond_decorators_end(); ++it) {
XYZDecorator d(it->get_particle());
Vector3D deriv;
// divide derivatives equally between endpoints
@@ -84,7 +87,7 @@
void CoverBondsScoreState::show(std::ostream &out) const
{
- out << "CoverBondsScoreState on " << bl_->get_number_of_bonds()
+ out << "CoverBondsScoreState on " << bl_->get_number_of_bond_decorators()
<< " bonds " << std::endl;
}
Index: kernel/src/score_states/BipartiteNonbondedListScoreState.cpp
===================================================================
--- kernel/src/score_states/BipartiteNonbondedListScoreState.cpp (revision 649)
+++ kernel/src/score_states/BipartiteNonbondedListScoreState.cpp (working copy)
@@ -15,8 +15,10 @@
namespace IMP
{
-//! Turn the default into an actual algorithm and work around missing algorithms
-/** This cannot be shared with AllNBL because that one has grid and this does
+
+
+// Turn the default into an actual algorithm and work around missing algorithms
+/* This cannot be shared with AllNBL because that one has grid and this does
not.
*/
static BipartiteNonbondedListScoreState::Algorithm
@@ -45,30 +47,21 @@
BipartiteNonbondedListScoreState
-::BipartiteNonbondedListScoreState(Float cut,
- FloatKey rk,
- Algorithm a):
- P(cut, rk), a_(translate_algorithm(a))
+::BipartiteNonbondedListScoreState(Float cut):
+ P(cut), a_(translate_algorithm(DEFAULT))
{
- FloatKeys ks;
- ks.push_back(rk);
mc0_= new MaxChangeScoreState(XYZDecorator::get_xyz_keys());
mc1_= new MaxChangeScoreState(XYZDecorator::get_xyz_keys());
- mcr_= new MaxChangeScoreState(ks);
}
BipartiteNonbondedListScoreState
-::BipartiteNonbondedListScoreState(Float cut, FloatKey rk,
+::BipartiteNonbondedListScoreState(Float cut,
const Particles &ps0,
- const Particles &ps1,
- Algorithm a):
- P(cut, rk), a_(translate_algorithm(a))
+ const Particles &ps1):
+ P(cut), a_(translate_algorithm(DEFAULT))
{
- FloatKeys ks;
- ks.push_back(rk);
mc0_= new MaxChangeScoreState(XYZDecorator::get_xyz_keys());
mc1_= new MaxChangeScoreState(XYZDecorator::get_xyz_keys());
- mcr_= new MaxChangeScoreState(ks);
set_particles(ps0, ps1);
}
@@ -79,7 +72,7 @@
{
mc0_->before_evaluate(ScoreState::get_before_evaluate_iteration());
mc0_->before_evaluate(ScoreState::get_before_evaluate_iteration());
- mcr_->before_evaluate(ScoreState::get_before_evaluate_iteration());
+ update_max_radius_change(ScoreState::get_before_evaluate_iteration());
Float mc= std::max(mc0_->get_max(), mc1_->get_max());
Float cost;
@@ -88,17 +81,19 @@
cost= 10*mc0_->get_number_of_particles()* mc1_->get_number_of_particles();
break;
case BBOX:
- cost= 1000 * mcr_->get_number_of_particles();
+ cost= 1000 * (mc0_->get_number_of_particles()
+ + mc1_->get_number_of_particles());
break;
default:
IMP_assert(0, "Bad algorithm");
- cost= 1000 * mcr_->get_number_of_particles();
+ cost= 1000 * (mc0_->get_number_of_particles()
+ + mc1_->get_number_of_particles());
}
- if (P::update(mc+ mcr_->get_max(), cost)) {
+ if (P::update(mc+ get_max_radius_change(), cost)) {
mc0_->reset();
mc1_->reset();
- mcr_->reset();
+ reset_max_radius_change();
}
IMP_IF_CHECK(EXPENSIVE) {
check_nbl();
@@ -150,9 +145,9 @@
mc0_->add_particles(ps0);
mc1_->clear_particles();
mc1_->add_particles(ps1);
- mcr_->clear_particles();
- mcr_->add_particles(P::particles_with_radius(ps0));
- mcr_->add_particles(P::particles_with_radius(ps1));
+ clear_particles_for_max_radius_change();
+ add_particles_for_max_radius_change(ps0);
+ add_particles_for_max_radius_change(ps1);
P::set_nbl_is_valid(false);
}
@@ -160,24 +155,36 @@
void BipartiteNonbondedListScoreState::add_particles_0(const Particles &ps)
{
if (P::get_nbl_is_valid()) process_sets(ps, mc1_->get_particles());
+ else P::set_nbl_is_valid(false);
mc0_->add_particles(ps);
- mcr_->add_particles(P::particles_with_radius(ps));
- P::set_nbl_is_valid(false);
+ add_particles_for_max_radius_change(ps);
}
void BipartiteNonbondedListScoreState::add_particles_1(const Particles &ps)
{
if (P::get_nbl_is_valid()) process_sets(ps, mc0_->get_particles());
+ else P::set_nbl_is_valid(false);
mc1_->add_particles(ps);
- mcr_->add_particles(P::particles_with_radius(ps));
- P::set_nbl_is_valid(false);
+ add_particles_for_max_radius_change(ps);
}
+
+void BipartiteNonbondedListScoreState::add_particle_0(Particle *ps)
+{
+ add_particles_0(Particles(1, ps));
+}
+
+void BipartiteNonbondedListScoreState::add_particle_1(Particle *ps)
+{
+ add_particles_1(Particles(1, ps));
+}
+
+
void BipartiteNonbondedListScoreState::clear_particles()
{
mc0_->clear_particles();
mc1_->clear_particles();
- mcr_->clear_particles();
+ clear_particles_for_max_radius_change();
P::set_nbl_is_valid(false);
P::set_nbl_is_valid(true);
}
Index: kernel/src/score_states/SConscript
===================================================================
--- kernel/src/score_states/SConscript (revision 649)
+++ kernel/src/score_states/SConscript (working copy)
@@ -1,10 +1,12 @@
Import('env')
-
-files = ['NonbondedListScoreState.cpp',
- 'MaxChangeScoreState.cpp', 'BondDecoratorListScoreState.cpp',
- 'AllNonbondedListScoreState.cpp',
- 'BipartiteNonbondedListScoreState.cpp',
- 'GravityCenterScoreState.cpp', 'CoverBondsScoreState.cpp']
-
-files = [File(x) for x in files]
+files=[
+ File( 'AllNonbondedListScoreState.cpp' ),
+ File( 'BipartiteNonbondedListScoreState.cpp' ),
+ File( 'BondDecoratorListScoreState.cpp' ),
+ File( 'CoverBondsScoreState.cpp' ),
+ File( 'GravityCenterScoreState.cpp' ),
+ File( 'ManualBondDecoratorListScoreState.cpp' ),
+ File( 'MaxChangeScoreState.cpp' ),
+ File( 'NonbondedListScoreState.cpp' ),
+ ]
Return('files')
Index: kernel/src/score_states/AllNonbondedListScoreState.cpp
===================================================================
--- kernel/src/score_states/AllNonbondedListScoreState.cpp (revision 649)
+++ kernel/src/score_states/AllNonbondedListScoreState.cpp (working copy)
@@ -18,7 +18,7 @@
{
-//! Turn the default into an actual algorithm and work around missing algorithms
+// Turn the default into an actual algorithm and work around missing algorithms
static AllNonbondedListScoreState::Algorithm
translate_algorithm(AllNonbondedListScoreState::Algorithm a)
{
@@ -44,29 +44,13 @@
AllNonbondedListScoreState::AllNonbondedListScoreState(Float cut,
- FloatKey rk,
- const Particles &ps,
- Algorithm a):
- P(cut, rk), a_(translate_algorithm(a))
+ const Particles &ps):
+ P(cut), a_(translate_algorithm(DEFAULT))
{
- FloatKeys ks;
- ks.push_back(rk);
mc_= new MaxChangeScoreState(XYZDecorator::get_xyz_keys());
- mcr_= new MaxChangeScoreState(ks);
add_particles(ps);
}
-AllNonbondedListScoreState::AllNonbondedListScoreState(Float cut,
- FloatKey rk,
- Algorithm a):
- P(cut, rk), a_(translate_algorithm(a))
-{
- FloatKeys ks;
- ks.push_back(rk);
- mc_= new MaxChangeScoreState(XYZDecorator::get_xyz_keys());
- mcr_= new MaxChangeScoreState(ks);
-}
-
AllNonbondedListScoreState::~AllNonbondedListScoreState()
{
}
@@ -74,8 +58,8 @@
void AllNonbondedListScoreState::do_before_evaluate()
{
mc_->before_evaluate(ScoreState::get_before_evaluate_iteration());
- mcr_->before_evaluate(ScoreState::get_before_evaluate_iteration());
- Float mc=mc_->get_max()+ mcr_->get_max();
+ update_max_radius_change(ScoreState::get_before_evaluate_iteration());
+ Float mc=mc_->get_max()+ get_max_radius_change();
Float cost;
switch (a_){
case QUADRATIC:
@@ -94,7 +78,7 @@
}
if (P::update(mc, cost)) {
mc_->reset();
- mcr_->reset();
+ reset_max_radius_change();
}
IMP_IF_CHECK(EXPENSIVE) {
check_nbl();
@@ -132,8 +116,8 @@
{
mc_->clear_particles();
mc_->add_particles(ps);
- mcr_->clear_particles();
- mcr_->add_particles(particles_with_radius(ps));
+ clear_particles_for_max_radius_change();
+ add_particles_for_max_radius_change(ps);
P::set_nbl_is_valid(false);
}
@@ -155,14 +139,14 @@
internal::NBLAddPairIfNonbonded(this));
}
}
- mcr_->add_particles(particles_with_radius(ps));
+ add_particles_for_max_radius_change(ps);
mc_->add_particles(ps);
}
void AllNonbondedListScoreState::clear_particles()
{
mc_->clear_particles();
- mcr_->clear_particles();
+ clear_particles_for_max_radius_change();
P::set_nbl_is_valid(false);
P::set_nbl_is_valid(true);
}
Index: kernel/src/score_states/MaxChangeScoreState.cpp
===================================================================
--- kernel/src/score_states/MaxChangeScoreState.cpp (revision 649)
+++ kernel/src/score_states/MaxChangeScoreState.cpp (working copy)
@@ -19,12 +19,6 @@
const Particles &ps): keys_(keys)
{
add_particles(ps);
- origkeys_.resize(keys_.size());
- std::ostringstream oss;
- oss << " MCSS base " << this;
- for (unsigned int i=0; i< keys_.size(); ++i) {
- origkeys_[i]= FloatKey((keys_[i].get_string()+oss.str()).c_str());
- }
}
@@ -35,12 +29,6 @@
<< obj,
ValueException);
};
- for (unsigned int i=0; i< origkeys_.size(); ++i) {
- if (!obj->has_attribute(origkeys_[i])) {
- obj->add_attribute(origkeys_[i],
- obj->get_value(keys_[i]), false);
- }
- }
}, {reset();});
void MaxChangeScoreState::do_before_evaluate()
@@ -48,14 +36,16 @@
max_change_=0;
// get rid of inactive particles and their stored values
internal::remove_inactive_particles(particle_vector_);
- for (ParticleIterator it= particles_begin(); it != particles_end(); ++it) {
- (*it)->assert_is_valid();
+ for (unsigned int i=0; i< get_number_of_particles(); ++i) {
+ get_particle(i)->assert_is_valid();
for (unsigned int j=0; j < keys_.size(); ++j) {
- Float v= (*it)->get_value(keys_[j]);
- Float ov= (*it)->get_value(origkeys_[j]);
- IMP_LOG(VERBOSE, "Particle " << (*it)->get_index()
- << " and attribute " << keys_[j]
- << " moved " << std::abs(v - ov) << std::endl);
+ Float v= get_particle(i)->get_value(keys_[j]);
+ Float ov= old_values_[i].get_value(keys_[j]);
+ if (std::abs(v - ov) != 0) {
+ IMP_LOG(VERBOSE, "Particle " << get_particle(i)->get_index()
+ << " and attribute " << keys_[j]
+ << " moved " << std::abs(v - ov) << std::endl);
+ }
max_change_= std::max(max_change_,
std::abs(v-ov));
}
@@ -66,9 +56,15 @@
void MaxChangeScoreState::reset()
{
- for (ParticleIterator it= particles_begin(); it != particles_end(); ++it) {
+ Table t;
+ for (unsigned int i=0; i < keys_.size(); ++i) {
+ t.insert(keys_[i], 0);
+ }
+ old_values_.resize(get_number_of_particles(), t);
+ for (unsigned int i=0; i< get_number_of_particles(); ++i) {
for (unsigned int j=0; j < keys_.size(); ++j) {
- (*it)->set_value(origkeys_[j], (*it)->get_value(keys_[j]));
+ old_values_[i].set_value(keys_[j],
+ get_particle(i)->get_value(keys_[j]));
}
}
max_change_=0;
Index: kernel/src/score_states/ManualBondDecoratorListScoreState.cpp
===================================================================
--- kernel/src/score_states/ManualBondDecoratorListScoreState.cpp (revision 0)
+++ kernel/src/score_states/ManualBondDecoratorListScoreState.cpp (revision 0)
@@ -0,0 +1,41 @@
+#include "IMP/score_states/ManualBondDecoratorListScoreState.h"
+
+namespace IMP
+{
+ ManualBondDecoratorListScoreState
+ ::ManualBondDecoratorListScoreState(const BondDecorators &ps) {
+ set_bond_decorators(ps);
+ }
+
+ IMP_LIST_IMPL(ManualBondDecoratorListScoreState, BondDecorator,
+ bond_decorator, BondDecorator,
+ ,);
+
+ bool ManualBondDecoratorListScoreState::are_bonded(Particle *a,
+ Particle *b) const
+ {
+ if (BondedDecorator::is_instance_of(a)
+ && BondedDecorator::is_instance_of(b)) {
+ BondedDecorator da(a);
+ BondedDecorator db(b);
+ return get_bond(da, db) != BondDecorator();
+ } else if (BondDecorator::is_instance_of(a)
+ && BondDecorator::is_instance_of(b)){
+ // check if they are bonds which share an endpoint
+ /**
+ \todo Decide if the check for adject bonds should go elsewhere
+ */
+ BondDecorator da(a);
+ BondDecorator db(b);
+ Particle *a0= da.get_bonded(0).get_particle();
+ Particle *a1= da.get_bonded(1).get_particle();
+ Particle *b0= db.get_bonded(0).get_particle();
+ Particle *b1= db.get_bonded(1).get_particle();
+
+ return (a0== b0 || a0 == b1 || a1== b0 || a1 == b1);
+ } else {
+ return false;
+ }
+ }
+
+}
Index: kernel/src/restraints/LowestNPairListRestraint.cpp
===================================================================
--- kernel/src/restraints/LowestNPairListRestraint.cpp (revision 0)
+++ kernel/src/restraints/LowestNPairListRestraint.cpp (revision 0)
@@ -0,0 +1,90 @@
+/**
+ * \file LowestNPairListRestraint.cpp
+ * \brief Apply a score function to a list of pairs of particles.
+ *
+ * Copyright 2007-8 Sali Lab. All rights reserved.
+ *
+ */
+
+#include "IMP/restraints/LowestNPairListRestraint.h"
+#include "IMP/PairScore.h"
+#include "IMP/log.h"
+#include "IMP/internal/MinimalSet.h"
+
+#include <cmath>
+#include <algorithm>
+
+namespace IMP
+{
+
+LowestNPairListRestraint::LowestNPairListRestraint(PairScore *s,
+ unsigned int n,
+ const ParticlePairs &ps) :
+ PairListRestraint(s, ps), n_(n)
+{
+}
+
+LowestNPairListRestraint::~LowestNPairListRestraint()
+{
+}
+
+namespace {
+
+struct CompareFirst {
+ template <class T>
+ bool operator()(const T &a, const T &b) const {
+ return a.first < b.first;
+ }
+};
+
+}
+
+Float LowestNPairListRestraint::evaluate(DerivativeAccumulator *accum)
+{
+ IMP_CHECK_OBJECT(ss_.get());
+ IMP_assert(get_number_of_particles()%2 == 0, "There should be an even number"
+ << " of particles");
+
+ internal::MinimalSet<float, unsigned int> bestn(n_);
+
+ IMP_LOG(TERSE, "Finding lowest " << n_ << " of "
+ << get_number_of_particles()/2 << std::endl);
+
+ for (unsigned int i=0; i< get_number_of_particles(); i+=2) {
+ float score= ss_->evaluate(get_particle(i), get_particle(i+1), NULL);
+ IMP_LOG(VERBOSE, "Pair " << get_particle(i)->get_index()
+ << " and " << get_particle(i+1)->get_index()
+ << " has score " << score << std::endl);
+
+ if (bestn.can_insert(score)) {
+ bestn.insert(score, i);
+ }
+ }
+
+ float score=0;
+ IMP_assert(bestn.size() == std::min(get_number_of_particles()/2, n_),
+ "Something wrong in finding best set of edges");
+ for (unsigned int i=0; i< bestn.size(); ++i) {
+ if (accum) {
+ // to set the derivs
+ ss_->evaluate(get_particle(bestn[i].second),
+ get_particle(bestn[i].second+1), accum);
+ }
+ score+= bestn[i].first;
+ IMP_LOG(VERBOSE, "Chose pair "
+ << get_particle(bestn[i].second)->get_index()
+ << " and " << get_particle(bestn[i].second+1)->get_index()
+ << " with score " << bestn[i].first << std::endl);
+ }
+
+ return score;
+}
+
+void LowestNPairListRestraint::show(std::ostream& out) const
+{
+ out << "Best N Pair list restraint with score function ";
+ ss_->show(out);
+ out << std::endl;
+}
+
+} // namespace IMP
Index: kernel/src/restraints/NonbondedRestraint.cpp
===================================================================
--- kernel/src/restraints/NonbondedRestraint.cpp (revision 649)
+++ kernel/src/restraints/NonbondedRestraint.cpp (working copy)
@@ -49,6 +49,12 @@
}
+ParticlesList NonbondedRestraint::get_interacting_particles() const {
+ IMP_failure("Not implemented yet because not needed", ErrorException);
+ return ParticlesList();
+}
+
+
void NonbondedRestraint::show(std::ostream& out) const
{
out << "Nonbonded restraint with score function ";
Index: kernel/src/restraints/SConscript
===================================================================
--- kernel/src/restraints/SConscript (revision 649)
+++ kernel/src/restraints/SConscript (working copy)
@@ -1,11 +1,18 @@
Import('env')
-
-files = ['ConnectivityRestraint.cpp', 'RestraintSet.cpp',
- 'DistanceRestraint.cpp', 'AngleRestraint.cpp', 'DihedralRestraint.cpp',
- 'NonbondedRestraint.cpp', 'BondDecoratorRestraint.cpp',
- 'SingletonListRestraint.cpp', 'PairListRestraint.cpp',
- 'TripletChainRestraint.cpp', 'PairChainRestraint.cpp',
- 'ConstantRestraint.cpp']
-
-files = [File(x) for x in files]
+files=[
+ File( 'AngleRestraint.cpp' ),
+ File( 'BondDecoratorRestraint.cpp' ),
+ File( 'ConnectivityRestraint.cpp' ),
+ File( 'ConstantRestraint.cpp' ),
+ File( 'DihedralRestraint.cpp' ),
+ File( 'DistanceRestraint.cpp' ),
+ File( 'LowestNPairListRestraint.cpp' ),
+ File( 'NonbondedRestraint.cpp' ),
+ File( 'PairChainRestraint.cpp' ),
+ File( 'PairListRestraint.cpp' ),
+ File( 'PairRestraint.cpp' ),
+ File( 'RestraintSet.cpp' ),
+ File( 'SingletonListRestraint.cpp' ),
+ File( 'TripletChainRestraint.cpp' ),
+ ]
Return('files')
Index: kernel/src/restraints/DihedralRestraint.cpp
===================================================================
--- kernel/src/restraints/DihedralRestraint.cpp (revision 649)
+++ 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
{
@@ -76,11 +78,10 @@
angle = -angle;
}
- Float score;
+ Float score, deriv;
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/restraints/ConnectivityRestraint.cpp
===================================================================
--- kernel/src/restraints/ConnectivityRestraint.cpp (revision 649)
+++ kernel/src/restraints/ConnectivityRestraint.cpp (working copy)
@@ -118,6 +118,8 @@
IMP_LOG(VERBOSE, "ConnectivityRestraint edge between "
<< get_particle(i)->get_index()
<< " and " << get_particle(j)->get_index() << std::endl);
+ IMP_LOG(VERBOSE, *get_particle(i) << std::endl
+ << *get_particle(j) << std::endl);
sum+= ps_->evaluate(get_particle(i),
get_particle(j),
accum);
Index: kernel/src/restraints/BondDecoratorRestraint.cpp
===================================================================
--- kernel/src/restraints/BondDecoratorRestraint.cpp (revision 649)
+++ kernel/src/restraints/BondDecoratorRestraint.cpp (working copy)
@@ -20,14 +20,15 @@
BondDecoratorRestraint
::BondDecoratorRestraint(UnaryFunction *f,
- BondDecoratorListScoreState *s): bl_(s),
- f_(f){}
+ ManualBondDecoratorListScoreState *s): bl_(s),
+ f_(f){}
Float BondDecoratorRestraint::evaluate(DerivativeAccumulator *accum)
{
Float sum=0;
- for (BondDecoratorListScoreState::BondIterator bi= bl_->bonds_begin();
- bi != bl_->bonds_end(); ++bi) {
+ for (BondDecoratorListScoreState::BondDecoratorIterator bi
+ = bl_->bond_decorators_begin();
+ bi != bl_->bond_decorators_end(); ++bi) {
BondDecorator bd= *bi;
Float l= bd.get_length();
Float s= bd.get_stiffness();
@@ -60,11 +61,29 @@
return sum;
}
+
+ParticlesList BondDecoratorRestraint::get_interacting_particles() const {
+ ParticlesList ret;
+ ret.reserve(std::distance(bl_->bond_decorators_begin(),
+ bl_->bond_decorators_end()));
+
+ for (BondDecoratorListScoreState::BondDecoratorConstIterator bi
+ = bl_->bond_decorators_begin();
+ bi != bl_->bond_decorators_end(); ++bi) {
+ ret.push_back(Particles(3));
+ ret.back()[0]= bi->get_particle();
+ ret.back()[1]= bi->get_bonded(0).get_particle();
+ ret.back()[2]= bi->get_bonded(1).get_particle();
+ }
+ return ret;
+}
+
+
void BondDecoratorRestraint::show(std::ostream& out) const
{
out << "Bond decorator restraint with unary function ";
f_->show(out);
- out << " on " << bl_->get_number_of_bonds() << " bonds";
+ out << " on " << bl_->get_number_of_bond_decorators() << " bonds";
out << std::endl;
}
Index: kernel/src/restraints/TripletChainRestraint.cpp
===================================================================
--- kernel/src/restraints/TripletChainRestraint.cpp (revision 649)
+++ kernel/src/restraints/TripletChainRestraint.cpp (working copy)
@@ -18,54 +18,34 @@
namespace IMP
{
-TripletChainRestraint::TripletChainRestraint(TripletScore* ts)
+TripletChainRestraint::TripletChainRestraint(TripletScore* ts,
+ const Particles &ps)
{
+ set_particles(ps);
ts_ = ts;
- clear_chains();
}
-void TripletChainRestraint::add_chain(const Particles &ps)
-{
- if (ps.size() <3) {
- IMP_WARN("Adding a chain of length 2 or less to the AnglesRestraint"
- << " doesn't accomplish anything."<< std::endl);
- } else {
- Restraint::add_particles(ps);
- chain_splits_.back()= Restraint::get_number_of_particles();
- chain_splits_.push_back(Restraint::get_number_of_particles());
- }
-}
-
Float TripletChainRestraint::evaluate(DerivativeAccumulator *accum)
{
- int cur_break=0;
- unsigned int i=2;
float score=0;
- while (i < Restraint::get_number_of_particles()) {
- /*IMP_LOG(VERBOSE, "Chain eval on "
- << Restraint::get_particle(i-2)->get_index()
- << Restraint::get_particle(i-1)->get_index()
- << Restraint::get_particle(i)->get_index()
- << " split is " << chain_splits_[cur_break]
- << std::endl);*/
+ for (unsigned int i=2; i< get_number_of_particles(); ++i) {
score += ts_->evaluate(Restraint::get_particle(i-2),
Restraint::get_particle(i-1),
Restraint::get_particle(i),
accum);
- if (chain_splits_[cur_break] == i) {
- i+=3;
- ++cur_break;
- } else {
- ++i;
- }
}
return score;
}
-void TripletChainRestraint::clear_chains() {
- Restraint::clear_particles();
- chain_splits_.clear();
- chain_splits_.push_back(0);
+ParticlesList TripletChainRestraint::get_interacting_particles() const {
+ ParticlesList ret(get_number_of_particles()-2);
+ for (unsigned int i=2; i< get_number_of_particles(); ++i) {
+ ret[i-2]= Particles(2);
+ ret[i-2][0]= get_particle(i);
+ ret[i-2][1]= get_particle(i-1);
+ ret[i-2][2]= get_particle(i-2);
+ }
+ return ret;
}
@@ -78,7 +58,6 @@
}
get_version_info().show(out);
- out << " " << chain_splits_.size()-1 << " chains" << std::endl;
ts_->show(out);
out << std::endl;
}
Index: kernel/src/restraints/PairListRestraint.cpp
===================================================================
--- kernel/src/restraints/PairListRestraint.cpp (revision 649)
+++ kernel/src/restraints/PairListRestraint.cpp (working copy)
@@ -59,7 +59,18 @@
}
}
+ParticlesList PairListRestraint::get_interacting_particles() const {
+ ParticlesList ret(get_number_of_particles()/2);
+ for (unsigned int i=0; 2*i< get_number_of_particles(); ++i) {
+ ret[i]= Particles(2);
+ ret[i][0]= get_particle(i*2);
+ ret[i][1]= get_particle(i*2+1);
+ }
+ return ret;
+}
+
+
void PairListRestraint::show(std::ostream& out) const
{
out << "Pair list restraint with score function ";
@@ -67,4 +78,5 @@
out << std::endl;
}
+
} // namespace IMP
Index: kernel/src/restraints/DistanceRestraint.cpp
===================================================================
--- kernel/src/restraints/DistanceRestraint.cpp (revision 649)
+++ kernel/src/restraints/DistanceRestraint.cpp (working copy)
@@ -17,23 +17,10 @@
DistanceRestraint::DistanceRestraint(UnaryFunction* score_func,
Particle* p1, Particle* p2) :
- dp_(score_func)
+ PairRestraint(new DistancePairScore(score_func), p1, p2)
{
- add_particle(p1);
- add_particle(p2);
}
-//! Calculate the score for this distance restraint.
-/** \param[in] accum If not NULL, use this object to accumulate partial first
- derivatives.
- \return Current score.
- */
-Float DistanceRestraint::evaluate(DerivativeAccumulator *accum)
-{
- return dp_.evaluate(get_particle(0), get_particle(1), accum);
-}
-
-
//! Show the current restraint.
/** \param[in] out Stream to send restraint description to.
*/
@@ -49,7 +36,7 @@
out << " particles: " << get_particle(0)->get_index();
out << " and " << get_particle(1)->get_index();
out << " ";
- dp_.show(out);
+ dp_->show(out);
out << std::endl;
}
Index: kernel/src/restraints/PairRestraint.cpp
===================================================================
--- kernel/src/restraints/PairRestraint.cpp (revision 0)
+++ kernel/src/restraints/PairRestraint.cpp (revision 0)
@@ -0,0 +1,51 @@
+/**
+ * \file PairRestraint.cpp \brief Pair restraint between two particles.
+ *
+ * Copyright 2007-8 Sali Lab. All rights reserved.
+ *
+ */
+
+#include "IMP/Particle.h"
+#include "IMP/Model.h"
+#include "IMP/log.h"
+#include "IMP/restraints/PairRestraint.h"
+#include "IMP/decorators/XYZDecorator.h"
+#include "IMP/PairScore.h"
+
+namespace IMP
+{
+
+PairRestraint::PairRestraint(PairScore* score_func,
+ Particle* p1, Particle* p2) :
+ dp_(score_func)
+{
+ add_particle(p1);
+ add_particle(p2);
+}
+
+Float PairRestraint::evaluate(DerivativeAccumulator *accum)
+{
+ return dp_->evaluate(get_particle(0), get_particle(1), accum);
+}
+
+
+//! Show the current restraint.
+/** \param[in] out Stream to send restraint description to.
+ */
+void PairRestraint::show(std::ostream& out) const
+{
+ if (get_is_active()) {
+ out << "pair restraint (active):" << std::endl;
+ } else {
+ out << "pair restraint (inactive):" << std::endl;
+ }
+
+ get_version_info().show(out);
+ out << " particles: " << get_particle(0)->get_index();
+ out << " and " << get_particle(1)->get_index();
+ out << " ";
+ dp_->show(out);
+ out << std::endl;
+}
+
+} // namespace IMP
Index: kernel/src/restraints/RestraintSet.cpp
===================================================================
--- kernel/src/restraints/RestraintSet.cpp (revision 649)
+++ kernel/src/restraints/RestraintSet.cpp (working copy)
@@ -65,7 +65,17 @@
return score * weight_;
}
+ParticlesList RestraintSet::get_interacting_particles() const {
+ ParticlesList ret;
+ for( RestraintConstIterator it= restraints_begin();
+ it != restraints_end(); ++it) {
+ ParticlesList c= (*it)->get_interacting_particles();
+ ret.insert(ret.end(), c.begin(), c.end());
+ }
+ return ret;
+}
+
//! Show the current restraint.
/** \param[in] out Stream to send restraint description to.
*/
Index: kernel/src/restraints/PairChainRestraint.cpp
===================================================================
--- kernel/src/restraints/PairChainRestraint.cpp (revision 649)
+++ kernel/src/restraints/PairChainRestraint.cpp (working copy)
@@ -21,53 +21,31 @@
PairChainRestraint::PairChainRestraint(PairScore* ts)
{
ts_ = ts;
- clear_chains();
}
-void PairChainRestraint::add_chain(const Particles &ps)
-{
- if (ps.size() <= 1) {
- IMP_WARN("Adding a chain of length 1 or less to the PairChainsRestraint"
- << " doesn't accomplish anything."<< std::endl);
- } else {
- Restraint::add_particles(ps);
- chain_splits_.back()= Restraint::get_number_of_particles();
- chain_splits_.push_back(Restraint::get_number_of_particles());
- }
-}
Float PairChainRestraint::evaluate(DerivativeAccumulator *accum)
{
- int cur_break=0;
- unsigned int i=1;
float score=0;
- while (i < Restraint::get_number_of_particles()) {
- /*IMP_LOG(VERBOSE, "Chain eval on "
- << Restraint::get_particle(i-2)->get_index()
- << Restraint::get_particle(i-1)->get_index()
- << Restraint::get_particle(i)->get_index()
- << " split is " << chain_splits_[cur_break]
- << std::endl);*/
+ for (unsigned int i=1; i< Restraint::get_number_of_particles(); ++i) {
score += ts_->evaluate(Restraint::get_particle(i-1),
Restraint::get_particle(i),
accum);
- if (chain_splits_[cur_break] == i) {
- i+=2;
- ++cur_break;
- } else {
- ++i;
- }
}
return score;
}
-void PairChainRestraint::clear_chains() {
- Restraint::clear_particles();
- chain_splits_.clear();
- chain_splits_.push_back(0);
+
+ParticlesList PairChainRestraint::get_interacting_particles() const {
+ ParticlesList ret(get_number_of_particles()-1);
+ for (unsigned int i=1; i< get_number_of_particles(); ++i) {
+ ret[i-1]= Particles(2);
+ ret[i-1][0]= get_particle(i);
+ ret[i-1][1]= get_particle(i-1);
+ }
+ return ret;
}
-
void PairChainRestraint::show(std::ostream& out) const
{
if (get_is_active()) {
@@ -77,7 +55,6 @@
}
get_version_info().show(out);
- out << " " << chain_splits_.size()-1 << " chains" << std::endl;
ts_->show(out);
out << std::endl;
}
Index: kernel/src/restraints/SingletonListRestraint.cpp
===================================================================
--- kernel/src/restraints/SingletonListRestraint.cpp (revision 649)
+++ kernel/src/restraints/SingletonListRestraint.cpp (working copy)
@@ -37,6 +37,13 @@
return score;
}
+ParticlesList SingletonListRestraint::get_interacting_particles() const {
+ ParticlesList ret(get_number_of_particles());
+ for( unsigned int i=0; i< get_number_of_particles(); ++i) {
+ ret[i]= Particles(1, get_particle(i));
+ }
+ return ret;
+}
void SingletonListRestraint::show(std::ostream& out) const
{
Index: kernel/src/particle_refiners/BondCoverParticleRefiner.cpp
===================================================================
--- kernel/src/particle_refiners/BondCoverParticleRefiner.cpp (revision 649)
+++ kernel/src/particle_refiners/BondCoverParticleRefiner.cpp (working copy)
@@ -16,12 +16,11 @@
{
BondCoverParticleRefiner::BondCoverParticleRefiner(FloatKey rk,
- FloatKey vk): rk_(rk),
- vk_(vk)
+ FloatKey vk,
+ IntKey tk): rk_(rk),
+ vk_(vk),
+ tk_(tk)
{
- if (0) {
- BondCoverParticleRefiner t(rk, vk);
- }
}
@@ -98,6 +97,10 @@
XYZDecorator d= XYZDecorator::create(np);
d.set_coordinates(vb+ (1+2*i)*r* ud);
np->add_attribute(rk_, r, false);
+ if (tk_ != IntKey() && p->has_attribute(tk_)) {
+ np->add_attribute(tk_,
+ p->get_value(tk_));
+ }
ret.push_back(np);
}
return ret;
Index: kernel/src/particle_refiners/SConscript
===================================================================
--- kernel/src/particle_refiners/SConscript (revision 649)
+++ kernel/src/particle_refiners/SConscript (working copy)
@@ -1,4 +1,6 @@
-files=['BondCoverParticleRefiner.cpp', 'ChildrenParticleRefiner.cpp']
-
-files = [File(x) for x in files]
+Import('env')
+files=[
+ File( 'BondCoverParticleRefiner.cpp' ),
+ File( 'ChildrenParticleRefiner.cpp' ),
+ ]
Return('files')
Index: kernel/src/triplet_scores/AngleTripletScore.cpp
===================================================================
--- kernel/src/triplet_scores/AngleTripletScore.cpp (revision 649)
+++ kernel/src/triplet_scores/AngleTripletScore.cpp (working copy)
@@ -8,6 +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 +47,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/src/triplet_scores/SConscript
===================================================================
--- kernel/src/triplet_scores/SConscript (revision 649)
+++ kernel/src/triplet_scores/SConscript (working copy)
@@ -1,6 +1,5 @@
Import('env')
-
-files = ['AngleTripletScore.cpp']
-
-files = [File(x) for x in files]
+files=[
+ File( 'AngleTripletScore.cpp' ),
+ ]
Return('files')
Index: kernel/src/Log.cpp
===================================================================
--- kernel/src/Log.cpp (revision 649)
+++ kernel/src/Log.cpp (working copy)
@@ -14,4 +14,13 @@
/* Initialize singleton pointer to NULL */
Log* Log::logpt_ = NULL;
+
+bool get_print_exception_messages() {
+ return internal::print_exceptions;
+}
+
+void set_print_exception_messages(bool tf) {
+ internal::print_exceptions=tf;
+}
+
} // namespace IMP
Index: kernel/src/OptimizerState.cpp
===================================================================
--- kernel/src/OptimizerState.cpp (revision 649)
+++ kernel/src/OptimizerState.cpp (working copy)
@@ -15,6 +15,7 @@
//! Constructor
OptimizerState::OptimizerState()
{
+ optimizer_=NULL;
IMP_LOG(VERBOSE, "OptimizerState constructed " << std::endl);
}
Index: kernel/src/decorators/XYZRDecorator.cpp
===================================================================
--- kernel/src/decorators/XYZRDecorator.cpp (revision 0)
+++ kernel/src/decorators/XYZRDecorator.cpp (revision 0)
@@ -0,0 +1,102 @@
+/**
+ * \file XYZRDecorator.cpp \brief Simple xyzr decorator.
+ *
+ * Copyright 2007-8 Sali Lab. All rights reserved.
+ *
+ */
+
+#include "IMP/decorators/XYZRDecorator.h"
+
+//#undef IMP_USE_CGAL
+
+#ifdef IMP_USE_CGAL
+#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
+#include <CGAL/Min_sphere_of_spheres_d.h>
+
+#endif
+
+namespace IMP
+{
+
+// These aren't statically initialized, as that way they may be initialized
+// before the table that caches them
+FloatKey XYZRDecorator::radius_;
+
+void XYZRDecorator::show(std::ostream &out, std::string prefix) const
+{
+ out << prefix << "(" << get_x()<< ", "
+ << get_y() << ", " << get_z() << ": " << get_radius() << ")";
+
+}
+
+
+
+IMP_DECORATOR_INITIALIZE(XYZRDecorator, XYZDecorator,
+ {
+ radius_= FloatKey("radius");
+ })
+
+
+Float distance(XYZRDecorator a, XYZRDecorator b)
+{
+ float d= distance(static_cast<XYZDecorator>(a),
+ static_cast<XYZDecorator>(b));
+ return d - a.get_radius() - b.get_radius();
+}
+
+
+void set_enclosing_sphere(const Particles &v,
+ XYZRDecorator out)
+{
+ IMP_check(!v.empty(), "Must pass some particles to have a bounding sphere",
+ ValueException);
+ for (unsigned int i=0; i< v.size(); ++i) {
+ XYZDecorator::cast(v[i]);
+ }
+ FloatKey rk= XYZRDecorator::get_radius_key();
+
+
+#ifdef IMP_USE_CGAL
+ typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
+ typedef CGAL::Min_sphere_of_spheres_d_traits_3<K, K::FT> Traits;
+ typedef CGAL::Min_sphere_of_spheres_d<Traits> Min_sphere;
+ typedef K::Point_3 Point;
+ typedef Traits::Sphere Sphere;
+
+ std::vector<Sphere> spheres;
+ for (unsigned int i=0; i< v.size(); ++i) {
+ // need cast to resolve ambiguity
+ XYZDecorator d(v[i]);
+ float r=0;
+ if (v[i]->has_attribute(rk)) r= v[i]->get_value(rk);
+
+ spheres.push_back(Sphere(Point(d.get_x(),
+ d.get_y(),
+ d.get_z()),
+ square(r)));
+ }
+ Min_sphere ms(spheres.begin(), spheres.end());
+
+ out.set_radius(ms.radius());
+ out.set_x(*ms.center_cartesian_begin());
+ out.set_y(*(ms.center_cartesian_begin()+1));
+ out.set_z(*(ms.center_cartesian_begin()+2));
+#else
+ Vector3D c(0,0,0);
+ for (unsigned int i=0; i< v.size(); ++i) {
+ XYZDecorator d(v[i]);
+ c+= d.get_vector();
+ }
+ c/=v.size();
+ out.set_coordinates(c);
+ float r=0;
+ for (unsigned int i=0; i< v.size(); ++i) {
+ float d= distance(XYZDecorator(v[i]), out);
+ if (v[i]->has_attribute(rk)) d+= v[i]->get_value(rk);
+ r= std::max(r, d);
+ }
+ out.set_radius(r);
+#endif
+}
+
+} // namespace IMP
Index: kernel/src/decorators/HierarchyDecorator.cpp
===================================================================
--- kernel/src/decorators/HierarchyDecorator.cpp (revision 649)
+++ kernel/src/decorators/HierarchyDecorator.cpp (working copy)
@@ -10,6 +10,7 @@
#include "IMP/decorators/NameDecorator.h"
#include <sstream>
+#include <memory>
namespace IMP
{
@@ -44,9 +45,9 @@
unsigned int count_hierarchy(HierarchyDecorator h)
{
- HierarchyCounter hc;
- depth_first_traversal(h,hc);
- return hc.get_count();
+ std::auto_ptr<HierarchyCounter> hc(new HierarchyCounter());
+ depth_first_traversal(h,hc.get());
+ return hc->get_count();
}
@@ -70,8 +71,8 @@
void HierarchyDecorator::validate() const
{
//std::cerr << "Checking hierarchy" << std::endl;
- internal::AssertHierarchy ah;
- depth_first_traversal(*this, ah);
+ std::auto_ptr<internal::AssertHierarchy> ah(new internal::AssertHierarchy());
+ depth_first_traversal(*this, ah.get());
}
@@ -95,7 +96,7 @@
})
-void breadth_first_traversal(HierarchyDecorator d, HierarchyVisitor &f)
+void breadth_first_traversal(HierarchyDecorator d, HierarchyVisitor *f)
{
std::deque<HierarchyDecorator> stack;
stack.push_back(d);
@@ -103,7 +104,7 @@
do {
HierarchyDecorator cur= stack.front();
stack.pop_front();
- if (f.visit(cur.get_particle())) {
+ if (f->visit(cur.get_particle())) {
//std::cerr << "Visiting particle " << cur.get_particle() << std::endl;
for (int i=cur.get_number_of_children()-1; i>=0; --i) {
stack.push_back(cur.get_child(i));
@@ -112,14 +113,14 @@
} while (!stack.empty());
}
-void depth_first_traversal(HierarchyDecorator d, HierarchyVisitor &f)
+void depth_first_traversal(HierarchyDecorator d, HierarchyVisitor *f)
{
std::vector<HierarchyDecorator> stack;
stack.push_back(d);
do {
HierarchyDecorator cur= stack.back();
stack.pop_back();
- if (f.visit(cur.get_particle())) {
+ if (f->visit(cur.get_particle())) {
for (int i=cur.get_number_of_children()-1; i>=0; --i) {
stack.push_back(cur.get_child(i));
}
@@ -144,27 +145,22 @@
Particles
-hierarchy_get_leaves(HierarchyDecorator mhd)
-{
- Particles out;
- hierarchy_gather(mhd, MHDMatchingLeaves(),
- std::back_inserter(out));
- return out;
+get_leaves(HierarchyDecorator mhd) {
+ return hierarchy_get(mhd, MHDMatchingLeaves());
}
-BondDecorators hierarchy_get_internal_bonds(HierarchyDecorator mhd)
-{
- Particles ps= hierarchy_get_all_descendants(mhd);
+BondDecorators get_internal_bonds(HierarchyDecorator mhd) {
+ Particles ps= get_all_descendants(mhd);
std::set<Particle*> sps(ps.begin(), ps.end());
BondDecorators ret;
for (unsigned int i=0; i< ps.size(); ++i) {
if (IMP::BondedDecorator::is_instance_of(ps[i])){
IMP::BondedDecorator b(ps[i]);
- for (unsigned int i=0; i< b.get_number_of_bonds(); ++i) {
- Particle *op= b.get_bonded(i).get_particle();
- if (op < ps[i]
+ for (unsigned int j=0; j< b.get_number_of_bonds(); ++j) {
+ Particle *op= b.get_bonded(j).get_particle();
+ if (op->get_index() < ps[i]->get_index()
&& sps.find(op) != sps.end()) {
- ret.push_back(b.get_bond(i));
+ ret.push_back(b.get_bond(j));
}
}
}
@@ -186,12 +182,8 @@
} // namespace
Particles
-hierarchy_get_all_descendants(HierarchyDecorator mhd)
-{
- Particles out;
- hierarchy_gather(mhd, MHDMatchingAll(),
- std::back_inserter(out));
- return out;
+get_all_descendants(HierarchyDecorator mhd) {
+ return hierarchy_get(mhd, MHDMatchingAll());
}
} // namespace IMP
Index: kernel/src/decorators/yaml.cpp
===================================================================
--- kernel/src/decorators/yaml.cpp (revision 0)
+++ kernel/src/decorators/yaml.cpp (revision 0)
@@ -0,0 +1,272 @@
+#include <IMP/decorators/yaml.h>
+
+namespace IMP
+{
+
+ namespace {
+ struct DefaultWrite {
+ template <class T>
+ std::ostream & operator()(const T t,
+ std::ostream &out) const {
+ return out << t;
+ }
+ };
+ struct ParticleWrite {
+ std::ostream & operator()(Particle *p,
+ std::ostream &out) const {
+ return out << p->get_index().get_index();
+ }
+ };
+ struct FloatWrite {
+ std::ostream & operator()(Float t,
+ std::ostream &out) const {
+ // set format better
+ return out << t;
+ }
+ };
+
+ template <class It, class Write>
+ void write_attributes(std::string indent,
+ Particle *p,
+ It b, It e, Write w,
+ std::ostream &out) {
+ for (It c= b; c != e; ++c) {
+ /** \todo should escape things properly */
+ out << indent << c->get_string() << ": ";
+ w(p->get_value(*c), out) << "\n";
+ }
+ }
+
+ // skips empty lines
+ struct LineStream{
+ std::istream ∈
+ std::vector<std::string> line;
+ bool has_indent(std::string str, unsigned int indent) {
+ for (unsigned int i=0; i< indent; ++i) {
+ if (i== str.size()) return false;
+ if (str[i] != ' ') return false;
+ }
+ return true;
+ }
+ bool not_white(char buf[]) const {
+ for (int i=0; buf[i] != '\0'; ++i) {
+ if (buf[i] != ' ') return true;
+ }
+ return false;
+ }
+
+ LineStream(std::istream &init): in(init){}
+ operator bool() const {return !line.empty() || static_cast<bool>(in);}
+ std::string get_line(unsigned int min_indent) {
+ while (line.empty()) {
+ char buf[3000];
+ in.getline(buf, 3000);
+ if (!in) return std::string();
+ if (buf[0] == '#') continue;
+ if (not_white(buf)) {
+ line.push_back(buf);
+ }
+ }
+ if (has_indent(line.back(), min_indent)) {
+ std::string ret(line.back(), min_indent);
+ line.pop_back();
+ return ret;
+ } else {
+ IMP_LOG(VERBOSE, "Line \"" << line.back() << "\" lacks "
+ << min_indent << " spaces" << std::endl);
+ return std::string();
+ }
+ }
+ void push_line(std::string s) {
+ if (s.empty()) return;
+ for (unsigned int i=0; i< s.size(); ++i) {
+ if (s[i] != ' ') {
+ line.push_back(s);
+ return;
+ }
+ }
+ }
+ };
+
+ template <class K, class V>
+ struct DefaultRead {
+ void operator()(Particle *p, std::string key, std::string value) const {
+ IMP_LOG(VERBOSE,
+ "Reading values from pair " << key << " "
+ << value << std::endl);
+ K k(key.c_str());
+ std::istringstream iss(value.c_str());
+ V v;
+ iss >> v;
+ IMP_check(iss, "Error reading value. Got " << v , ValueException);
+ p->set_value(k, v);
+ }
+ };
+
+ struct ParticleRead {
+ void operator()(Particle *p, std::string key, std::string value) const {
+ IMP_LOG(VERBOSE,
+ "Reading values from pair " << key << " "
+ << value << std::endl);
+ ParticleKey k(key.c_str());
+ std::istringstream iss(value.c_str());
+ int i;
+ iss >> i;
+ IMP_check(iss, "Error reading value" , ValueException);
+ Particle *op= p->get_model()->get_particle(ParticleIndex(i));
+ p->set_value(k, op);
+ }
+ };
+
+ int get_next_indent(LineStream &in) {
+ std::string buf= in.get_line(0);
+ if (buf.empty()) return 0;
+ unsigned int i=0;
+ for (; i < buf.size() && buf[i] == ' '; ++i) {
+ }
+ in.push_line(buf);
+ return i;
+ }
+
+
+ template <class Read>
+ void read_attributes(Particle *p, LineStream &in,
+ int indent,
+ Read read) {
+ IMP_LOG(VERBOSE, "Reading attributes " << indent << std::endl);
+ int nindent= get_next_indent(in);
+ if (nindent <= indent) return;
+ indent=nindent;
+ IMP_LOG(VERBOSE, "Required indent is " << indent<< std::endl);
+ do {
+ std::string buf = in.get_line(indent);
+ if (buf.empty()) {
+ IMP_LOG(VERBOSE, "Done reading attributes" << std::endl);
+ return;
+ }
+ IMP_check(buf[0] != ' ', "Extra white space on line "
+ << buf, InvalidStateException);
+ std::istringstream iss(buf.c_str());
+ char key[2000];
+ iss.get(key, 2000, ':');
+ IMP_check(iss, "no : found in line " << buf,
+ ValueException);
+ char colon;
+ iss >> colon;
+ IMP_check(colon == ':', "No colon found" << buf,
+ ValueException);
+
+ char value[2000];
+ iss.getline(value, 2000);
+ IMP_check(iss, "Error reading line " << buf,
+ ValueException);
+ read(p, key, value);
+
+ } while (true);
+ }
+
+ void read_yaml(Model *m,
+ LineStream &in,
+ unsigned int indent) {
+ std::string buf=in.get_line(indent);
+ if (buf.empty()) return;
+ //IMP_LOG(VERBOSE, "Got line " << buf << std::endl);
+ //IMP_check(in, "Error reading particle line from yaml", ValueException);
+ int id;
+ int nread=sscanf(buf.c_str(), "particle: %d", &id);
+ IMP_check(nread==1, "Couldn't read id", InvalidStateException);
+ Particle *p= m->get_particle(id);
+ IMP_LOG(VERBOSE, "Reading particle " << id << std::endl);
+ unsigned int nindent= get_next_indent(in);
+ if (nindent <= indent) return;
+ indent=nindent;
+ while (in) {
+ std::string buf=in.get_line(indent);
+ if (buf.empty()) break;
+ IMP_check(buf[0] != ' ', "Indent error" << buf, InvalidStateException);
+
+ IMP_LOG(VERBOSE, "Looking for attributes in line " << buf << std::endl);
+ std::istringstream iss(buf);
+ std::string type;
+ iss >> type;
+ if (type.compare("float-attributes:")==0) {
+ read_attributes(p, in, indent, DefaultRead<FloatKey, Float>());
+ } else if (type.compare("int-attributes:")==0) {
+ read_attributes(p, in, indent, DefaultRead<IntKey, Int>());
+ } else if (type.compare("string-attributes:")==0) {
+ read_attributes(p, in, indent, DefaultRead<StringKey, String>());
+ } else if (type.compare("particle-attributes:")==0) {
+ read_attributes(p, in, indent, ParticleRead());
+ } else {
+ break;
+ }
+ }
+ IMP_LOG(VERBOSE, "Done reading particle " << id << std::endl);
+ }
+ }
+
+ static std::string indent_level=" ";
+
+ void write_yaml(Particle *p,
+ std::ostream &out,
+ std::string indent) {
+ out << indent << "particle: " << p->get_index().get_index() << "\n";
+ out << indent << indent_level << "float-attributes:\n";
+ write_attributes(indent+indent_level+" ",
+ p,
+ p->float_keys_begin(),
+ p->float_keys_end(),
+ FloatWrite(),
+ out);
+ out << indent << indent_level << "int-attributes:\n";
+ write_attributes(indent+indent_level+" ",
+ p,
+ p->int_keys_begin(),
+ p->int_keys_end(),
+ DefaultWrite(),
+ out);
+ out << indent << indent_level << "string-attributes:\n";
+ write_attributes(indent+indent_level+" ",
+ p,
+ p->string_keys_begin(),
+ p->string_keys_end(),
+ DefaultWrite(),
+ out);
+ out << indent << indent_level << "particle-attributes:\n";
+ write_attributes(indent+indent_level+" ",
+ p,
+ p->particle_keys_begin(),
+ p->particle_keys_end(),
+ ParticleWrite(),
+ out);
+ }
+
+ void write_yaml(Model *m,
+ std::ostream &out,
+ std::string indent) {
+ for (Model::ParticleIterator pit= m->particles_begin();
+ pit != m->particles_end(); ++pit) {
+ write_yaml(*pit, out, indent);
+ }
+ }
+
+ void read_yaml(std::istream &in,
+ Model *m) {
+ LineStream r(in);
+ unsigned int nread=0;
+ do {
+ read_yaml(m, r, get_next_indent(r));
+ ++nread;
+ } while (r);
+ IMP_check(nread== m->get_number_of_particles(),
+ "Read wrong number of particles. Model is corrupted. Bye.",
+ ErrorException);
+ }
+
+ void read_yaml(std::string contents,
+ Model *m) {
+ std::istringstream iss(contents.c_str());
+ IMP_LOG(VERBOSE, "Got string:\n" << contents << std::endl);
+ read_yaml(iss, m);
+ }
+}
Index: kernel/src/decorators/AtomDecorator.cpp
===================================================================
--- kernel/src/decorators/AtomDecorator.cpp (revision 649)
+++ kernel/src/decorators/AtomDecorator.cpp (working copy)
@@ -5,13 +5,13 @@
*
*/
+#include "IMP/decorators/AtomDecorator.h"
+#include "IMP/log.h"
+
#include <sstream>
#include <vector>
#include <limits>
-#include "IMP/decorators/AtomDecorator.h"
-#include "IMP/log.h"
-
namespace IMP
{
Index: kernel/src/decorators/ResidueDecorator.cpp
===================================================================
--- kernel/src/decorators/ResidueDecorator.cpp (revision 649)
+++ kernel/src/decorators/ResidueDecorator.cpp (working copy)
@@ -12,6 +12,7 @@
#include <sstream>
#include <vector>
+#include <iterator>
namespace IMP
{
Index: kernel/src/decorators/MolecularHierarchyDecorator.cpp
===================================================================
--- kernel/src/decorators/MolecularHierarchyDecorator.cpp (revision 649)
+++ kernel/src/decorators/MolecularHierarchyDecorator.cpp (working copy)
@@ -54,34 +54,12 @@
})
-namespace
+Particles get_by_type(MolecularHierarchyDecorator mhd,
+ MolecularHierarchyDecorator::Type t)
{
-
-struct MHDMatchingType
-{
- MHDMatchingType(MolecularHierarchyDecorator::Type t): t_(t){}
-
- bool operator()(Particle *p) const {
- MolecularHierarchyDecorator mhd= MolecularHierarchyDecorator::cast(p);
- if (mhd== MolecularHierarchyDecorator()) {
- return false;
- } else {
- return mhd.get_type()==t_;
- }
- }
-
- MolecularHierarchyDecorator::Type t_;
-};
-
-} // namespace
-
-Particles molecular_hierarchy_get_by_type(MolecularHierarchyDecorator mhd,
- MolecularHierarchyDecorator::Type t)
-{
- Particles out;
- hierarchy_gather(mhd, MHDMatchingType(t),
- std::back_inserter(out));
- return out;
+ return hierarchy_get_by_attribute(mhd,
+ MolecularHierarchyDecorator::get_type_key(),
+ t);
}
@@ -108,7 +86,7 @@
ResidueDecorator
-molecular_hierarchy_get_residue(MolecularHierarchyDecorator mhd,
+get_residue(MolecularHierarchyDecorator mhd,
unsigned int index)
{
IMP_check(mhd.get_type() == MolecularHierarchyDecorator::PROTEIN
@@ -156,5 +134,4 @@
return fd;
}
-
} // namespace IMP
Index: kernel/src/decorators/SConscript
===================================================================
--- kernel/src/decorators/SConscript (revision 649)
+++ kernel/src/decorators/SConscript (working copy)
@@ -1,8 +1,13 @@
Import('env')
-
-files = ['AtomDecorator.cpp', 'bond_decorators.cpp',
- 'HierarchyDecorator.cpp', 'NameDecorator.cpp', 'ResidueDecorator.cpp',
- 'MolecularHierarchyDecorator.cpp', 'XYZDecorator.cpp']
-
-files = [File(x) for x in files]
+files=[
+ File( 'AtomDecorator.cpp' ),
+ File( 'bond_decorators.cpp' ),
+ File( 'HierarchyDecorator.cpp' ),
+ File( 'MolecularHierarchyDecorator.cpp' ),
+ File( 'NameDecorator.cpp' ),
+ File( 'ResidueDecorator.cpp' ),
+ File( 'XYZDecorator.cpp' ),
+ File( 'XYZRDecorator.cpp' ),
+ File( 'yaml.cpp' ),
+ ]
Return('files')
Index: kernel/src/decorators/bond_decorators.cpp
===================================================================
--- kernel/src/decorators/bond_decorators.cpp (revision 649)
+++ kernel/src/decorators/bond_decorators.cpp (working copy)
@@ -75,6 +75,11 @@
IMP_DECORATOR_INITIALIZE(BondedDecorator, DecoratorBase,
bond_initialize_static_data());
+ParticleKey BondDecorator::get_bonded_key(unsigned int i) {
+ decorator_initialize_static_data();
+ IMP_check(i < 2, "Invalid bond key requested", IndexException);
+ return internal::bond_graph_data_.get_node_key(i);
+}
BondDecorator bond(BondedDecorator a, BondedDecorator b, Int t)
{
Index: kernel/src/decorators/XYZDecorator.cpp
===================================================================
--- kernel/src/decorators/XYZDecorator.cpp (revision 649)
+++ 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/optimizers/ConjugateGradients.cpp
===================================================================
--- kernel/src/optimizers/ConjugateGradients.cpp (revision 649)
+++ kernel/src/optimizers/ConjugateGradients.cpp (working copy)
@@ -16,7 +16,7 @@
{
//! Estimate of limit of machine precision
-static const float eps = 1.2e-7;
+static const Float eps = 1.2e-7;
//! Get the score for a given model state.
/** \param[in] model The model to score.
@@ -26,14 +26,15 @@
\param[out] dscore First derivatives for current state.
\return The model score.
*/
-Float ConjugateGradients::get_score(std::vector<FloatIndex> float_indices,
- std::vector<Float> &x,
- std::vector<Float> &dscore)
+Float ConjugateGradients
+::get_score(const std::vector<FloatIndex>& float_indices,
+ std::vector<Float> &x,
+ std::vector<Float> &dscore)
{
int i, opt_var_cnt = float_indices.size();
/* set model state */
for (i = 0; i < opt_var_cnt; i++) {
- IMP_check(x[i] == x[i], "Got NaN in CG",
+ IMP_check(!is_nan(x[i]), "Got NaN in CG",
ValueException);
if (std::abs(x[i] - get_value(float_indices[i])) > max_change_) {
if (x[i] < get_value(float_indices[i])) {
@@ -76,16 +77,16 @@
*/
bool ConjugateGradients::line_search(std::vector<Float> &x,
std::vector<Float> &dx,
- float &alpha,
+ Float &alpha,
const std::vector<FloatIndex>
&float_indices,
- int &ifun, float &f,
- float &dg, float &dg1,
+ int &ifun, Float &f,
+ Float &dg, Float &dg1,
int max_steps,
const std::vector<Float> &search,
const std::vector<Float> &estimate)
{
- float ap, fp, dp, step, minf, u1, u2;
+ Float ap, fp, dp, step, minf, u1, u2;
int i, n, ncalls = ifun;
n = float_indices.size();
@@ -109,7 +110,6 @@
/* BEGIN THE LINEAR SEARCH ITERATION. */
while (true) {
- float dal, at;
/* TEST FOR FAILURE OF THE LINEAR SEARCH. */
if (alpha * step <= eps) {
@@ -130,7 +130,7 @@
}
/* COMPUTE THE DERIVATIVE OF F AT ALPHA. */
- dal = 0.0;
+ Float dal = 0.0;
for (i = 0; i < n; i++) {
dal += dx[i] * search[i];
}
@@ -167,8 +167,12 @@
u2 = 0.;
}
u2 = sqrt(u2);
- at = alpha - (alpha - ap) * (dal + u2 - u1) / (dal - dp + 2. * u2);
+ Float denom=(dal - dp + 2. * u2);
+ IMP_assert(denom != 0, "Error in ConjugateGradients. Divide by zero.");
+
+ Float at = alpha - (alpha - ap) * (dal + u2 - u1) / denom;
+
/* TEST WHETHER THE LINE MINIMUM HAS BEEN BRACKETED. */
if (dal / dp <= 0.) {
@@ -252,8 +256,8 @@
// Initialize optimization variables
int ifun = 0;
int nrst, nflag = 0;
- float dg1, xsq, dxsq, alpha, step, u1, u2, u3, u4;
- float f = 0., dg = 1., w1 = 0., w2 = 0., rtst, bestf;
+ Float dg1, xsq, dxsq, alpha, step, u1, u2, u3, u4;
+ Float f = 0., dg = 1., w1 = 0., w2 = 0., rtst, bestf;
bool gradient_direction;
// dx holds the gradient at x
@@ -293,7 +297,7 @@
dxsq = -dg1;
/* Test if the initial point is the minimizer. */
- if (dxsq <= eps * eps * std::max(1.0f, xsq)) {
+ if (dxsq <= eps * eps * std::max(Float(1.0), xsq)) {
goto end;
}
@@ -354,6 +358,8 @@
for (i = 0; i < n; i++) {
rtst += dx[i] * destimate[i];
}
+ IMP_check(dxsq != 0, "No gradient found in conjugate gradients",
+ ValueException);
if (fabs(rtst / dxsq) > 0.2) {
nrst = n;
}
@@ -371,6 +377,9 @@
}
}
+ IMP_check(w1 != 0 && w2 != 0, "Error in conjugate gradients. w1 or w2 is 0.",
+ ValueException);
+
/* CALCULATE THE RESTART HESSIAN TIMES THE CURRENT GRADIENT. */
u1 = u2 = 0.0;
for (i = 0; i < n; i++) {
Index: kernel/src/optimizers/states/SConscript
===================================================================
--- kernel/src/optimizers/states/SConscript (revision 649)
+++ kernel/src/optimizers/states/SConscript (working copy)
@@ -1,7 +1,7 @@
Import('env')
-
-files = ['VRMLLogOptimizerState.cpp', 'CMMLogOptimizerState.cpp',
- 'VelocityScalingOptimizerState.cpp']
-
-files = [File(x) for x in files]
+files=[
+ File( 'CMMLogOptimizerState.cpp' ),
+ File( 'VelocityScalingOptimizerState.cpp' ),
+ File( 'VRMLLogOptimizerState.cpp' ),
+ ]
Return('files')
Index: kernel/src/optimizers/states/VRMLLogOptimizerState.cpp
===================================================================
--- kernel/src/optimizers/states/VRMLLogOptimizerState.cpp (revision 649)
+++ kernel/src/optimizers/states/VRMLLogOptimizerState.cpp (working copy)
@@ -10,15 +10,43 @@
#include <sstream>
#include "IMP/optimizers/states/VRMLLogOptimizerState.h"
-#include "IMP/decorators/XYZDecorator.h"
+#include "IMP/decorators/XYZRDecorator.h"
namespace IMP
{
+//! A class for managing VRML IO
+class VRMLFile {
+ std::ostream &out_;
+ IntKey color_;
+ FloatKey radius_;
+ typedef std::map<int, Vector3D > Colors;
+ Colors colors_;
+
+ void write_color(std::string prefix, Particle *p) const;
+ void write_colors(std::string prefix) const;
+ void write_ball(std::string prefix, Particle *p) const;
+ public:
+ VRMLFile(std::ostream &out);
+
+ void add_balls(const Particles &ps) const;
+ void add_derivatives(const Particles &ps) const;
+
+ void set_color_key(IntKey ik) { color_=ik;}
+ void set_radius_key(FloatKey ik) { radius_=ik;}
+ template <class It>
+ void add_colors(It b, It e) {
+ colors_.insert(b,e);
+ }
+};
+
+
+
+
VRMLLogOptimizerState::VRMLLogOptimizerState(std::string filename,
const Particles &pis) :
filename_(filename), file_number_(0), call_number_(0),
- skip_steps_(0)
+ skip_steps_(0), radius_(XYZRDecorator::get_radius_key())
{
set_particles(pis);
}
@@ -47,14 +75,30 @@
} else {
IMP_LOG(VERBOSE, "Writing " << get_number_of_particles()
<< " particles to file " << buf << "..." << std::flush);
- write(out, get_particles());
+ VRMLFile file(out);
+ file.set_color_key(color_);
+ file.set_radius_key(radius_);
+ file.add_colors(colors_.begin(), colors_.end());
+ Particles ps;
+ std::vector<ParticleRefinerGuard > pgs(get_number_of_particles());
+ if (pr_) {
+ for (unsigned int i=0; i< get_number_of_particles(); ++i) {
+ pgs[i].refine(pr_, get_particle(i));
+ ps.insert(ps.end(),
+ pgs[i].get_particles().begin(),
+ pgs[i].get_particles().end());
+ }
+ } else {
+ ps= get_particles();
+ }
+ file.add_balls(ps);
+ file.add_derivatives(get_derivatives());
//IMP_LOG(TERSE, "done" << std::endl);
}
}
IMP_LIST_IMPL(VRMLLogOptimizerState, Particle, particle, Particle*, ,);
-IMP_CONTAINER_IMPL(VRMLLogOptimizerState, ParticleRefiner, particle_refiner,
- ParticleRefinerIndex ,,,);
+IMP_LIST_IMPL(VRMLLogOptimizerState, Derivative, derivative, Particle*, ,);
static Float snap(Float f)
{
@@ -69,82 +113,147 @@
snap(v[2]));
}
-void VRMLLogOptimizerState::write(std::ostream &out, const Particles &ps) const
+
+
+void VRMLLogOptimizerState::show(std::ostream &out) const
{
- out << "#VRML V2.0 utf8\n";
- out << "Group {\n";
- out << "children [\n";
+ out << "Writing VRML files " << filename_ << std::endl;
+}
- for (Particles::const_iterator it = ps.begin(); it != ps.end(); ++it) {
- Particle *p = *it;
- bool wasrefined=false;
- for (ParticleRefinerConstIterator prit= particle_refiners_begin();
- prit != particle_refiners_end(); ++prit) {
- if ((*prit)->get_can_refine(p)) {
- Particles refined= (*prit)->get_refined(p);
- write(out, refined);
- (*prit)->cleanup_refined(p, refined, NULL);
- wasrefined=true;
- break;
- }
+
+
+VRMLFile::VRMLFile(std::ostream &out): out_(out){
+ out_ << "#VRML V2.0 utf8\n";
+}
+
+
+
+void VRMLFile::write_color(std::string prefix,
+ Particle *p) const {
+ if (color_ != IntKey()
+ && p->has_attribute(color_)) {
+ int cv = p->get_value(color_);
+ if (colors_.find(cv) != colors_.end()) {
+ Float rv= colors_.find(cv)->second[0];
+ Float gv= colors_.find(cv)->second[1];
+ Float bv= colors_.find(cv)->second[2];
+ out_ << prefix << "appearance Appearance {\n";
+ out_ << prefix << " material Material {\n";
+ out_ << prefix << " diffuseColor " << rv << " " << gv
+ << " " << bv << "\n";
+ out_ << prefix << " }\n";
+ out_ << prefix << "}\n";
}
- if (wasrefined) continue;
- try {
- XYZDecorator xyz = XYZDecorator::cast(p);
- float x = xyz.get_x();
- float y = xyz.get_y();
- float z = xyz.get_z();
- Float rv = -1, gv = -1, bv = -1;
- if (color_ != IntKey()
- && p->has_attribute(color_)) {
- int cv = p->get_value(color_);
- if (colors_.find(cv) != colors_.end()) {
- rv= colors_.find(cv)->second[0];
- gv= colors_.find(cv)->second[1];
- bv= colors_.find(cv)->second[2];
- }
- }
- Float radius = .1;
- if (radius_ != FloatKey() && p->has_attribute(radius_)) {
- radius = p->get_value(radius_);
- //oss << ".sphere " << x << " " << y << " " << z << " " << r << "\n";
- }
+ }
+}
- out << "Transform {\n";
- out << " translation " << x << " " << y << " " << z << std::endl;
- out << " children [\n";
- out << " Shape {\n";
- if (rv >= 0 && gv >= 0 && bv >= 0) {
- /* appearance Appearance {
- material Material {
- diffuseColor 1 0 0
- }
- }
- */
- out << " appearance Appearance {\n";
- out << " material Material {\n";
- out << " diffuseColor " << rv << " " << gv;
- out << " " << bv << "\n";
- out << " }\n";
- out << " }\n";
- }
- out << " geometry Sphere { radius " << radius << "}\n";
- out << " }\n";
- out << " ]\n";
- out << "}\n";
- } catch (InvalidStateException &e) {
- IMP_WARN("Particle " << p << " does not have "
- << " cartesian coordinates");
+void VRMLFile::write_colors(std::string prefix) const {
+ int max=0;
+ for (std::map<int, Vector3D >::const_iterator it= colors_.begin();
+ it != colors_.end(); ++it) {
+ max= std::max(max, it->first);
+ }
+ out_ << prefix << "color Color {\n";
+ out_ << prefix << " color [\n";
+ for (int i=0; i<= max; ++i) {
+ Colors::const_iterator cit= colors_.find(i);
+ if (cit != colors_.end()) {
+ out_ << prefix << " " << cit->second[0] << " "
+ << cit->second[1] << " "
+ << cit->second[2] << ",\n";
+ } else {
+ out_ << prefix << " .5 .5 .5,\n";
}
}
- out << "]\n";
- out << "}\n";
+ out_ << prefix << " ]\n";
+ out_ << prefix << "}\n";
}
-void VRMLLogOptimizerState::show(std::ostream &out) const
+void VRMLFile::write_ball(std::string prefix, Particle *p) const
{
- out << "Writing VRML files " << filename_ << std::endl;
+ try {
+ XYZDecorator xyz = XYZDecorator::cast(p);
+ float x = xyz.get_x();
+ float y = xyz.get_y();
+ float z = xyz.get_z();
+
+ Float radius = .1;
+ if (radius_ != FloatKey() && p->has_attribute(radius_)) {
+ radius = p->get_value(radius_);
+ //oss << ".sphere " << x << " " << y << " " << z << " " << r << "\n";
+ }
+
+ out_ << "Transform {\n";
+ out_ << " translation " << x << " " << y << " " << z << std::endl;
+ out_ << " children [\n";
+ out_ << " Shape {\n";
+ write_color(" ", p);
+ out_ << " geometry Sphere { radius " << radius << "}\n";
+ out_ << " }\n";
+ out_ << " ]\n";
+ out_ << "}\n";
+
+ } catch (InvalidStateException &e) {
+ IMP_WARN("Particle " << p << " does not have "
+ << " cartesian coordinates");
+ }
}
+
+void VRMLFile::add_balls(const Particles &ps) const
+{
+ out_ << "Group {\n";
+ out_ << "children [\n";
+
+ for (Particles::const_iterator it = ps.begin(); it != ps.end(); ++it) {
+ Particle *p = *it;
+ write_ball(" ", p);
+ }
+ out_ << "]\n";
+ out_ << "}\n";
+}
+
+void VRMLFile::add_derivatives(const Particles &ps) const
+{
+ out_ << "Shape {\n";
+ out_ << "geometry IndexedLineSet {\n";
+ write_colors(" ");
+ out_ << " coord Coordinate { point[\n";
+ for (Particles::const_iterator it = ps.begin(); it != ps.end(); ++it) {
+ Particle *p = *it;
+ IMP_LOG(TERSE, "Writing derivatives for particle " << p->get_index()
+ << std::endl);
+ XYZDecorator d= XYZDecorator::cast(p);
+ Vector3D vc=d.get_vector();
+ Vector3D vd=d.get_derivative_vector();
+ Vector3D pd= vc+vd;
+ out_ << " " << vc[0] << ", " << vc[1] << ", " << vc[2]
+ << ", " << pd[0] << ", " << pd[1] << ", " << pd[2] << ",\n";
+ }
+ out_ << " ]}\n";
+ out_ << " coordIndex [\n";
+ unsigned int i=0;
+ for (Particles::const_iterator it = ps.begin(); it != ps.end(); ++it) {
+ out_ << " " << i << ", " << i+1 << ", -1,\n";
+ i+=2;
+ }
+ out_ << " ]\n";
+ if (color_ != IntKey()) {
+ out_ << " colorPerVertex FALSE\n";
+ out_ << " colorIndex [\n";
+ for (Particles::const_iterator it = ps.begin(); it != ps.end(); ++it) {
+ if ((*it)->has_attribute(color_)) {
+ int c= (*it)->get_value(color_);
+ out_ << " " << c << ",\n";
+ } else {
+ out_ << " " << 0 << ",\n";
+ }
+ }
+ out_ << " ]\n";
+ }
+ out_ << "}\n";
+ out_ << "}\n";
+}
+
} // namespace IMP
Index: kernel/src/optimizers/MonteCarlo.cpp
===================================================================
--- kernel/src/optimizers/MonteCarlo.cpp (revision 649)
+++ 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 "\
@@ -53,12 +53,13 @@
for (unsigned int i=0; i< max_steps; ++i) {
//make it a parameter
for (MoverIterator it = movers_begin(); it != movers_end(); ++it) {
- IMP_LOG(VERBOSE, "MC Trying move " << **it << std::endl);
+ IMP_LOG(VERBOSE, "MC Trying move "
+ << std::distance(movers_begin(), it) << std::endl);
IMP_CHECK_OBJECT(*it);
(*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());
@@ -69,7 +70,7 @@
}
bool accept= (next_energy < prior_energy_);
- if (!accept) {
+ if (!accept && temp_ != 0) {
Float diff= next_energy- prior_energy_;
Float e= std::exp(-diff/temp_);
Float r= rand(random_number_generator);
Index: kernel/src/optimizers/SConscript
===================================================================
--- kernel/src/optimizers/SConscript (revision 649)
+++ kernel/src/optimizers/SConscript (working copy)
@@ -1,11 +1,12 @@
Import('env')
-
-states_files = SConscript('states/SConscript')
-movers_files = SConscript('movers/SConscript')
-
-files = ['SteepestDescent.cpp', 'ConjugateGradients.cpp',
- 'MolecularDynamics.cpp', 'MonteCarlo.cpp', 'MoverBase.cpp',
- 'BrownianDynamics.cpp'] + states_files + movers_files
-
-files = [File(x) for x in files]
+files=[
+ File( 'BrownianDynamics.cpp' ),
+ File( 'ConjugateGradients.cpp' ),
+ File( 'MolecularDynamics.cpp' ),
+ File( 'MonteCarlo.cpp' ),
+ File( 'MoverBase.cpp' ),
+ File( 'SteepestDescent.cpp' ),
+ SConscript( 'movers/SConscript' ),
+ SConscript( 'states/SConscript' ),
+ ]
Return('files')
Index: kernel/src/optimizers/movers/BallMover.cpp
===================================================================
--- kernel/src/optimizers/movers/BallMover.cpp (revision 649)
+++ kernel/src/optimizers/movers/BallMover.cpp (working copy)
@@ -14,20 +14,21 @@
// These functions probably should be exposed at some point
-static void random_point_in_sphere(unsigned int D,
- Float radius,
- std::vector<Float> v)
+static std::vector<Float> random_point_in_sphere(unsigned int D,
+ Float radius)
{
IMP_assert(radius > 0, "No volume there");
::boost::uniform_real<> rand(-radius, radius);
Float norm;
+ std::vector<Float> ret(D);
do {
norm = 0;
for (unsigned int i = 0; i < D; ++i) {
- v[i] = rand(random_number_generator);
- norm += v[i] * v[i];
+ ret[i] = rand(random_number_generator);
+ norm += ret[i] * ret[i];
}
} while (norm > radius*radius);
+ return ret;
}
static std::vector<Float>
@@ -35,8 +36,7 @@
Float radius)
{
IMP_assert(radius > 0, "No volume there");
- std::vector<Float> v(center.size());
- random_point_in_sphere(center.size(), radius, v);
+ std::vector<Float> v= random_point_in_sphere(center.size(), radius);
std::vector<Float> r(center.size());
for (unsigned int i = 0; i < center.size(); ++i) {
r[i] = center[i] + v[i];
@@ -66,7 +66,20 @@
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/optimizers/movers/SConscript
===================================================================
--- kernel/src/optimizers/movers/SConscript (revision 649)
+++ kernel/src/optimizers/movers/SConscript (working copy)
@@ -1,6 +1,6 @@
Import('env')
-
-files = [ 'BallMover.cpp', 'NormalMover.cpp']
-
-files = [File(x) for x in files]
+files=[
+ File( 'BallMover.cpp' ),
+ File( 'NormalMover.cpp' ),
+ ]
Return('files')
Index: kernel/src/pair_scores/RefineHighPairScore.cpp
===================================================================
--- kernel/src/pair_scores/RefineHighPairScore.cpp (revision 0)
+++ kernel/src/pair_scores/RefineHighPairScore.cpp (revision 0)
@@ -0,0 +1,71 @@
+/**
+ * \file RefineHighPairScore.cpp
+ * \brief Refine particles at most high with a ParticleRefiner.
+ *
+ *
+ * Copyright 2007-8 Sali Lab. All rights reserved.
+ */
+
+#include "IMP/pair_scores/RefineHighPairScore.h"
+#include "IMP/decorators/bond_decorators.h"
+#include "IMP/decorators/XYZDecorator.h"
+#include "IMP/internal/constants.h"
+
+#include <cmath>
+
+namespace IMP
+{
+
+RefineHighPairScore::RefineHighPairScore(PairScore *f,
+ ParticleRefiner *pr,
+ Float threshold): f_(f),
+ pr_(pr),
+ threshold_(threshold){}
+
+
+Float RefineHighPairScore::evaluate(Particle *a, Particle *b,
+ DerivativeAccumulator *da) const
+{
+ Float score=f_->evaluate(a,b, NULL);
+ if (score > threshold_) {
+ ParticleRefinerGuard prg[2];
+ prg[0].refine(pr_, a);
+ prg[1].refine(pr_, b);
+ if (!prg[0].get_was_refined() && !prg[1].get_was_refined()) {
+ return da? f_->evaluate(a, b, da): score;
+ } else {
+ Float ret=0;
+ for (unsigned int i=0; i< prg[0].get_number_of_particles(); ++i) {
+ for (unsigned int j=0; j< prg[1].get_number_of_particles(); ++j) {
+ ret+=evaluate(prg[0].get_particle(i), prg[1].get_particle(j), da);
+ }
+ }
+ return ret;
+ }
+ } else {
+ IMP_IF_CHECK(EXPENSIVE) {
+ ParticleRefinerGuard prg[2];
+ prg[0].refine(pr_, a);
+ prg[1].refine(pr_, b);
+ for (unsigned int i=0; i< prg[0].get_number_of_particles(); ++i) {
+ for (unsigned int j=0; j< prg[1].get_number_of_particles(); ++j) {
+ IMP_assert(f_->evaluate(prg[0].get_particle(i),
+ prg[1].get_particle(j), NULL)
+ <= threshold_,
+ "Children are not under threshold even though parent is");
+ }
+ }
+ }
+
+ return 0;
+ }
+}
+
+void RefineHighPairScore::show(std::ostream &out) const
+{
+ out << "RefineHighPairScore using ";
+ f_->show(out);
+ out << threshold_ << std::endl;
+}
+
+} // namespace IMP
Index: kernel/src/pair_scores/DistancePairScore.cpp
===================================================================
--- kernel/src/pair_scores/DistancePairScore.cpp (revision 649)
+++ 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/RefineOncePairScore.cpp
===================================================================
--- kernel/src/pair_scores/RefineOncePairScore.cpp (revision 649)
+++ kernel/src/pair_scores/RefineOncePairScore.cpp (working copy)
@@ -23,33 +23,20 @@
DerivativeAccumulator *da) const
{
Particle* p[2]={a,b};
- Particles ps[2];
+ ParticleRefinerGuard ps[2];
for (unsigned int i=0; i< 2; ++i) {
- if (r_->get_can_refine(p[i])) {
- ps[i]= r_->get_refined(p[i]);
- } else {
- ps[i].push_back(p[i]);
- }
+ ps[i].refine(r_, p[i]);
IMP_LOG(VERBOSE, "Refining " << p[i]->get_index()
- << " resulted in " << ps[i].size() << " particles"
+ << " resulted in " << ps[i].get_particles().size() << " particles"
<< std::endl);
}
Float ret=0;
- for (unsigned int i=0; i< ps[0].size(); ++i) {
- for (unsigned int j=0; j< ps[1].size(); ++j) {
- ret+=f_->evaluate(ps[0][i], ps[1][j], da);
+ for (unsigned int i=0; i< ps[0].get_number_of_particles(); ++i) {
+ for (unsigned int j=0; j< ps[1].get_number_of_particles(); ++j) {
+ ret+=f_->evaluate(ps[0].get_particle(i), ps[1].get_particle(j), da);
}
}
-
- for (unsigned int i=0; i< 2; ++i) {
- if (ps[i].size() != 1 || (ps[i].size()==1 && ps[i][0] != p[i])) {
- IMP_LOG(VERBOSE, "Refining " << p[i]->get_index()
- << " resulted in " << ps[i].size() << " particles"
- << std::endl);
- r_->cleanup_refined(p[i], ps[i], da);
- }
- }
return ret;
}
Index: kernel/src/pair_scores/SConscript
===================================================================
--- kernel/src/pair_scores/SConscript (revision 649)
+++ kernel/src/pair_scores/SConscript (working copy)
@@ -1,8 +1,10 @@
Import('env')
-
-files = ['DistancePairScore.cpp', 'SphereDistancePairScore.cpp',
- 'RefineOncePairScore.cpp', 'TypedPairScore.cpp',
- 'TransformedDistancePairScore.cpp']
-
-files = [File(x) for x in files]
+files=[
+ File( 'DistancePairScore.cpp' ),
+ File( 'RefineHighPairScore.cpp' ),
+ File( 'RefineOncePairScore.cpp' ),
+ File( 'SphereDistancePairScore.cpp' ),
+ File( 'TransformedDistancePairScore.cpp' ),
+ File( 'TypedPairScore.cpp' ),
+ ]
Return('files')
Index: kernel/src/pair_scores/TransformedDistancePairScore.cpp
===================================================================
--- kernel/src/pair_scores/TransformedDistancePairScore.cpp (revision 649)
+++ kernel/src/pair_scores/TransformedDistancePairScore.cpp (working copy)
@@ -52,43 +52,36 @@
}
+
/*
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 +111,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,62 @@
+/**
+ * \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 <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/src/internal/SConscript
===================================================================
--- kernel/src/internal/SConscript (revision 649)
+++ kernel/src/internal/SConscript (working copy)
@@ -1,7 +1,9 @@
Import('env')
-
-files = ['graph_base.cpp', 'ParticleGrid.cpp',
- 'kernel_version_info.cpp', 'constants.cpp', 'bbox_nbl_helpers.cpp']
-
-files = [File(x) for x in files]
+files=[
+ File( 'bbox_nbl_helpers.cpp' ),
+ File( 'constants.cpp' ),
+ File( 'graph_base.cpp' ),
+ File( 'kernel_version_info.cpp' ),
+ File( 'ParticleGrid.cpp' ),
+ ]
Return('files')
Index: kernel/src/internal/bbox_nbl_helpers.cpp
===================================================================
--- kernel/src/internal/bbox_nbl_helpers.cpp (revision 649)
+++ kernel/src/internal/bbox_nbl_helpers.cpp (working copy)
@@ -12,7 +12,9 @@
/* compile the CGAL code with NDEBUG since it doesn't have the
same level of control over errors as IMP
*/
+#ifndef NDEBUG
#define NDEBUG
+#endif
#ifdef IMP_USE_CGAL
#include <CGAL/box_intersection_d.h>
#include <vector>
@@ -21,7 +23,7 @@
namespace IMP
{
-namespace internal
+namespace
{
struct NBLBbox
@@ -49,7 +51,7 @@
#ifdef IMP_USE_CGAL
static void copy_particles_to_boxes(const Particles &ps,
FloatKey rk, Float slack, Float cutoff,
- std::vector<internal::NBLBbox> &boxes)
+ std::vector<NBLBbox> &boxes)
{
boxes.resize(ps.size());
for (unsigned int i=0; i< ps.size(); ++i) {
@@ -59,17 +61,22 @@
if (rk != FloatKey() && p->has_attribute(rk)) {
r+= p->get_value(rk);
}
- boxes[i]=internal::NBLBbox(p, r);
+ boxes[i]=NBLBbox(p, r);
}
}
#endif
+} // anon
+
+namespace internal
+{
+
void bipartite_bbox_scan(const Particles &ps0, const Particles &ps1,
FloatKey rk, Float slack, Float cutoff,
const NBLAddPairIfNonbonded &ap)
{
#ifdef IMP_USE_CGAL
- std::vector<internal::NBLBbox> boxes0, boxes1;
+ std::vector<NBLBbox> boxes0, boxes1;
copy_particles_to_boxes(ps0, rk, slack, cutoff, boxes0);
copy_particles_to_boxes(ps1, rk, slack, cutoff, boxes1);
@@ -85,7 +92,7 @@
const NBLAddPairIfNonbonded &ap)
{
#ifdef IMP_USE_CGAL
- std::vector<internal::NBLBbox> boxes;
+ std::vector<NBLBbox> boxes;
copy_particles_to_boxes(ps, rk, slack, cutoff, boxes);
Index: kernel/pyext/IMP/test.py
===================================================================
--- kernel/pyext/IMP/test.py (revision 649)
+++ kernel/pyext/IMP/test.py (working copy)
@@ -3,6 +3,40 @@
import IMP
+def create_particles_in_box( model, num=10,
+ lb= [0,0,0],
+ ub= [10,10,10]):
+ """Create a bunch of particles in a box"""
+ lbv=IMP.Vector3D(lb[0],lb[1],lb[2])
+ ubv=IMP.Vector3D(ub[0],ub[1],ub[2])
+ ps= IMP.Particles()
+ for i in range(0,num):
+ p= IMP.Particle()
+ model.add_particle(p)
+ d= IMP.XYZDecorator.create(p)
+ d.set_coordinates(IMP.random_vector_in_box(lbv, ubv))
+ ps.append(p)
+ d.set_coordinates_are_optimized(True)
+ return ps
+
+def randomize_particles(particles, deviation):
+ """Randomize the xyz coordinates of a list of particles"""
+ for p in particles:
+ d= IMP.XYZDecorator.cast(p)
+ d.set_x(random.uniform(-deviation, deviation))
+ d.set_y(random.uniform(-deviation, deviation))
+ d.set_z(random.uniform(-deviation, deviation))
+
+def create_chain(ps, length=1, stiffness=1):
+ bds= []
+ IMP.BondedDecorator.create(ps[0])
+ for i in range(1,len(ps)):
+ ba= IMP.BondedDecorator.cast(ps[i-1])
+ bb= IMP.BondedDecorator.create(ps[i])
+ bds.append(IMP.custom_bond(ba, bb, length, stiffness))
+ return bds
+
+
class TestCase(unittest.TestCase):
"""Super class for IMP test cases"""
@@ -10,6 +44,8 @@
self.__check_level = IMP.get_check_level()
# Turn on expensive runtime checks while running the test suite:
IMP.set_check_level(IMP.EXPENSIVE)
+ # should go in library init code
+ IMP.set_print_exception_messages(False)
def tearDown(self):
# Restore original check level
@@ -33,14 +69,6 @@
p.add_attribute(IMP.FloatKey("z"), z, True)
return p
- def randomize_particles(self, particles, deviation):
- """Randomize the xyz coordinates of a list of particles"""
- for p in particles:
- d= IMP.XYZDecorator.cast(p)
- d.set_x(random.uniform(-deviation, deviation))
- d.set_y(random.uniform(-deviation, deviation))
- d.set_z(random.uniform(-deviation, deviation))
-
def particle_distance(self, p1, p2):
"""Return distance between two given particles"""
xkey = IMP.FloatKey("x")
@@ -73,30 +101,7 @@
fmin, vmin = f, v
self.assertInTolerance(fmin, expected_fmin, step)
- def create_particles_in_box(self, model, num=10,
- lb= [0,0,0],
- ub= [10,10,10]):
- """Create a bunch of particles in a box"""
- lbv=IMP.Vector3D(lb[0],lb[1],lb[2])
- ubv=IMP.Vector3D(ub[0],ub[1],ub[2])
- ps= IMP.Particles()
- for i in range(0,num):
- p= IMP.Particle()
- model.add_particle(p)
- d= IMP.XYZDecorator.create(p)
- d.randomize_in_box(lbv, ubv)
- ps.append(p)
- d.set_coordinates_are_optimized(True)
- return ps
- def create_chain(self, ps, length=1, stiffness=1):
- bds= []
- IMP.BondedDecorator.create(ps[0])
- for i in range(1,len(ps)):
- ba= IMP.BondedDecorator.cast(ps[i-1])
- bb= IMP.BondedDecorator.create(ps[i])
- bds.append(IMP.custom_bond(ba, bb, length, stiffness))
- return bds
class ConstPairScore(IMP.PairScore):
Index: kernel/pyext/IMP.i
===================================================================
--- kernel/pyext/IMP.i (revision 649)
+++ 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()):
@@ -26,6 +35,7 @@
raise IndexError("Particle does not have attribute")
%}
+
namespace IMP {
%pythonprepend Model::add_restraint %{
args[1].thisown=0
@@ -36,6 +46,9 @@
%pythonprepend Optimizer::add_optimizer_state %{
args[1].thisown=0
%}
+ %pythonprepend Optimizer::set_model %{
+ args[1].thisown=0
+ %}
%pythonprepend RestraintSet::add_restraint %{
args[1].thisown=0
%}
@@ -45,6 +58,9 @@
%pythonprepend DistanceRestraint::DistanceRestraint %{
args[0].thisown=0
%}
+ %pythonprepend DistanceRestraint::PairRestraint %{
+ args[0].thisown=0
+ %}
%pythonprepend AngleRestraint::AngleRestraint %{
args[0].thisown=0
%}
@@ -72,6 +88,9 @@
%pythonprepend PairChainRestraint::PairChainRestraint %{
args[0].thisown=0
%}
+ %pythonprepend PairRestraint::PairRestraint %{
+ args[0].thisown=0
+ %}
%pythonprepend ConnectivityRestraint::ConnectivityRestraint %{
args[0].thisown=0
%}
@@ -91,6 +110,10 @@
args[0].thisown=0
args[1].thisown=0
%}
+ %pythonprepend RefineHighPairScore::RefineHighPairScore %{
+ args[0].thisown=0
+ args[1].thisown=0
+ %}
%pythonprepend DistanceToSingletonScore::DistanceToSingletonScore %{
args[0].thisown=0
%}
@@ -152,7 +175,10 @@
%feature("ref") Particle "$this->ref();"
%feature("unref") Particle "$this->unref(); if (! $this->get_has_ref()) delete $this;"
+/*%feature("ref") Model "$this->ref();"
+%feature("unref") Model "$this->unref(); if (! $this->get_has_ref()) delete $this;"*/
+
/* Don't wrap internal functions */
%ignore IMP::internal::evaluate_distance_pair_score;
%ignore IMP::internal::check_particles_active;
@@ -194,6 +220,7 @@
%include "IMP/SingletonScore.h"
%include "IMP/TripletScore.h"
%include "IMP/Particle.h"
+%include "IMP/ParticleRefiner.h"
%include "Vector3D.i"
%include "IMP/DecoratorBase.h"
%include "IMP/decorators/bond_decorators.h"
@@ -202,8 +229,9 @@
%include "IMP/decorators/NameDecorator.h"
%include "IMP/decorators/ResidueDecorator.h"
%include "IMP/decorators/XYZDecorator.h"
+%include "IMP/decorators/XYZRDecorator.h"
%include "IMP/decorators/AtomDecorator.h"
-%include "IMP/ParticleRefiner.h"
+%include "IMP/decorators/yaml.h"
%include "IMP/particle_refiners/BondCoverParticleRefiner.h"
%include "IMP/particle_refiners/ChildrenParticleRefiner.h"
%include "IMP/Optimizer.h"
@@ -221,6 +249,7 @@
%include "IMP/optimizers/states/VelocityScalingOptimizerState.h"
%include "IMP/pair_scores/DistancePairScore.h"
%include "IMP/pair_scores/RefineOncePairScore.h"
+%include "IMP/pair_scores/RefineHighPairScore.h"
%include "IMP/pair_scores/SphereDistancePairScore.h"
%include "IMP/pair_scores/TypedPairScore.h"
%include "IMP/pair_scores/TransformedDistancePairScore.h"
@@ -234,11 +263,13 @@
%include "IMP/score_states/MaxChangeScoreState.h"
%include "IMP/score_states/NonbondedListScoreState.h"
%include "IMP/score_states/AllNonbondedListScoreState.h"
+%include "IMP/score_states/ManualBondDecoratorListScoreState.h"
%include "IMP/score_states/BondDecoratorListScoreState.h"
%include "IMP/score_states/BipartiteNonbondedListScoreState.h"
%include "IMP/score_states/GravityCenterScoreState.h"
%include "IMP/score_states/CoverBondsScoreState.h"
%include "IMP/restraints/AngleRestraint.h"
+%include "IMP/restraints/PairRestraint.h"
%include "IMP/restraints/BondDecoratorRestraint.h"
%include "IMP/restraints/ConnectivityRestraint.h"
%include "IMP/restraints/ConstantRestraint.h"
@@ -247,6 +278,7 @@
%include "IMP/restraints/NonbondedRestraint.h"
%include "IMP/restraints/PairChainRestraint.h"
%include "IMP/restraints/PairListRestraint.h"
+%include "IMP/restraints/LowestNPairListRestraint.h"
%include "IMP/restraints/RestraintSet.h"
%include "IMP/restraints/SingletonListRestraint.h"
%include "IMP/restraints/TripletChainRestraint.h"
@@ -275,6 +307,7 @@
%template(ParticleIndexes) ::std::vector<ParticleIndex>;
%template(BondDecorators) ::std::vector<BondDecorator>;
%template(MolecularHiearchyDecorators) ::std::vector<MolecularHierarchyDecorator>;
+ %template(XYZDecorators) ::std::vector<XYZDecorator>;
%template(FloatKeys) ::std::vector<FloatKey>;
%template(StringKeys) ::std::vector<StringKey>;
%template(IntKeys) ::std::vector<IntKey>;
@@ -282,4 +315,5 @@
%template(Floats) ::std::vector<Float>;
%template(Strings) ::std::vector<String>;
%template(Ints) ::std::vector<Int>;
+ %template(Segment) ::std::pair< IMP::Vector3D, IMP::Vector3D >;
}
Index: SConstruct
===================================================================
--- SConstruct (revision 649)
+++ SConstruct (working copy)
@@ -43,7 +43,7 @@
#SConscript('domino/SConscript')
# bin script first requires kernel libraries to be built:
-env.Depends(bin, [src, pyext, pymod])
+env.Depends(bin, [pyext, src, pymod])
# Build the binaries by default:
env.Default(bin)