cycles.py 40 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660661662663664665666667668669670671672673674675676677678679680681682683684685686687688689690691692693694695696697698699700701702703704705706707708709710711712713714715716717718719720721722723724725726727728729730731732733734735736737738739740741742743744745746747748749750751752753754755756757758759760761762763764765766767768769770771772773774775776777778779780781782783784785786787788789790791792793794795796797798799800801802803804805806807808809810811812813814815816817818819820821822823824825826827828829830831832833834835836837838839840841842843844845846847848849850851852853854855856857858859860861862863864865866867868869870871872873874875876877878879880881882883884885886887888889890891892893894895896897898899900901902903904905906907908909910911912913914915916917918919920921922923924925926927928929930931932933934935936937938939940941942943944945946947948949950951952953954955956957958959960961962963964965966967968969970971972973974975976977978979980981982983984985986987988989990991992993994995996997998999100010011002100310041005100610071008100910101011101210131014101510161017101810191020102110221023102410251026102710281029103010311032103310341035103610371038103910401041104210431044104510461047104810491050105110521053105410551056105710581059106010611062106310641065106610671068106910701071107210731074107510761077107810791080108110821083108410851086108710881089109010911092109310941095109610971098109911001101110211031104110511061107110811091110111111121113111411151116111711181119112011211122112311241125112611271128112911301131113211331134113511361137113811391140114111421143114411451146114711481149115011511152
  1. """
  2. ========================
  3. Cycle finding algorithms
  4. ========================
  5. """
  6. from collections import defaultdict
  7. from itertools import combinations, product
  8. import networkx as nx
  9. from networkx.utils import not_implemented_for, pairwise
  10. __all__ = [
  11. "cycle_basis",
  12. "simple_cycles",
  13. "recursive_simple_cycles",
  14. "find_cycle",
  15. "minimum_cycle_basis",
  16. "chordless_cycles",
  17. ]
  18. @not_implemented_for("directed")
  19. @not_implemented_for("multigraph")
  20. def cycle_basis(G, root=None):
  21. """Returns a list of cycles which form a basis for cycles of G.
  22. A basis for cycles of a network is a minimal collection of
  23. cycles such that any cycle in the network can be written
  24. as a sum of cycles in the basis. Here summation of cycles
  25. is defined as "exclusive or" of the edges. Cycle bases are
  26. useful, e.g. when deriving equations for electric circuits
  27. using Kirchhoff's Laws.
  28. Parameters
  29. ----------
  30. G : NetworkX Graph
  31. root : node, optional
  32. Specify starting node for basis.
  33. Returns
  34. -------
  35. A list of cycle lists. Each cycle list is a list of nodes
  36. which forms a cycle (loop) in G.
  37. Examples
  38. --------
  39. >>> G = nx.Graph()
  40. >>> nx.add_cycle(G, [0, 1, 2, 3])
  41. >>> nx.add_cycle(G, [0, 3, 4, 5])
  42. >>> print(nx.cycle_basis(G, 0))
  43. [[3, 4, 5, 0], [1, 2, 3, 0]]
  44. Notes
  45. -----
  46. This is adapted from algorithm CACM 491 [1]_.
  47. References
  48. ----------
  49. .. [1] Paton, K. An algorithm for finding a fundamental set of
  50. cycles of a graph. Comm. ACM 12, 9 (Sept 1969), 514-518.
  51. See Also
  52. --------
  53. simple_cycles
  54. """
  55. gnodes = set(G.nodes())
  56. cycles = []
  57. while gnodes: # loop over connected components
  58. if root is None:
  59. root = gnodes.pop()
  60. stack = [root]
  61. pred = {root: root}
  62. used = {root: set()}
  63. while stack: # walk the spanning tree finding cycles
  64. z = stack.pop() # use last-in so cycles easier to find
  65. zused = used[z]
  66. for nbr in G[z]:
  67. if nbr not in used: # new node
  68. pred[nbr] = z
  69. stack.append(nbr)
  70. used[nbr] = {z}
  71. elif nbr == z: # self loops
  72. cycles.append([z])
  73. elif nbr not in zused: # found a cycle
  74. pn = used[nbr]
  75. cycle = [nbr, z]
  76. p = pred[z]
  77. while p not in pn:
  78. cycle.append(p)
  79. p = pred[p]
  80. cycle.append(p)
  81. cycles.append(cycle)
  82. used[nbr].add(z)
  83. gnodes -= set(pred)
  84. root = None
  85. return cycles
  86. def simple_cycles(G, length_bound=None):
  87. """Find simple cycles (elementary circuits) of a graph.
  88. A `simple cycle`, or `elementary circuit`, is a closed path where
  89. no node appears twice. In a directed graph, two simple cycles are distinct
  90. if they are not cyclic permutations of each other. In an undirected graph,
  91. two simple cycles are distinct if they are not cyclic permutations of each
  92. other nor of the other's reversal.
  93. Optionally, the cycles are bounded in length. In the unbounded case, we use
  94. a nonrecursive, iterator/generator version of Johnson's algorithm [1]_. In
  95. the bounded case, we use a version of the algorithm of Gupta and
  96. Suzumura[2]_. There may be better algorithms for some cases [3]_ [4]_ [5]_.
  97. The algorithms of Johnson, and Gupta and Suzumura, are enhanced by some
  98. well-known preprocessing techniques. When G is directed, we restrict our
  99. attention to strongly connected components of G, generate all simple cycles
  100. containing a certain node, remove that node, and further decompose the
  101. remainder into strongly connected components. When G is undirected, we
  102. restrict our attention to biconnected components, generate all simple cycles
  103. containing a particular edge, remove that edge, and further decompose the
  104. remainder into biconnected components.
  105. Note that multigraphs are supported by this function -- and in undirected
  106. multigraphs, a pair of parallel edges is considered a cycle of length 2.
  107. Likewise, self-loops are considered to be cycles of length 1. We define
  108. cycles as sequences of nodes; so the presence of loops and parallel edges
  109. does not change the number of simple cycles in a graph.
  110. Parameters
  111. ----------
  112. G : NetworkX DiGraph
  113. A directed graph
  114. length_bound : int or None, optional (default=None)
  115. If length_bound is an int, generate all simple cycles of G with length at
  116. most length_bound. Otherwise, generate all simple cycles of G.
  117. Yields
  118. ------
  119. list of nodes
  120. Each cycle is represented by a list of nodes along the cycle.
  121. Examples
  122. --------
  123. >>> edges = [(0, 0), (0, 1), (0, 2), (1, 2), (2, 0), (2, 1), (2, 2)]
  124. >>> G = nx.DiGraph(edges)
  125. >>> sorted(nx.simple_cycles(G))
  126. [[0], [0, 1, 2], [0, 2], [1, 2], [2]]
  127. To filter the cycles so that they don't include certain nodes or edges,
  128. copy your graph and eliminate those nodes or edges before calling.
  129. For example, to exclude self-loops from the above example:
  130. >>> H = G.copy()
  131. >>> H.remove_edges_from(nx.selfloop_edges(G))
  132. >>> sorted(nx.simple_cycles(H))
  133. [[0, 1, 2], [0, 2], [1, 2]]
  134. Notes
  135. -----
  136. When length_bound is None, the time complexity is $O((n+e)(c+1))$ for $n$
  137. nodes, $e$ edges and $c$ simple circuits. Otherwise, when length_bound > 1,
  138. the time complexity is $O((c+n)(k-1)d^k)$ where $d$ is the average degree of
  139. the nodes of G and $k$ = length_bound.
  140. Raises
  141. ------
  142. ValueError
  143. when length_bound < 0.
  144. References
  145. ----------
  146. .. [1] Finding all the elementary circuits of a directed graph.
  147. D. B. Johnson, SIAM Journal on Computing 4, no. 1, 77-84, 1975.
  148. https://doi.org/10.1137/0204007
  149. .. [2] Finding All Bounded-Length Simple Cycles in a Directed Graph
  150. A. Gupta and T. Suzumura https://arxiv.org/abs/2105.10094
  151. .. [3] Enumerating the cycles of a digraph: a new preprocessing strategy.
  152. G. Loizou and P. Thanish, Information Sciences, v. 27, 163-182, 1982.
  153. .. [4] A search strategy for the elementary cycles of a directed graph.
  154. J.L. Szwarcfiter and P.E. Lauer, BIT NUMERICAL MATHEMATICS,
  155. v. 16, no. 2, 192-204, 1976.
  156. .. [5] Optimal Listing of Cycles and st-Paths in Undirected Graphs
  157. R. Ferreira and R. Grossi and A. Marino and N. Pisanti and R. Rizzi and
  158. G. Sacomoto https://arxiv.org/abs/1205.2766
  159. See Also
  160. --------
  161. cycle_basis
  162. chordless_cycles
  163. """
  164. if length_bound is not None:
  165. if length_bound == 0:
  166. return
  167. elif length_bound < 0:
  168. raise ValueError("length bound must be non-negative")
  169. directed = G.is_directed()
  170. yield from ([v] for v, Gv in G.adj.items() if v in Gv)
  171. if length_bound is not None and length_bound == 1:
  172. return
  173. if G.is_multigraph() and not directed:
  174. visited = set()
  175. for u, Gu in G.adj.items():
  176. multiplicity = ((v, len(Guv)) for v, Guv in Gu.items() if v in visited)
  177. yield from ([u, v] for v, m in multiplicity if m > 1)
  178. visited.add(u)
  179. # explicitly filter out loops; implicitly filter out parallel edges
  180. if directed:
  181. G = nx.DiGraph((u, v) for u, Gu in G.adj.items() for v in Gu if v != u)
  182. else:
  183. G = nx.Graph((u, v) for u, Gu in G.adj.items() for v in Gu if v != u)
  184. # this case is not strictly necessary but improves performance
  185. if length_bound is not None and length_bound == 2:
  186. if directed:
  187. visited = set()
  188. for u, Gu in G.adj.items():
  189. yield from (
  190. [v, u] for v in visited.intersection(Gu) if G.has_edge(v, u)
  191. )
  192. visited.add(u)
  193. return
  194. if directed:
  195. yield from _directed_cycle_search(G, length_bound)
  196. else:
  197. yield from _undirected_cycle_search(G, length_bound)
  198. def _directed_cycle_search(G, length_bound):
  199. """A dispatch function for `simple_cycles` for directed graphs.
  200. We generate all cycles of G through binary partition.
  201. 1. Pick a node v in G which belongs to at least one cycle
  202. a. Generate all cycles of G which contain the node v.
  203. b. Recursively generate all cycles of G \\ v.
  204. This is accomplished through the following:
  205. 1. Compute the strongly connected components SCC of G.
  206. 2. Select and remove a biconnected component C from BCC. Select a
  207. non-tree edge (u, v) of a depth-first search of G[C].
  208. 3. For each simple cycle P containing v in G[C], yield P.
  209. 4. Add the biconnected components of G[C \\ v] to BCC.
  210. If the parameter length_bound is not None, then step 3 will be limited to
  211. simple cycles of length at most length_bound.
  212. Parameters
  213. ----------
  214. G : NetworkX DiGraph
  215. A directed graph
  216. length_bound : int or None
  217. If length_bound is an int, generate all simple cycles of G with length at most length_bound.
  218. Otherwise, generate all simple cycles of G.
  219. Yields
  220. ------
  221. list of nodes
  222. Each cycle is represented by a list of nodes along the cycle.
  223. """
  224. scc = nx.strongly_connected_components
  225. components = [c for c in scc(G) if len(c) >= 2]
  226. while components:
  227. c = components.pop()
  228. Gc = G.subgraph(c)
  229. v = next(iter(c))
  230. if length_bound is None:
  231. yield from _johnson_cycle_search(Gc, [v])
  232. else:
  233. yield from _bounded_cycle_search(Gc, [v], length_bound)
  234. # delete v after searching G, to make sure we can find v
  235. G.remove_node(v)
  236. components.extend(c for c in scc(Gc) if len(c) >= 2)
  237. def _undirected_cycle_search(G, length_bound):
  238. """A dispatch function for `simple_cycles` for undirected graphs.
  239. We generate all cycles of G through binary partition.
  240. 1. Pick an edge (u, v) in G which belongs to at least one cycle
  241. a. Generate all cycles of G which contain the edge (u, v)
  242. b. Recursively generate all cycles of G \\ (u, v)
  243. This is accomplished through the following:
  244. 1. Compute the biconnected components BCC of G.
  245. 2. Select and remove a biconnected component C from BCC. Select a
  246. non-tree edge (u, v) of a depth-first search of G[C].
  247. 3. For each (v -> u) path P remaining in G[C] \\ (u, v), yield P.
  248. 4. Add the biconnected components of G[C] \\ (u, v) to BCC.
  249. If the parameter length_bound is not None, then step 3 will be limited to simple paths
  250. of length at most length_bound.
  251. Parameters
  252. ----------
  253. G : NetworkX Graph
  254. An undirected graph
  255. length_bound : int or None
  256. If length_bound is an int, generate all simple cycles of G with length at most length_bound.
  257. Otherwise, generate all simple cycles of G.
  258. Yields
  259. ------
  260. list of nodes
  261. Each cycle is represented by a list of nodes along the cycle.
  262. """
  263. bcc = nx.biconnected_components
  264. components = [c for c in bcc(G) if len(c) >= 3]
  265. while components:
  266. c = components.pop()
  267. Gc = G.subgraph(c)
  268. uv = list(next(iter(Gc.edges)))
  269. G.remove_edge(*uv)
  270. # delete (u, v) before searching G, to avoid fake 3-cycles [u, v, u]
  271. if length_bound is None:
  272. yield from _johnson_cycle_search(Gc, uv)
  273. else:
  274. yield from _bounded_cycle_search(Gc, uv, length_bound)
  275. components.extend(c for c in bcc(Gc) if len(c) >= 3)
  276. class _NeighborhoodCache(dict):
  277. """Very lightweight graph wrapper which caches neighborhoods as list.
  278. This dict subclass uses the __missing__ functionality to query graphs for
  279. their neighborhoods, and store the result as a list. This is used to avoid
  280. the performance penalty incurred by subgraph views.
  281. """
  282. def __init__(self, G):
  283. self.G = G
  284. def __missing__(self, v):
  285. Gv = self[v] = list(self.G[v])
  286. return Gv
  287. def _johnson_cycle_search(G, path):
  288. """The main loop of the cycle-enumeration algorithm of Johnson.
  289. Parameters
  290. ----------
  291. G : NetworkX Graph or DiGraph
  292. A graph
  293. path : list
  294. A cycle prefix. All cycles generated will begin with this prefix.
  295. Yields
  296. ------
  297. list of nodes
  298. Each cycle is represented by a list of nodes along the cycle.
  299. References
  300. ----------
  301. .. [1] Finding all the elementary circuits of a directed graph.
  302. D. B. Johnson, SIAM Journal on Computing 4, no. 1, 77-84, 1975.
  303. https://doi.org/10.1137/0204007
  304. """
  305. G = _NeighborhoodCache(G)
  306. blocked = set(path)
  307. B = defaultdict(set) # graph portions that yield no elementary circuit
  308. start = path[0]
  309. stack = [iter(G[path[-1]])]
  310. closed = [False]
  311. while stack:
  312. nbrs = stack[-1]
  313. for w in nbrs:
  314. if w == start:
  315. yield path[:]
  316. closed[-1] = True
  317. elif w not in blocked:
  318. path.append(w)
  319. closed.append(False)
  320. stack.append(iter(G[w]))
  321. blocked.add(w)
  322. break
  323. else: # no more nbrs
  324. stack.pop()
  325. v = path.pop()
  326. if closed.pop():
  327. if closed:
  328. closed[-1] = True
  329. unblock_stack = {v}
  330. while unblock_stack:
  331. u = unblock_stack.pop()
  332. if u in blocked:
  333. blocked.remove(u)
  334. unblock_stack.update(B[u])
  335. B[u].clear()
  336. else:
  337. for w in G[v]:
  338. B[w].add(v)
  339. def _bounded_cycle_search(G, path, length_bound):
  340. """The main loop of the cycle-enumeration algorithm of Gupta and Suzumura.
  341. Parameters
  342. ----------
  343. G : NetworkX Graph or DiGraph
  344. A graph
  345. path : list
  346. A cycle prefix. All cycles generated will begin with this prefix.
  347. length_bound: int
  348. A length bound. All cycles generated will have length at most length_bound.
  349. Yields
  350. ------
  351. list of nodes
  352. Each cycle is represented by a list of nodes along the cycle.
  353. References
  354. ----------
  355. .. [1] Finding All Bounded-Length Simple Cycles in a Directed Graph
  356. A. Gupta and T. Suzumura https://arxiv.org/abs/2105.10094
  357. """
  358. G = _NeighborhoodCache(G)
  359. lock = {v: 0 for v in path}
  360. B = defaultdict(set)
  361. start = path[0]
  362. stack = [iter(G[path[-1]])]
  363. blen = [length_bound]
  364. while stack:
  365. nbrs = stack[-1]
  366. for w in nbrs:
  367. if w == start:
  368. yield path[:]
  369. blen[-1] = 1
  370. elif len(path) < lock.get(w, length_bound):
  371. path.append(w)
  372. blen.append(length_bound)
  373. lock[w] = len(path)
  374. stack.append(iter(G[w]))
  375. break
  376. else:
  377. stack.pop()
  378. v = path.pop()
  379. bl = blen.pop()
  380. if blen:
  381. blen[-1] = min(blen[-1], bl)
  382. if bl < length_bound:
  383. relax_stack = [(bl, v)]
  384. while relax_stack:
  385. bl, u = relax_stack.pop()
  386. if lock.get(u, length_bound) < length_bound - bl + 1:
  387. lock[u] = length_bound - bl + 1
  388. relax_stack.extend((bl + 1, w) for w in B[u].difference(path))
  389. else:
  390. for w in G[v]:
  391. B[w].add(v)
  392. def chordless_cycles(G, length_bound=None):
  393. """Find simple chordless cycles of a graph.
  394. A `simple cycle` is a closed path where no node appears twice. In a simple
  395. cycle, a `chord` is an additional edge between two nodes in the cycle. A
  396. `chordless cycle` is a simple cycle without chords. Said differently, a
  397. chordless cycle is a cycle C in a graph G where the number of edges in the
  398. induced graph G[C] is equal to the length of `C`.
  399. Note that some care must be taken in the case that G is not a simple graph
  400. nor a simple digraph. Some authors limit the definition of chordless cycles
  401. to have a prescribed minimum length; we do not.
  402. 1. We interpret self-loops to be chordless cycles, except in multigraphs
  403. with multiple loops in parallel. Likewise, in a chordless cycle of
  404. length greater than 1, there can be no nodes with self-loops.
  405. 2. We interpret directed two-cycles to be chordless cycles, except in
  406. multi-digraphs when any edge in a two-cycle has a parallel copy.
  407. 3. We interpret parallel pairs of undirected edges as two-cycles, except
  408. when a third (or more) parallel edge exists between the two nodes.
  409. 4. Generalizing the above, edges with parallel clones may not occur in
  410. chordless cycles.
  411. In a directed graph, two chordless cycles are distinct if they are not
  412. cyclic permutations of each other. In an undirected graph, two chordless
  413. cycles are distinct if they are not cyclic permutations of each other nor of
  414. the other's reversal.
  415. Optionally, the cycles are bounded in length.
  416. We use an algorithm strongly inspired by that of Dias et al [1]_. It has
  417. been modified in the following ways:
  418. 1. Recursion is avoided, per Python's limitations
  419. 2. The labeling function is not necessary, because the starting paths
  420. are chosen (and deleted from the host graph) to prevent multiple
  421. occurrences of the same path
  422. 3. The search is optionally bounded at a specified length
  423. 4. Support for directed graphs is provided by extending cycles along
  424. forward edges, and blocking nodes along forward and reverse edges
  425. 5. Support for multigraphs is provided by omitting digons from the set
  426. of forward edges
  427. Parameters
  428. ----------
  429. G : NetworkX DiGraph
  430. A directed graph
  431. length_bound : int or None, optional (default=None)
  432. If length_bound is an int, generate all simple cycles of G with length at
  433. most length_bound. Otherwise, generate all simple cycles of G.
  434. Yields
  435. ------
  436. list of nodes
  437. Each cycle is represented by a list of nodes along the cycle.
  438. Examples
  439. --------
  440. >>> sorted(list(nx.chordless_cycles(nx.complete_graph(4))))
  441. [[1, 0, 2], [1, 0, 3], [2, 0, 3], [2, 1, 3]]
  442. Notes
  443. -----
  444. When length_bound is None, and the graph is simple, the time complexity is
  445. $O((n+e)(c+1))$ for $n$ nodes, $e$ edges and $c$ chordless cycles.
  446. Raises
  447. ------
  448. ValueError
  449. when length_bound < 0.
  450. References
  451. ----------
  452. .. [1] Efficient enumeration of chordless cycles
  453. E. Dias and D. Castonguay and H. Longo and W.A.R. Jradi
  454. https://arxiv.org/abs/1309.1051
  455. See Also
  456. --------
  457. simple_cycles
  458. """
  459. if length_bound is not None:
  460. if length_bound == 0:
  461. return
  462. elif length_bound < 0:
  463. raise ValueError("length bound must be non-negative")
  464. directed = G.is_directed()
  465. multigraph = G.is_multigraph()
  466. if multigraph:
  467. yield from ([v] for v, Gv in G.adj.items() if len(Gv.get(v, ())) == 1)
  468. else:
  469. yield from ([v] for v, Gv in G.adj.items() if v in Gv)
  470. if length_bound is not None and length_bound == 1:
  471. return
  472. # Nodes with loops cannot belong to longer cycles. Let's delete them here.
  473. # also, we implicitly reduce the multiplicity of edges down to 1 in the case
  474. # of multiedges.
  475. if directed:
  476. F = nx.DiGraph((u, v) for u, Gu in G.adj.items() if u not in Gu for v in Gu)
  477. B = F.to_undirected(as_view=False)
  478. else:
  479. F = nx.Graph((u, v) for u, Gu in G.adj.items() if u not in Gu for v in Gu)
  480. B = None
  481. # If we're given a multigraph, we have a few cases to consider with parallel
  482. # edges.
  483. #
  484. # 1. If we have 2 or more edges in parallel between the nodes (u, v), we
  485. # must not construct longer cycles along (u, v).
  486. # 2. If G is not directed, then a pair of parallel edges between (u, v) is a
  487. # chordless cycle unless there exists a third (or more) parallel edge.
  488. # 3. If G is directed, then parallel edges do not form cyles, but do
  489. # preclude back-edges from forming cycles (handled in the next section),
  490. # Thus, if an edge (u, v) is duplicated and the reverse (v, u) is also
  491. # present, then we remove both from F.
  492. #
  493. # In directed graphs, we need to consider both directions that edges can
  494. # take, so iterate over all edges (u, v) and possibly (v, u). In undirected
  495. # graphs, we need to be a little careful to only consider every edge once,
  496. # so we use a "visited" set to emulate node-order comparisons.
  497. if multigraph:
  498. if not directed:
  499. B = F.copy()
  500. visited = set()
  501. for u, Gu in G.adj.items():
  502. if directed:
  503. multiplicity = ((v, len(Guv)) for v, Guv in Gu.items())
  504. for v, m in multiplicity:
  505. if m > 1:
  506. F.remove_edges_from(((u, v), (v, u)))
  507. else:
  508. multiplicity = ((v, len(Guv)) for v, Guv in Gu.items() if v in visited)
  509. for v, m in multiplicity:
  510. if m == 2:
  511. yield [u, v]
  512. if m > 1:
  513. F.remove_edge(u, v)
  514. visited.add(u)
  515. # If we're given a directed graphs, we need to think about digons. If we
  516. # have two edges (u, v) and (v, u), then that's a two-cycle. If either edge
  517. # was duplicated above, then we removed both from F. So, any digons we find
  518. # here are chordless. After finding digons, we remove their edges from F
  519. # to avoid traversing them in the search for chordless cycles.
  520. if directed:
  521. for u, Fu in F.adj.items():
  522. digons = [[u, v] for v in Fu if F.has_edge(v, u)]
  523. yield from digons
  524. F.remove_edges_from(digons)
  525. F.remove_edges_from(e[::-1] for e in digons)
  526. if length_bound is not None and length_bound == 2:
  527. return
  528. # Now, we prepare to search for cycles. We have removed all cycles of
  529. # lengths 1 and 2, so F is a simple graph or simple digraph. We repeatedly
  530. # separate digraphs into their strongly connected components, and undirected
  531. # graphs into their biconnected components. For each component, we pick a
  532. # node v, search for chordless cycles based at each "stem" (u, v, w), and
  533. # then remove v from that component before separating the graph again.
  534. if directed:
  535. separate = nx.strongly_connected_components
  536. # Directed stems look like (u -> v -> w), so we use the product of
  537. # predecessors of v with successors of v.
  538. def stems(C, v):
  539. for u, w in product(C.pred[v], C.succ[v]):
  540. if not G.has_edge(u, w): # omit stems with acyclic chords
  541. yield [u, v, w], F.has_edge(w, u)
  542. else:
  543. separate = nx.biconnected_components
  544. # Undirected stems look like (u ~ v ~ w), but we must not also search
  545. # (w ~ v ~ u), so we use combinations of v's neighbors of length 2.
  546. def stems(C, v):
  547. yield from (([u, v, w], F.has_edge(w, u)) for u, w in combinations(C[v], 2))
  548. components = [c for c in separate(F) if len(c) > 2]
  549. while components:
  550. c = components.pop()
  551. v = next(iter(c))
  552. Fc = F.subgraph(c)
  553. Fcc = Bcc = None
  554. for S, is_triangle in stems(Fc, v):
  555. if is_triangle:
  556. yield S
  557. else:
  558. if Fcc is None:
  559. Fcc = _NeighborhoodCache(Fc)
  560. Bcc = Fcc if B is None else _NeighborhoodCache(B.subgraph(c))
  561. yield from _chordless_cycle_search(Fcc, Bcc, S, length_bound)
  562. components.extend(c for c in separate(F.subgraph(c - {v})) if len(c) > 2)
  563. def _chordless_cycle_search(F, B, path, length_bound):
  564. """The main loop for chordless cycle enumeration.
  565. This algorithm is strongly inspired by that of Dias et al [1]_. It has been
  566. modified in the following ways:
  567. 1. Recursion is avoided, per Python's limitations
  568. 2. The labeling function is not necessary, because the starting paths
  569. are chosen (and deleted from the host graph) to prevent multiple
  570. occurrences of the same path
  571. 3. The search is optionally bounded at a specified length
  572. 4. Support for directed graphs is provided by extending cycles along
  573. forward edges, and blocking nodes along forward and reverse edges
  574. 5. Support for multigraphs is provided by omitting digons from the set
  575. of forward edges
  576. Parameters
  577. ----------
  578. F : _NeighborhoodCache
  579. A graph of forward edges to follow in constructing cycles
  580. B : _NeighborhoodCache
  581. A graph of blocking edges to prevent the production of chordless cycles
  582. path : list
  583. A cycle prefix. All cycles generated will begin with this prefix.
  584. length_bound : int
  585. A length bound. All cycles generated will have length at most length_bound.
  586. Yields
  587. ------
  588. list of nodes
  589. Each cycle is represented by a list of nodes along the cycle.
  590. References
  591. ----------
  592. .. [1] Efficient enumeration of chordless cycles
  593. E. Dias and D. Castonguay and H. Longo and W.A.R. Jradi
  594. https://arxiv.org/abs/1309.1051
  595. """
  596. blocked = defaultdict(int)
  597. target = path[0]
  598. blocked[path[1]] = 1
  599. for w in path[1:]:
  600. for v in B[w]:
  601. blocked[v] += 1
  602. stack = [iter(F[path[2]])]
  603. while stack:
  604. nbrs = stack[-1]
  605. for w in nbrs:
  606. if blocked[w] == 1 and (length_bound is None or len(path) < length_bound):
  607. Fw = F[w]
  608. if target in Fw:
  609. yield path + [w]
  610. else:
  611. Bw = B[w]
  612. if target in Bw:
  613. continue
  614. for v in Bw:
  615. blocked[v] += 1
  616. path.append(w)
  617. stack.append(iter(Fw))
  618. break
  619. else:
  620. stack.pop()
  621. for v in B[path.pop()]:
  622. blocked[v] -= 1
  623. @not_implemented_for("undirected")
  624. def recursive_simple_cycles(G):
  625. """Find simple cycles (elementary circuits) of a directed graph.
  626. A `simple cycle`, or `elementary circuit`, is a closed path where
  627. no node appears twice. Two elementary circuits are distinct if they
  628. are not cyclic permutations of each other.
  629. This version uses a recursive algorithm to build a list of cycles.
  630. You should probably use the iterator version called simple_cycles().
  631. Warning: This recursive version uses lots of RAM!
  632. It appears in NetworkX for pedagogical value.
  633. Parameters
  634. ----------
  635. G : NetworkX DiGraph
  636. A directed graph
  637. Returns
  638. -------
  639. A list of cycles, where each cycle is represented by a list of nodes
  640. along the cycle.
  641. Example:
  642. >>> edges = [(0, 0), (0, 1), (0, 2), (1, 2), (2, 0), (2, 1), (2, 2)]
  643. >>> G = nx.DiGraph(edges)
  644. >>> nx.recursive_simple_cycles(G)
  645. [[0], [2], [0, 1, 2], [0, 2], [1, 2]]
  646. Notes
  647. -----
  648. The implementation follows pp. 79-80 in [1]_.
  649. The time complexity is $O((n+e)(c+1))$ for $n$ nodes, $e$ edges and $c$
  650. elementary circuits.
  651. References
  652. ----------
  653. .. [1] Finding all the elementary circuits of a directed graph.
  654. D. B. Johnson, SIAM Journal on Computing 4, no. 1, 77-84, 1975.
  655. https://doi.org/10.1137/0204007
  656. See Also
  657. --------
  658. simple_cycles, cycle_basis
  659. """
  660. # Jon Olav Vik, 2010-08-09
  661. def _unblock(thisnode):
  662. """Recursively unblock and remove nodes from B[thisnode]."""
  663. if blocked[thisnode]:
  664. blocked[thisnode] = False
  665. while B[thisnode]:
  666. _unblock(B[thisnode].pop())
  667. def circuit(thisnode, startnode, component):
  668. closed = False # set to True if elementary path is closed
  669. path.append(thisnode)
  670. blocked[thisnode] = True
  671. for nextnode in component[thisnode]: # direct successors of thisnode
  672. if nextnode == startnode:
  673. result.append(path[:])
  674. closed = True
  675. elif not blocked[nextnode]:
  676. if circuit(nextnode, startnode, component):
  677. closed = True
  678. if closed:
  679. _unblock(thisnode)
  680. else:
  681. for nextnode in component[thisnode]:
  682. if thisnode not in B[nextnode]: # TODO: use set for speedup?
  683. B[nextnode].append(thisnode)
  684. path.pop() # remove thisnode from path
  685. return closed
  686. path = [] # stack of nodes in current path
  687. blocked = defaultdict(bool) # vertex: blocked from search?
  688. B = defaultdict(list) # graph portions that yield no elementary circuit
  689. result = [] # list to accumulate the circuits found
  690. # Johnson's algorithm exclude self cycle edges like (v, v)
  691. # To be backward compatible, we record those cycles in advance
  692. # and then remove from subG
  693. for v in G:
  694. if G.has_edge(v, v):
  695. result.append([v])
  696. G.remove_edge(v, v)
  697. # Johnson's algorithm requires some ordering of the nodes.
  698. # They might not be sortable so we assign an arbitrary ordering.
  699. ordering = dict(zip(G, range(len(G))))
  700. for s in ordering:
  701. # Build the subgraph induced by s and following nodes in the ordering
  702. subgraph = G.subgraph(node for node in G if ordering[node] >= ordering[s])
  703. # Find the strongly connected component in the subgraph
  704. # that contains the least node according to the ordering
  705. strongcomp = nx.strongly_connected_components(subgraph)
  706. mincomp = min(strongcomp, key=lambda ns: min(ordering[n] for n in ns))
  707. component = G.subgraph(mincomp)
  708. if len(component) > 1:
  709. # smallest node in the component according to the ordering
  710. startnode = min(component, key=ordering.__getitem__)
  711. for node in component:
  712. blocked[node] = False
  713. B[node][:] = []
  714. dummy = circuit(startnode, startnode, component)
  715. return result
  716. def find_cycle(G, source=None, orientation=None):
  717. """Returns a cycle found via depth-first traversal.
  718. The cycle is a list of edges indicating the cyclic path.
  719. Orientation of directed edges is controlled by `orientation`.
  720. Parameters
  721. ----------
  722. G : graph
  723. A directed/undirected graph/multigraph.
  724. source : node, list of nodes
  725. The node from which the traversal begins. If None, then a source
  726. is chosen arbitrarily and repeatedly until all edges from each node in
  727. the graph are searched.
  728. orientation : None | 'original' | 'reverse' | 'ignore' (default: None)
  729. For directed graphs and directed multigraphs, edge traversals need not
  730. respect the original orientation of the edges.
  731. When set to 'reverse' every edge is traversed in the reverse direction.
  732. When set to 'ignore', every edge is treated as undirected.
  733. When set to 'original', every edge is treated as directed.
  734. In all three cases, the yielded edge tuples add a last entry to
  735. indicate the direction in which that edge was traversed.
  736. If orientation is None, the yielded edge has no direction indicated.
  737. The direction is respected, but not reported.
  738. Returns
  739. -------
  740. edges : directed edges
  741. A list of directed edges indicating the path taken for the loop.
  742. If no cycle is found, then an exception is raised.
  743. For graphs, an edge is of the form `(u, v)` where `u` and `v`
  744. are the tail and head of the edge as determined by the traversal.
  745. For multigraphs, an edge is of the form `(u, v, key)`, where `key` is
  746. the key of the edge. When the graph is directed, then `u` and `v`
  747. are always in the order of the actual directed edge.
  748. If orientation is not None then the edge tuple is extended to include
  749. the direction of traversal ('forward' or 'reverse') on that edge.
  750. Raises
  751. ------
  752. NetworkXNoCycle
  753. If no cycle was found.
  754. Examples
  755. --------
  756. In this example, we construct a DAG and find, in the first call, that there
  757. are no directed cycles, and so an exception is raised. In the second call,
  758. we ignore edge orientations and find that there is an undirected cycle.
  759. Note that the second call finds a directed cycle while effectively
  760. traversing an undirected graph, and so, we found an "undirected cycle".
  761. This means that this DAG structure does not form a directed tree (which
  762. is also known as a polytree).
  763. >>> G = nx.DiGraph([(0, 1), (0, 2), (1, 2)])
  764. >>> nx.find_cycle(G, orientation="original")
  765. Traceback (most recent call last):
  766. ...
  767. networkx.exception.NetworkXNoCycle: No cycle found.
  768. >>> list(nx.find_cycle(G, orientation="ignore"))
  769. [(0, 1, 'forward'), (1, 2, 'forward'), (0, 2, 'reverse')]
  770. See Also
  771. --------
  772. simple_cycles
  773. """
  774. if not G.is_directed() or orientation in (None, "original"):
  775. def tailhead(edge):
  776. return edge[:2]
  777. elif orientation == "reverse":
  778. def tailhead(edge):
  779. return edge[1], edge[0]
  780. elif orientation == "ignore":
  781. def tailhead(edge):
  782. if edge[-1] == "reverse":
  783. return edge[1], edge[0]
  784. return edge[:2]
  785. explored = set()
  786. cycle = []
  787. final_node = None
  788. for start_node in G.nbunch_iter(source):
  789. if start_node in explored:
  790. # No loop is possible.
  791. continue
  792. edges = []
  793. # All nodes seen in this iteration of edge_dfs
  794. seen = {start_node}
  795. # Nodes in active path.
  796. active_nodes = {start_node}
  797. previous_head = None
  798. for edge in nx.edge_dfs(G, start_node, orientation):
  799. # Determine if this edge is a continuation of the active path.
  800. tail, head = tailhead(edge)
  801. if head in explored:
  802. # Then we've already explored it. No loop is possible.
  803. continue
  804. if previous_head is not None and tail != previous_head:
  805. # This edge results from backtracking.
  806. # Pop until we get a node whose head equals the current tail.
  807. # So for example, we might have:
  808. # (0, 1), (1, 2), (2, 3), (1, 4)
  809. # which must become:
  810. # (0, 1), (1, 4)
  811. while True:
  812. try:
  813. popped_edge = edges.pop()
  814. except IndexError:
  815. edges = []
  816. active_nodes = {tail}
  817. break
  818. else:
  819. popped_head = tailhead(popped_edge)[1]
  820. active_nodes.remove(popped_head)
  821. if edges:
  822. last_head = tailhead(edges[-1])[1]
  823. if tail == last_head:
  824. break
  825. edges.append(edge)
  826. if head in active_nodes:
  827. # We have a loop!
  828. cycle.extend(edges)
  829. final_node = head
  830. break
  831. else:
  832. seen.add(head)
  833. active_nodes.add(head)
  834. previous_head = head
  835. if cycle:
  836. break
  837. else:
  838. explored.update(seen)
  839. else:
  840. assert len(cycle) == 0
  841. raise nx.exception.NetworkXNoCycle("No cycle found.")
  842. # We now have a list of edges which ends on a cycle.
  843. # So we need to remove from the beginning edges that are not relevant.
  844. for i, edge in enumerate(cycle):
  845. tail, head = tailhead(edge)
  846. if tail == final_node:
  847. break
  848. return cycle[i:]
  849. @not_implemented_for("directed")
  850. @not_implemented_for("multigraph")
  851. def minimum_cycle_basis(G, weight=None):
  852. """Returns a minimum weight cycle basis for G
  853. Minimum weight means a cycle basis for which the total weight
  854. (length for unweighted graphs) of all the cycles is minimum.
  855. Parameters
  856. ----------
  857. G : NetworkX Graph
  858. weight: string
  859. name of the edge attribute to use for edge weights
  860. Returns
  861. -------
  862. A list of cycle lists. Each cycle list is a list of nodes
  863. which forms a cycle (loop) in G. Note that the nodes are not
  864. necessarily returned in a order by which they appear in the cycle
  865. Examples
  866. --------
  867. >>> G = nx.Graph()
  868. >>> nx.add_cycle(G, [0, 1, 2, 3])
  869. >>> nx.add_cycle(G, [0, 3, 4, 5])
  870. >>> print([sorted(c) for c in nx.minimum_cycle_basis(G)])
  871. [[0, 1, 2, 3], [0, 3, 4, 5]]
  872. References:
  873. [1] Kavitha, Telikepalli, et al. "An O(m^2n) Algorithm for
  874. Minimum Cycle Basis of Graphs."
  875. http://link.springer.com/article/10.1007/s00453-007-9064-z
  876. [2] de Pina, J. 1995. Applications of shortest path methods.
  877. Ph.D. thesis, University of Amsterdam, Netherlands
  878. See Also
  879. --------
  880. simple_cycles, cycle_basis
  881. """
  882. # We first split the graph in connected subgraphs
  883. return sum(
  884. (_min_cycle_basis(G.subgraph(c), weight) for c in nx.connected_components(G)),
  885. [],
  886. )
  887. def _min_cycle_basis(comp, weight):
  888. cb = []
  889. # We extract the edges not in a spanning tree. We do not really need a
  890. # *minimum* spanning tree. That is why we call the next function with
  891. # weight=None. Depending on implementation, it may be faster as well
  892. spanning_tree_edges = list(nx.minimum_spanning_edges(comp, weight=None, data=False))
  893. edges_excl = [frozenset(e) for e in comp.edges() if e not in spanning_tree_edges]
  894. N = len(edges_excl)
  895. # We maintain a set of vectors orthogonal to sofar found cycles
  896. set_orth = [{edge} for edge in edges_excl]
  897. for k in range(N):
  898. # kth cycle is "parallel" to kth vector in set_orth
  899. new_cycle = _min_cycle(comp, set_orth[k], weight=weight)
  900. cb.append(list(set().union(*new_cycle)))
  901. # now update set_orth so that k+1,k+2... th elements are
  902. # orthogonal to the newly found cycle, as per [p. 336, 1]
  903. base = set_orth[k]
  904. set_orth[k + 1 :] = [
  905. orth ^ base if len(orth & new_cycle) % 2 else orth
  906. for orth in set_orth[k + 1 :]
  907. ]
  908. return cb
  909. def _min_cycle(G, orth, weight=None):
  910. """
  911. Computes the minimum weight cycle in G,
  912. orthogonal to the vector orth as per [p. 338, 1]
  913. """
  914. T = nx.Graph()
  915. nodes_idx = {node: idx for idx, node in enumerate(G.nodes())}
  916. idx_nodes = {idx: node for node, idx in nodes_idx.items()}
  917. nnodes = len(nodes_idx)
  918. # Add 2 copies of each edge in G to T. If edge is in orth, add cross edge;
  919. # otherwise in-plane edge
  920. for u, v, data in G.edges(data=True):
  921. uidx, vidx = nodes_idx[u], nodes_idx[v]
  922. edge_w = data.get(weight, 1)
  923. if frozenset((u, v)) in orth:
  924. T.add_edges_from(
  925. [(uidx, nnodes + vidx), (nnodes + uidx, vidx)], weight=edge_w
  926. )
  927. else:
  928. T.add_edges_from(
  929. [(uidx, vidx), (nnodes + uidx, nnodes + vidx)], weight=edge_w
  930. )
  931. all_shortest_pathlens = dict(nx.shortest_path_length(T, weight=weight))
  932. cross_paths_w_lens = {
  933. n: all_shortest_pathlens[n][nnodes + n] for n in range(nnodes)
  934. }
  935. # Now compute shortest paths in T, which translates to cyles in G
  936. start = min(cross_paths_w_lens, key=cross_paths_w_lens.get)
  937. end = nnodes + start
  938. min_path = nx.shortest_path(T, source=start, target=end, weight="weight")
  939. # Now we obtain the actual path, re-map nodes in T to those in G
  940. min_path_nodes = [node if node < nnodes else node - nnodes for node in min_path]
  941. # Now remove the edges that occur two times
  942. mcycle_pruned = _path_to_cycle(min_path_nodes)
  943. return {frozenset((idx_nodes[u], idx_nodes[v])) for u, v in mcycle_pruned}
  944. def _path_to_cycle(path):
  945. """
  946. Removes the edges from path that occur even number of times.
  947. Returns a set of edges
  948. """
  949. edges = set()
  950. for edge in pairwise(path):
  951. # Toggle whether to keep the current edge.
  952. edges ^= {edge}
  953. return edges