Shortcuts on this page
r m x toggle line displays
j k next/prev highlighted chunk
0 (zero) top of page
1 (one) first highlighted chunk
1"""
2ISMAGS Algorithm
3================
5Provides a Python implementation of the ISMAGS algorithm. [1]_
7It is capable of finding (subgraph) isomorphisms between two graphs, taking the
8symmetry of the subgraph into account. In most cases the VF2 algorithm is
9faster (at least on small graphs) than this implementation, but in some cases
10there is an exponential number of isomorphisms that are symmetrically
11equivalent. In that case, the ISMAGS algorithm will provide only one solution
12per symmetry group.
14>>> petersen = nx.petersen_graph()
15>>> ismags = nx.isomorphism.ISMAGS(petersen, petersen)
16>>> isomorphisms = list(ismags.isomorphisms_iter(symmetry=False))
17>>> len(isomorphisms)
18120
19>>> isomorphisms = list(ismags.isomorphisms_iter(symmetry=True))
20>>> answer = [{0: 0, 1: 1, 2: 2, 3: 3, 4: 4, 5: 5, 6: 6, 7: 7, 8: 8, 9: 9}]
21>>> answer == isomorphisms
22True
24In addition, this implementation also provides an interface to find the
25largest common induced subgraph [2]_ between any two graphs, again taking
26symmetry into account. Given `graph` and `subgraph` the algorithm will remove
27nodes from the `subgraph` until `subgraph` is isomorphic to a subgraph of
28`graph`. Since only the symmetry of `subgraph` is taken into account it is
29worth thinking about how you provide your graphs:
31>>> graph1 = nx.path_graph(4)
32>>> graph2 = nx.star_graph(3)
33>>> ismags = nx.isomorphism.ISMAGS(graph1, graph2)
34>>> ismags.is_isomorphic()
35False
36>>> largest_common_subgraph = list(ismags.largest_common_subgraph())
37>>> answer = [{1: 0, 0: 1, 2: 2}, {2: 0, 1: 1, 3: 2}]
38>>> answer == largest_common_subgraph
39True
40>>> ismags2 = nx.isomorphism.ISMAGS(graph2, graph1)
41>>> largest_common_subgraph = list(ismags2.largest_common_subgraph())
42>>> answer = [
43... {1: 0, 0: 1, 2: 2},
44... {1: 0, 0: 1, 3: 2},
45... {2: 0, 0: 1, 1: 2},
46... {2: 0, 0: 1, 3: 2},
47... {3: 0, 0: 1, 1: 2},
48... {3: 0, 0: 1, 2: 2},
49... ]
50>>> answer == largest_common_subgraph
51True
53However, when not taking symmetry into account, it doesn't matter:
55>>> largest_common_subgraph = list(ismags.largest_common_subgraph(symmetry=False))
56>>> answer = [
57... {1: 0, 0: 1, 2: 2},
58... {1: 0, 2: 1, 0: 2},
59... {2: 0, 1: 1, 3: 2},
60... {2: 0, 3: 1, 1: 2},
61... {1: 0, 0: 1, 2: 3},
62... {1: 0, 2: 1, 0: 3},
63... {2: 0, 1: 1, 3: 3},
64... {2: 0, 3: 1, 1: 3},
65... {1: 0, 0: 2, 2: 3},
66... {1: 0, 2: 2, 0: 3},
67... {2: 0, 1: 2, 3: 3},
68... {2: 0, 3: 2, 1: 3},
69... ]
70>>> answer == largest_common_subgraph
71True
72>>> largest_common_subgraph = list(ismags2.largest_common_subgraph(symmetry=False))
73>>> answer = [
74... {1: 0, 0: 1, 2: 2},
75... {1: 0, 0: 1, 3: 2},
76... {2: 0, 0: 1, 1: 2},
77... {2: 0, 0: 1, 3: 2},
78... {3: 0, 0: 1, 1: 2},
79... {3: 0, 0: 1, 2: 2},
80... {1: 1, 0: 2, 2: 3},
81... {1: 1, 0: 2, 3: 3},
82... {2: 1, 0: 2, 1: 3},
83... {2: 1, 0: 2, 3: 3},
84... {3: 1, 0: 2, 1: 3},
85... {3: 1, 0: 2, 2: 3},
86... ]
87>>> answer == largest_common_subgraph
88True
90Notes
91-----
92- The current implementation works for undirected graphs only. The algorithm
93 in general should work for directed graphs as well though.
94- Node keys for both provided graphs need to be fully orderable as well as
95 hashable.
96- Node and edge equality is assumed to be transitive: if A is equal to B, and
97 B is equal to C, then A is equal to C.
99References
100----------
101.. [1] M. Houbraken, S. Demeyer, T. Michoel, P. Audenaert, D. Colle,
102 M. Pickavet, "The Index-Based Subgraph Matching Algorithm with General
103 Symmetries (ISMAGS): Exploiting Symmetry for Faster Subgraph
104 Enumeration", PLoS One 9(5): e97896, 2014.
105 https://doi.org/10.1371/journal.pone.0097896
106.. [2] https://en.wikipedia.org/wiki/Maximum_common_induced_subgraph
107"""
109__all__ = ["ISMAGS"]
111import itertools
112from collections import Counter, defaultdict
113from functools import reduce, wraps
116def are_all_equal(iterable):
117 """
118 Returns ``True`` if and only if all elements in `iterable` are equal; and
119 ``False`` otherwise.
121 Parameters
122 ----------
123 iterable: collections.abc.Iterable
124 The container whose elements will be checked.
126 Returns
127 -------
128 bool
129 ``True`` iff all elements in `iterable` compare equal, ``False``
130 otherwise.
131 """
132 try:
133 shape = iterable.shape
134 except AttributeError:
135 pass
136 else:
137 if len(shape) > 1:
138 message = "The function does not works on multidimensional arrays."
139 raise NotImplementedError(message) from None
141 iterator = iter(iterable)
142 first = next(iterator, None)
143 return all(item == first for item in iterator)
146def make_partitions(items, test):
147 """
148 Partitions items into sets based on the outcome of ``test(item1, item2)``.
149 Pairs of items for which `test` returns `True` end up in the same set.
151 Parameters
152 ----------
153 items : collections.abc.Iterable[collections.abc.Hashable]
154 Items to partition
155 test : collections.abc.Callable[collections.abc.Hashable, collections.abc.Hashable]
156 A function that will be called with 2 arguments, taken from items.
157 Should return `True` if those 2 items need to end up in the same
158 partition, and `False` otherwise.
160 Returns
161 -------
162 list[set]
163 A list of sets, with each set containing part of the items in `items`,
164 such that ``all(test(*pair) for pair in itertools.combinations(set, 2))
165 == True``
167 Notes
168 -----
169 The function `test` is assumed to be transitive: if ``test(a, b)`` and
170 ``test(b, c)`` return ``True``, then ``test(a, c)`` must also be ``True``.
171 """
172 partitions = []
173 for item in items:
174 for partition in partitions:
175 p_item = next(iter(partition))
176 if test(item, p_item):
177 partition.add(item)
178 break
179 else: # No break
180 partitions.append({item})
181 return partitions
184def partition_to_color(partitions):
185 """
186 Creates a dictionary that maps each item in each partition to the index of
187 the partition to which it belongs.
189 Parameters
190 ----------
191 partitions: collections.abc.Sequence[collections.abc.Iterable]
192 As returned by :func:`make_partitions`.
194 Returns
195 -------
196 dict
197 """
198 colors = {}
199 for color, keys in enumerate(partitions):
200 for key in keys:
201 colors[key] = color
202 return colors
205def intersect(collection_of_sets):
206 """
207 Given an collection of sets, returns the intersection of those sets.
209 Parameters
210 ----------
211 collection_of_sets: collections.abc.Collection[set]
212 A collection of sets.
214 Returns
215 -------
216 set
217 An intersection of all sets in `collection_of_sets`. Will have the same
218 type as the item initially taken from `collection_of_sets`.
219 """
220 collection_of_sets = list(collection_of_sets)
221 first = collection_of_sets.pop()
222 out = reduce(set.intersection, collection_of_sets, set(first))
223 return type(first)(out)
226class ISMAGS:
227 """
228 Implements the ISMAGS subgraph matching algorithm. [1]_ ISMAGS stands for
229 "Index-based Subgraph Matching Algorithm with General Symmetries". As the
230 name implies, it is symmetry aware and will only generate non-symmetric
231 isomorphisms.
233 Notes
234 -----
235 The implementation imposes additional conditions compared to the VF2
236 algorithm on the graphs provided and the comparison functions
237 (:attr:`node_equality` and :attr:`edge_equality`):
239 - Node keys in both graphs must be orderable as well as hashable.
240 - Equality must be transitive: if A is equal to B, and B is equal to C,
241 then A must be equal to C.
243 Attributes
244 ----------
245 graph: networkx.Graph
246 subgraph: networkx.Graph
247 node_equality: collections.abc.Callable
248 The function called to see if two nodes should be considered equal.
249 It's signature looks like this:
250 ``f(graph1: networkx.Graph, node1, graph2: networkx.Graph, node2) -> bool``.
251 `node1` is a node in `graph1`, and `node2` a node in `graph2`.
252 Constructed from the argument `node_match`.
253 edge_equality: collections.abc.Callable
254 The function called to see if two edges should be considered equal.
255 It's signature looks like this:
256 ``f(graph1: networkx.Graph, edge1, graph2: networkx.Graph, edge2) -> bool``.
257 `edge1` is an edge in `graph1`, and `edge2` an edge in `graph2`.
258 Constructed from the argument `edge_match`.
260 References
261 ----------
262 .. [1] M. Houbraken, S. Demeyer, T. Michoel, P. Audenaert, D. Colle,
263 M. Pickavet, "The Index-Based Subgraph Matching Algorithm with General
264 Symmetries (ISMAGS): Exploiting Symmetry for Faster Subgraph
265 Enumeration", PLoS One 9(5): e97896, 2014.
266 https://doi.org/10.1371/journal.pone.0097896
267 """
269 def __init__(self, graph, subgraph, node_match=None, edge_match=None, cache=None):
270 """
271 Parameters
272 ----------
273 graph: networkx.Graph
274 subgraph: networkx.Graph
275 node_match: collections.abc.Callable or None
276 Function used to determine whether two nodes are equivalent. Its
277 signature should look like ``f(n1: dict, n2: dict) -> bool``, with
278 `n1` and `n2` node property dicts. See also
279 :func:`~networkx.algorithms.isomorphism.categorical_node_match` and
280 friends.
281 If `None`, all nodes are considered equal.
282 edge_match: collections.abc.Callable or None
283 Function used to determine whether two edges are equivalent. Its
284 signature should look like ``f(e1: dict, e2: dict) -> bool``, with
285 `e1` and `e2` edge property dicts. See also
286 :func:`~networkx.algorithms.isomorphism.categorical_edge_match` and
287 friends.
288 If `None`, all edges are considered equal.
289 cache: collections.abc.Mapping
290 A cache used for caching graph symmetries.
291 """
292 # TODO: graph and subgraph setter methods that invalidate the caches.
293 # TODO: allow for precomputed partitions and colors
294 self.graph = graph
295 self.subgraph = subgraph
296 self._symmetry_cache = cache
297 # Naming conventions are taken from the original paper. For your
298 # sanity:
299 # sg: subgraph
300 # g: graph
301 # e: edge(s)
302 # n: node(s)
303 # So: sgn means "subgraph nodes".
304 self._sgn_partitions_ = None
305 self._sge_partitions_ = None
307 self._sgn_colors_ = None
308 self._sge_colors_ = None
310 self._gn_partitions_ = None
311 self._ge_partitions_ = None
313 self._gn_colors_ = None
314 self._ge_colors_ = None
316 self._node_compat_ = None
317 self._edge_compat_ = None
319 if node_match is None:
320 self.node_equality = self._node_match_maker(lambda n1, n2: True)
321 self._sgn_partitions_ = [set(self.subgraph.nodes)]
322 self._gn_partitions_ = [set(self.graph.nodes)]
323 self._node_compat_ = {0: 0}
324 else:
325 self.node_equality = self._node_match_maker(node_match)
326 if edge_match is None:
327 self.edge_equality = self._edge_match_maker(lambda e1, e2: True)
328 self._sge_partitions_ = [set(self.subgraph.edges)]
329 self._ge_partitions_ = [set(self.graph.edges)]
330 self._edge_compat_ = {0: 0}
331 else:
332 self.edge_equality = self._edge_match_maker(edge_match)
334 @property
335 def _sgn_partitions(self):
336 if self._sgn_partitions_ is None:
338 def nodematch(node1, node2):
339 return self.node_equality(self.subgraph, node1, self.subgraph, node2)
341 self._sgn_partitions_ = make_partitions(self.subgraph.nodes, nodematch)
342 return self._sgn_partitions_
344 @property
345 def _sge_partitions(self):
346 if self._sge_partitions_ is None:
348 def edgematch(edge1, edge2):
349 return self.edge_equality(self.subgraph, edge1, self.subgraph, edge2)
351 self._sge_partitions_ = make_partitions(self.subgraph.edges, edgematch)
352 return self._sge_partitions_
354 @property
355 def _gn_partitions(self):
356 if self._gn_partitions_ is None:
358 def nodematch(node1, node2):
359 return self.node_equality(self.graph, node1, self.graph, node2)
361 self._gn_partitions_ = make_partitions(self.graph.nodes, nodematch)
362 return self._gn_partitions_
364 @property
365 def _ge_partitions(self):
366 if self._ge_partitions_ is None:
368 def edgematch(edge1, edge2):
369 return self.edge_equality(self.graph, edge1, self.graph, edge2)
371 self._ge_partitions_ = make_partitions(self.graph.edges, edgematch)
372 return self._ge_partitions_
374 @property
375 def _sgn_colors(self):
376 if self._sgn_colors_ is None:
377 self._sgn_colors_ = partition_to_color(self._sgn_partitions)
378 return self._sgn_colors_
380 @property
381 def _sge_colors(self):
382 if self._sge_colors_ is None:
383 self._sge_colors_ = partition_to_color(self._sge_partitions)
384 return self._sge_colors_
386 @property
387 def _gn_colors(self):
388 if self._gn_colors_ is None:
389 self._gn_colors_ = partition_to_color(self._gn_partitions)
390 return self._gn_colors_
392 @property
393 def _ge_colors(self):
394 if self._ge_colors_ is None:
395 self._ge_colors_ = partition_to_color(self._ge_partitions)
396 return self._ge_colors_
398 @property
399 def _node_compatibility(self):
400 if self._node_compat_ is not None:
401 return self._node_compat_
402 self._node_compat_ = {}
403 for sgn_part_color, gn_part_color in itertools.product(
404 range(len(self._sgn_partitions)), range(len(self._gn_partitions))
405 ):
406 sgn = next(iter(self._sgn_partitions[sgn_part_color]))
407 gn = next(iter(self._gn_partitions[gn_part_color]))
408 if self.node_equality(self.subgraph, sgn, self.graph, gn):
409 self._node_compat_[sgn_part_color] = gn_part_color
410 return self._node_compat_
412 @property
413 def _edge_compatibility(self):
414 if self._edge_compat_ is not None:
415 return self._edge_compat_
416 self._edge_compat_ = {}
417 for sge_part_color, ge_part_color in itertools.product(
418 range(len(self._sge_partitions)), range(len(self._ge_partitions))
419 ):
420 sge = next(iter(self._sge_partitions[sge_part_color]))
421 ge = next(iter(self._ge_partitions[ge_part_color]))
422 if self.edge_equality(self.subgraph, sge, self.graph, ge):
423 self._edge_compat_[sge_part_color] = ge_part_color
424 return self._edge_compat_
426 @staticmethod
427 def _node_match_maker(cmp):
428 @wraps(cmp)
429 def comparer(graph1, node1, graph2, node2):
430 return cmp(graph1.nodes[node1], graph2.nodes[node2])
432 return comparer
434 @staticmethod
435 def _edge_match_maker(cmp):
436 @wraps(cmp)
437 def comparer(graph1, edge1, graph2, edge2):
438 return cmp(graph1.edges[edge1], graph2.edges[edge2])
440 return comparer
442 def find_isomorphisms(self, symmetry=True):
443 """Find all subgraph isomorphisms between subgraph and graph
445 Finds isomorphisms where :attr:`subgraph` <= :attr:`graph`.
447 Parameters
448 ----------
449 symmetry: bool
450 Whether symmetry should be taken into account. If False, found
451 isomorphisms may be symmetrically equivalent.
453 Yields
454 ------
455 dict
456 The found isomorphism mappings of {graph_node: subgraph_node}.
457 """
458 # The networkx VF2 algorithm is slightly funny in when it yields an
459 # empty dict and when not.
460 if not self.subgraph:
461 yield {}
462 return
463 elif not self.graph:
464 return
465 elif len(self.graph) < len(self.subgraph):
466 return
468 if symmetry:
469 _, cosets = self.analyze_symmetry(
470 self.subgraph, self._sgn_partitions, self._sge_colors
471 )
472 constraints = self._make_constraints(cosets)
473 else:
474 constraints = []
476 candidates = self._find_nodecolor_candidates()
477 la_candidates = self._get_lookahead_candidates()
478 for sgn in self.subgraph:
479 extra_candidates = la_candidates[sgn]
480 if extra_candidates:
481 candidates[sgn] = candidates[sgn] | {frozenset(extra_candidates)}
483 if any(candidates.values()):
484 start_sgn = min(candidates, key=lambda n: min(candidates[n], key=len))
485 candidates[start_sgn] = (intersect(candidates[start_sgn]),)
486 yield from self._map_nodes(start_sgn, candidates, constraints)
487 else:
488 return
490 @staticmethod
491 def _find_neighbor_color_count(graph, node, node_color, edge_color):
492 """
493 For `node` in `graph`, count the number of edges of a specific color
494 it has to nodes of a specific color.
495 """
496 counts = Counter()
497 neighbors = graph[node]
498 for neighbor in neighbors:
499 n_color = node_color[neighbor]
500 if (node, neighbor) in edge_color:
501 e_color = edge_color[node, neighbor]
502 else:
503 e_color = edge_color[neighbor, node]
504 counts[e_color, n_color] += 1
505 return counts
507 def _get_lookahead_candidates(self):
508 """
509 Returns a mapping of {subgraph node: collection of graph nodes} for
510 which the graph nodes are feasible candidates for the subgraph node, as
511 determined by looking ahead one edge.
512 """
513 g_counts = {}
514 for gn in self.graph:
515 g_counts[gn] = self._find_neighbor_color_count(
516 self.graph, gn, self._gn_colors, self._ge_colors
517 )
518 candidates = defaultdict(set)
519 for sgn in self.subgraph:
520 sg_count = self._find_neighbor_color_count(
521 self.subgraph, sgn, self._sgn_colors, self._sge_colors
522 )
523 new_sg_count = Counter()
524 for (sge_color, sgn_color), count in sg_count.items():
525 try:
526 ge_color = self._edge_compatibility[sge_color]
527 gn_color = self._node_compatibility[sgn_color]
528 except KeyError:
529 pass
530 else:
531 new_sg_count[ge_color, gn_color] = count
533 for gn, g_count in g_counts.items():
534 if all(new_sg_count[x] <= g_count[x] for x in new_sg_count):
535 # Valid candidate
536 candidates[sgn].add(gn)
537 return candidates
539 def largest_common_subgraph(self, symmetry=True):
540 """
541 Find the largest common induced subgraphs between :attr:`subgraph` and
542 :attr:`graph`.
544 Parameters
545 ----------
546 symmetry: bool
547 Whether symmetry should be taken into account. If False, found
548 largest common subgraphs may be symmetrically equivalent.
550 Yields
551 ------
552 dict
553 The found isomorphism mappings of {graph_node: subgraph_node}.
554 """
555 # The networkx VF2 algorithm is slightly funny in when it yields an
556 # empty dict and when not.
557 if not self.subgraph:
558 yield {}
559 return
560 elif not self.graph:
561 return
563 if symmetry:
564 _, cosets = self.analyze_symmetry(
565 self.subgraph, self._sgn_partitions, self._sge_colors
566 )
567 constraints = self._make_constraints(cosets)
568 else:
569 constraints = []
571 candidates = self._find_nodecolor_candidates()
573 if any(candidates.values()):
574 yield from self._largest_common_subgraph(candidates, constraints)
575 else:
576 return
578 def analyze_symmetry(self, graph, node_partitions, edge_colors):
579 """
580 Find a minimal set of permutations and corresponding co-sets that
581 describe the symmetry of `graph`, given the node and edge equalities
582 given by `node_partitions` and `edge_colors`, respectively.
584 Parameters
585 ----------
586 graph : networkx.Graph
587 The graph whose symmetry should be analyzed.
588 node_partitions : list of sets
589 A list of sets containing node keys. Node keys in the same set
590 are considered equivalent. Every node key in `graph` should be in
591 exactly one of the sets. If all nodes are equivalent, this should
592 be ``[set(graph.nodes)]``.
593 edge_colors : dict mapping edges to their colors
594 A dict mapping every edge in `graph` to its corresponding color.
595 Edges with the same color are considered equivalent. If all edges
596 are equivalent, this should be ``{e: 0 for e in graph.edges}``.
599 Returns
600 -------
601 set[frozenset]
602 The found permutations. This is a set of frozensets of pairs of node
603 keys which can be exchanged without changing :attr:`subgraph`.
604 dict[collections.abc.Hashable, set[collections.abc.Hashable]]
605 The found co-sets. The co-sets is a dictionary of
606 ``{node key: set of node keys}``.
607 Every key-value pair describes which ``values`` can be interchanged
608 without changing nodes less than ``key``.
609 """
610 if self._symmetry_cache is not None:
611 key = hash(
612 (
613 tuple(graph.nodes),
614 tuple(graph.edges),
615 tuple(map(tuple, node_partitions)),
616 tuple(edge_colors.items()),
617 )
618 )
619 if key in self._symmetry_cache:
620 return self._symmetry_cache[key]
621 node_partitions = list(
622 self._refine_node_partitions(graph, node_partitions, edge_colors)
623 )
624 assert len(node_partitions) == 1
625 node_partitions = node_partitions[0]
626 permutations, cosets = self._process_ordered_pair_partitions(
627 graph, node_partitions, node_partitions, edge_colors
628 )
629 if self._symmetry_cache is not None:
630 self._symmetry_cache[key] = permutations, cosets
631 return permutations, cosets
633 def is_isomorphic(self, symmetry=False):
634 """
635 Returns True if :attr:`graph` is isomorphic to :attr:`subgraph` and
636 False otherwise.
638 Returns
639 -------
640 bool
641 """
642 return len(self.subgraph) == len(self.graph) and self.subgraph_is_isomorphic(
643 symmetry
644 )
646 def subgraph_is_isomorphic(self, symmetry=False):
647 """
648 Returns True if a subgraph of :attr:`graph` is isomorphic to
649 :attr:`subgraph` and False otherwise.
651 Returns
652 -------
653 bool
654 """
655 # symmetry=False, since we only need to know whether there is any
656 # example; figuring out all symmetry elements probably costs more time
657 # than it gains.
658 isom = next(self.subgraph_isomorphisms_iter(symmetry=symmetry), None)
659 return isom is not None
661 def isomorphisms_iter(self, symmetry=True):
662 """
663 Does the same as :meth:`find_isomorphisms` if :attr:`graph` and
664 :attr:`subgraph` have the same number of nodes.
665 """
666 if len(self.graph) == len(self.subgraph):
667 yield from self.subgraph_isomorphisms_iter(symmetry=symmetry)
669 def subgraph_isomorphisms_iter(self, symmetry=True):
670 """Alternative name for :meth:`find_isomorphisms`."""
671 return self.find_isomorphisms(symmetry)
673 def _find_nodecolor_candidates(self):
674 """
675 Per node in subgraph find all nodes in graph that have the same color.
676 """
677 candidates = defaultdict(set)
678 for sgn in self.subgraph.nodes:
679 sgn_color = self._sgn_colors[sgn]
680 if sgn_color in self._node_compatibility:
681 gn_color = self._node_compatibility[sgn_color]
682 candidates[sgn].add(frozenset(self._gn_partitions[gn_color]))
683 else:
684 candidates[sgn].add(frozenset())
685 candidates = dict(candidates)
686 for sgn, options in candidates.items():
687 candidates[sgn] = frozenset(options)
688 return candidates
690 @staticmethod
691 def _make_constraints(cosets):
692 """
693 Turn cosets into constraints.
694 """
695 constraints = []
696 for node_i, node_ts in cosets.items():
697 for node_t in node_ts:
698 if node_i != node_t:
699 # Node i must be smaller than node t.
700 constraints.append((node_i, node_t))
701 return constraints
703 @staticmethod
704 def _find_node_edge_color(graph, node_colors, edge_colors):
705 """
706 For every node in graph, come up with a color that combines 1) the
707 color of the node, and 2) the number of edges of a color to each type
708 of node.
709 """
710 counts = defaultdict(lambda: defaultdict(int))
711 for node1, node2 in graph.edges:
712 if (node1, node2) in edge_colors:
713 # FIXME directed graphs
714 ecolor = edge_colors[node1, node2]
715 else:
716 ecolor = edge_colors[node2, node1]
717 # Count per node how many edges it has of what color to nodes of
718 # what color
719 counts[node1][ecolor, node_colors[node2]] += 1
720 counts[node2][ecolor, node_colors[node1]] += 1
722 node_edge_colors = {}
723 for node in graph.nodes:
724 node_edge_colors[node] = node_colors[node], set(counts[node].items())
726 return node_edge_colors
728 @staticmethod
729 def _get_permutations_by_length(items):
730 """
731 Get all permutations of items, but only permute items with the same
732 length.
734 >>> found = list(ISMAGS._get_permutations_by_length([[1], [2], [3, 4], [4, 5]]))
735 >>> answer = [
736 ... (([1], [2]), ([3, 4], [4, 5])),
737 ... (([1], [2]), ([4, 5], [3, 4])),
738 ... (([2], [1]), ([3, 4], [4, 5])),
739 ... (([2], [1]), ([4, 5], [3, 4])),
740 ... ]
741 >>> found == answer
742 True
743 """
744 by_len = defaultdict(list)
745 for item in items:
746 by_len[len(item)].append(item)
748 yield from itertools.product(
749 *(itertools.permutations(by_len[l]) for l in sorted(by_len))
750 )
752 @classmethod
753 def _refine_node_partitions(cls, graph, node_partitions, edge_colors, branch=False):
754 """
755 Given a partition of nodes in graph, make the partitions smaller such
756 that all nodes in a partition have 1) the same color, and 2) the same
757 number of edges to specific other partitions.
758 """
760 def equal_color(node1, node2):
761 return node_edge_colors[node1] == node_edge_colors[node2]
763 node_partitions = list(node_partitions)
764 node_colors = partition_to_color(node_partitions)
765 node_edge_colors = cls._find_node_edge_color(graph, node_colors, edge_colors)
766 if all(
767 are_all_equal(node_edge_colors[node] for node in partition)
768 for partition in node_partitions
769 ):
770 yield node_partitions
771 return
773 new_partitions = []
774 output = [new_partitions]
775 for partition in node_partitions:
776 if not are_all_equal(node_edge_colors[node] for node in partition):
777 refined = make_partitions(partition, equal_color)
778 if (
779 branch
780 and len(refined) != 1
781 and len({len(r) for r in refined}) != len([len(r) for r in refined])
782 ):
783 # This is where it breaks. There are multiple new cells
784 # in refined with the same length, and their order
785 # matters.
786 # So option 1) Hit it with a big hammer and simply make all
787 # orderings.
788 permutations = cls._get_permutations_by_length(refined)
789 new_output = []
790 for n_p in output:
791 for permutation in permutations:
792 new_output.append(n_p + list(permutation[0]))
793 output = new_output
794 else:
795 for n_p in output:
796 n_p.extend(sorted(refined, key=len))
797 else:
798 for n_p in output:
799 n_p.append(partition)
800 for n_p in output:
801 yield from cls._refine_node_partitions(graph, n_p, edge_colors, branch)
803 def _edges_of_same_color(self, sgn1, sgn2):
804 """
805 Returns all edges in :attr:`graph` that have the same colour as the
806 edge between sgn1 and sgn2 in :attr:`subgraph`.
807 """
808 if (sgn1, sgn2) in self._sge_colors:
809 # FIXME directed graphs
810 sge_color = self._sge_colors[sgn1, sgn2]
811 else:
812 sge_color = self._sge_colors[sgn2, sgn1]
813 if sge_color in self._edge_compatibility:
814 ge_color = self._edge_compatibility[sge_color]
815 g_edges = self._ge_partitions[ge_color]
816 else:
817 g_edges = []
818 return g_edges
820 def _map_nodes(self, sgn, candidates, constraints, mapping=None, to_be_mapped=None):
821 """
822 Find all subgraph isomorphisms honoring constraints.
823 """
824 if mapping is None:
825 mapping = {}
826 else:
827 mapping = mapping.copy()
828 if to_be_mapped is None:
829 to_be_mapped = set(self.subgraph.nodes)
831 # Note, we modify candidates here. Doesn't seem to affect results, but
832 # remember this.
833 # candidates = candidates.copy()
834 sgn_candidates = intersect(candidates[sgn])
835 candidates[sgn] = frozenset([sgn_candidates])
836 for gn in sgn_candidates:
837 # We're going to try to map sgn to gn.
838 if gn in mapping.values() or sgn not in to_be_mapped:
839 # gn is already mapped to something
840 continue # pragma: no cover
842 # REDUCTION and COMBINATION
843 mapping[sgn] = gn
844 # BASECASE
845 if to_be_mapped == set(mapping.keys()):
846 yield {v: k for k, v in mapping.items()}
847 continue
848 left_to_map = to_be_mapped - set(mapping.keys())
850 new_candidates = candidates.copy()
851 sgn_nbrs = set(self.subgraph[sgn])
852 not_gn_nbrs = set(self.graph.nodes) - set(self.graph[gn])
853 for sgn2 in left_to_map:
854 if sgn2 not in sgn_nbrs:
855 gn2_options = not_gn_nbrs
856 else:
857 # Get all edges to gn of the right color:
858 g_edges = self._edges_of_same_color(sgn, sgn2)
859 # FIXME directed graphs
860 # And all nodes involved in those which are connected to gn
861 gn2_options = {n for e in g_edges for n in e if gn in e}
862 # Node color compatibility should be taken care of by the
863 # initial candidate lists made by find_subgraphs
865 # Add gn2_options to the right collection. Since new_candidates
866 # is a dict of frozensets of frozensets of node indices it's
867 # a bit clunky. We can't do .add, and + also doesn't work. We
868 # could do |, but I deem union to be clearer.
869 new_candidates[sgn2] = new_candidates[sgn2].union(
870 [frozenset(gn2_options)]
871 )
873 if (sgn, sgn2) in constraints:
874 gn2_options = {gn2 for gn2 in self.graph if gn2 > gn}
875 elif (sgn2, sgn) in constraints:
876 gn2_options = {gn2 for gn2 in self.graph if gn2 < gn}
877 else:
878 continue # pragma: no cover
879 new_candidates[sgn2] = new_candidates[sgn2].union(
880 [frozenset(gn2_options)]
881 )
883 # The next node is the one that is unmapped and has fewest
884 # candidates
885 next_sgn = min(left_to_map, key=lambda n: min(new_candidates[n], key=len))
886 yield from self._map_nodes(
887 next_sgn,
888 new_candidates,
889 constraints,
890 mapping=mapping,
891 to_be_mapped=to_be_mapped,
892 )
893 # Unmap sgn-gn. Strictly not necessary since it'd get overwritten
894 # when making a new mapping for sgn.
895 # del mapping[sgn]
897 def _largest_common_subgraph(self, candidates, constraints, to_be_mapped=None):
898 """
899 Find all largest common subgraphs honoring constraints.
900 """
901 if to_be_mapped is None:
902 to_be_mapped = {frozenset(self.subgraph.nodes)}
904 # The LCS problem is basically a repeated subgraph isomorphism problem
905 # with smaller and smaller subgraphs. We store the nodes that are
906 # "part of" the subgraph in to_be_mapped, and we make it a little
907 # smaller every iteration.
909 current_size = len(next(iter(to_be_mapped), []))
911 found_iso = False
912 if current_size <= len(self.graph):
913 # There's no point in trying to find isomorphisms of
914 # graph >= subgraph if subgraph has more nodes than graph.
916 # Try the isomorphism first with the nodes with lowest ID. So sort
917 # them. Those are more likely to be part of the final
918 # correspondence. This makes finding the first answer(s) faster. In
919 # theory.
920 for nodes in sorted(to_be_mapped, key=sorted):
921 # Find the isomorphism between subgraph[to_be_mapped] <= graph
922 next_sgn = min(nodes, key=lambda n: min(candidates[n], key=len))
923 isomorphs = self._map_nodes(
924 next_sgn, candidates, constraints, to_be_mapped=nodes
925 )
927 # This is effectively `yield from isomorphs`, except that we look
928 # whether an item was yielded.
929 try:
930 item = next(isomorphs)
931 except StopIteration:
932 pass
933 else:
934 yield item
935 yield from isomorphs
936 found_iso = True
938 # BASECASE
939 if found_iso or current_size == 1:
940 # Shrinking has no point because either 1) we end up with a smaller
941 # common subgraph (and we want the largest), or 2) there'll be no
942 # more subgraph.
943 return
945 left_to_be_mapped = set()
946 for nodes in to_be_mapped:
947 for sgn in nodes:
948 # We're going to remove sgn from to_be_mapped, but subject to
949 # symmetry constraints. We know that for every constraint we
950 # have those subgraph nodes are equal. So whenever we would
951 # remove the lower part of a constraint, remove the higher
952 # instead. This is all dealth with by _remove_node. And because
953 # left_to_be_mapped is a set, we don't do double work.
955 # And finally, make the subgraph one node smaller.
956 # REDUCTION
957 new_nodes = self._remove_node(sgn, nodes, constraints)
958 left_to_be_mapped.add(new_nodes)
959 # COMBINATION
960 yield from self._largest_common_subgraph(
961 candidates, constraints, to_be_mapped=left_to_be_mapped
962 )
964 @staticmethod
965 def _remove_node(node, nodes, constraints):
966 """
967 Returns a new set where node has been removed from nodes, subject to
968 symmetry constraints. We know, that for every constraint we have
969 those subgraph nodes are equal. So whenever we would remove the
970 lower part of a constraint, remove the higher instead.
971 """
972 while True:
973 for low, high in constraints:
974 if low == node and high in nodes:
975 node = high
976 break
977 else: # no break, couldn't find node in constraints
978 break
979 return frozenset(nodes - {node})
981 @staticmethod
982 def _find_permutations(top_partitions, bottom_partitions):
983 """
984 Return the pairs of top/bottom partitions where the partitions are
985 different. Ensures that all partitions in both top and bottom
986 partitions have size 1.
987 """
988 # Find permutations
989 permutations = set()
990 for top, bot in zip(top_partitions, bottom_partitions):
991 # top and bot have only one element
992 if len(top) != 1 or len(bot) != 1:
993 raise IndexError(
994 "Not all nodes are coupled. This is"
995 f" impossible: {top_partitions}, {bottom_partitions}"
996 )
997 if top != bot:
998 permutations.add(frozenset((next(iter(top)), next(iter(bot)))))
999 return permutations
1001 @staticmethod
1002 def _update_orbits(orbits, permutations):
1003 """
1004 Update orbits based on permutations. Orbits is modified in place.
1005 For every pair of items in permutations their respective orbits are
1006 merged.
1007 """
1008 for permutation in permutations:
1009 node, node2 = permutation
1010 # Find the orbits that contain node and node2, and replace the
1011 # orbit containing node with the union
1012 first = second = None
1013 for idx, orbit in enumerate(orbits):
1014 if first is not None and second is not None:
1015 break
1016 if node in orbit:
1017 first = idx
1018 if node2 in orbit:
1019 second = idx
1020 if first != second:
1021 orbits[first].update(orbits[second])
1022 del orbits[second]
1024 def _couple_nodes(
1025 self,
1026 top_partitions,
1027 bottom_partitions,
1028 pair_idx,
1029 t_node,
1030 b_node,
1031 graph,
1032 edge_colors,
1033 ):
1034 """
1035 Generate new partitions from top and bottom_partitions where t_node is
1036 coupled to b_node. pair_idx is the index of the partitions where t_ and
1037 b_node can be found.
1038 """
1039 t_partition = top_partitions[pair_idx]
1040 b_partition = bottom_partitions[pair_idx]
1041 assert t_node in t_partition and b_node in b_partition
1042 # Couple node to node2. This means they get their own partition
1043 new_top_partitions = [top.copy() for top in top_partitions]
1044 new_bottom_partitions = [bot.copy() for bot in bottom_partitions]
1045 new_t_groups = {t_node}, t_partition - {t_node}
1046 new_b_groups = {b_node}, b_partition - {b_node}
1047 # Replace the old partitions with the coupled ones
1048 del new_top_partitions[pair_idx]
1049 del new_bottom_partitions[pair_idx]
1050 new_top_partitions[pair_idx:pair_idx] = new_t_groups
1051 new_bottom_partitions[pair_idx:pair_idx] = new_b_groups
1053 new_top_partitions = self._refine_node_partitions(
1054 graph, new_top_partitions, edge_colors
1055 )
1056 new_bottom_partitions = self._refine_node_partitions(
1057 graph, new_bottom_partitions, edge_colors, branch=True
1058 )
1059 new_top_partitions = list(new_top_partitions)
1060 assert len(new_top_partitions) == 1
1061 new_top_partitions = new_top_partitions[0]
1062 for bot in new_bottom_partitions:
1063 yield list(new_top_partitions), bot
1065 def _process_ordered_pair_partitions(
1066 self,
1067 graph,
1068 top_partitions,
1069 bottom_partitions,
1070 edge_colors,
1071 orbits=None,
1072 cosets=None,
1073 ):
1074 """
1075 Processes ordered pair partitions as per the reference paper. Finds and
1076 returns all permutations and cosets that leave the graph unchanged.
1077 """
1078 if orbits is None:
1079 orbits = [{node} for node in graph.nodes]
1080 else:
1081 # Note that we don't copy orbits when we are given one. This means
1082 # we leak information between the recursive branches. This is
1083 # intentional!
1084 orbits = orbits
1085 if cosets is None:
1086 cosets = {}
1087 else:
1088 cosets = cosets.copy()
1090 if not all(
1091 len(t_p) == len(b_p) for t_p, b_p in zip(top_partitions, bottom_partitions)
1092 ):
1093 # This used to be an assertion, but it gets tripped in rare cases:
1094 # 5 - 4 \ / 12 - 13
1095 # 0 - 3
1096 # 9 - 8 / \ 16 - 17
1097 # Assume 0 and 3 are coupled and no longer equivalent. At that point
1098 # {4, 8} and {12, 16} are no longer equivalent, and neither are
1099 # {5, 9} and {13, 17}. Coupling 4 and refinement results in 5 and 9
1100 # getting their own partitions, *but not 13 and 17*. Further
1101 # iterations will attempt to couple 5 to {13, 17}, which cannot
1102 # result in more symmetries?
1103 return [], cosets
1105 # BASECASE
1106 if all(len(top) == 1 for top in top_partitions):
1107 # All nodes are mapped
1108 permutations = self._find_permutations(top_partitions, bottom_partitions)
1109 self._update_orbits(orbits, permutations)
1110 if permutations:
1111 return [permutations], cosets
1112 else:
1113 return [], cosets
1115 permutations = []
1116 unmapped_nodes = {
1117 (node, idx)
1118 for idx, t_partition in enumerate(top_partitions)
1119 for node in t_partition
1120 if len(t_partition) > 1
1121 }
1122 node, pair_idx = min(unmapped_nodes)
1123 b_partition = bottom_partitions[pair_idx]
1125 for node2 in sorted(b_partition):
1126 if len(b_partition) == 1:
1127 # Can never result in symmetry
1128 continue
1129 if node != node2 and any(
1130 node in orbit and node2 in orbit for orbit in orbits
1131 ):
1132 # Orbit prune branch
1133 continue
1134 # REDUCTION
1135 # Couple node to node2
1136 partitions = self._couple_nodes(
1137 top_partitions,
1138 bottom_partitions,
1139 pair_idx,
1140 node,
1141 node2,
1142 graph,
1143 edge_colors,
1144 )
1145 for opp in partitions:
1146 new_top_partitions, new_bottom_partitions = opp
1148 new_perms, new_cosets = self._process_ordered_pair_partitions(
1149 graph,
1150 new_top_partitions,
1151 new_bottom_partitions,
1152 edge_colors,
1153 orbits,
1154 cosets,
1155 )
1156 # COMBINATION
1157 permutations += new_perms
1158 cosets.update(new_cosets)
1160 mapped = {
1161 k
1162 for top, bottom in zip(top_partitions, bottom_partitions)
1163 for k in top
1164 if len(top) == 1 and top == bottom
1165 }
1166 ks = {k for k in graph.nodes if k < node}
1167 # Have all nodes with ID < node been mapped?
1168 find_coset = ks <= mapped and node not in cosets
1169 if find_coset:
1170 # Find the orbit that contains node
1171 for orbit in orbits:
1172 if node in orbit:
1173 cosets[node] = orbit.copy()
1174 return permutations, cosets