diff --git a/src/doc/en/reference/graphs/index.rst b/src/doc/en/reference/graphs/index.rst index 60044beb6e6..0485f3a6565 100644 --- a/src/doc/en/reference/graphs/index.rst +++ b/src/doc/en/reference/graphs/index.rst @@ -61,6 +61,7 @@ Libraries of algorithms sage/graphs/graph_coloring sage/graphs/cliquer + sage/graphs/independent_sets sage/graphs/comparability sage/graphs/line_graph sage/graphs/spanning_tree diff --git a/src/module_list.py b/src/module_list.py index ecffd2608d5..b4a18427f4a 100755 --- a/src/module_list.py +++ b/src/module_list.py @@ -372,6 +372,9 @@ def uname_specific(name, value, alternative): sources = ['sage/graphs/cliquer.pyx', 'sage/graphs/cliquer/cl.c'], libraries = ['cliquer']), + Extension('sage.graphs.independent_sets', + sources = ['sage/graphs/independent_sets.pyx']), + Extension('sage.graphs.graph_decompositions.vertex_separation', sources = ['sage/graphs/graph_decompositions/vertex_separation.pyx']), diff --git a/src/sage/graphs/base/static_dense_graph.pxd b/src/sage/graphs/base/static_dense_graph.pxd new file mode 100644 index 00000000000..aeec83609a1 --- /dev/null +++ b/src/sage/graphs/base/static_dense_graph.pxd @@ -0,0 +1,3 @@ +include "sage/misc/binary_matrix.pxi" + +cdef dict dense_graph_init(binary_matrix_t m, g, translation = ?) diff --git a/src/sage/graphs/base/static_dense_graph.pyx b/src/sage/graphs/base/static_dense_graph.pyx index 2ff09a3b415..6e6867583aa 100644 --- a/src/sage/graphs/base/static_dense_graph.pyx +++ b/src/sage/graphs/base/static_dense_graph.pyx @@ -40,7 +40,6 @@ Index Functions --------- """ -include "sage/misc/binary_matrix.pxi" cdef dict dense_graph_init(binary_matrix_t m, g, translation = False): r""" diff --git a/src/sage/graphs/graph.py b/src/sage/graphs/graph.py index ddfcc5c2e3c..73c31cae0a8 100644 --- a/src/sage/graphs/graph.py +++ b/src/sage/graphs/graph.py @@ -4830,7 +4830,6 @@ def topological_minor(self, H, vertices = False, paths = False, solver=None, ver v_repr = p.get_values(v_repr) flow = p.get_values(flow) - for u,v in minor.edges(labels = False): used = False for C in H.edges(labels = False): @@ -4856,16 +4855,33 @@ def topological_minor(self, H, vertices = False, paths = False, solver=None, ver ### Cliques - def cliques_maximal(self): + def def cliques_maximal(self, algorithm = "native"): """ Returns the list of all maximal cliques, with each clique represented by a list of vertices. A clique is an induced complete subgraph, and a maximal clique is one not contained in a larger one. + INPUT: + + - ``algorithm`` -- can be set to ``"native"`` (default) to use Sage's + own implementation, or to ``"NetworkX"`` to use NetworkX' + implementation of the Bron and Kerbosch Algorithm [BroKer1973]_. + + .. NOTE:: - Currently only implemented for undirected graphs. Use to_undirected - to convert a digraph to an undirected graph. + This method sorts its output before returning it. If you prefer to + save the extra time, you can call + :class:`sage.graphs.independent_sets.IndependentSets` directly. + + .. NOTE:: + + Sage's implementation of the enumeration of *maximal* independent + sets is not much faster than NetworkX' (expect a 2x speedup), which + is surprising as it is written in Cython. This being said, the + algorithm from NetworkX appears to be sligthly different from this + one, and that would be a good thing to explore if one wants to + improve the implementation. ALGORITHM: @@ -4882,21 +4898,37 @@ def cliques_maximal(self): EXAMPLES:: sage: graphs.ChvatalGraph().cliques_maximal() - [[0, 1], [0, 4], [0, 6], [0, 9], [2, 1], [2, 3], [2, 6], [2, 8], [3, 4], [3, 7], [3, 9], [5, 1], [5, 4], [5, 10], [5, 11], [7, 1], [7, 8], [7, 11], [8, 4], [8, 10], [10, 6], [10, 9], [11, 6], [11, 9]] + [[0, 1], [0, 4], [0, 6], [0, 9], [1, 2], [1, 5], [1, 7], [2, 3], + [2, 6], [2, 8], [3, 4], [3, 7], [3, 9], [4, 5], [4, 8], [5, 10], + [5, 11], [6, 10], [6, 11], [7, 8], [7, 11], [8, 10], [9, 10], [9, 11]] sage: G = Graph({0:[1,2,3], 1:[2], 3:[0,1]}) sage: G.show(figsize=[2,2]) sage: G.cliques_maximal() [[0, 1, 2], [0, 1, 3]] sage: C=graphs.PetersenGraph() sage: C.cliques_maximal() - [[0, 1], [0, 4], [0, 5], [2, 1], [2, 3], [2, 7], [3, 4], [3, 8], [6, 1], [6, 8], [6, 9], [7, 5], [7, 9], [8, 5], [9, 4]] + [[0, 1], [0, 4], [0, 5], [1, 2], [1, 6], [2, 3], [2, 7], [3, 4], + [3, 8], [4, 9], [5, 7], [5, 8], [6, 8], [6, 9], [7, 9]] sage: C = Graph('DJ{') sage: C.cliques_maximal() - [[4, 0], [4, 1, 2, 3]] + [[0, 4], [4, 1, 2, 3]] + Comparing the two implementations:: + + sage: g = graphs.RandomGNP(20,.7) + sage: s1 = Set(map(Set, g.cliques_maximal(algorithm="NetworkX"))) + sage: s2 = Set(map(Set, g.cliques_maximal(algorithm="native"))) + sage: s1 == s2 + True """ - import networkx - return sorted(networkx.find_cliques(self.networkx_graph(copy=False))) + if algorithm == "native": + from sage.graphs.independent_sets import IndependentSets + return sorted(IndependentSets(self, maximal = True, complement = True)) + elif algorithm == "NetworkX": + import networkx + return sorted(networkx.find_cliques(self.networkx_graph(copy=False))) + else: + raise ValueError("Algorithm must be equal to 'native' or to 'NetworkX'.") cliques = deprecated_function_alias(5739, cliques_maximal) @@ -5060,7 +5092,7 @@ def cliques_number_of(self, vertices=None, cliques=None): {0: 1, 1: 1, 2: 1, 3: 1, 4: 2} sage: E = C.cliques_maximal() sage: E - [[4, 0], [4, 1, 2, 3]] + [[0, 4], [1, 2, 3, 4]] sage: C.cliques_number_of(cliques=E) {0: 1, 1: 1, 2: 1, 3: 1, 4: 2} sage: F = graphs.Grid2dGraph(2,3) @@ -5089,6 +5121,8 @@ def cliques_get_max_clique_graph(self, name=''): edges between maximal cliques with common members in the original graph. + For more information, see the :wikipedia:`Clique_graph`. + .. NOTE:: Currently only implemented for undirected graphs. Use to_undirected @@ -5150,6 +5184,10 @@ def independent_set(self, algorithm = "Cliquer", value_only = False, reduction_r Equivalently, an independent set is defined as the complement of a vertex cover. + For more information, see the + :wikipedia:`Independent_set_(graph_theory)` and the + :wikipedia:`Vertex_cover`. + INPUT: - ``algorithm`` -- the algorithm to be used @@ -5523,7 +5561,7 @@ def cliques_vertex_clique_number(self, algorithm="cliquer", vertices=None, {0: 2, 1: 4, 2: 4, 3: 4, 4: 4} sage: E = C.cliques_maximal() sage: E - [[4, 0], [4, 1, 2, 3]] + [[0, 4], [1, 2, 3, 4]] sage: C.cliques_vertex_clique_number(cliques=E,algorithm="networkx") {0: 2, 1: 4, 2: 4, 3: 4, 4: 4} sage: F = graphs.Grid2dGraph(2,3) @@ -5586,9 +5624,9 @@ def cliques_containing_vertex(self, vertices=None, cliques=None): {0: [[4, 0]], 1: [[4, 1, 2, 3]], 2: [[4, 1, 2, 3]], 3: [[4, 1, 2, 3]], 4: [[4, 0], [4, 1, 2, 3]]} sage: E = C.cliques_maximal() sage: E - [[4, 0], [4, 1, 2, 3]] + [[0, 4], [1, 2, 3, 4]] sage: C.cliques_containing_vertex(cliques=E) - {0: [[4, 0]], 1: [[4, 1, 2, 3]], 2: [[4, 1, 2, 3]], 3: [[4, 1, 2, 3]], 4: [[4, 0], [4, 1, 2, 3]]} + {0: [[0, 4]], 1: [[1, 2, 3, 4]], 2: [[1, 2, 3, 4]], 3: [[1, 2, 3, 4]], 4: [[0, 4], [1, 2, 3, 4]]} sage: F = graphs.Grid2dGraph(2,3) sage: X = F.cliques_containing_vertex() sage: for v in sorted(X.iterkeys()): diff --git a/src/sage/graphs/independent_sets.pxd b/src/sage/graphs/independent_sets.pxd new file mode 100644 index 00000000000..e2a63de61c5 --- /dev/null +++ b/src/sage/graphs/independent_sets.pxd @@ -0,0 +1,14 @@ +from libc.stdint cimport uint32_t, uint64_t +include "sage/misc/binary_matrix.pxi" +include "sage/misc/bitset_pxd.pxi" +#include "sage/misc/bitset_pxd.pxi" + +cdef class IndependentSets: + cdef binary_matrix_t g + cdef list vertices + cdef dict vertex_to_int + cdef int n + cdef int i + cdef int count_only + cdef int maximal + diff --git a/src/sage/graphs/independent_sets.pyx b/src/sage/graphs/independent_sets.pyx new file mode 100644 index 00000000000..f04bc1488af --- /dev/null +++ b/src/sage/graphs/independent_sets.pyx @@ -0,0 +1,353 @@ +r""" +Independent sets + +This module implements the :class:`IndependentSets` class which can be used to : + +- List the independent sets (or cliques) of a graph +- Count them (which is obviously faster) +- Test whether a set of vertices is an independent set + +It can also be restricted to focus on (inclusionwise) maximal independent +sets. See the documentation of :class:`IndependentSets` for actual examples. + +Classes and methods +------------------- + +""" + +include "sage/misc/bitset.pxi" +from sage.misc.cachefunc import cached_method +from sage.graphs.base.static_dense_graph cimport dense_graph_init + +cdef inline int bitset_are_disjoint(unsigned long * b1, unsigned long * b2, int width): + r""" + Tests whether two bitsets of length width*sizeof(unsigned int) have an empty + intersection. + """ + cdef int i + for i in range(width): + if b1[i]&b2[i]: + return False + return True + +cdef inline int ismaximal(binary_matrix_t g, int n, bitset_t s): + cdef int i + for i in range(n): + if (not bitset_in(s,i)) and bitset_are_disjoint(g.rows[i],s.bits,g.width): + return False + + return True + +cdef class IndependentSets: + r""" + An instance of this class represents the set of independent sets of a graph. + + For more information on independent sets, see the + :wikipedia:`Independent_set_(graph_theory)`. + + INPUT: + + - ``G`` -- a graph + + - ``maximal`` (boolean) -- whether to only consider (inclusionwise) maximal + independent sets. Set to ``False`` by default. + + - ``complement`` (boolean) -- whether to consider the graph's complement + (i.e. cliques instead of independent sets). Set to ``False`` by default. + + ALGORITHM: + + The enumeration of independent sets is done naively : given an independent + sets, this implementation considers all ways to add a new vertex to it + (while keeping it an independent sets), and then creates new independent + sets from all those that were created this way. + + The implementation, however, is not recursive. + + .. NOTE:: + + This implementation of the enumeration of *maximal* independent sets is + not much faster than NetworkX', which is surprising as it is written in + Cython. This being said, the algorithm from NetworkX appears to be + sligthly different from this one, and that would be a good thing to + explore if one wants to improve the implementation. + + EXAMPLES: + + Listing all independent sets of the Claw graph:: + + sage: from sage.graphs.independent_sets import IndependentSets + sage: g = graphs.ClawGraph() + sage: list(IndependentSets(g)) + [[0], [1], [1, 2], [1, 2, 3], [1, 3], [2], [2, 3], [3], []] + + Count them:: + + sage: IndependentSets(g).cardinality() + 9 + + List only the maximal independent sets:: + + sage: list(IndependentSets(g, maximal = True)) + [[0], [1, 2, 3]] + + And count them:: + + sage: IndependentSets(g, maximal = True).cardinality() + 2 + + One can easily obtain an iterator over all independent sets of given + cardinality:: + + sage: g = graphs.PetersenGraph() + sage: is4 = (x for x in IndependentSets(g) if len(x) == 4) + sage: list(is4) + [[0, 2, 8, 9], [0, 3, 6, 7], [1, 3, 5, 9], [1, 4, 7, 8], [2, 4, 5, 6]] + + Similarly, one can easily count the number of independent sets of each + cardinality:: + + sage: number_of_is = [0] * g.order() + sage: for x in IndependentSets(g): + ... number_of_is[len(x)] += 1 + sage: print number_of_is + [1, 10, 30, 30, 5, 0, 0, 0, 0, 0] + """ + def __init__(self, G, maximal = False, complement = False): + r""" + Constructor for this class + + TESTS:: + + sage: from sage.graphs.independent_sets import IndependentSets + sage: IndependentSets(graphs.PetersenGraph()) + i are zero. + + while True: + + # If i is in current_set + if bitset_in(current_set,i): + + # We have found an independent set ! + if bitset_are_disjoint(self.g.rows[i],current_set.bits,self.g.width): + + # Saving that set + bitset_copy(tmp, current_set) + + # Preparing for the next set, except if we set the last bit. + if i < self.n-1: + + # Adding (i+1)th bit + bitset_add(current_set,i+1) + i += 1 + else: + bitset_discard(current_set,i) + + # Returning the result if necessary ... + + if self.maximal and not ismaximal(self.g,self.n, tmp): + continue + + count += 1 + + if not self.count_only: + yield [self.vertices[j] for j in range(i+1) if bitset_in(tmp,j)] + continue + + else: + # Removing the ith bit + bitset_discard(current_set, i) + + # Preparing for the next set ! + if i < self.n-1: + bitset_add(current_set, i+1) + i += 1 + + # Not already included in the set + else: + if i == 0: + if not self.maximal: + count+=1 + yield [] + + break + + # Going backward, we explored all we could there ! + if bitset_in(current_set,i-1): + bitset_discard(current_set, i-1) + bitset_add(current_set,i) + else: + i -= 1 + if i == -1: + break + + + if self.count_only: + yield count + + bitset_free(current_set) + bitset_free(tmp) + + + def __dealloc__(self): + r""" + Frees everything we ever allocated + """ + binary_matrix_free(self.g) + + @cached_method + def cardinality(self): + r""" + Computes and returns the number of independent sets + + TESTS:: + + sage: from sage.graphs.independent_sets import IndependentSets + sage: IndependentSets(graphs.PetersenGraph()).cardinality() + 76 + + Only maximal ones:: + + sage: from sage.graphs.independent_sets import IndependentSets + sage: IndependentSets(graphs.PetersenGraph(), maximal = True).cardinality() + 15 + """ + self.count_only = 1 + + for i in self: + pass + + self.count_only = 0 + + from sage.rings.integer import Integer + return Integer(i) + + def __contains__(self, S): + r""" + Checks whether the set is an independent set (possibly maximal) + + - ``S`` -- a set of vertices to be tested. + + EXAMPLE:: + + All independent sets of PetersenGraph are ... independent sets :: + + sage: from sage.graphs.independent_sets import IndependentSets + sage: IS = IndependentSets(graphs.PetersenGraph()) + sage: all( s in IS for s in IS) + True + + Same with maximal independent sets:: + + sage: IS = IndependentSets(graphs.PetersenGraph(), maximal = True) + sage: all( s in IS for s in IS) + True + """ + # Set of vertices as a bitset + cdef bitset_t s + + bitset_init(s,self.n) + bitset_set_first_n(s,0) + + cdef int i + for I in S: + try: + i = self.vertex_to_int[I] + except KeyError: + raise ValueError("An element of the set being tested does not belong to ") + + # Adding the new vertex to s + bitset_add(s, i) + + # Checking that the set s is independent + if not bitset_are_disjoint(self.g.rows[i],s.bits,self.g.width): + return False + + if self.maximal and not ismaximal(self.g, self.n,s): + return False + + return True diff --git a/src/sage/misc/binary_matrix.pxi b/src/sage/misc/binary_matrix.pxi index 60eacb03500..af145f27ab0 100644 --- a/src/sage/misc/binary_matrix.pxi +++ b/src/sage/misc/binary_matrix.pxi @@ -22,6 +22,7 @@ a ``binary_matrix_t`` structure contains : width``. """ include "binary_matrix_pxd.pxi" +include "sage/ext/stdsage.pxi" cdef inline binary_matrix_init(binary_matrix_t m, long n_rows, long n_cols): r""" @@ -56,7 +57,20 @@ cdef inline binary_matrix_fill(binary_matrix_t m, bint bit): r""" Fill the whole matrix with a bit """ - memset(m.rows[0],-( bit),m.width * m.n_rows * sizeof(unsigned long)) + memset(m.rows[0],-( bit),m.width * m.n_rows * sizeof(unsigned long)) + +cdef inline binary_matrix_complement(binary_matrix_t m): + r""" + Complements all of the matrix' bits. + """ + cdef int i,j + for i in range(m.n_rows): + for j in range(m.width): + m.rows[i][j] = ~m.rows[i][j] + + # Set the "useless" bits to 0 too. + if ((m.n_cols) & offset_mask): + m.rows[i][m.width-1] &= ( 1<<((m.n_cols) & offset_mask))-1 cdef inline binary_matrix_set1(binary_matrix_t m, long row, long col): r""" @@ -90,6 +104,10 @@ cdef inline binary_matrix_print(binary_matrix_t m): cdef int i,j import sys for i in range(m.n_rows): + # If you want to print the *whole* matrix, including the useless bits, + # use the following line instead + # + # for j in (m.width*8*sizeof(unsigned long)): for j in range(m.n_cols): sys.stdout.write("1" if binary_matrix_get(m, i, j) else ".",) print "" diff --git a/src/sage/misc/binary_matrix_pxd.pxi b/src/sage/misc/binary_matrix_pxd.pxi index 99ab2d0680d..ac106fa0561 100644 --- a/src/sage/misc/binary_matrix_pxd.pxi +++ b/src/sage/misc/binary_matrix_pxd.pxi @@ -3,6 +3,7 @@ include 'sage/ext/cdefs.pxi' cdef extern from *: int __builtin_popcountl(unsigned long) + void *memset(void *, int, size_t) # Constants from bitset.pxd cdef extern from *: