IMP-dev

imp-dev@salilab.org

869 discussions

Nonbonded list
by Daniel Russel 15 Jan '08

15 Jan '08

This patch adds a non-bonded list. The implementation is dumb, but functional. The code more or less follows what was discussed earlier on the list. There is also a NonbondedRestraint which applies a PairScore to each nonbonded pair. There are PairScores for distance between the centers of the two particles and distance between the two spheres. The DistanceRestraint now uses the DistancePairScore internally to remove redundant code. There is some trivial testing.

1 0

Minor cleanup to bonds
by Daniel Russel 11 Jan '08

11 Jan '08

I cleaned up the bond decorator code. Among other things I removed the length parameter with the expectation that I will add another similar set of decorators for custom connectors (not quite sure what to call them). These will represent arbitrary length attachments between particles. I also added some test code for the bond decorators. The connectors will have a length and stiffness parameter and will be handled all together by a restraint (along the lines of what was discussed with the non-bonded list). Normal inter-atomic bonds could be handled by generating the appropriate lengths from the names of the atoms involved in each atomic bond and adding the connector information.

1 0

Re: [IMP-dev] minor quibble
by Ben Webb 11 Jan '08

11 Jan '08

Daniel Russel wrote: > It probably better to call it UnaryFunction rather than UnaryFunctor as > unary functors already mean something more general in C++ (any functor > that takes one variable). I guess strictly speaking it should be > UnaryRealFunction, but we can get elide the Real. Fair enough; as of r319 they are UnaryFunctions. Ben -- ben(a)salilab.org http://salilab.org/~ben/ "It is a capital mistake to theorize before one has data." - Sir Arthur Conan Doyle

1 0

Proposal: replace BasicScoreFuncParams with PairPotential
by Ben Webb 04 Jan '08

04 Jan '08

I propose to replace BasicScoreFuncParams with classes derived from PairPotential (of course the name is also up for discussion if you hate that one). This would then be used by restraints that need a ScoreFunc generator (Connectivity, ExclusionVolume, PairConnectivity, Proximity). PairPotential classes are simple generators that create new ScoreFunc objects, given two Particles: class PairPotential { public: virtual ScoreFunc* create_score_func(Particle *p1, Particle *p2) = 0; }; class UpperBoundPairPotential : public PairPotential { public: UpperBoundPairPotential(Float mean, Float stdev); UpperBoundPairPotential(FloatKey radius, Float stdev); ScoreFunc* create_score_func(Particle*, Particle*) { return new HarmonicUpperBound(some_mean, some_stdev); } }; Restraints such as ProximityRestraint then take a pointer to a PairPotential, and own it. The advantage is that ProximityRestraint's two constructors (one using a fixed mean, and the other a fixed particle radius attribute) can now be collapsed into one; the 'add two radii together to get a mean distance' logic is now in the PairPotential. The other advantage is that pair potentials which do not use mean/stdev at all are straightforward. For example, a DOPE potential might look something like: class DOPEPotential : public PairPotential { public: DOPEPotential(); ScoreFunc* create_score_func(Particle* p1, Particle* p2) { return new CubicSpline(p1->get_value("atom_type"), p2->get_value("atom_type")); } }; Any thoughts/problems? I can of course clarify if any of this is unclear. Ben -- ben(a)salilab.org http://salilab.org/~ben/ "It is a capital mistake to theorize before one has data." - Sir Arthur Conan Doyle

2 11

SWIG %import
by Daniel Russel 03 Jan '08

03 Jan '08

