IMP  2.0.1
The Integrative Modeling Platform
domino/six_particles_optimization.py

Optimize six particles on a 2D unit grid. In order to remove translation degrees of freedom, the 0th particle is pinned at the origin by allowing it only a single conformation. To remove flips, the first particle is restrained to have a positive x coordinate.

1 ## \example domino/six_particles_optimization.py
2 ## Optimize six particles on a 2D unit grid. In order to remove translation degrees
3 ## of freedom, the 0th particle is pinned at the origin by allowing it only a
4 ## single conformation. To remove flips, the first particle is restrained to
5 ## have a positive x coordinate.
6 
7 import IMP
8 import IMP.domino
9 import IMP.core
10 import IMP.container
11 
12 #set restraints
13 def create_scoring(m, ps):
14  pairs=[[0,1],[0,2],[1,3],[2,3],[3,4],[4,5],[1,5]]
15  # we will restrain various pairs to be 1 apart
17  # the restraint will be broken apart during optimization
18  # map the indices above to actual particles
19  pc= IMP.container.ListPairContainer([(ps[p[0]], ps[p[1]]) for p in pairs],
20  "Restrained pairs")
21  pr= IMP.container.PairsRestraint(score, pc)
22  pr.set_maximum_score(.01)
23  pr.set_model(m)
25  IMP.algebra.Vector3D(2,0,0))
26  # force ps[1] to be on the positive side to remove flip degree of freedom
27  dr= IMP.core.SingletonRestraint(d, ps[1])
28  dr.set_model(m)
29  # we are not interested in conformations which don't fit the distances
30  # exactly, but using 0 is tricky
31  dr.set_maximum_score(.01)
32  print m.get_root_restraint_set()
33  return [pr, dr]
34 
35 def create_representation(m):
36  ps=[]
37  # create size particles, initial coordinates don't matter.
38  for i in range(0,6):
39  p=IMP.Particle(m)
41  ps.append(p)
42  return ps
43 
44 def create_discrete_states(ps):
46  vs=[IMP.algebra.Vector3D(1,0,0),
47  IMP.algebra.Vector3D(0,1,0),
48  IMP.algebra.Vector3D(1,1,0),
49  IMP.algebra.Vector3D(2,1,0),
50  IMP.algebra.Vector3D(2,0,0)]
51  vs= vs+[-v for v in vs]
52  print len(vs), "states for each particle"
53  states= IMP.domino.XYZStates(vs)
54  # special case ps[0] to remove a sliding degree of freedom
55  # all other particles are given the same set of states
56  for p in ps[1:]:
57  pst.set_particle_states(p, states)
58  return pst
59 
60 def create_sampler(m, r, pst):
61  # create the sampler and pass it the states for each patricle
62  s=IMP.domino.DominoSampler(m, pst)
63  # the following lines recreate the defaults and so are optional
64  filters=[]
65  # create a restraint cache to avoid re-evaluating restraints
67  # add the list of restraints we want to use
68  rc.add_restraints(r)
69  # do not allow particles with the same ParticleStates object
70  # to have the same state index
71  filters.append(IMP.domino.ExclusionSubsetFilterTable(pst))
72  # filter states that score worse than the cutoffs in the Model
74  filters[-1].set_log_level(IMP.base.SILENT)
75  # try to be intelligent about enumerating the states in each subset
76  states= IMP.domino.BranchAndBoundAssignmentsTable(pst, filters);
77  states.set_log_level(IMP.base.SILENT);
78  s.set_assignments_table(states)
79  s.set_subset_filter_tables(filters)
80 
81  return s
82 
83 IMP.base.set_log_level(IMP.base.TERSE)
84 m=IMP.Model()
85 # don't print information during Model.evaluate
86 m.set_log_level(IMP.base.SILENT)
87 
88 print "creating representation"
89 ps=create_representation(m)
90 print "creating discrete states"
91 pst=create_discrete_states(ps)
92 print "creating score function"
93 rs=create_scoring(m, ps)
94 print "creating sampler"
95 s=create_sampler(m, rs, pst)
96 
97 print "sampling"
98 # get an IMP.ConfigurationSet with the sampled states. If there are very
99 # many, it might be better to use s.get_sample_states() and then
100 # IMP.domino.load_particle_states() to handle the states as that takes
101 # much less memory, and time.
102 cs=s.get_sample()
103 
104 print "found ", cs.get_number_of_configurations(), "solutions"
105 for i in range(cs.get_number_of_configurations()):
106  cs.load_configuration(i)
107  print "solution number:",i," is:", m.evaluate(False)
108  for p in ps:
109  print IMP.core.XYZ(p).get_x(), IMP.core.XYZ(p).get_y()