Finding all cycles in undirected graphs

Question

I need a working algorithm for finding all simple cycles in an undirected graph. I know the cost can be exponential and the problem is NP-complete, but I am going to use it in a small graph (up to 20-30 vertices) and the cycles are small in number.

After a long research (mainly here) I still don't have a working approach. Here is a summary of my search:

Finding all cycles in an undirected graph

Cycles in an Undirected Graph -> detects only whether there is a cycle or not

Finding polygons within an undirected Graph -> very nice description, but no solution

Finding all cycles in a directed graph -> finds cycles only in directed graphs

Detect cycles in undirected graph using boost graph library

The only answer I found, which approaches my problem, is this one:

Find all cycles in graph, redux

It seems that finding a basic set of cycles and XOR-ing them could do the trick. Finding a basic set of cycles is easy, but I don't understand how to combine them in order to obtain all cycles in the graph...

Answer 1

For an undirected graph the standard approach is to look for a so called cycle base : a set of simple cycles from which one can generate through combinations all other cycles. These are not necessarily all simple cycles in the graph. Consider for example the following graph:

    A 
  /   \
B ----- C
  \   /
    D

There are 3 simple cycles here : A-B-C-A, B-C-D-B and A-B-D-C-A. You can however take each 2 of these as a basis and obtain the 3rd as a combination of the 2. This is a substantial difference from directed graphs where one can not combine so freely cycles due to the need to observe edge direction.

The standard baseline algorithm for finding a cycle base for an undirected graph is this : Build a spanning tree and then for each edge which is not part of the tree build a cycle from that edge and some edges on the tree. Such cycle must exist because otherwise the edge would be part of the tree.

For example one of the possible spanning trees for the sample graph above is this:

    A 
  /   \
B      C
  \ 
    D

The 2 edges not in the tree are B-C and C-D. And the corresponding simple cycles are A-B-C-A and A-B-D-C-A.

You can also build the following spanning tree:

    A 
  /   
B ----- C
  \   
    D

And for this spanning tree the simple cycles would be A-B-C-A and B-C-D-B.

The baseline algorithm can be refined in different ways. To the best of my knowledge the best refinement belongs to Paton (K. Paton, An algorithm for finding a fundamental set of cycles for an undirected linear graph, Comm. ACM 12 (1969), pp. 514-518.). An open source implementation in Java is available here : http://code.google.com/p/niographs/ .

I should have mentioned how you combine simple cycles from the cycle base to form new simple cycles. You start off by listing in any (but fixed hereafter) order all edges of the graph. Then you represent cycles by sequences of zeros and ones by placing ones in the positions of edges which belong to the cycle and zeros in the positions of edges which are not part of the cycle. Then you do bitwise exclusive OR (XOR) of the sequences. The reason you do XOR is that you want to exclude edges which belong to both cycles and thus make the combined cycle non-simple. You need to check also that the 2 cycles have SOME common edges by checking that the bitwise AND of the sequences is not all zeros. Otherwise the result of XOR will be 2 disjoint cycles rather than a new simple cycle.

Here is an example for the sample graph above:

We start by listing the edges : ((AB), (AC), (BC), (BD), (CD)). Then the simple cycles A-B-C-A, B-D-C-B and A-B-D-C-A are represented as (1, 1, 1, 0, 0), (0, 0, 1, 1, 1) and (1, 1, 0, 1, 1). Now we can for example XOR A-B-C-A with B-D-C-B and the result is (1, 1, 0, 1, 1) which is exactly A-B-D-C-A. Or we can XOR A-B-C-A and A-B-D-C-A with the result being (0, 0, 1, 1, 1). Which is exactly B-D-C-B.

Given a cycle base you can discover all simple cycles by examining all possible combinations of 2 or more distinct base cycles. The procedure is described in more detail here : http://dspace.mit.edu/bitstream/handle/1721.1/68106/FTL_R_1982_07.pdf on page 14.

For the sake of completeness, I would notice that it seems possible (and inefficient) to use algorithms for finding all simple cycles of a directed graph. Every edge of the undirected graph can be replaced by 2 directed edges going in opposite directions. Then algorithms for directed graphs should work. There will be 1 "false" 2-node cycle for every edge of the undirected graph which will have to be ignored and there will be a clockwise and a counterclockwise version of every simple cycle of the undirected graph. Open source implementation in Java of algorithms for finding all cycles in a directed graph can be found at the link I already quoted.

Answer 2

Here is a vb .net version of the python code above:

