IMP  2.4.0
The Integrative Modeling Platform
domino/multiscale.py

We are interested in applying domino to problems systematically in a multiscale manner. This script experiments with those approaches.

1 ## \example domino/multiscale.py
2 # We are interested in applying domino to problems systematically in a
3 # multiscale manner. This script experiments with those approaches.
4 
5 import IMP.domino
6 import IMP.core
7 
8 m = IMP.kernel.Model()
9 m.set_log_level(IMP.base.SILENT)
10 ds = [IMP.core.XYZR.setup_particle(IMP.kernel.Particle(m))
11  for i in range(0, 3)]
12 for i, d in enumerate(ds):
13  d.set_radius(1)
14  IMP.display.Colored.setup_particle(d, IMP.display.get_display_color(i))
15 
16 k = 1
17 h = IMP.core.Harmonic(0, k)
18 r0 = IMP.core.SingletonRestraint(
19  IMP.core.DistanceToSingletonScore(h, IMP.algebra.Vector3D(0, 0, 0)),
20  ds[0], "0 at origin")
21 
22 r1 = IMP.core.SingletonRestraint(IMP.core.AttributeSingletonScore(h,
23  IMP.core.XYZ.get_xyz_keys(
24  )[0]),
25  ds[1], "1 on axis")
26 
27 rs = [r0, r1]
28 for pr in [(0, 1), (1, 2), (0, 2)]:
29  r = IMP.core.PairRestraint(IMP.core.HarmonicSphereDistancePairScore(0, k),
30  (ds[pr[0]], ds[pr[1]]),
31  "R for " + str(pr))
32  rs.append(r)
33 
34 
35 bb = IMP.algebra.BoundingBox2D(IMP.algebra.Vector2D(0, 0),
36  IMP.algebra.Vector2D(4, 4))
37 
38 covers = []
39 for i in range(0, 6):
40  cur = IMP.algebra.get_grid_interior_cover_by_spacing(bb, 4.0 / 2 ** i)
41  print(cur)
42  covers.append([IMP.algebra.Vector3D(x[0], x[1], 0) for x in cur])
43 
44 
45 def setup(cover, scale):
46  pst = IMP.domino.ParticleStatesTable()
47  st = IMP.domino.XYZStates(cover)
48  for p in ds:
49  pst.set_particle_states(p, st)
50  for r in rs:
51  r.set_maximum_score(.5 * scale ** 2)
52  lf = IMP.domino.ListSubsetFilterTable(pst)
53  rc = IMP.domino.RestraintCache(pst)
54  rc.add_restraints(rs)
55  fs = [IMP.domino.RestraintScoreSubsetFilterTable(rc),
56  lf]
57  sampler = IMP.domino.DominoSampler(m, pst)
58  sampler.set_subset_filter_tables(fs)
59  sampler.set_log_level(IMP.base.SILENT)
60  return (sampler, lf, pst)
61 
62 (sampler, lf, pst) = setup(covers[0], 4.0)
63 
64 subset = IMP.domino.Subset(ds)
65 ac = sampler.get_sample_assignments(subset)
66 
67 print(ac)
68 
69 
70 def get_mapping(cover0, cover1):
71  nn = IMP.algebra.NearestNeighbor3D(cover0)
72  ret = [[] for c in cover0]
73  for i, p in enumerate(cover1):
74  nns = nn.get_nearest_neighbor(p)
75  ret[nns].append(i)
76  return ret
77 
78 mw = IMP.display.PymolWriter("mapping.pym")
79 
80 
81 def display_mapping(index, cover0, cover1, mapping):
82  mw.set_frame(index)
83  for i, c in enumerate(mapping):
84  for p in c:
85  g = IMP.display.PointGeometry(cover1[p])
86  g.set_color(IMP.display.get_display_color(i))
87  g.set_name("fine")
88  mw.add_geometry(g)
89  for i, c in enumerate(cover0):
90  g = IMP.display.PointGeometry(c)
91  g.set_color(IMP.display.get_display_color(i))
92  g.set_name("coarse")
93  mw.add_geometry(g)
94 
95 
96 for curi in range(1, len(covers)):
97  scale = 4.0 / 2 ** curi
98  print(scale)
99  mapping = get_mapping(covers[curi - 1], covers[curi])
100  print(mapping)
101  display_mapping(curi - 1, covers[curi - 1], covers[curi], mapping)
102  (sampler, lf, pst) = setup(covers[curi], scale)
103  lac = ac
104  cac = []
105  for a in lac:
106  for i, p in enumerate(subset):
107  s = a[i]
108  allowed = mapping[s]
109  lf.set_allowed_states(p, allowed)
110  ccac = sampler.get_sample_assignments(subset)
111  print(a, ccac)
112  cac = cac + ccac
113  ac = list(set(cac))
114  print("for scale", scale, "got", ac)
115  sw = IMP.display.PymolWriter("solutions." + str(curi) + ".pym")
116  for i, a in enumerate(ac):
117  IMP.domino.load_particle_states(subset, a, pst)
118  sw.set_frame(i)
119  for p in ds:
120  g = IMP.core.XYZRGeometry(p)
121  sw.add_geometry(g)
122  for c in covers[curi]:
123  g = IMP.display.PointGeometry(c)
124  g.set_color(IMP.display.Color(1, 1, 1))
125  g.set_name("grid")
126  sw.add_geometry(g)