Apparently the %import directive is meant to be used with .i files rather than .h files. impEM isn't doing the right thing (nor was mine). I haven't yet checked if that solves my problems linking multiple swig libs as I have only had my mac which doesn't exhibit the problems. In addition, we are not properly sharing the swig runtime info between modules. Apparently we need to add " -DSWIG_TYPE_TABLE=IMP" to each compilation of a _wrap.cc so that the information is properly shared. Unfortunately, this includes EMbed (although EMbed doesn't have to actually see any of the IMP code as the IMP there is just used as a string). I would also propose changing IMP_exceptions.i as follows in order to allow other libraries to reuse the exception handling function by %including IMP_exceptions.i rather than redefine a new one: IMP_exceptions.i should contain: %{ void IMP_swig_handle_exception(void); %} %exception { try { $action } catch (...) { IMP_swig_handle_exception(); /* This should be unnecessary, since handle_imp_exception cannot return; here just to quell lots of warnings about the 'result' variable not being initialized. */ SWIG_fail; } } and putting the following in IMP.i to define the function uniquely. %{ /* Code to convert C++ exceptions into scripting language errors. Saves having lots of catch statements in every single wrapper. */ void IMP_swig_handle_exception(void) { try { throw; } catch (std::out_of_range &e) { SWIG_exception(SWIG_IndexError, e.what()); } catch (IMP::IndexException &e) { SWIG_exception(SWIG_IndexError, e.what()); } catch (IMP::InvalidStateException &e) { SWIG_exception(SWIG_ValueError, e.what()); } catch (IMP::ErrorException &e) { SWIG_exception(SWIG_RuntimeError, e.what()); } /* SWIG_exception contains "goto fail" so make sure the label is defined */ fail: return; } %} Here is a patch for everything in IMP (not for EMbed).

2 4

Removing names from restraints
by Daniel Russel 21 Dec '07

21 Dec '07

