trac #18864: Merged with 6.8.beta8

src/sage/graphs/distances_all_pairs.pxd

@@ -1,8 +1,10 @@
+from libc.stdint cimport uint32_t
+
 cdef unsigned short * c_shortest_path_all_pairs(G) except NULL
 cdef unsigned short * c_distances_all_pairs(G)
 cdef all_pairs_shortest_path_BFS(gg,
                                  unsigned short * predecessors,
                                  unsigned short * distances,
-                                 int            * eccentricity)
+                                 uint32_t       * eccentricity)
 
-cdef int * c_eccentricity(G) except NULL
+cdef uint32_t * c_eccentricity(G) except NULL

src/sage/graphs/distances_all_pairs.pyx

@@ -131,6 +131,9 @@ REFERENCE:
   tight bounds for the diameter of massive graphs. *ACM Journal of Experimental
   Algorithms* 13 (2008) http://dx.doi.org/10.1145/1412228.1455266
 
+.. [FK13] F. W. Takes and W. A. Kosters. Computing the eccentricity distribution
+  of large graphs. *Algorithms* 6:100-118 (2013)
+  http://dx.doi.org/10.3390/a6010100
 
 Functions
 ---------
@@ -152,13 +155,13 @@ from libc.stdint cimport uint64_t, uint32_t, INT32_MAX, UINT32_MAX
 from sage.graphs.base.c_graph cimport CGraphBackend
 from sage.graphs.base.c_graph cimport CGraph
 
-from sage.graphs.base.static_sparse_graph cimport short_digraph, init_short_digraph, free_short_digraph
+from sage.graphs.base.static_sparse_graph cimport short_digraph, init_short_digraph, free_short_digraph, out_degree
 
 
 cdef inline all_pairs_shortest_path_BFS(gg,
                                         unsigned short * predecessors,
                                         unsigned short * distances,
-                                        int            * eccentricity):
+                                        uint32_t       * eccentricity):
     """
     See the module's documentation.
     """
