I would like to modify the networkx implementation of Johnson's algorithm for finding all elementary cycles in a graph (also copied below) so that is does not search for cycles larger than some maximum length.
def simple_cycles(G):
def _unblock(thisnode,blocked,B):
stack=set([thisnode])
while stack:
node=stack.pop()
if node in blocked:
blocked.remove(node)
stack.update(B[node])
B[node].clear()
# Johnson's algorithm requires some ordering of the nodes.
# We assign the arbitrary ordering given by the strongly connected comps
# There is no need to track the ordering as each node removed as processed.
subG = type(G)(G.edges_iter()) # save the actual graph so we can mutate it here
# We only take the edges because we do not want to
# copy edge and node attributes here.
sccs = list(nx.strongly_connected_components(subG))
while sccs:
scc=sccs.pop()
# order of scc determines ordering of nodes
startnode = scc.pop()
# Processing node runs "circuit" routine from recursive version
path=[startnode]
blocked = set() # vertex: blocked from search?
closed = set() # nodes involved in a cycle
blocked.add(startnode)
B=defaultdict(set) # graph portions that yield no elementary circuit
stack=[ (startnode,list(subG[startnode])) ] # subG gives component nbrs
while stack:
thisnode,nbrs = stack[-1]
if nbrs:
nextnode = nbrs.pop()
# print thisnode,nbrs,":",nextnode,blocked,B,path,stack,startnode
# f=raw_input("pause")
if nextnode == startnode:
yield path[:]
closed.update(path)
# print "Found a cycle",path,closed
elif nextnode not in blocked:
path.append(nextnode)
stack.append( (nextnode,list(subG[nextnode])) )
closed.discard(nextnode)
blocked.add(nextnode)
continue
# done with nextnode... look for more neighbors
if not nbrs: # no more nbrs
if thisnode in closed:
_unblock(thisnode,blocked,B)
else:
for nbr in subG[thisnode]:
if thisnode not in B[nbr]:
B[nbr].add(thisnode)
stack.pop()
assert path[-1]==thisnode
path.pop()
# done processing this node
subG.remove_node(startnode)
H=subG.subgraph(scc) # make smaller to avoid work in SCC routine
sccs.extend(list(nx.strongly_connected_components(H)))
Of course, I'd also accept a suggestion that differs from the implementation above but runs in similar time. Also, my project uses networkx, so feel free to use any other function from that library, such as shortest_path.
(Note: not homework!)
Edit
Dorijan Cirkveni suggested (if I understood correctly):
if len(blocked) >= limit + 1:
continue
elif nextnode == startnode:
yield path[:]
However, that doesn't work. Here's a counterexample:
G = nx.DiGraph()
G.add_edge(1, 2)
G.add_edge(2, 3)
G.add_edge(3, 1)
G.add_edge(3, 2)
G.add_edge(3, 4)
my_cycles = list(simple_cycles(G, limit = 3)) # Modification
nx_cycles = list(nx.simple_cycles(G)) # Original networkx code
print("MY:", my_cycles)
print("NX:", nx_cycles)
Will output
MY: [[2, 3]]
NX: [[1, 2, 3], [2, 3]]
Also, if we substitute blocked by stack or path, the result will be correct for this example, but will give the wrong answer for other graphs.