Drastic performance improvement of computing the chromatic polynomial

[sagetrac-azi]

Hello!

I am not sure I have the energy to do that at the moment but its a nice task that we might want to address at some point. Hence I am leaving some stuff here for further reference.

The current implementation of the chromatic_polynomial uses the deletion-contraction recurrence directly. This gives us an exponential branching tree. The neat thing is that at some point some of the branches compute the chromatic polynomial for the same graphs (multiple - exponential) times!

A very dumb implementation in Python (see below) that is faster than the C code. It takes 700us to compute the chromatic polynomial of a random dense graph of order 15 and 74 minutes (!) with the C code currently used by Sage.

Of course this code needs more work and can be optimized even further (perhaps by calling the C implementation under a threshold) but for now here it is.

global cache
cache = {}
def chrompoly(G):
    
    global cache 

    s = tuple(G.canonical_label().edges(labels=False))

    if s in cache:
        return cache[s]
    if not G.is_connected():
        return prod([chrompoly(C) for C in G.connected_components_subgraphs()])

    R = ZZ['x']
    x = R.gen()
    
    n = G.order()
    m = G.size()

    #  Since we know that G is connected it is enough to check the number of edges
    #  to test if it is a tree
    if m == n-1:
        return x*(x-1)**(n-1)

    u,v = G.edge_iterator(labels=False).next()

    G.delete_edge(u,v)
    p = chrompoly(G)
    G.add_edge(u,v)

    G2 = G.copy()
    G2.add_edges(((u,x) for x in G.neighbor_iterator(v)))
    G2.delete_vertex(v)

    ret = p - chrompoly(G2)
    
    cache[s] = ret 

    return ret 

CC: @nathanncohen @Stefan @lobiCode @rbeezer @orlitzky

Component: graph theory

Author: Jernej Azarija, David Coudert

Branch/Commit: 66708cc

Reviewer: Michael Orlitzky

Issue created by migration from https://trac.sagemath.org/ticket/14529

[sagetrac-azi]

comment:1

Hm... Maybe I am being too fast with the time results. I should really check what happens if I delete the cache after each run

[sagetrac-azi]

comment:2

Looks fine!!

sage: G = graphs.RandomGNP(14,0.8)
sage: %timeit G.chromatic_polynomial()
1 loops, best of 3: 50.7 s per loop

sage: %timeit chrompoly(G)
1 loops, best of 3: 790 us per loop

[sagetrac-azi]

Description changed:

--- 
+++ 
@@ -4,7 +4,7 @@
 
 The current implementation of the chromatic_polynomial uses the deletion-contraction recurrence directly. This gives us an exponential branching tree. The neat thing is that at some point some of the branches compute the chromatic polynomial for the same graphs (multiple - exponential) times!
 
-A very dumb implementation in Python (see below) is already faster on graphs of order 13 than the C code. It takes 2ms to compute the chromatic polynomial with this code and 1 minute with the one we have now!
+A very dumb implementation in Python (see below) that is faster than the C code. It takes 700us to compute the chromatic polynomial of a random dense graph of order 15 and 74 minutes (!) with the C code currently used by Sage.
 
 Of course this code needs more work and can be optimized even further  (perhaps by calling the C implementation under a threshold) but for now here it is.
 
@@ -13,37 +13,40 @@
 global cache
 cache = {}
 def chrompoly(G):
-        
+    
     global cache 
 
-    s = G.canonical_label().sparse6_string()
+    s = tuple(G.canonical_label().edges(labels=False))
 
     if s in cache:
         return cache[s]
     if not G.is_connected():
         return prod([chrompoly(C) for C in G.connected_components_subgraphs()])
+
     R = ZZ['x']
     x = R.gen()
-    if G.is_tree():
-        return x*(x-1)**(G.order()-1)
+    
+    n = G.order()
+    m = G.size()
 
-    u,v = G.edges(labels=False)[0]    
+    #  Since we know that G is connected it is enough to check the number of edges
+    #  to test if it is a tree
+    if m == n-1:
+        return x*(x-1)**(n-1)
+
+    u,v = G.edge_iterator(labels=False).next()
 
     G.delete_edge(u,v)
-
-    p1 = chrompoly(G)
-
+    p = chrompoly(G)
     G.add_edge(u,v)
 
-    ret = p1 - chrompoly(contract_edge(G, u,v))
-        
+    G2 = G.copy()
+    G2.add_edges(((u,x) for x in G.neighbor_iterator(v)))
+    G2.delete_vertex(v)
+
+    ret = p - chrompoly(G2)
+    
     cache[s] = ret 
 
     return ret 
-
-def contract_edge(G, u,v):
-    Gret = G.copy()
-    Gret.add_edges([(u,x) for x in G.neighbors(v)])
-    Gret.delete_vertex(v)
-    return Gret
 ```

[sagetrac-azi]

comment:4

Does anyone of you guys sees a reason not to remove the current Cython implementation and adding this code

[nathanncohen]

comment:5

Yep : the authors of the Cython code may be so impressed by your caching that they may want to adapt their code to it. Worth sending an email. Otherwise, if the results are better and the code clearer.... :-P

Nathann

[sagetrac-azi]

comment:6

Good!

How do you think we should handle the global cache if this is to be made into a proper patch?

What I would like to do is have a global ____sage_very_very_private__cache__chrompoly___dict = {} variable outside the function and delete the cache each time the computation is over. In addition I would have a parameter keep_cache that keeps the cache alive (convenient when you're computing the polynomial in a loop)

What do you think?

Also whom would you send an email to?

[nathanncohen]

comment:7

Helloooooooo !

How do you think we should handle the global cache if this is to be made into a proper patch?

Hmmm... Well, there is the @cached_method decorator, that automatically manages the cache of any function. But you can't empty it automatically when the computations are done though.

Technically, what you can do is write a function A (without the decorator) which call function B (the real one, with the cached_method decorator). Before returning B's result, A would clear the cache of B.

I wouldn't go to such lengths before there is a real need, actually... If somebody computes one chromatic polynomial perhaps he will compute two, and ... well. It's up to you. To me caching the method with the cached_method decorator would do the job. and you can add in its documentation an explanation of how the cache can be cleared manually, with the clear_cache method.

http://www.sagemath.org/doc/reference/misc/sage/misc/cachefunc.html

Also whom would you send an email to?

The file seems to have been written by Robert Miller and Gordon Royle (just reading the copyright in the file itself). Robert Miller does not work on Sage anymore but should be interested (he's the one who initially wrote all the graph-related stuff in Sage. Backends included). And Gordon Royle may be interested too, though I probably never wrote to him :-)

Nathann

[sagetrac-azi]

comment:8

Great thanks for the info! I'll read through that stuff and produce a proper patch.

As for the "compute one or two,.." I would for example like to test the Birkhoff-Lewis conjecture for all planar graphs of some order and having such a cache alive through all the computation would be very nice.