API changes: Restraint.set_name/Restraint.get_name went away because no one was using them. There are a couple of minor changes: - RestraintSet uses the IMP_CONTAINER macros. - RestraintSet uses an auto_ptr instead of catching all exceptions Index: test/xml/test_xml.py =================================================================== --- test/xml/test_xml.py (revision 664) +++ test/xml/test_xml.py (working copy) @@ -35,18 +35,6 @@ #restraint_sets = self.imp_model.get_restraints self.assertEqual(self.imp_model.number_of_restraints(), 3, "xml file contains three restraint sets") - self.assertEqual(self.imp_model.get_restraint(IMP.RestraintIndex(0)).get_name(), - 'exclusion_volumes', - "not expecting restraint set name: " - + self.imp_model.get_restraint(IMP.RestraintIndex(0)).get_name()) - self.assertEqual(self.imp_model.get_restraint(IMP.RestraintIndex(1)).get_name(), - 'torus', - "not expecting restraint set name: " - + self.imp_model.get_restraint(IMP.RestraintIndex(1)).get_name()) - self.assertEqual(self.imp_model.get_restraint(IMP.RestraintIndex(2)).get_name(), - 'connections', - "not expecting restraint set name: " - + self.imp_model.get_restraint(IMP.RestraintIndex(2)).get_name()) # test restraints score = self.imp_model.get_restraint(IMP.RestraintIndex(0)).evaluate(None) Index: include/IMP/utility.h =================================================================== --- include/IMP/utility.h (revision 664) +++ include/IMP/utility.h (working copy) @@ -164,7 +164,7 @@ \param[in] init Code to modify the passed in object. The object is obj its index index. */ -#define IMP_CONTAINER_IMPL(Class, Ucname, lcname, IndexType, init) \ +#define IMP_CONTAINER_IMPL(Class, Ucname, lcname, IndexType, init) \ IndexType Class::add_##lcname(Ucname *obj) { \ IndexType index(lcname##_vector_.size()); \ lcname##_vector_.push_back(obj); \ Index: include/IMP/restraints/RestraintSet.h =================================================================== --- include/IMP/restraints/RestraintSet.h (revision 664) +++ include/IMP/restraints/RestraintSet.h (working copy) @@ -25,32 +25,20 @@ IMP_RESTRAINT("0.5", "Daniel Russel") - //! The type to use to retrieve a restraint - typedef int RestraintIndex; + IMP_CONTAINER(Restraint, restraint, RestraintIndex); + public: + //! Set weight for all restraints contained by this set. - /** \param[in] weight The new value of the weight. + /** Setting the weight to 0 disables the restraints in the set. + + \param[in] weight The new value of the weight. */ void set_weight(Float weight) { weight_ = weight; } - //! Add restraint to the restraint set. - /** \param[in] restraint The restraint to add to the restraint set. - \return the index of the newly-added restraint in the restraint set. - */ - RestraintIndex add_restraint(Restraint *restraint); + //! Get weight for all restraints contained by this set. + Float get_weight() { return weight_; } - //! Access a restraint in the restraint set. - /** \param[in] i The RestraintIndex of the restraint to retrieve. - \exception std::out_of_range the index is out of range. - \return Pointer to the Restraint. - */ - Restraint *get_restraint(RestraintIndex i) const; - - //! Return the total number of restraints - unsigned int number_of_restraints() const { - return restraints_.size(); - } - //! Called when at least one particle has been inactivated. /** Check each restraint to see if it changes its active status. */ @@ -58,15 +46,10 @@ protected: - //! Restraints to evaluate. - /** These can be accessed with an iterator by a filter. - In the case where the restraint is a single simple restraint, - this vector contains a pointer to itself. - */ - std::vector<Restraint *> restraints_; - //! Weight for all restraints. Float weight_; + + std::string name_; }; } // namespace IMP Index: include/IMP/Restraint.h =================================================================== --- include/IMP/Restraint.h (revision 664) +++ include/IMP/Restraint.h (working copy) @@ -29,8 +29,8 @@ class IMPDLLEXPORT Restraint : public Object { public: - //! Initialize the Restraint and its model pointer - Restraint(std::string name=std::string()); + //! Initialize the Restraint + Restraint(); virtual ~Restraint(); //! Return the score for this restraint for the current state of the model. @@ -65,16 +65,6 @@ virtual std::string last_modified_by() const = 0; - //! Get the name of the restraint - const std::string& get_name() const { - return name_; - } - - //! Set the name of the restraint - void set_name(const std::string &name) { - name_=name; - } - //! The model the restraint is part of. /** \param[in] model The model. */ @@ -134,8 +124,6 @@ bool are_particles_active_; std::vector<Particle*> particles_; - - std::string name_; }; IMP_OUTPUT_OPERATOR(Restraint); Index: impEM/pyext/IMPEM.i =================================================================== --- impEM/pyext/IMPEM.i (revision 664) +++ impEM/pyext/IMPEM.i (working copy) @@ -28,6 +28,7 @@ } /* Get definitions of IMP base classes (but do not wrap; that is done by IMP) */ +%import "IMP/Object.h" %import "IMP/Restraint.h" /* Get definitions of EMLIB base classes (but do not wrap) */ Index: src/Model.cpp =================================================================== --- src/Model.cpp (revision 664) +++ src/Model.cpp (working copy) @@ -105,7 +105,7 @@ for (RestraintIterator it = restraints_begin(); it != restraints_end(); ++it) { IMP_CHECK_OBJECT(*it); - IMP_LOG(VERBOSE, (*it)->get_name() << ": " << std::flush); + IMP_LOG(VERBOSE, **it); Float tscore=0; if ((*it)->get_is_active()) { tscore = (*it)->evaluate(accpt); @@ -203,12 +203,6 @@ out << "version: " << version() << " "; out << "last_modified_by: " << last_modified_by() << std::endl; out << number_of_particles() << " particles" << std::endl; - out << "Restraints:" << std::endl; - for (RestraintConstIterator it = restraints_begin(); - it != restraints_end(); ++it) { - IMP_CHECK_OBJECT(*it); - out << (*it)->get_name() << std::endl; - } internal::show_attributes(out); } Index: src/Restraint.cpp =================================================================== --- src/Restraint.cpp (revision 664) +++ src/Restraint.cpp (working copy) @@ -15,7 +15,7 @@ namespace IMP { -Restraint::Restraint(std::string name): name_(name) +Restraint::Restraint() { model_ = NULL; is_active_ = true; // active by default @@ -57,8 +57,8 @@ out << "unknown restraint (inactive):" << std::endl; } - out << "version: " << version() << std::endl; - out << "last_modified_by: " << last_modified_by() << std::endl; + out << " version: " << version() << std::endl; + out << " last_modified_by: " << last_modified_by() << std::endl; } Index: src/restraints/RestraintSet.cpp =================================================================== --- src/restraints/RestraintSet.cpp (revision 664) +++ src/restraints/RestraintSet.cpp (working copy) @@ -10,6 +10,7 @@ #include "IMP/log.h" #include "IMP/restraints/RestraintSet.h" #include "../mystdexcept.h" +#include <utility> namespace IMP { @@ -17,12 +18,14 @@ //! Constructor RestraintSet::RestraintSet(const std::string& name) - : Restraint(name), weight_(1.0) + : weight_(1.0) { IMP_LOG(VERBOSE, "Restraint set constructed"); + name_=name; } + //! Destructor /** \todo Should reference count restraints correctly, to avoid deleting restraints here which live in two or more RestraintSets. @@ -30,39 +33,15 @@ RestraintSet::~RestraintSet() { IMP_LOG(VERBOSE, "Delete RestraintSet"); - for (unsigned int i = 0; i < restraints_.size(); ++i) { - restraints_[i]->set_model(NULL); - delete restraints_[i]; - } + IMP_CONTAINER_DELETE(Restraint, restraint); } +IMP_CONTAINER_IMPL(RestraintSet, Restraint, restraint, RestraintIndex, + obj->set_model(get_model());); -//! Add restraint to the restraint set. -/** \param[in] restraint The restraint to add to the restraint set. - \return the index of the newly-added restraint in the restraint set. - */ -RestraintSet::RestraintIndex RestraintSet::add_restraint(Restraint* restraint) -{ - restraints_.push_back(restraint); - restraint->set_model(get_model()); - return restraints_.size() - 1; -} -//! Access a restraint in the restraint set. -/** \param[in] i The RestraintIndex of the restraint to retrieve. - \exception std::out_of_range the index is out of range. - \return Pointer to the Restraint. - */ -Restraint* RestraintSet::get_restraint(RestraintIndex i) const -{ - IMP_check(static_cast<unsigned int>(i) < restraints_.size(), - "Invalid restraint requested", - std::out_of_range("Invalid restraint")); - return restraints_[i]; -} - //! Calculate the score for this restraint for the current model state. /** \param[in] accum If not NULL, use this object to accumulate partial first derivatives. @@ -70,25 +49,21 @@ */ Float RestraintSet::evaluate(DerivativeAccumulator *accum) { + if (get_weight() == 0) return 0; Float score; - + typedef std::auto_ptr<DerivativeAccumulator> DAP; // Use a local copy of the accumulator for our sub-restraints - DerivativeAccumulator *ouracc = NULL; + DAP ouracc; if (accum) { - ouracc = new DerivativeAccumulator(*accum, weight_); + ouracc = DAP(new DerivativeAccumulator(*accum, weight_)); } score = (Float) 0.0; - try { - for (size_t i = 0; i < restraints_.size(); i++) { - if (restraints_[i]->get_is_active()) - score += restraints_[i]->evaluate(ouracc); + for (RestraintIterator it= restraints_begin(); it != restraints_end(); ++it) { + if ((*it)->get_is_active()) { + score += (*it)->evaluate(ouracc.get()); } - } catch (...) { - delete ouracc; - throw; } - delete ouracc; return score * weight_; } @@ -99,8 +74,8 @@ */ void RestraintSet::check_particles_active() { - for (size_t i = 0; i < restraints_.size(); i++) { - restraints_[i]->check_particles_active(); + for (RestraintIterator it= restraints_begin(); it != restraints_end(); ++it) { + (*it)->check_particles_active(); } } @@ -111,11 +86,12 @@ */ void RestraintSet::show(std::ostream& out) const { - out << "restraint set " << get_name() << ":" << std::endl; - for (size_t i = 0; i < restraints_.size(); i++) { - restraints_[i]->show(out); + out << "restraint set " << name_ << ":..." << std::endl; + for (RestraintConstIterator it= restraints_begin(); + it != restraints_end(); ++it) { + (*it)->show(out); } - out << "... end restraint set " << get_name() << std::endl; + out << "... end restraint set " << name_ << std::endl; }

