Index: kernel/include/IMP/score_states/QuadraticNonbondedListScoreState.h
===================================================================
--- kernel/include/IMP/score_states/QuadraticNonbondedListScoreState.h (revision 0)
+++ kernel/include/IMP/score_states/QuadraticNonbondedListScoreState.h (revision 0)
@@ -0,0 +1,53 @@
+/**
+ * \file QuadraticNonbondedListScoreState.h
+ * \brief Allow iteration through pairs of a set of spheres.
+ *
+ * Copyright 2007-8 Sali Lab. All rights reserved.
+ */
+
+#ifndef __IMP_QUADRATIC_NONBONDED_LIST_SCORE_STATE_H
+#define __IMP_QUADRATIC_NONBONDED_LIST_SCORE_STATE_H
+
+#include "NonbondedListScoreState.h"
+#include "../internal/kernel_version_info.h"
+#include "MaxChangeScoreState.h"
+
+#include <vector>
+#include <limits>
+
+namespace IMP
+{
+
+//! A base class for nonbonded lists which test all pairs of particles
+/** This should not be used by end users. But since it needs to be in the
+ inheritance hierarchy, it should be in the IMP namespace.
+ */
+class IMPDLLEXPORT QuadraticNonbondedListScoreState:
+ public IMP::NonbondedListScoreState
+{
+ typedef NonbondedListScoreState P;
+ internal::ObjectPointer<MaxChangeScoreState, true> mc_;
+ float slack_;
+ protected:
+ void handle_nbl_pair( Particle *a, Particle *b, float d);
+ const Particles &get_particles() const {
+ return mc_->get_particles();
+ }
+ void set_particles(const Particles &ps) {
+ P::clear_nbl();
+ mc_->clear_particles();
+ mc_->add_particles(ps);
+ }
+
+ QuadraticNonbondedListScoreState(FloatKey radius);
+ ~QuadraticNonbondedListScoreState();
+
+ public:
+ void update();
+};
+
+//! This class maintains a list of non-bonded pairs of spheres
+
+} // namespace IMP
+
+#endif /* __IMP_QUADRATIC_NONBONDED_LIST_SCORE_STATE_H */
Index: kernel/include/IMP/score_states/BondDecoratorListScoreState.h
===================================================================
--- kernel/include/IMP/score_states/BondDecoratorListScoreState.h (revision 489)
+++ kernel/include/IMP/score_states/BondDecoratorListScoreState.h (working copy)
@@ -31,7 +31,7 @@
//! Find bonds amongst the following points.
/** \param [in] ps The set of particles to use.
*/
- BondDecoratorListScoreState(const Particles &ps);
+ BondDecoratorListScoreState(const Particles &ps= Particles());
virtual ~BondDecoratorListScoreState(){}
virtual void set_particles(const Particles &ps);
Index: kernel/include/IMP/score_states/NonbondedListScoreState.h
===================================================================
--- kernel/include/IMP/score_states/NonbondedListScoreState.h (revision 489)
+++ kernel/include/IMP/score_states/NonbondedListScoreState.h (working copy)
@@ -22,18 +22,19 @@
typedef std::vector<BondedListScoreState*> BondedListScoreStates;
//! A base class for classes that maintain a list of non-bonded pairs.
-/** \ingroup restraint
+/**
*/
class IMPDLLEXPORT NonbondedListScoreState: public ScoreState
{
private:
+ FloatKey rk_;
+
protected:
// made protected for debugging code, do not use otherwise
typedef std::vector<std::pair<Particle*, Particle*> > NBL;
NBL nbl_;
float last_cutoff_;
-
unsigned int size_nbl() const {return nbl_.size();}
//! rebuild the nonbonded list
@@ -46,14 +47,25 @@
nbl_.clear();
}
+ bool are_bonded(Particle *a, Particle *b) const {
+ for (BondedListConstIterator bli= bonded_lists_begin();
+ bli != bonded_lists_end(); ++bli) {
+ if ((*bli)->are_bonded(a, b)) {
+ return true;
+ }
+ }
+ return false;
+ }
+
struct AddToNBL;
friend struct AddToNBL;
+ // these should be two separate classes, but they are not
struct AddToNBL{
NonbondedListScoreState *state_;
Particle *p_;
AddToNBL(NonbondedListScoreState *s, Particle *p): state_(s),
- p_(p){}
+ p_(p){}
void operator()(Particle *p) {
operator()(p_, p);
}
@@ -62,16 +74,6 @@
}
};
- bool are_bonded(Particle *a, Particle *b) const {
- for (BondedListConstIterator bli= bonded_lists_begin();
- bli != bonded_lists_end(); ++bli) {
- if ((*bli)->are_bonded(a, b)) {
- return true;
- }
- }
- return false;
- }
-
//! tell the bonded lists what particles to pay attention to
void propagate_particles(const Particles &ps);
@@ -82,14 +84,29 @@
IMP_LOG(VERBOSE, "Found pair " << a->get_index()
<< " " << b->get_index() << std::endl);
nbl_.push_back(std::make_pair(a, b));
+ } else {
+ IMP_LOG(VERBOSE, "Pair " << a->get_index()
+ << " and " << b->get_index() << " rejected on bond"
+ <<std::endl);
}
}
+ Float get_radius(Particle *a) const {
+ if (rk_ != FloatKey() && a->has_attribute(rk_)) {
+ return a->get_value(rk_);
+ } else {
+ return 0;
+ }
+ }
+
public:
/**
*/
- NonbondedListScoreState();
+ NonbondedListScoreState(FloatKey rk=FloatKey());
+ FloatKey get_radius_key() const {return rk_;}
+ void set_radius_key(FloatKey rk) {rk_=rk;}
+
IMP_CONTAINER(BondedListScoreState, bonded_list,
BondedListIndex);
// kind of evil hack to make the names better
Index: kernel/include/IMP/score_states/AllSphereNonbondedListScoreState.h
===================================================================
--- kernel/include/IMP/score_states/AllSphereNonbondedListScoreState.h (revision 489)
+++ kernel/include/IMP/score_states/AllSphereNonbondedListScoreState.h (working copy)
@@ -1,76 +0,0 @@
-/**
- * \file AllSphereNonbondedListScoreState.h
- * \brief Allow iteration through pairs of a set of spheres.
- *
- * Copyright 2007-8 Sali Lab. All rights reserved.
- */
-
-#ifndef __IMP_ALL_SPHERE_NONBONDED_LIST_SCORE_STATE_H
-#define __IMP_ALL_SPHERE_NONBONDED_LIST_SCORE_STATE_H
-
-#include "NonbondedListScoreState.h"
-#include "../internal/ParticleGrid.h"
-#include "../internal/kernel_version_info.h"
-#include "MaxChangeScoreState.h"
-
-#include <vector>
-#include <limits>
-
-namespace IMP
-{
-
-//! This class maintains a list of non-bonded pairs of spheres
-/** \note The radius is currently assumed not to change. This could
- be fixed later.
- \note Points whose coordinates are not optimized are assumed to
- stay that way (and vice
- versa, although that direction doesn't matter so much).
- \note The structure is slightly dumb about rebuilding and will
- rebuild the whole list of any of the grids become invalidated.
- This could be improved as each piece is computed separately (so
- they could be cached).
-
- \ingroup restraint
- */
-class IMPDLLEXPORT AllSphereNonbondedListScoreState:
- public NonbondedListScoreState
-{
-protected:
- //! \internal
- struct Bin
- {
- internal::ParticleGrid *grid;
- Float rmax;
- Bin(): grid(NULL), rmax(-1){}
- Bin(const Bin &o): grid(o.grid), rmax(o.rmax){}
- };
- std::vector<Bin> bins_;
- std::vector<Bin> fixed_bins_;
- FloatKey rk_;
-
- //! \internal
- void rebuild_nbl(Float cut);
-
- void repartition_points(const Particles &ps, std::vector<Bin> &out);
-
- float side_from_r(float r) const {return r*1.6;}
-
- void generate_nbl(const Bin &particle_bin, const Bin &grid_bin, float cut);
-
- void cleanup(std::vector<Bin> &bins);
-
-public:
- /**
- \param[in] ps A list of particles to use.
- \param[in] radius The key to use to get the radius
- */
- AllSphereNonbondedListScoreState(const Particles &ps, FloatKey radius);
- ~AllSphereNonbondedListScoreState();
- IMP_SCORE_STATE(internal::kernel_version_info)
-
- void set_particles(const Particles &ps);
-};
-
-} // namespace IMP
-
-#endif /* __IMP_ALL_SPHERE_NONBONDED_LIST_SCORE_STATE_H */
Index: kernel/include/IMP/score_states/SConscript
===================================================================
--- kernel/include/IMP/score_states/SConscript (revision 489)
+++ kernel/include/IMP/score_states/SConscript (working copy)
@@ -1,13 +1,15 @@
+Import('env')
import os.path
-Import('env')
-
-files = ['BondedListScoreState.h', 'NonbondedListScoreState.h',
- 'BipartiteNonbondedListScoreState.h', 'MaxChangeScoreState.h',
- 'BondDecoratorListScoreState.h', 'AllNonbondedListScoreState.h',
- 'AllSphereNonbondedListScoreState.h']
-
-
-# Install the include files:
-includedir = os.path.join(env['includedir'], 'IMP', 'score_states')
+files=[]
+files.append( 'AllNonbondedListScoreState.h' )
+files.append( 'BipartiteNonbondedListScoreState.h' )
+files.append( 'BondDecoratorListScoreState.h' )
+files.append( 'BondedListScoreState.h' )
+files.append( 'MaxChangeScoreState.h' )
+files.append( 'NonbondedListScoreState.h' )
+files.append( 'QuadraticAllNonbondedListScoreState.h' )
+files.append( 'QuadraticBipartiteNonbondedListScoreState.h' )
+files.append( 'QuadraticNonbondedListScoreState.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/QuadraticBipartiteNonbondedListScoreState.h
===================================================================
--- kernel/include/IMP/score_states/QuadraticBipartiteNonbondedListScoreState.h (revision 0)
+++ kernel/include/IMP/score_states/QuadraticBipartiteNonbondedListScoreState.h (revision 0)
@@ -0,0 +1,51 @@
+/**
+ * \file QuadraticBipartiteNonbondedListScoreState.h
+ * \brief Allow iteration through pairs of a set of atoms.
+ *
+ * Copyright 2007-8 Sali Lab. All rights reserved.
+ */
+
+#ifndef __IMP_QUADRATIC_BIPARTITE_NONBONDED_LIST_SCORE_STATE_H
+#define __IMP_QUADRATIC_BIPARTITE_NONBONDED_LIST_SCORE_STATE_H
+
+#include "QuadraticNonbondedListScoreState.h"
+#include "../internal/ParticleGrid.h"
+#include "../internal/kernel_version_info.h"
+
+#include <vector>
+#include <limits>
+
+namespace IMP
+{
+
+//! This class maintains a list of non-bonded pairs between two sets.
+/** The class works roughly like the NonbondedListScoreState except
+ only pairs where one particle is taken from each set are returned.
+
+ \note QuadraticBipartiteNonbondedListScoreState is basically an
+ implementation detail for performance analysis and should not
+ be used by end users.
+ */
+class IMPDLLEXPORT QuadraticBipartiteNonbondedListScoreState:
+ public QuadraticNonbondedListScoreState
+{
+ typedef QuadraticNonbondedListScoreState P;
+ unsigned int na_;
+
+ virtual void rebuild_nbl(float cut);
+public:
+ QuadraticBipartiteNonbondedListScoreState(FloatKey rk,
+ const Particles &ps0,
+ const Particles &ps1);
+ QuadraticBipartiteNonbondedListScoreState(const Particles &ps0,
+ const Particles &ps1);
+ QuadraticBipartiteNonbondedListScoreState(FloatKey rk=FloatKey());
+
+ IMP_SCORE_STATE(internal::kernel_version_info)
+
+ void set_particles(const Particles &ps0, const Particles &ps1);
+};
+
+} // namespace IMP
+
+#endif /* __IMP_QUADRATIC_BIPARTITE_NONBONDED_LIST_SCORE_STATE_H */
Index: kernel/include/IMP/score_states/BipartiteNonbondedListScoreState.h
===================================================================
--- kernel/include/IMP/score_states/BipartiteNonbondedListScoreState.h (revision 489)
+++ kernel/include/IMP/score_states/BipartiteNonbondedListScoreState.h (working copy)
@@ -8,39 +8,27 @@
#ifndef __IMP_BIPARTITE_NONBONDED_LIST_SCORE_STATE_H
#define __IMP_BIPARTITE_NONBONDED_LIST_SCORE_STATE_H
-#include "NonbondedListScoreState.h"
-#include "../internal/ParticleGrid.h"
-#include "../internal/kernel_version_info.h"
+#include "QuadraticBipartiteNonbondedListScoreState.h"
-#include <vector>
-#include <limits>
-
namespace IMP
{
+//! Maintain a nonbonded list between two disjoint sets.
+/**
+ \note If no value for the radius key is specified, all radii are
+ considered to be zero.
-class BondedListScoreState;
-
-//! This class maintains a list of non-bonded pairs between two sets.
-/** The class works roughly like the NonbondedListScoreState except
- only pairs where one particle is taken from each set are returned.
- \ingroup restraint
+ \ingroup restraint
*/
class IMPDLLEXPORT BipartiteNonbondedListScoreState:
- public NonbondedListScoreState
-{
- typedef NonbondedListScoreState P;
- internal::ObjectPointer<MaxChangeScoreState, true> mc_;
- internal::ParticleGrid grid_;
-
- virtual void rebuild_nbl(float cut);
-public:
+ public QuadraticBipartiteNonbondedListScoreState {
+ typedef QuadraticBipartiteNonbondedListScoreState P;
+ public:
+ BipartiteNonbondedListScoreState(FloatKey rk,
+ const Particles &ps0,
+ const Particles &ps1): P(rk, ps0, ps1){}
BipartiteNonbondedListScoreState(const Particles &ps0,
- const Particles &ps1,
- float target_side);
-
- IMP_SCORE_STATE(internal::kernel_version_info)
-
- void set_particles(const Particles &ps0, const Particles &ps1);
+ const Particles &ps1): P(ps0, ps1){}
+ BipartiteNonbondedListScoreState(FloatKey rk=FloatKey()): P(rk){}
};
} // namespace IMP
Index: kernel/include/IMP/score_states/QuadraticAllNonbondedListScoreState.h
===================================================================
--- kernel/include/IMP/score_states/QuadraticAllNonbondedListScoreState.h (revision 0)
+++ kernel/include/IMP/score_states/QuadraticAllNonbondedListScoreState.h (revision 0)
@@ -0,0 +1,54 @@
+/**
+ * \file QuadraticAllNonbondedListScoreState.h
+ * \brief Allow iteration through pairs of a set of spheres.
+ *
+ * Copyright 2007-8 Sali Lab. All rights reserved.
+ */
+
+#ifndef __IMP_QUADRATIC_ALL_SPHERE_NONBONDED_LIST_SCORE_STATE_H
+#define __IMP_QUADRATIC_ALL_SPHERE_NONBONDED_LIST_SCORE_STATE_H
+
+#include "NonbondedListScoreState.h"
+#include "../internal/kernel_version_info.h"
+#include "QuadraticNonbondedListScoreState.h"
+
+#include <vector>
+#include <limits>
+
+namespace IMP
+{
+
+//! This class maintains a list of non-bonded pairs of spheres
+/** \note Points whose coordinates are not optimized are assumed to
+ stay that way and are pairs involving only fixed points are skipped
+
+ \note QuadraticBipartiteNonbondedListScoreState is basically an
+ implementation detail for performance analysis and should not
+ be used by end users.
+ */
+class IMPDLLEXPORT QuadraticAllNonbondedListScoreState:
+ public QuadraticNonbondedListScoreState
+{
+ typedef QuadraticNonbondedListScoreState P;
+ Particles fixed_;
+
+ //! \internal
+ void rebuild_nbl(Float cut);
+
+public:
+ /**
+ \param[in] ps A list of particles to use.
+ \param[in] radius The key to use to get the radius
+ */
+ QuadraticAllNonbondedListScoreState(FloatKey radius= FloatKey(),
+ const Particles &ps= Particles());
+ QuadraticAllNonbondedListScoreState(const Particles &ps);
+ ~QuadraticAllNonbondedListScoreState();
+ IMP_SCORE_STATE(internal::kernel_version_info)
+
+ void set_particles(const Particles &ps);
+};
+
+} // namespace IMP
+
+#endif /* __IMP_QUADRATIC_ALL_NONBONDED_LIST_SCORE_STATE_H */
Index: kernel/include/IMP/score_states/MaxChangeScoreState.h
===================================================================
--- kernel/include/IMP/score_states/MaxChangeScoreState.h (revision 489)
+++ kernel/include/IMP/score_states/MaxChangeScoreState.h (working copy)
@@ -18,20 +18,31 @@
{
//! 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.
+
+ */
class IMPDLLEXPORT MaxChangeScoreState: public ScoreState
{
FloatKeys keys_;
- std::vector<std::vector<float> > orig_;
+ std::vector<Floats > orig_;
Particles ps_;
float max_change_;
public:
- MaxChangeScoreState(const FloatKeys &keys):keys_(keys){}
+ //! Track the changes with the specified keys.
+ MaxChangeScoreState(const FloatKeys &keys,
+ const Particles &ps= Particles()):keys_(keys){
+ add_particles(ps);
+ }
virtual ~MaxChangeScoreState(){}
virtual void update();
+ //! Measure differences from the current value.
void reset();
+ //! Return the maximum amount any attribute has changed.
float get_max() const {
return max_change_;
}
Index: kernel/include/IMP/score_states/AllNonbondedListScoreState.h
===================================================================
--- kernel/include/IMP/score_states/AllNonbondedListScoreState.h (revision 489)
+++ kernel/include/IMP/score_states/AllNonbondedListScoreState.h (working copy)
@@ -1,6 +1,6 @@
/**
* \file AllNonbondedListScoreState.h
- * \brief Allow iteration through pairs of a set of particles.
+ * \brief Allow iteration through pairs of a set of s.
*
* Copyright 2007-8 Sali Lab. All rights reserved.
*/
@@ -9,7 +9,9 @@
#define __IMP_ALL_NONBONDED_LIST_SCORE_STATE_H
#include "NonbondedListScoreState.h"
+#include "../internal/ParticleGrid.h"
#include "../internal/kernel_version_info.h"
+#include "MaxChangeScoreState.h"
#include <vector>
#include <limits>
@@ -17,24 +19,59 @@
namespace IMP
{
-//! This class maintains a list of non-bonded pairs.
-/** The nonbonded list returns all nonboned pairs in the set of points.
+//! This class maintains a list of non-bonded pairs of particles
+/** \note If no value for the radius key is specified, all radii are
+ considered to be zero.
+
+ \note The radius is currently assumed not to change. This could
+ be fixed later.
+
+ \note Points whose coordinates are not optimized are assumed to
+ stay that way (and vice versa, although that direction doesn't
+ matter so much).
+ \todo The structure is slightly dumb about rebuilding and will
+ rebuild the whole list of any of the grids become invalidated.
+ This could be improved as each piece is computed separately (so
+ they could be cached).
+
\ingroup restraint
*/
-class IMPDLLEXPORT AllNonbondedListScoreState: public NonbondedListScoreState
+class IMPDLLEXPORT AllNonbondedListScoreState:
+ public NonbondedListScoreState
{
+ typedef NonbondedListScoreState P;
protected:
- internal::ParticleGrid grid_;
+ //! \internal
+ struct Bin
+ {
+ internal::ParticleGrid *grid;
+ Float rmax;
+ Bin(): grid(NULL), rmax(-1){}
+ Bin(const Bin &o): grid(o.grid), rmax(o.rmax){}
+ };
+ std::vector<Bin> bins_;
+ std::vector<Bin> fixed_bins_;
- virtual void rebuild_nbl(float cut);
+ //! \internal
+ void rebuild_nbl(Float cut);
+ void repartition_points(const Particles &ps, std::vector<Bin> &out);
+
+ float side_from_r(float r) const;
+
+ void generate_nbl(const Bin &particle_bin, const Bin &grid_bin, float cut);
+
+ void cleanup(std::vector<Bin> &bins);
+
public:
/**
\param[in] ps A list of particles to use.
- \param[in] tvs A suggested size for the voxel side.
- */
- AllNonbondedListScoreState(const Particles &ps, float tvs=1);
-
+ \param[in] radius The key to use to get the radius
+ */
+ AllNonbondedListScoreState(FloatKey radius= FloatKey(),
+ const Particles &ps=Particles());
+ AllNonbondedListScoreState(const Particles &ps);
+ ~AllNonbondedListScoreState();
IMP_SCORE_STATE(internal::kernel_version_info)
void set_particles(const Particles &ps);
Index: kernel/test/states/test_nonbonded_list.py
===================================================================
--- kernel/test/states/test_nonbonded_list.py (revision 489)
+++ kernel/test/states/test_nonbonded_list.py (working copy)
@@ -33,14 +33,17 @@
p= IMP.Particle()
m.add_particle(p)
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));
- s= IMP.AllNonbondedListScoreState(m.get_particles())
+ s= IMP.AllNonbondedListScoreState(IMP.FloatKey(),
+ m.get_particles())
m.add_score_state(s)
o= OnePair()
r= IMP.NonbondedRestraint(o, s, 15)
m.add_restraint(r)
score= m.evaluate(False)
+ print score
self.assertEqual(score, 4950, "Wrong score")
def test_bl(self):
@@ -51,6 +54,7 @@
p= IMP.Particle()
m.add_particle(p)
d=IMP.XYZDecorator.create(p)
+ d.set_coordinates_are_optimized(True)
d.set_x(0)
d.set_y(10)
d.set_z(10)
@@ -60,7 +64,8 @@
pts.append(p)
for i in range(1,100):
IMP.custom_bond(bds[i-1], bds[i], 1, .1)
- s= IMP.AllNonbondedListScoreState(pts, 1)
+ s= IMP.AllNonbondedListScoreState(IMP.FloatKey(),
+ pts)
b= IMP.BondDecoratorListScoreState(pts)
s.add_bonded_list(b)
m.add_score_state(s)
@@ -84,19 +89,21 @@
m.add_particle(p)
d=IMP.XYZDecorator.create(p)
ps.append(p)
+ d.set_coordinates_are_optimized(True)
if (i < 25):
d.randomize_in_box(IMP.Vector3D(0,0,0),
IMP.Vector3D(10,10,10));
else:
d.randomize_in_box(IMP.Vector3D(60,60,60),
IMP.Vector3D(70,70,70));
- s= IMP.AllNonbondedListScoreState(ps)
+ s= IMP.AllNonbondedListScoreState(IMP.FloatKey(), ps)
m.add_score_state(s)
o= OnePair()
r= IMP.NonbondedRestraint(o, s, 15)
m.add_restraint(r)
score= m.evaluate(False)
self.assertEqual(score, 1225, "Wrong score")
+
def test_bi(self):
"""Test the bipartite nonbonded list and restraint which uses it"""
m= IMP.Model()
@@ -106,13 +113,15 @@
p= IMP.Particle()
m.add_particle(p)
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));
if (i < 5):
ps0.append(p)
else:
ps1.append(p)
- s= IMP.BipartiteNonbondedListScoreState(ps0, ps1, 1000)
+ s= IMP.QuadraticBipartiteNonbondedListScoreState(IMP.FloatKey(),
+ ps0, ps1)
m.add_score_state(s)
o= OnePair()
r= IMP.NonbondedRestraint(o, s, 15)
@@ -120,7 +129,7 @@
score= m.evaluate(False)
self.assertEqual(score, 25, "Wrong score")
def test_spheres2(self):
- """Test the nonbonded list of spheres (num pairs)"""
+ """Test the (quadratic and optimized) nonbonded lists of spheres (num pairs)"""
m= IMP.Model()
rk= IMP.FloatKey("radius")
for i in range(0,100):
@@ -131,12 +140,37 @@
IMP.Vector3D(10,10,10));
p.add_attribute(rk, i, False)
d.set_coordinates_are_optimized(True)
- s= IMP.AllSphereNonbondedListScoreState(m.get_particles(), rk)
+ s= IMP.AllNonbondedListScoreState(rk,
+ m.get_particles())
+ s2= IMP.QuadraticAllNonbondedListScoreState(rk, m.get_particles())
m.add_score_state(s)
+ m.add_score_state(s2)
o= OnePair()
+ o2= OnePair()
r= IMP.NonbondedRestraint(o, s, 15)
+ r2= IMP.NonbondedRestraint(o2, s2, 15)
m.add_restraint(r)
+ m.add_restraint(r2)
score= m.evaluate(False)
+ self.assertEqual(score, 2*4950, "Wrong score")
+ def test_spheres3(self):
+ """Test the quadratic nonbonded list of spheres (num pairs)"""
+ m= IMP.Model()
+ rk= IMP.FloatKey("radius")
+ for i in range(0,100):
+ 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));
+ p.add_attribute(rk, i, False)
+ d.set_coordinates_are_optimized(True)
+ s= IMP.QuadraticAllNonbondedListScoreState(rk, m.get_particles())
+ m.add_score_state(s)
+ o= OnePair()
+ r= IMP.NonbondedRestraint(o, s, 15)
+ m.add_restraint(r)
+ score= m.evaluate(False)
self.assertEqual(score, 4950, "Wrong score")
def test_frido_spheres(self):
"""Test the nonbonded list with frido's spheres"""
@@ -262,7 +296,7 @@
print ".sphere "+str(d.get_x())+ " " + str(d.get_y())\
+ " " + str(d.get_z()) + " " + str(p.get_value(rk))
- s= IMP.AllSphereNonbondedListScoreState(m.get_particles(), rk)
+ s= IMP.AllNonbondedListScoreState(rk, m.get_particles())
m.add_score_state(s)
sd= IMP.SphereDistancePairScore(IMP.HarmonicLowerBound(0,1),
rk)
@@ -283,7 +317,7 @@
IMP.Vector3D(20,20,20));
p.add_attribute(rk, random.uniform(0,100), False)
d.set_coordinates_are_optimized(True)
- s= IMP.AllSphereNonbondedListScoreState(m.get_particles(), rk)
+ s= IMP.AllNonbondedListScoreState(rk, m.get_particles())
m.add_score_state(s)
sd= IMP.SphereDistancePairScore(IMP.HarmonicLowerBound(0,1),
rk)
Index: kernel/src/score_states/BondDecoratorListScoreState.cpp
===================================================================
--- kernel/src/score_states/BondDecoratorListScoreState.cpp (revision 489)
+++ kernel/src/score_states/BondDecoratorListScoreState.cpp (working copy)
@@ -19,7 +19,7 @@
void BondDecoratorListScoreState::update()
{
- IMP_LOG(VERBOSE, "Updating BondDecoratorList for "
+ IMP_LOG(TERSE, "Updating BondDecoratorList for "
<< ps_.size() << " particles" << std::endl);
bonds_.clear();
for (unsigned int i=0; i< ps_.size(); ++i) {
@@ -35,11 +35,12 @@
continue;
}
if (di < dj) {
+ IMP_LOG(VERBOSE, "Found bond " << di.get_bond(j) << std::endl);
bonds_.push_back(di.get_bond(j));
}
}
}
- IMP_LOG(VERBOSE, "Found " << bonds_.size() << " bonds"<< std::endl);
+ IMP_LOG(TERSE, "Found " << bonds_.size() << " bonds"<< std::endl);
}
void BondDecoratorListScoreState::set_particles(const Particles &ps)
Index: kernel/src/score_states/NonbondedListScoreState.cpp
===================================================================
--- kernel/src/score_states/NonbondedListScoreState.cpp (revision 489)
+++ kernel/src/score_states/NonbondedListScoreState.cpp (working copy)
@@ -13,7 +13,7 @@
namespace IMP
{
-NonbondedListScoreState::NonbondedListScoreState()
+NonbondedListScoreState::NonbondedListScoreState(FloatKey rk): rk_(rk)
{
last_cutoff_=-1;
}
Index: kernel/src/score_states/AllSphereNonbondedListScoreState.cpp
===================================================================
--- kernel/src/score_states/AllSphereNonbondedListScoreState.cpp (revision 489)
+++ kernel/src/score_states/AllSphereNonbondedListScoreState.cpp (working copy)
@@ -1,242 +0,0 @@
-/**
- * \file AllSphereNonbondedListScoreState.cpp
- * \brief Allow iteration through pairs of a set of spheres.
- *
- * Copyright 2007-8 Sali Lab. All rights reserved.
- */
-
-#include "IMP/score_states/AllSphereNonbondedListScoreState.h"
-#include "IMP/decorators/XYZDecorator.h"
-
-#include <algorithm>
-
-namespace IMP
-{
-
-static unsigned int min_grid_size=20;
-
-AllSphereNonbondedListScoreState
-::AllSphereNonbondedListScoreState(const Particles &ps, FloatKey radius):
- rk_(radius)
-{
- set_particles(ps);
-}
-
-AllSphereNonbondedListScoreState::~AllSphereNonbondedListScoreState()
-{
- cleanup(bins_);
- cleanup(fixed_bins_);
-}
-
-void AllSphereNonbondedListScoreState::set_particles(const Particles &ps)
-{
- NonbondedListScoreState::clear_nbl();
- // split into moving and fixed and call repartition twice
- Particles moving, fixed;
- for (unsigned int i=0; i< ps.size(); ++i) {
- XYZDecorator d= XYZDecorator::cast(ps[i]);
- if (d.get_coordinates_are_optimized()) {
- moving.push_back(ps[i]);
- } else {
- fixed.push_back(ps[i]);
- }
- }
- cleanup(bins_);
- cleanup(fixed_bins_);
- repartition_points(moving, bins_);
- repartition_points(fixed, fixed_bins_);
- for (unsigned int i=0; i< fixed_bins_.size(); ++i) {
- fixed_bins_[i].grid->update();
- }
-}
-
-void AllSphereNonbondedListScoreState::repartition_points(const Particles &ps,
- std::vector<Bin> &out)
-{
- cleanup(out);
- if (ps.empty()) return;
- float minr=std::numeric_limits<float>::max(), maxr=0;
- for (unsigned int i=0; i< ps.size(); ++i) {
- float r= ps[i]->get_value(rk_);
- if ( r > maxr) maxr=r;
- if ( r > 0 && r < minr) minr=r;
- }
- float curr=minr*2;
- Floats cuts;
- do {
- cuts.push_back(curr);
- curr *= 2;
- } while (curr < maxr);
- cuts.push_back(2*maxr);
-
- std::vector<Particles> ops(cuts.size());
- for (unsigned int i=0; i< ps.size(); ++i) {
- float r= ps[i]->get_value(rk_);
- for (unsigned int j=0; ; ++j) {
- IMP_assert(j< cuts.size(), "Internal error in ASNBLSS");
- if (cuts[j] > r) {
- ops[j].push_back(ps[i]);
- break;
- }
- }
- }
- // consolidate
- for (unsigned int i=1; i< ops.size(); ++i) {
- if (ops[i-1].size() + ops[i].size() < min_grid_size) {
- ops[i].insert(ops[i].end(), ops[i-1].begin(), ops[i-1].end());
- ops[i-1].clear();
- }
- }
- for (unsigned int i=0; i< cuts.size(); ++i) {
- if (ops[i].empty()) continue;
- out.push_back(Bin());
- float rmax=0;
- for (unsigned int j=0; j< ops[i].size(); ++j) {
- rmax= std::max(rmax, ops[i][j]->get_value(rk_));
- }
- out.back().rmax= rmax;
- internal::ParticleGrid *pg
- = new internal::ParticleGrid(side_from_r(out.back().rmax));
- out.back().grid= pg;
- out.back().grid->add_particles(ops[i]);
- }
- IMP_LOG(VERBOSE, "Created " << out.size() << " grids" << std::endl);
- for (unsigned int i=0; i< out.size(); ++i) {
- IMP_LOG(VERBOSE, out[i].rmax
- << ": " << *out[i].grid << std::endl);
- }
-
-#ifndef NDEBUG
- Particles all;
- for (unsigned int i=0; i< out.size(); ++i) {
- all.insert(all.end(), out[i].grid->get_particles().begin(),
- out[i].grid->get_particles().end());
- }
- std::sort(all.begin(), all.end());
- Particles::iterator it = std::unique(all.begin(), all.end());
- IMP_assert(it == all.end(), "Duplicate particles " << all.size());
- IMP_assert(all.size() == ps.size(), "Wrong number of particles at end "
- << all.size() << " " << ps.size());
-#endif
-}
-
-void AllSphereNonbondedListScoreState::cleanup(std::vector<Bin> &bins)
-{
- for (unsigned int i=0; i< bins.size(); ++i) {
- delete bins[i].grid;
- }
- bins.clear();
-}
-
-void AllSphereNonbondedListScoreState::update()
-{
- NonbondedListScoreState::update();
- bool bad=false;
- for (unsigned int i=0; i< bins_.size(); ++i) {
- if (bins_[i].grid->update()) bad=true;
- IMP_LOG(VERBOSE, bins_[i].rmax << std::endl << *bins_[i].grid << std::endl);
- }
- if (bad) {
- IMP_LOG(VERBOSE, "Destroying nbl in Sphere list"<< std::endl);
- NonbondedListScoreState::clear_nbl();
- }
-}
-
-void AllSphereNonbondedListScoreState::generate_nbl(const Bin &particle_bin,
- const Bin &grid_bin,
- float cut)
-{
- IMP_CHECK_OBJECT(particle_bin.grid);
- IMP_CHECK_OBJECT(grid_bin.grid);
- IMP_LOG(VERBOSE, "Generate nbl for " << particle_bin.rmax
- << " and " << grid_bin.rmax << std::endl);
- for (unsigned int k=0; k< particle_bin.grid->get_particles().size(); ++k) {
- Particle *p= particle_bin.grid->get_particles()[k];
- NonbondedListScoreState::AddToNBL f(this, p);
- XYZDecorator d= XYZDecorator::cast(p);
- internal::ParticleGrid::VirtualIndex index
- = grid_bin.grid->get_virtual_index(Vector3D(d.get_x(),
- d.get_y(),
- d.get_z()));
- IMP_LOG(VERBOSE, "Searching for " << p->get_index()
- << " from " << index
- << " of bin " << grid_bin.rmax << std::endl);
- grid_bin.grid->apply_to_nearby(f, index,
- cut+2*particle_bin.rmax + grid_bin.rmax,
- false);
- }
-
-}
-
-void AllSphereNonbondedListScoreState::rebuild_nbl(Float cut)
-{
- IMP_LOG(VERBOSE, "Rebuilding NBL with " << bins_.size() << " dynamic and "
- << fixed_bins_.size() << " fixed" << std::endl);
- for (unsigned int i=0; i< bins_.size(); ++i) {
- for (unsigned int j=i+1; j< bins_.size(); ++j) {
- generate_nbl(bins_[i], bins_[j], cut);
- }
-
- for (unsigned int j=0; j< fixed_bins_.size(); ++j) {
- generate_nbl(bins_[i], fixed_bins_[j], cut);
- }
-
- // same code as AllNonbonded. Would be nice to consolidate
- internal::ParticleGrid::Index last_index;
- for (internal::ParticleGrid::ParticleVoxelIterator it
- = bins_[i].grid->particle_voxels_begin();
- it != bins_[i].grid->particle_voxels_end(); ++it) {
- IMP_LOG(VERBOSE, "Searching with particle " << it->first->get_index()
- << std::endl);
- NonbondedListScoreState::AddToNBL f(this, it->first);
- bins_[i].grid->apply_to_nearby(f, it->second,
- cut+3*bins_[i].rmax,
- true);
- if (it->second != last_index) {
- IMP_LOG(VERBOSE, "Searching in " << it->second
- << std::endl);
- bins_[i].grid->apply_to_cell_pairs(f, it->second);
- }
- last_index=it->second;
- }
- }
-
-#ifndef NDEBUG
- Particles ps;
- for (unsigned int i=0; i< bins_.size(); ++i) {
- if (bins_[i].grid) {
- ps.insert(ps.end(), bins_[i].grid->get_particles().begin(),
- bins_[i].grid->get_particles().end());
- }
- }
- for (unsigned int i=0; i< ps.size(); ++i) {
- XYZDecorator di= XYZDecorator::cast(ps[i]);
- for (unsigned int j=0; j< i; ++j) {
- XYZDecorator dj= XYZDecorator::cast(ps[j]);
- if (distance(di, dj) - ps[i]->get_value(rk_) - ps[j]->get_value(rk_)
- <= cut && !are_bonded(ps[i], ps[j])) {
- bool found=false;
- for (NonbondedIterator nit= nbl_.begin();
- nit != nbl_.end(); ++nit) {
- if (nit->first == ps[i] && nit->second == ps[j]
- || nit->first == ps[j] && nit->second == ps[i]) {
- IMP_assert(!found, "Entry is in list twice");
- found=true;
- }
- }
- IMP_assert(found, "Nonbonded list is missing "
- << ps[i]->get_index() << " " << di
- << " and " << ps[j]->get_index() << " " << dj << std::endl);
- }
- }
- }
-#endif
-}
-
-
-void AllSphereNonbondedListScoreState::show(std::ostream &out) const
-{
- out << "AllSphereNonbondedListScoreState" << std::endl;
-}
-
-} // namespace IMP
Index: kernel/src/score_states/QuadraticBipartiteNonbondedListScoreState.cpp
===================================================================
--- kernel/src/score_states/QuadraticBipartiteNonbondedListScoreState.cpp (revision 0)
+++ kernel/src/score_states/QuadraticBipartiteNonbondedListScoreState.cpp (revision 0)
@@ -0,0 +1,74 @@
+/**
+ * \file QuadraticBipartiteNonbondedListScoreState.cpp
+ * \brief Allow iteration through pairs of a set of atoms.
+ *
+ * Copyright 2007-8 Sali Lab. All rights reserved.
+ */
+
+#include "IMP/score_states/QuadraticBipartiteNonbondedListScoreState.h"
+#include "IMP/decorators/XYZDecorator.h"
+#include "IMP/score_states/MaxChangeScoreState.h"
+
+namespace IMP
+{
+
+QuadraticBipartiteNonbondedListScoreState
+::QuadraticBipartiteNonbondedListScoreState(FloatKey rk,
+ const Particles &ps0,
+ const Particles &ps1): P(rk)
+{
+ set_particles(ps0, ps1);
+}
+
+
+QuadraticBipartiteNonbondedListScoreState
+::QuadraticBipartiteNonbondedListScoreState(const Particles &ps0,
+ const Particles &ps1): P(FloatKey())
+{
+ set_particles(ps0, ps1);
+}
+
+
+QuadraticBipartiteNonbondedListScoreState
+::QuadraticBipartiteNonbondedListScoreState(FloatKey rk): P(rk)
+{
+}
+
+void QuadraticBipartiteNonbondedListScoreState::rebuild_nbl(float cut)
+{
+ IMP_LOG(TERSE, "Rebuilding QBNBL with lists of size " << na_
+ << " and " << P::get_particles().size() -na_
+ << " and cutoff " << cut << std::endl);
+ for (unsigned int i=0; i< na_; ++i) {
+ for (unsigned int j=na_; j < P::get_particles().size(); ++j) {
+ P::handle_nbl_pair(P::get_particles()[i],
+ P::get_particles()[j],
+ cut);
+ }
+ }
+ IMP_LOG(TERSE, "NBL has " << P::nbl_.size() << " pairs" << std::endl);
+}
+
+void QuadraticBipartiteNonbondedListScoreState
+::set_particles(const Particles &ps0,
+ const Particles &ps1)
+{
+ IMP_LOG(TERSE, "Setting QBNBLSS particles " << ps0.size()
+ << " and " << ps1.size() << std::endl);
+ Particles all(ps0);
+ all.insert(all.end(), ps1.begin(), ps1.end());
+ P::set_particles(all);
+ na_= ps0.size();
+}
+
+void QuadraticBipartiteNonbondedListScoreState::update()
+{
+ P::update();
+}
+
+void QuadraticBipartiteNonbondedListScoreState::show(std::ostream &out) const
+{
+ out << "QuadraticBipartiteNonbondedList" << std::endl;
+}
+
+} // namespace IMP
Index: kernel/src/score_states/BipartiteNonbondedListScoreState.cpp
===================================================================
--- kernel/src/score_states/BipartiteNonbondedListScoreState.cpp (revision 489)
+++ kernel/src/score_states/BipartiteNonbondedListScoreState.cpp (working copy)
@@ -1,94 +0,0 @@
-/**
- * \file BipartiteNonbondedListScoreState.cpp
- * \brief Allow iteration through pairs of a set of atoms.
- *
- * Copyright 2007-8 Sali Lab. All rights reserved.
- */
-
-#include "IMP/score_states/BipartiteNonbondedListScoreState.h"
-#include "IMP/decorators/XYZDecorator.h"
-#include "IMP/score_states/MaxChangeScoreState.h"
-
-namespace IMP
-{
-
- /*
- Iterator needs:
- iterator into ps_, end iterator into ps_
- cell for ps_[i]
- -> actually has value current iterator around cell
- -put in logic for current
- number of bins to search
-
-
- nonbipartite:
- iterator through all cells
- index in current cell
- -> actual value in iterator around cell
- num bins
- */
-
-BipartiteNonbondedListScoreState
-::BipartiteNonbondedListScoreState(const Particles &ps0,
- const Particles &ps1, float target):
- grid_(target)
-{
- FloatKeys fks;
- fks.push_back(FloatKey("x"));
- fks.push_back(FloatKey("y"));
- fks.push_back(FloatKey("z"));
- mc_= new MaxChangeScoreState(fks);
- set_particles(ps0, ps1);
-}
-
-void BipartiteNonbondedListScoreState::rebuild_nbl(float cut)
-{
- for (unsigned int i=0; i< mc_->get_particles().size(); ++i) {
- Particle *p= mc_->get_particles()[i];
- NonbondedListScoreState::AddToNBL f(this, p);
- XYZDecorator d= XYZDecorator::cast(p);
- internal::ParticleGrid::VirtualIndex index
- = grid_.get_virtual_index(Vector3D(d.get_x(), d.get_y(), d.get_z()));
- IMP_LOG(VERBOSE, "Searching for " << p->get_index()
- << " from cell " << index << std::endl);
- grid_.apply_to_nearby(f, index, cut, false);
- }
- IMP_LOG(VERBOSE, "Found " << P::size_nbl()
- << " bipartite nonbonded pairs"
- << " for " << mc_->get_particles().size()
- << " nongridded particles" << std::endl);
-}
-
-void BipartiteNonbondedListScoreState::set_particles(const Particles &ps0,
- const Particles &ps1)
-{
- IMP_LOG(VERBOSE, "Setting bipartite particles " << ps0.size()
- << " and " << ps1.size() << std::endl);
- mc_->clear_particles();
- mc_->add_particles(ps1);
- grid_.clear_particles();
- grid_.add_particles(ps0);
- Particles all(ps0);
- all.insert(all.end(), ps1.begin(), ps1.end());
- NonbondedListScoreState::propagate_particles(all);
- /*IMP_assert(ps0.size() == P::get_particles().size(),
- "Where did they go?");*/
- IMP_assert(ps1.size() == mc_->get_particles().size(),
- "Where did we go?");
-}
-
-void BipartiteNonbondedListScoreState::update()
-{
- NonbondedListScoreState::update();
- if (mc_->get_max() > grid_.get_voxel_size()
- || grid_.update()) {
- NonbondedListScoreState::clear_nbl();
- }
-}
-
-void BipartiteNonbondedListScoreState::show(std::ostream &out) const
-{
- out << "BipartiteNonbondedList" << std::endl;
-}
-
-} // namespace IMP
Index: kernel/src/score_states/SConscript
===================================================================
--- kernel/src/score_states/SConscript (revision 489)
+++ kernel/src/score_states/SConscript (working copy)
@@ -1,9 +1,10 @@
Import('env')
-
-files = ['NonbondedListScoreState.cpp', 'BipartiteNonbondedListScoreState.cpp',
- 'MaxChangeScoreState.cpp', 'BondDecoratorListScoreState.cpp',
- 'AllNonbondedListScoreState.cpp',
- 'AllSphereNonbondedListScoreState.cpp']
-
-files = [File(x) for x in files]
+files=[]
+files.append(File( 'AllNonbondedListScoreState.cpp' ))
+files.append(File( 'BondDecoratorListScoreState.cpp' ))
+files.append(File( 'MaxChangeScoreState.cpp' ))
+files.append(File( 'NonbondedListScoreState.cpp' ))
+files.append(File( 'QuadraticAllNonbondedListScoreState.cpp' ))
+files.append(File( 'QuadraticBipartiteNonbondedListScoreState.cpp' ))
+files.append(File( 'QuadraticNonbondedListScoreState.cpp' ))
Return('files')
Index: kernel/src/score_states/QuadraticAllNonbondedListScoreState.cpp
===================================================================
--- kernel/src/score_states/QuadraticAllNonbondedListScoreState.cpp (revision 0)
+++ kernel/src/score_states/QuadraticAllNonbondedListScoreState.cpp (revision 0)
@@ -0,0 +1,77 @@
+/**
+ * \file QuadraticAllNonbondedListScoreState.cpp
+ * \brief Allow iteration through pairs of a set of s.
+ *
+ * Copyright 2007-8 Sali Lab. All rights reserved.
+ */
+
+#include "IMP/score_states/QuadraticAllNonbondedListScoreState.h"
+#include "IMP/decorators/XYZDecorator.h"
+
+#include <algorithm>
+
+namespace IMP
+{
+
+
+QuadraticAllNonbondedListScoreState
+::QuadraticAllNonbondedListScoreState(FloatKey radius,
+ const Particles &ps): P(radius)
+{
+ set_particles(ps);
+}
+
+QuadraticAllNonbondedListScoreState
+::QuadraticAllNonbondedListScoreState(const Particles &ps): P(FloatKey())
+{
+ set_particles(ps);
+}
+
+QuadraticAllNonbondedListScoreState::~QuadraticAllNonbondedListScoreState()
+{
+}
+
+void QuadraticAllNonbondedListScoreState::set_particles(const Particles &ps)
+{
+ Particles moving;
+ fixed_.clear();
+ for (unsigned int i=0; i< ps.size(); ++i) {
+ XYZDecorator d= XYZDecorator::cast(ps[i]);
+ if (d.get_coordinates_are_optimized()) {
+ moving.push_back(ps[i]);
+ } else {
+ fixed_.push_back(ps[i]);
+ }
+ }
+ P::set_particles(moving);
+}
+
+
+void QuadraticAllNonbondedListScoreState::update()
+{
+ P::update();
+}
+
+
+void QuadraticAllNonbondedListScoreState::rebuild_nbl(Float cut)
+{
+ IMP_LOG(TERSE, "Rebuilding QNBL with cutoff " << cut << std::endl);
+ const Particles &moving= P::get_particles();
+ for (unsigned int j=0; j< moving.size(); ++j) {
+ for (unsigned int i=0; i< fixed_.size(); ++i) {
+ P::handle_nbl_pair(fixed_[i], moving[j], cut);
+ }
+ for (unsigned int i=0; i< j; ++i) {
+ P::handle_nbl_pair(moving[i], moving[j], cut);
+ }
+ }
+ IMP_LOG(TERSE, "NBL has " << P::nbl_.size() << " pairs" << std::endl);
+}
+
+
+void QuadraticAllNonbondedListScoreState::show(std::ostream &out) const
+{
+ out << "QuadraticAllNonbondedListScoreState" << std::endl;
+}
+
+} // namespace IMP
Index: kernel/src/score_states/MaxChangeScoreState.cpp
===================================================================
--- kernel/src/score_states/MaxChangeScoreState.cpp (revision 489)
+++ kernel/src/score_states/MaxChangeScoreState.cpp (working copy)
@@ -19,11 +19,16 @@
max_change_=0;
for (unsigned int i=0; i < orig_.size(); ++i) {
for (unsigned int j=0; j < keys_.size(); ++j) {
+ IMP_LOG(VERBOSE, "Particle " << ps_[i]->get_index()
+ << " and attribute " << keys_[j]
+ << " moved " << std::abs(ps_[i]->get_value(keys_[j])
+ - orig_[i][j]) << std::endl);
max_change_= std::max(max_change_,
std::abs(ps_[i]->get_value(keys_[j])
- orig_[i][j]));
}
}
+ IMP_LOG(TERSE, "MaxChange update got " << max_change_ << std::endl);
}
Index: kernel/src/score_states/AllNonbondedListScoreState.cpp
===================================================================
--- kernel/src/score_states/AllNonbondedListScoreState.cpp (revision 489)
+++ kernel/src/score_states/AllNonbondedListScoreState.cpp (working copy)
@@ -1,54 +1,243 @@
/**
* \file AllNonbondedListScoreState.cpp
- * \brief Allow iteration through pairs of a set of atoms.
+ * \brief Allow iteration through pairs of a set of s.
*
* Copyright 2007-8 Sali Lab. All rights reserved.
*/
#include "IMP/score_states/AllNonbondedListScoreState.h"
#include "IMP/decorators/XYZDecorator.h"
-#include "IMP/internal/Grid3D.h"
-#include "IMP/score_states/MaxChangeScoreState.h"
+#include <algorithm>
+
namespace IMP
{
-AllNonbondedListScoreState::AllNonbondedListScoreState(const Particles &ps,
- float tvs):
- grid_(tvs)
+static unsigned int min_grid_size=20;
+
+AllNonbondedListScoreState
+::AllNonbondedListScoreState(FloatKey radius, const Particles &ps):
+ P(radius)
{
set_particles(ps);
}
-void AllNonbondedListScoreState::rebuild_nbl(float cut)
+AllNonbondedListScoreState
+::AllNonbondedListScoreState(const Particles &ps):
+ P(FloatKey())
{
- IMP_LOG(VERBOSE, "Scanning for nonbonded pairs..." << std::flush);
- internal::ParticleGrid::Index last_index;
- for (internal::ParticleGrid::ParticleVoxelIterator it
- = grid_.particle_voxels_begin();
- it != grid_.particle_voxels_end(); ++it) {
- IMP_LOG(VERBOSE, "Searching with particle " << it->first->get_index()
- << std::endl);
- NonbondedListScoreState::AddToNBL f(this, it->first);
- grid_.apply_to_nearby(f, it->second, cut, true);
- if (it->second != last_index) {
- grid_.apply_to_cell_pairs(f, it->second);
+ set_particles(ps);
+}
+
+
+AllNonbondedListScoreState::~AllNonbondedListScoreState()
+{
+ cleanup(bins_);
+ cleanup(fixed_bins_);
+}
+
+float AllNonbondedListScoreState::side_from_r(float r) const {
+ if (r==0) return 1;
+ else return r*1.6;
+}
+
+
+
+void AllNonbondedListScoreState::set_particles(const Particles &ps)
+{
+ NonbondedListScoreState::clear_nbl();
+ // split into moving and fixed and call repartition twice
+ Particles moving, fixed;
+ for (unsigned int i=0; i< ps.size(); ++i) {
+ XYZDecorator d= XYZDecorator::cast(ps[i]);
+ if (d.get_coordinates_are_optimized()) {
+ moving.push_back(ps[i]);
+ } else {
+ fixed.push_back(ps[i]);
}
- last_index=it->second;
}
+ cleanup(bins_);
+ cleanup(fixed_bins_);
+ repartition_points(moving, bins_);
+ repartition_points(fixed, fixed_bins_);
+ for (unsigned int i=0; i< fixed_bins_.size(); ++i) {
+ fixed_bins_[i].grid->update();
+ }
+}
- IMP_LOG(VERBOSE, "done" << std::endl);
- IMP_LOG(VERBOSE, "Found " << NonbondedListScoreState::size_nbl()
- << " nonbonded pairs"
- << std::endl);
+void AllNonbondedListScoreState::repartition_points(const Particles &ps,
+ std::vector<Bin> &out)
+{
+ cleanup(out);
+ if (ps.empty()) return;
+ float minr=std::numeric_limits<float>::max(), maxr=0;
+ for (unsigned int i=0; i< ps.size(); ++i) {
+ float r= P::get_radius(ps[i]);
+ if ( r > maxr) maxr=r;
+ if ( r > 0 && r < minr) minr=r;
+ }
+ minr= std::min(maxr, minr);
+ float curr=minr*2;
+ Floats cuts;
+ cuts.push_back(0);
+ do {
+ cuts.push_back(curr);
+ curr *= 2;
+ } while (curr < maxr);
+ cuts.push_back(2*maxr);
+ std::vector<Particles> ops(cuts.size());
+ for (unsigned int i=0; i< ps.size(); ++i) {
+ float r= P::get_radius(ps[i]);
+ bool found=false;
+ for (unsigned int j=0; ; ++j) {
+ IMP_assert(j< cuts.size(), "Internal error in ASNBLSS");
+ if (cuts[j] >= r) {
+ ops[j].push_back(ps[i]);
+ found=true;
+ break;
+ }
+ }
+ IMP_assert(found, "Didn't put particle anywhere");
+ }
+ // consolidate
+ for (unsigned int i=1; i< ops.size(); ++i) {
+ if (ops[i-1].size() + ops[i].size() < min_grid_size) {
+ ops[i].insert(ops[i].end(), ops[i-1].begin(), ops[i-1].end());
+ ops[i-1].clear();
+ }
+ }
+ for (unsigned int i=0; i< cuts.size(); ++i) {
+ if (ops[i].empty()) continue;
+ out.push_back(Bin());
+ float rmax=0;
+ for (unsigned int j=0; j< ops[i].size(); ++j) {
+ rmax= std::max(rmax, P::get_radius(ops[i][j]));
+ }
+ out.back().rmax= rmax;
+ internal::ParticleGrid *pg
+ = new internal::ParticleGrid(side_from_r(out.back().rmax));
+ out.back().grid= pg;
+ out.back().grid->add_particles(ops[i]);
+ }
+ IMP_LOG(VERBOSE, "Created " << out.size() << " grids" << std::endl);
+ for (unsigned int i=0; i< out.size(); ++i) {
+ IMP_LOG(VERBOSE, out[i].rmax
+ << ": " << *out[i].grid << std::endl);
+ }
+
#ifndef NDEBUG
- Particles ps= grid_.get_particles();
+ Particles all;
+ for (unsigned int i=0; i< out.size(); ++i) {
+ all.insert(all.end(), out[i].grid->get_particles().begin(),
+ out[i].grid->get_particles().end());
+ }
+ std::sort(all.begin(), all.end());
+ Particles::iterator it = std::unique(all.begin(), all.end());
+ IMP_assert(it == all.end(), "Duplicate particles " << all.size());
+ IMP_assert(all.size() == ps.size(), "Wrong number of particles at end "
+ << all.size() << " " << ps.size());
+#endif
+}
+
+void AllNonbondedListScoreState::cleanup(std::vector<Bin> &bins)
+{
+ for (unsigned int i=0; i< bins.size(); ++i) {
+ delete bins[i].grid;
+ }
+ bins.clear();
+}
+
+void AllNonbondedListScoreState::update()
+{
+ NonbondedListScoreState::update();
+ bool bad=false;
+ for (unsigned int i=0; i< bins_.size(); ++i) {
+ if (bins_[i].grid->update()) bad=true;
+ IMP_LOG(VERBOSE, bins_[i].rmax << std::endl << *bins_[i].grid << std::endl);
+ }
+ if (bad) {
+ IMP_LOG(VERBOSE, "Destroying nbl in list"<< std::endl);
+ NonbondedListScoreState::clear_nbl();
+ }
+}
+
+void AllNonbondedListScoreState::generate_nbl(const Bin &particle_bin,
+ const Bin &grid_bin,
+ float cut)
+{
+ IMP_CHECK_OBJECT(particle_bin.grid);
+ IMP_CHECK_OBJECT(grid_bin.grid);
+ IMP_LOG(VERBOSE, "Generate nbl for " << particle_bin.rmax
+ << " and " << grid_bin.rmax << std::endl);
+ for (unsigned int k=0; k< particle_bin.grid->get_particles().size(); ++k) {
+ Particle *p= particle_bin.grid->get_particles()[k];
+ NonbondedListScoreState::AddToNBL f(this, p);
+ XYZDecorator d= XYZDecorator::cast(p);
+ internal::ParticleGrid::VirtualIndex index
+ = grid_bin.grid->get_virtual_index(Vector3D(d.get_x(),
+ d.get_y(),
+ d.get_z()));
+ IMP_LOG(VERBOSE, "Searching for " << p->get_index()
+ << " from " << index
+ << " of bin " << grid_bin.rmax << std::endl);
+ grid_bin.grid->apply_to_nearby(f, index,
+ cut+2*particle_bin.rmax + grid_bin.rmax,
+ false);
+ }
+
+}
+
+void AllNonbondedListScoreState::rebuild_nbl(Float cut)
+{
+ IMP_LOG(TERSE, "Rebuilding NBL with " << bins_.size() << " dynamic and "
+ << fixed_bins_.size() << " fixed "
+ << " and cutoff " << cut << std::endl);
+ for (unsigned int i=0; i< bins_.size(); ++i) {
+ for (unsigned int j=i+1; j< bins_.size(); ++j) {
+ generate_nbl(bins_[i], bins_[j], cut);
+ }
+
+ for (unsigned int j=0; j< fixed_bins_.size(); ++j) {
+ generate_nbl(bins_[i], fixed_bins_[j], cut);
+ }
+
+ // same code as AllNonbonded. Would be nice to consolidate
+ internal::ParticleGrid::Index last_index;
+ for (internal::ParticleGrid::ParticleVoxelIterator it
+ = bins_[i].grid->particle_voxels_begin();
+ it != bins_[i].grid->particle_voxels_end(); ++it) {
+ IMP_LOG(VERBOSE, "Searching with particle " << it->first->get_index()
+ << std::endl);
+ NonbondedListScoreState::AddToNBL f(this, it->first);
+ bins_[i].grid->apply_to_nearby(f, it->second,
+ cut+3*bins_[i].rmax,
+ true);
+ if (it->second != last_index) {
+ IMP_LOG(VERBOSE, "Searching in " << it->second
+ << std::endl);
+ bins_[i].grid->apply_to_cell_pairs(f, it->second);
+ }
+ last_index=it->second;
+ }
+ }
+
+ IMP_LOG(TERSE, "NBL has " << P::nbl_.size() << " pairs" << std::endl);
+
+#ifndef NDEBUG
+ Particles ps;
+ for (unsigned int i=0; i< bins_.size(); ++i) {
+ if (bins_[i].grid) {
+ ps.insert(ps.end(), bins_[i].grid->get_particles().begin(),
+ bins_[i].grid->get_particles().end());
+ }
+ }
for (unsigned int i=0; i< ps.size(); ++i) {
XYZDecorator di= XYZDecorator::cast(ps[i]);
for (unsigned int j=0; j< i; ++j) {
XYZDecorator dj= XYZDecorator::cast(ps[j]);
- if (distance(di, dj) <= cut && !are_bonded(ps[i], ps[j]) ) {
+ if (distance(di, dj) - P::get_radius(ps[i]) - P::get_radius(ps[j])
+ <= cut && !are_bonded(ps[i], ps[j])) {
bool found=false;
for (NonbondedIterator nit= nbl_.begin();
nit != nbl_.end(); ++nit) {
@@ -67,25 +256,10 @@
#endif
}
-void AllNonbondedListScoreState::set_particles(const Particles &ps)
-{
- grid_.clear_particles();
- grid_.add_particles(ps);
- NonbondedListScoreState::propagate_particles(ps);
-}
-void AllNonbondedListScoreState::update()
-{
- NonbondedListScoreState::update();
- if (grid_.update()) {
- NonbondedListScoreState::clear_nbl();
- }
- IMP_LOG(VERBOSE, grid_);
-}
-
void AllNonbondedListScoreState::show(std::ostream &out) const
{
- out << "NonbondedList" << std::endl;
+ out << "AllNonbondedListScoreState" << std::endl;
}
} // namespace IMP
Index: kernel/src/score_states/QuadraticNonbondedListScoreState.cpp
===================================================================
--- kernel/src/score_states/QuadraticNonbondedListScoreState.cpp (revision 0)
+++ kernel/src/score_states/QuadraticNonbondedListScoreState.cpp (revision 0)
@@ -0,0 +1,55 @@
+/**
+ * \file QuadraticNonbondedListScoreState.cpp
+ * \brief Allow iteration through pairs of a set of s.
+ *
+ * Copyright 2007-8 Sali Lab. All rights reserved.
+ */
+
+#include "IMP/score_states/QuadraticNonbondedListScoreState.h"
+#include "IMP/decorators/XYZDecorator.h"
+
+#include <algorithm>
+
+namespace IMP
+{
+
+QuadraticNonbondedListScoreState
+::QuadraticNonbondedListScoreState(FloatKey radius):
+ P(radius), slack_(.1){
+ mc_= new MaxChangeScoreState(XYZDecorator::get_xyz_keys());
+}
+
+QuadraticNonbondedListScoreState::~QuadraticNonbondedListScoreState(){}
+
+
+void QuadraticNonbondedListScoreState::update() {
+ // placeholder to do tuning of slack
+ NonbondedListScoreState::update();
+ if (mc_->get_max() > slack_) {
+ NonbondedListScoreState::clear_nbl();
+ mc_->reset();
+ }
+}
+
+void QuadraticNonbondedListScoreState
+::handle_nbl_pair( Particle *a, Particle *b,
+ float d) {
+ XYZDecorator da= XYZDecorator::cast(a);
+ XYZDecorator db= XYZDecorator::cast(b);
+ float ra= P::get_radius(a);
+ float rb= P::get_radius(b);
+ for (unsigned int i=0; i< 3; ++i) {
+ float delta=std::abs(da.get_coordinate(i) - db.get_coordinate(i));
+ if (delta -ra -rb > d-slack_) {
+ IMP_LOG(VERBOSE, "Pair " << a->get_index()
+ << " and " << b->get_index() << " rejected on coordinate "
+ << i << std::endl);
+ return ;
+ }
+ }
+ IMP_LOG(VERBOSE, "Adding pair " << a->get_index()
+ << " and " << b->get_index() << std::endl);
+ P::add_to_nbl(a, b);
+}
+
+} // namespace IMP