@@ -196,7 +199,6 @@ cdef inline all_pairs_shortest_path_BFS(gg,
     cdef uint32_t * end
 
     cdef unsigned short * c_predecessors = predecessors
-    cdef int * c_eccentricity = eccentricity
     cdef int * c_distances = <int *> sage_malloc( n * sizeof(int))
     if c_distances == NULL:
         sage_free(waiting_list)
@@ -269,20 +271,25 @@ cdef inline all_pairs_shortest_path_BFS(gg,
             waiting_beginning += 1
 
         # If not all the vertices have been met
-        for v in range(n):
-            if not bitset_in(seen, v):
+        if bitset_len(seen) < n:
+            bitset_complement(seen, seen)
+            v = bitset_next(seen, 0)
+            while v >= 0:
                 c_distances[v] = INT32_MAX
                 if predecessors != NULL:
                     c_predecessors[v] = -1
+                v = bitset_next(seen, v+1)
+
+            if eccentricity != NULL:
+                eccentricity[source] = UINT32_MAX
+
+        elif eccentricity != NULL:
+            eccentricity[source] = c_distances[waiting_list[n-1]]
+
 
         if predecessors != NULL:
             c_predecessors += n
 
-        if eccentricity != NULL:
-            c_eccentricity[source] = 0
-            for i in range(n):
-                c_eccentricity[source] = c_eccentricity[source] if c_eccentricity[source] >= c_distances[i] else c_distances[i]
-
         if distances != NULL:
             for i in range(n):
                 distances[i] = <unsigned short> c_distances[i]
@@ -710,25 +717,112 @@ def distances_and_predecessors_all_pairs(G):
 # Eccentricity #
 ################
 
-cdef int * c_eccentricity(G) except NULL:
+cdef uint32_t * c_eccentricity(G) except NULL:
     r"""
-    Returns the vector of eccentricities in G.
+    Return the vector of eccentricities in G.
 
     The array returned is of length n, and its ith component is the eccentricity
     of the ith vertex in ``G.vertices()``.
     """
     cdef unsigned int n = G.order()
 
-    cdef int * ecc = <int *> sage_malloc(n*sizeof(int))
+    cdef uint32_t * ecc = <uint32_t *> sage_calloc(n, sizeof(uint32_t))
     if ecc == NULL:
         raise MemoryError()
     all_pairs_shortest_path_BFS(G, NULL, NULL, ecc)
 
     return ecc
 
-def eccentricity(G):
+cdef uint32_t * c_eccentricity_bounding(G) except NULL:
+    r"""
+    Return the vector of eccentricities in G using the algorithm of [FK13]_.
+
+    The array returned is of length n, and its ith component is the eccentricity
+    of the ith vertex in ``G.vertices()``.
+
+    The algorithm proposed in [FK13]_ is based on the observation that for all
+    nodes `v,w\in V`, we have `\max(ecc[v]-d(v,w), d(v,w))\leq ecc[w] \leq
+    ecc[v] + d(v,w)`.Also the algorithms iteratively improves upper and lower
+    bounds on the eccentricity of each node until no further improvements can be
+    done. This algorithm offers good running time reduction on scale-free graphs.
+    """
+    # Copying the whole graph to obtain the list of neighbors quicker than by
+    # calling out_neighbors.  This data structure is well documented in the
+    # module sage.graphs.base.static_sparse_graph
+    cdef unsigned int n = G.order()
+    cdef short_digraph sd
+    init_short_digraph(sd, G)
+
+    # allocated some data structures
+    cdef bitset_t seen
+    bitset_init(seen, n)
+    cdef uint32_t * distances = <uint32_t *>sage_malloc(2 * n * sizeof(uint32_t))
+    if distances==NULL:
+        bitset_free(seen)
+        raise MemoryError()
+    cdef uint32_t * waiting_list = distances + n
+
+    cdef uint32_t * LB = <uint32_t *>sage_calloc(n, sizeof(uint32_t))
+    cdef uint32_t * UB = <uint32_t *>sage_malloc(n * sizeof(uint32_t))
+    if LB==NULL or UB==NULL:
+        bitset_free(seen)
+        sage_free(LB)
+        sage_free(UB)
+        sage_free(distances)
+        raise MemoryError()
+
+    memset(UB, -1, n * sizeof(uint32_t))
+
+    cdef list W = range(n)
+    cdef uint32_t v, w, next_v, tmp, cpt = 0
+
+    # The first vertex is the one with largest degree
+    next_v = max((out_degree(sd, v), v) for v in range(n))[1]
+
+    sig_on()
+    while W:
+
+        v = next_v
+        W.remove(v)
+        cpt += 1
+
+        # Compute the exact eccentricity of v
+        LB[v] = simple_BFS(n, sd.neighbors, v, distances, NULL, waiting_list, seen)
+
+        if LB[v]==UINT32_MAX:
+            # The (di) graph is not (strongly) connected. We set maximum value
+            # and exit.
+            for w in range(n):
+                LB[w] = UINT32_MAX
+            break
+
+        # Improve the bounds on the eccentricity of other vertices and select
+        # source of the next BFS
+        if not W:
+            break
+        next_v = W[0]
+        for w in W:
+            LB[w] = max(LB[w], max(LB[v] - distances[w], distances[w]))
+            UB[w] = min(UB[w], LB[v] + distances[w])
+            if LB[w]==UB[w]:
+                W.remove(w)
+            elif (cpt%2==0 and LB[w]<LB[next_v]) or (cpt%2==1 and UB[w]>UB[next_v]):
+                # The next vertex is either the vertex with largest upper bound
+                # or smallest lower bound
+                next_v = w
+
+    sig_off()
+
+    sage_free(distances)
+    sage_free(UB)
+    bitset_free(seen)
+    free_short_digraph(sd)
+
+    return LB
+
+def eccentricity(G, method="bounds"):
     r"""
-    Returns the vector of eccentricities in G.
+    Return the vector of eccentricities in G.
 
     The array returned is of length n, and its ith component is the eccentricity
     of the ith vertex in ``G.vertices()``.
@@ -739,22 +833,56 @@ def eccentricity(G):
         sage: g = graphs.PetersenGraph()
         sage: eccentricity(g)
         [2, 2, 2, 2, 2, 2, 2, 2, 2, 2]
+
+    TEST:
+
+    All methods are valid::
+
+        sage: from sage.graphs.distances_all_pairs import eccentricity
+        sage: g = graphs.RandomGNP(50,.1)
+        sage: eccentricity(g, method='standard')==eccentricity(g, method='bounds')
+        True
+
+    Case of not (strongly) connected (directed) graph::
+
+        sage: from sage.graphs.distances_all_pairs import eccentricity
+        sage: g = 2*graphs.PathGraph(2)
+        sage: eccentricity(g, method='bounds')
+        [+Infinity, +Infinity, +Infinity, +Infinity]
+        sage: g = digraphs.Path(3)
+        sage: eccentricity(g, method='bounds')
+        [+Infinity, +Infinity, +Infinity]
+
+    Asking for unknown method::
+
+        sage: from sage.graphs.distances_all_pairs import eccentricity
+        sage: g = graphs.PathGraph(2)
+        sage: eccentricity(g, method='Nice Jazz Festival')
+        Traceback (most recent call last):
+        ...
+        ValueError: Unknown method 'Nice Jazz Festival'. Please contribute.
     """
     from sage.rings.infinity import Infinity
     cdef int n = G.order()
 
+    # Trivial cases
     if n == 0:
         return []
+    elif G.is_directed():
+        if not G.is_strongly_connected():
+            return [Infinity]*n
+    elif not G.is_connected():
+        return [Infinity]*n
+
+    cdef uint32_t * ecc
+    if method=="bounds":
+        ecc = c_eccentricity_bounding(G)
+    elif method=="standard":
+        ecc = c_eccentricity(G)
+    else:
+        raise ValueError("Unknown method '{}'. Please contribute.".format(method))
 
-    cdef int * ecc = c_eccentricity(G)
-
-    cdef list l_ecc = []
-    cdef int i
-    for i in range(n):
-        if ecc[i] == INT32_MAX:
-            l_ecc.append(Infinity)
-        else:
-            l_ecc.append(ecc[i])
+    cdef list l_ecc = [ecc[i] for i in range(n)]
 
     sage_free(ecc)