1 0

Several small patches
by Daniel Russel 20 Dec '07

20 Dec '07

First I added a pair of functions internal::assert_failure and internal::check_failure which get called when assertions and checks fail. You can have the debugger break in these to find failures. The derivative accumulator checks for NaN. And finally ConjugateGradients wasn't initializing its threshold_ member variable.

1 0

Restraint niggles
by Daniel Russel 20 Dec '07

20 Dec '07

No one seems to ever use the name on Restraints as opposed to restraint sets. I propose we get rid of the name field to make restraint lighter. I also propose we get rid of the active flag as it is not really used and can also be pushed into a restraint set if the functionality is desired. Alternatively, if we like names, it should be made non-optional.

2 1

Proposal: new OptimizerState class
by Ben Webb 20 Dec '07

20 Dec '07

We currently have a State class, which you can use to update some kind of shared state every time the score is calculated. You could also in principle use State to generate a trajectory during an optimization (e.g. see test/md_optimizer/test_md_optimizer.py). This works fine for optimizers such as molecular dynamics, where each step of the optimization has a single evaluation of the score function. But it's no good for optimizers like CG and SD, which may do multiple trial evaluations per step - you only want to include the accepted points in the trajectory. You could also imagine doing other things during an optimization - for example, periodic recentering of the system, resetting of thermostats, or rescaling of particle velocities. Thus, I propose an OptimizerState class, which is similar to State but is added to an optimizer rather than a model. Its update method would then be called at every successful optimization step, rather than at every score evaluation. To prevent people adding OptimizerStates to the model, the current State should be renamed ScoreState and both OptimizerState and ScoreState should derive from a common State base class. A derived MDOptimizerState or similar class could be used in turn as a base class for MD-specific stuff such as velocity scaling. Ben -- ben(a)salilab.org http://salilab.org/~ben/ "It is a capital mistake to theorize before one has data." - Sir Arthur Conan Doyle