Module Module1
'  Graph modelled as list of edges
Public graph As Integer(,) = {{{1, 2}, {1, 3}, {1, 4}, {2, 3},
        {3, 4}, {2, 6}, {4, 6}, {7, 8},
        {8, 9}, {9, 7}}

Public cycles As New List(Of Integer())()
Sub Main()
    For i As Integer = 0 To graph.GetLength(0) - 1
        For j As Integer = 0 To graph.GetLength(1) - 1
            findNewCycles(New Integer() {graph(i, j)})
        Next
    Next

    For Each cy As Integer() In cycles
        Dim s As String
        s = cy(0)
        For i As Integer = 1 To cy.Length - 1
            s = s & "," & cy(i)
        Next

        Console.WriteLine(s)
        Debug.Print(s)
    Next

End Sub
Private Sub findNewCycles(path As Integer())
    Dim n As Integer = path(0)
    Dim x As Integer
    Dim [sub] As Integer() = New Integer(path.Length) {}

    For i As Integer = 0 To graph.GetLength(0) - 1
        For y As Integer = 0 To 1
            If graph(i, y) = n Then
                '  edge referes to our current node
                x = graph(i, (y + 1) Mod 2)
                If Not visited(x, path) Then
                    '  neighbor node not on path yet
                    [sub](0) = x
                    Array.Copy(path, 0, [sub], 1, path.Length)
                    '  explore extended path
                    findNewCycles([sub])
                ElseIf (path.Length > 2) AndAlso (x = path(path.Length - 1)) Then
                    '  cycle found
                    Dim p As Integer() = normalize(path)
                    Dim inv As Integer() = invert(p)
                    If isNew(p) AndAlso isNew(inv) Then
                        cycles.Add(p)
                    End If
                End If
            End If
        Next
    Next
End Sub

Private Function equals(a As Integer(), b As Integer()) As Boolean
    Dim ret As Boolean = (a(0) = b(0)) AndAlso (a.Length = b.Length)

    Dim i As Integer = 1
    While ret AndAlso (i < a.Length)
        If a(i) <> b(i) Then
            ret = False
        End If
        i += 1
    End While

    Return ret
End Function

Private Function invert(path As Integer()) As Integer()
    Dim p As Integer() = New Integer(path.Length - 1) {}

    For i As Integer = 0 To path.Length - 1
        p(i) = path(path.Length - 1 - i)
    Next

    Return normalize(p)
End Function

'  rotate cycle path such that it begins with the smallest node
Private Function normalize(path As Integer()) As Integer()
    Dim p As Integer() = New Integer(path.Length - 1) {}
    Dim x As Integer = smallest(path)
    Dim n As Integer

    Array.Copy(path, 0, p, 0, path.Length)

    While p(0) <> x
        n = p(0)
        Array.Copy(p, 1, p, 0, p.Length - 1)
        p(p.Length - 1) = n
    End While

    Return p
End Function

Private Function isNew(path As Integer()) As Boolean
    Dim ret As Boolean = True

    For Each p As Integer() In cycles
        If equals(p, path) Then
            ret = False
            Exit For
        End If
    Next

    Return ret
End Function

Private Function smallest(path As Integer()) As Integer
    Dim min As Integer = path(0)

    For Each p As Integer In path
        If p < min Then
            min = p
        End If
    Next

    Return min
End Function

Private Function visited(n As Integer, path As Integer()) As Boolean
    Dim ret As Boolean = False

    For Each p As Integer In path
        If p = n Then
            ret = True
            Exit For
        End If
    Next

    Return ret
End Function

End Module

Answer 3

The Matlab version missed something, function findNewCycles(path) should be:

function findNewCycles(path)

global graph cycles numCycles;
startNode = path(1);
nextNode = nan;
sub = [];

% visit each edge and each node of each edge
for i = 1:size(graph,1)
    node1 = graph(i,1);
    node2 = graph(i,2);
    if (node1 == startNode) || (node2==startNode) %% this if is required
        if node1 == startNode
            nextNode = node2;
        elseif node2 == startNode
            nextNode = node1;
        end
        if ~(visited(nextNode, path))
            % neighbor node not on path yet
            sub = nextNode;
            sub = [sub path];
            % explore extended path
            findNewCycles(sub);
        elseif size(path,2) > 2 && nextNode == path(end)
            % cycle found
            p = rotate_to_smallest(path);
            inv = invert(p);
            if isNew(p) && isNew(inv)
                numCycles = numCycles + 1;
                cycles{numCycles} = p;
            end
        end
    end
end

Answer 4

Axel, I've translated your code to python. About 1/4th the lines of code and clearer to read.

graph = [[1, 2], [1, 3], [1, 4], [2, 3], [3, 4], [2, 6], [4, 6], [8, 7], [8, 9], [9, 7]]
cycles = []

def main():
    global graph
    global cycles
    for edge in graph:
        for node in edge:
            findNewCycles([node])
    for cy in cycles:
        path = [str(node) for node in cy]
        s = ",".join(path)
        print(s)

def findNewCycles(path):
    start_node = path[0]
    next_node= None
    sub = []

    #visit each edge and each node of each edge
    for edge in graph:
        node1, node2 = edge
        if start_node in edge:
                if node1 == start_node:
                    next_node = node2
                else:
                    next_node = node1
        if not visited(next_node, path):
                # neighbor node not on path yet
                sub = [next_node]
                sub.extend(path)
                # explore extended path
                findNewCycles(sub);
        elif len(path) > 2  and next_node == path[-1]:
                # cycle found
                p = rotate_to_smallest(path);
                inv = invert(p)
                if isNew(p) and isNew(inv):
                    cycles.append(p)

def invert(path):
    return rotate_to_smallest(path[::-1])

#  rotate cycle path such that it begins with the smallest node
def rotate_to_smallest(path):
    n = path.index(min(path))
    return path[n:]+path[:n]

def isNew(path):
    return not path in cycles

def visited(node, path):
    return node in path

main()

Answer 5

Here's just a very lame MATLAB version of this algorithm adapted from the python code above, for anyone who might need it as well.

function cycleList = searchCycles(edgeMap)

    tic
    global graph cycles numCycles;
    graph = edgeMap;
    numCycles = 0;
    cycles = {};
    for i = 1:size(graph,1)
        for j = 1:2
            findNewCycles(graph(i,j))
        end
    end
    % print out all found cycles
    for i = 1:size(cycles,2)
        cycles{i}
    end
    % return the result
    cycleList = cycles;
    toc

function findNewCycles(path)

    global graph cycles numCycles;
    startNode = path(1);
    nextNode = nan;
    sub = [];

    % visit each edge and each node of each edge
    for i = 1:size(graph,1)
        node1 = graph(i,1);
        node2 = graph(i,2);
        if node1 == startNode
            nextNode = node2;
        elseif node2 == startNode
            nextNode = node1;
        end
        if ~(visited(nextNode, path))
            % neighbor node not on path yet
            sub = nextNode;
            sub = [sub path];
            % explore extended path
            findNewCycles(sub);
        elseif size(path,2) > 2 && nextNode == path(end)
            % cycle found
            p = rotate_to_smallest(path);
            inv = invert(p);
            if isNew(p) && isNew(inv)
                numCycles = numCycles + 1;
                cycles{numCycles} = p;
            end
        end
    end

function inv = invert(path)
    inv = rotate_to_smallest(path(end:-1:1));

% rotate cycle path such that it begins with the smallest node
function new_path = rotate_to_smallest(path)
    [~,n] = min(path);
    new_path = [path(n:end), path(1:n-1)];

function result = isNew(path)
    global cycles
    result = 1;
    for i = 1:size(cycles,2)
        if size(path,2) == size(cycles{i},2) && all(path == cycles{i})
            result = 0;
            break;
        end
    end

function result = visited(node,path)
    result = 0;
    if isnan(node) && any(isnan(path))
        result = 1;
        return
    end
    for i = 1:size(path,2)
        if node == path(i)
            result = 1;
            break
        end
    end

Answer 6

Here is a C++ version of the python code above:

std::vector< std::vector<vertex_t> > Graph::findAllCycles()
{
    std::vector< std::vector<vertex_t> > cycles;

    std::function<void(std::vector<vertex_t>)> findNewCycles = [&]( std::vector<vertex_t> sub_path )
    {
        auto visisted = []( vertex_t v, const std::vector<vertex_t> & path ){
            return std::find(path.begin(),path.end(),v) != path.end();
        };

        auto rotate_to_smallest = []( std::vector<vertex_t> path ){
            std::rotate(path.begin(), std::min_element(path.begin(), path.end()), path.end());
            return path;
        };

        auto invert = [&]( std::vector<vertex_t> path ){
            std::reverse(path.begin(),path.end());
            return rotate_to_smallest(path);
        };

        auto isNew = [&cycles]( const std::vector<vertex_t> & path ){
            return std::find(cycles.begin(), cycles.end(), path) == cycles.end();
        };

        vertex_t start_node = sub_path[0];
        vertex_t next_node;

        // visit each edge and each node of each edge
        for(auto edge : edges)
        {
            if( edge.has(start_node) )
            {
                vertex_t node1 = edge.v1, node2 = edge.v2;

                if(node1 == start_node)
                    next_node = node2;
                else
                    next_node = node1;

                if( !visisted(next_node, sub_path) )
                {
                    // neighbor node not on path yet
                    std::vector<vertex_t> sub;
                    sub.push_back(next_node);
                    sub.insert(sub.end(), sub_path.begin(), sub_path.end());
                    findNewCycles( sub );
                } 
                else if( sub_path.size() > 2 && next_node == sub_path.back() )
                {
                    // cycle found
                    auto p = rotate_to_smallest(sub_path);
                    auto inv = invert(p);

                    if( isNew(p) && isNew(inv) )
                        cycles.push_back( p );
                }
            }
        }
    };

    for(auto edge : edges)
    {
        findNewCycles( std::vector<vertex_t>(1,edge.v1) );
        findNewCycles( std::vector<vertex_t>(1,edge.v2) );
    }
}