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