2 2

Proposed bonded/nonbonded list interfaces
by Daniel Russel 20 Dec '07

20 Dec '07

Since I needed caching of my bonds so I decided to refactor things to have bonded and nonbonded lists. Here are my proposed interfaces for people's perusal. I then have a DynamicExclusionVolumeRestraint and DynamicBondRestraint which take pointers to these states and get their particle pairs from them. I don't know the best way to make things accessible from python (I don't think the C++ iterators will work from python, but I haven't tried). For reference, a BondDecorator is a decorator wrapping a particle which contains information (length, stiffness, type) for a bond. //! A pair of Particle * which can be used in maps struct ParticlePair { ParticlePair(Particle *a, Particle *b); Particle *get_particle(unsigned int i) const ; IMP_COMPARISONS_2(p_[0], p_[1]); }; //! A state that maintains a list of bonded particles. /** The state uses the BondedDecorator to find bonded pairs in its list of particles and builds a fast lookup structure on them. */ class IMPDLLEXPORT BondedListState : public State { public: BondedListState(const Particles &pis); virtual ~BondedListState(); IMP_STATE("0.5", "Daniel Russel"); void set_particles(const Particles &pis) ; typedef std::map<ParticlePair, BondDecorator>::const_iterator BondsIterator; BondsIterator bonds_begin() const ; BondsIterator bonds_end() const ; //! Get a bonds between two particles /** \returns BondDecorator() if there is not bond or the appropriate one if there is a bond. */ BondDecorator get_bond(Particle *a, Particle *b) const ; }; //! A state that maintains a list of nonbonded particles. /** */ class IMPDLLEXPORT NonbondedListState : public State { public: /** \note The distance_threshold parameter is just a guideline. Particle pairs with larger distance can be returned. But all with smaller distance will be. BondedListState will be updated when this one is and deleted when this State is and so should not be added seperately to the Model. */ NonbondedListState(const Particles &pis, Float distance_threshold, BondedList *bl=NULL); virtual ~NonbondedListState(); IMP_STATE("0.5", "Daniel Russel"); void set_particles(const Particles &pis); /*Eventually these iterators could generate the pairs on the fly from the voxel grid*/ typedef std::vector<ParticlePair>::const_iterator NonbondsIterator; NonbondsIterator nonbonds_begin() const ; NonbondsIterator nonbonds_end() const ; };

2 4

← Newer
1
...
80
81
82
83
84
85
86
87
Older →

Jump to page:

2024

2023

2022

2021

2020

2019

2018

2017

2016

2015

2014

2013

2012

2011

2010

2009

2008

2007

IMP-dev