IMP  2.2.1
The Integrative Modeling Platform
core/optimize_balls.py

This example optimizes a set of a balls to form 100 chains packed into a box. It illustrates using Monte Carlo (incremental) and conjugate gradients in conjunction in a non-trivial optimization.

1 ## \example core/optimize_balls.py
2 # This example optimizes a set of a balls to form 100 chains packed into a
3 # box. It illustrates using Monte Carlo (incremental) and conjugate
4 # gradients in conjunction in a non-trivial optimization.
5 
6 import IMP.core
7 import IMP.display
8 import IMP.container
9 import IMP.base
10 import sys
11 
12 IMP.base.setup_from_argv(sys.argv, "Optimize balls example")
13 
15  IMP.algebra.Vector3D(10, 10, 10))
17  ni = 2
18  nj = 2
19  np = 2
20  radius = .45
21  k = 100
22  ncg = 10
23  nmc = 1
24  ninner = 1
25  nouter = 1
26 else:
27  ni = 10
28  nj = 10
29  np = 10
30  radius = .45
31  k = 100
32  ncg = 1000
33  nmc = ni * nj * np * 100
34  ninner = 5
35  nouter = 11
36 
37 print IMP.base.get_is_quick_test(), ni, nj, np, ninner, nouter
38 # using a HarmonicDistancePairScore for fixed length links is more
39 # efficient than using a HarmonicSphereDistnacePairScore and works
40 # better with the optimizer
41 lps = IMP.core.HarmonicDistancePairScore(1.5 * radius, k)
43 
44 m = IMP.kernel.Model()
45 # IMP.base.set_log_level(IMP.base.SILENT)
46 aps = []
47 filters = []
48 movers = []
49 restraints = []
50 for i in range(0, ni):
51  for j in range(0, nj):
52  base = IMP.algebra.Vector3D(i, j, 0)
53  chain = []
54  for k in range(0, np):
55  p = IMP.kernel.Particle(m)
56  p.set_name("P" + str(i) + " " + str(j) + " " + str(k))
60  d.set_coordinates_are_optimized(True)
61  movers.append(IMP.core.BallMover([p], radius * 2))
62  movers[-1].set_was_used(True)
64  p, IMP.display.get_display_color(i * nj + j))
65  if k == 0:
66  d.set_coordinates(base)
67  else:
68  d.set_coordinates_are_optimized(True)
69  chain.append(p)
70  aps.append(p)
71  # set up a chain of bonds
73  r = IMP.container.PairsRestraint(lps, cpc)
74  restraints.append(r)
75 
76 # don't apply excluded volume to consecutive particles
80 bbr = IMP.container.SingletonsRestraint(ibss, aps)
81 restraints.append(bbr)
82 
84 mc = IMP.core.MonteCarlo(m)
85 mc.set_name("MC")
86 sm = IMP.core.SerialMover(movers)
87 mc.add_mover(sm)
88 # we are special casing the nbl term
89 isf = IMP.core.IncrementalScoringFunction(aps, restraints)
90 isf.set_name("I")
91 # use special incremental support for the non-bonded part
92 # apply the pair score sps to all touching ball pairs from the list of particles
93 # aps, using the filters to remove undersired pairs
94 # this is equivalent to the nbl construction above but optimized for
95 # incremental
96 isf.add_close_pair_score(sps, 0, aps, filters)
97 
98 # create a scoring function for conjugate gradients that includes the
99 # ExcludedVolumeRestraint
100 nbl = IMP.core.ExcludedVolumeRestraint(aps, k, 1)
101 nbl.set_pair_filters(filters)
102 sf = IMP.core.RestraintsScoringFunction(restraints + [nbl], "RSF")
103 
104 if True:
105  mc.set_incremental_scoring_function(isf)
106 else:
107  # we could, instead do non-incremental scoring
108  mc.set_scoring_function(sf)
109 
110 # first relax the bonds a bit
111 rs = []
112 for p in aps:
114  0))
115 cg.set_scoring_function(sf)
116 cg.optimize(ncg)
117 for r in restraints:
118  print r.get_name(), r.evaluate(False)
119 
120 # shrink each of the particles, relax the configuration, repeat
121 for i in range(1, nouter):
122  rs = []
123  factor = .1 * i
124  for p in aps:
125  rs.append(
127  IMP.core.XYZR(p).get_radius() * factor))
128  # move each particle nmc times
129  print factor
130  for j in range(0, ninner):
131  print "stage", j
132  mc.set_kt(100.0 / (3 * j + 1))
133  print "mc", mc.optimize((j + 1) * nmc), cg.optimize(nmc)
134  del rs
135  for r in restraints:
136  print r.get_name(), r.evaluate(False)
137 
138 w = IMP.display.PymolWriter("final.pym")
139 for p in aps:
140  g = IMP.core.XYZRGeometry(p)
141  w.add_geometry(g)
143 g.set_name("bb")
144 w.add_geometry(g)