Ticket #12052: trac_12052.patch

File trac_12052.patch, 52.8 KB (added by ncohen, 9 years ago)
  • doc/en/reference/graphs.rst

    # HG changeset patch
    # User Nathann Cohen <nathann.cohen@gmail.com>
    # Date 1321640996 -3600
    # Node ID ee1b7de8245fa974a6f8d0d39ce0f6cf58c41c31
    # Parent  9276337e5cd9643ceb185388014bd03858b4cb86
    trac #12052 -- Some distance computations remained *slow*
    
    diff --git a/doc/en/reference/graphs.rst b/doc/en/reference/graphs.rst
    a b  
    4949   sage/graphs/matchpoly
    5050   sage/graphs/graph_decompositions/vertex_separation
    5151   sage/graphs/convexity_properties
     52   sage/graphs/distances_all_pairs
    5253   sage/graphs/graph_latex
    5354   sage/graphs/graph_list
    5455
  • module_list.py

    diff --git a/module_list.py b/module_list.py
    a b  
    329329              sources = ['sage/graphs/generic_graph_pyx.pyx'],
    330330              libraries = ['gmp']),
    331331
     332    Extension('sage.graphs.distances_all_pairs',
     333              sources = ['sage/graphs/distances_all_pairs.pyx']),
     334
    332335    Extension('sage.graphs.modular_decomposition.modular_decomposition',
    333336              sources = ['sage/graphs/modular_decomposition/modular_decomposition.pyx'],
    334337              depends = ['sage/graphs/modular_decomposition/src/dm.c',
  • sage/graphs/all.py

    diff --git a/sage/graphs/all.py b/sage/graphs/all.py
    a b  
    1515from sage.graphs.cliquer import *
    1616from graph_database import graph_db_info
    1717from graph_editor import graph_editor
     18import sage.graphs.distances_all_pairs
  • sage/graphs/base/c_graph.pyx

    diff --git a/sage/graphs/base/c_graph.pyx b/sage/graphs/base/c_graph.pyx
    a b  
    25602560        r"""
    25612561        Returns whether the graph is connected.
    25622562
    2563         INPUT:
    2564 
    2565         - None.
    2566 
    2567         OUTPUT:
    2568 
    2569         - ``True`` if this graph is connected; ``False`` otherwise.
    2570 
    25712563        EXAMPLES:
    25722564
    25732565        Petersen's graph is connected::
     
    25962588        r"""
    25972589        Returns whether the graph is strongly connected.
    25982590
    2599         INPUT:
    2600 
    2601         - None.
    2602 
    2603         OUTPUT:
    2604 
    2605         - ``True`` if this graph is strongly connected; ``False`` otherwise.
    2606 
    26072591        EXAMPLES:
    26082592
    26092593        The circuit on 3 vertices is obviously strongly connected::
     
    26682652        Returns whether the graph is both directed and acylic (possibly with a
    26692653        certificate)
    26702654
    2671 
    26722655        INPUT:
    26732656
    26742657        - ``certificate`` -- whether to return a certificate (``False`` by
     
    28512834        else:
    28522835            return True
    28532836
    2854 def floyd_warshall(gg, paths = True, distances = False):
    2855     r"""
    2856     Computes the shortest path/distances between all pairs of vertices.
    2857 
    2858     For more information on the Floyd-Warshall algorithm, see the `Wikipedia
    2859     article on Floyd-Warshall
    2860     <http://en.wikipedia.org/wiki/Floyd%E2%80%93Warshall_algorithm>`_.
    2861 
    2862     INPUT:
    2863 
    2864         - ``gg`` -- the graph on which to work.
    2865 
    2866         - ``paths`` (boolean) -- whether to return the dictionary of shortest
    2867           paths. Set to ``True`` by default.
    2868 
    2869         - ``distances`` (boolean) -- whether to return the dictionary of
    2870           distances. Set to ``False`` by default.
    2871 
    2872     OUTPUT:
    2873 
    2874         Depending on the input, this function return the dictionary of paths,
    2875         the dictionary of distances, or a pair of dictionaries
    2876         ``(distances, paths)`` where ``distance[u][v]`` denotes the distance of a
    2877         shortest path from `u` to `v` and ``paths[u][v]`` denotes an inneighbor
    2878         `w` of `v` such that `dist(u,v)= 1 + dist(u,w)`.
    2879 
    2880     .. WARNING::
    2881 
    2882         Because this function works on matrices whose size is quadratic compared
    2883         to the number of vertices, it uses short variables instead of long ones
    2884         to divide by 2 the size in memory. This means that the current
    2885         implementation does not run on a graph of more than 65536 nodes (this
    2886         can be easily changed if necessary, but would require much more
    2887         memory. It may be worth writing two versions). For information, the
    2888         current version of the algorithm on a graph with `65536=2^{16}` nodes
    2889         creates in memory `2` tables on `2^{32}` short elements (2bytes each),
    2890         for a total of `2^{34}` bytes or `16` gigabytes. Let us also remember
    2891         that if the memory size is quadratic, the algorithm runs in cubic time.
    2892 
    2893     .. NOTE::
    2894 
    2895         When ``paths = False`` the algorithm saves roughly half of the memory as
    2896         it does not have to maintain the matrix of predecessors. However,
    2897         setting ``distances=False`` produces no such effect as the algorithm can
    2898         not run without computing them. They will not be returned, but they will
    2899         be stored while the method is running.
    2900 
    2901     EXAMPLES:
    2902 
    2903     Shortest paths in a small grid ::
    2904 
    2905         sage: g = graphs.Grid2dGraph(2,2)
    2906         sage: from sage.graphs.base.c_graph import floyd_warshall
    2907         sage: print floyd_warshall(g)
    2908         {(0, 1): {(0, 1): None, (1, 0): (0, 0), (0, 0): (0, 1), (1, 1): (0, 1)},
    2909         (1, 0): {(0, 1): (0, 0), (1, 0): None, (0, 0): (1, 0), (1, 1): (1, 0)},
    2910         (0, 0): {(0, 1): (0, 0), (1, 0): (0, 0), (0, 0): None, (1, 1): (0, 1)},
    2911         (1, 1): {(0, 1): (1, 1), (1, 0): (1, 1), (0, 0): (0, 1), (1, 1): None}}
    2912 
    2913     Checking the distances are correct ::
    2914 
    2915         sage: g = graphs.Grid2dGraph(5,5)
    2916         sage: dist,path = floyd_warshall(g, distances = True)
    2917         sage: all( dist[u][v] == g.distance(u,v) for u in g for v in g )
    2918         True
    2919 
    2920     Checking a random path is valid ::
    2921 
    2922         sage: u,v = g.random_vertex(), g.random_vertex()
    2923         sage: p = [v]
    2924         sage: while p[0] != None:
    2925         ...     p.insert(0,path[u][p[0]])
    2926         sage: len(p) == dist[u][v] + 2
    2927         True
    2928 
    2929     Distances for all pairs of vertices in a diamond::
    2930 
    2931         sage: g = graphs.DiamondGraph()
    2932         sage: floyd_warshall(g, paths = False, distances = True)
    2933         {0: {0: 0, 1: 1, 2: 1, 3: 2},
    2934          1: {0: 1, 1: 0, 2: 1, 3: 1},
    2935          2: {0: 1, 1: 1, 2: 0, 3: 1},
    2936          3: {0: 2, 1: 1, 2: 1, 3: 0}}
    2937 
    2938     TESTS:
    2939 
    2940     Too large graphs::
    2941 
    2942         sage: from sage.graphs.base.c_graph import floyd_warshall
    2943         sage: floyd_warshall(Graph(65536))
    2944         Traceback (most recent call last):
    2945         ...
    2946         ValueError: The graph backend contains more than 65535 nodes
    2947     """
    2948 
    2949     from sage.rings.infinity import Infinity
    2950     cdef CGraph g = <CGraph> gg._backend._cg
    2951 
    2952     cdef list gverts = g.verts()
    2953 
    2954     cdef int n = max(gverts) + 1
    2955 
    2956     if n >= <unsigned short> -1:
    2957         raise ValueError("The graph backend contains more than "+str(<unsigned short> -1)+" nodes")
    2958 
    2959     # All this just creates two tables prec[n][n] and dist[n][n]
    2960     cdef unsigned short * t_prec
    2961     cdef unsigned short * t_dist
    2962     cdef unsigned short ** prec
    2963     cdef unsigned short ** dist
    2964 
    2965     cdef int i
    2966     cdef int v_int
    2967     cdef int u_int
    2968     cdef int w_int
    2969 
    2970     # init dist
    2971     t_dist = <unsigned short *> sage_malloc(n*n*sizeof(short))
    2972     dist = <unsigned short **> sage_malloc(n*sizeof(short *))
    2973     dist[0] = t_dist
    2974     for 1 <= i< n:
    2975         dist[i] = dist[i-1] + n
    2976     memset(t_dist, -1, n*n*sizeof(short))
    2977     # Copying the adjacency matrix (vertices at distance 1)
    2978     for v_int in gverts:
    2979         dist[v_int][v_int] =  0
    2980         for u_int in g.out_neighbors(v_int):
    2981             dist[v_int][u_int] = 1
    2982 
    2983     if paths:
    2984         # init prec
    2985         t_prec = <unsigned short *> sage_malloc(n*n*sizeof(short))
    2986         prec = <unsigned short **> sage_malloc(n*sizeof(short *))
    2987         prec[0] = t_prec
    2988         for 1 <= i< n:
    2989             prec[i] = prec[i-1] + n
    2990         memset(t_prec, 0, n*n*sizeof(short))
    2991         # Copying the adjacency matrix (vertices at distance 1)
    2992         for v_int in gverts:
    2993             prec[v_int][v_int] = v_int
    2994             for u_int in g.out_neighbors(v_int):
    2995                 prec[v_int][u_int] = v_int
    2996 
    2997     # The algorithm itself.
    2998     cdef unsigned short *dv, *dw
    2999     cdef int dvw
    3000     cdef int val
    3001 
    3002     for w_int in gverts:
    3003         dw = dist[w_int]
    3004         for v_int in gverts:
    3005             dv = dist[v_int]
    3006             dvw = dv[w_int]
    3007             for u_int in gverts:
    3008                 val = dvw + dw[u_int]
    3009                 # If it is shorter to go from u to v through w, do it
    3010                 if dv[u_int] > val:
    3011                     dv[u_int] = val
    3012                     if paths:
    3013                         prec[v_int][u_int] = prec[w_int][u_int]
    3014 
    3015     # Dictionaries of distance, precedent element, and integers
    3016     cdef dict d_prec
    3017     cdef dict d_dist
    3018     cdef dict tmp_prec
    3019     cdef dict tmp_dist
    3020 
    3021     cdef dict ggbvi = gg._backend.vertex_ints
    3022     cdef dict ggbvl = gg._backend.vertex_labels
    3023 
    3024     if paths: d_prec = {}
    3025     if distances: d_dist = {}
    3026     for v_int in gverts:
    3027         if paths: tmp_prec = {}
    3028         if distances: tmp_dist = {}
    3029         v = vertex_label(v_int, ggbvi, ggbvl, g)
    3030         dv = dist[v_int]
    3031         for u_int in gverts:
    3032             u = vertex_label(u_int, ggbvi, ggbvl, g)
    3033             if paths:
    3034                 tmp_prec[u] = (None if v == u
    3035                                else vertex_label(prec[v_int][u_int], ggbvi, ggbvl, g))
    3036             if distances:
    3037                 tmp_dist[u] = (dv[u_int] if (dv[u_int] != <unsigned short> -1)
    3038                                else Infinity)
    3039         if paths: d_prec[v] = tmp_prec
    3040         if distances: d_dist[v] = tmp_dist
    3041 
    3042     if paths:
    3043         sage_free(t_prec)
    3044         sage_free(prec)
    3045 
    3046     sage_free(t_dist)
    3047     sage_free(dist)
    3048 
    3049     if distances and paths:
    3050         return d_dist, d_prec
    3051     if paths:
    3052         return d_prec
    3053     if distances:
    3054         return d_dist
    3055 
    3056 cpdef all_pairs_shortest_path_BFS(gg):
    3057     """
    3058     Computes the shortest path and the distances between each pair of vertices
    3059     through successive breadth-first-searches
    3060 
    3061     OUTPUT:
    3062 
    3063         This function return a pair of dictionaries ``(distances,paths)`` where
    3064         ``distance[u][v]`` denotes the distance of a shortest path from `u` to
    3065         `v` and ``paths[u][v]`` denotes an inneighbor `w` of `v` such that
    3066         `dist(u,v)= 1 + dist(u,w)`.
    3067 
    3068     .. WARNING::
    3069 
    3070         Because this function works on matrices whose size is quadratic compared
    3071         to the number of vertices, it uses short variables instead of long ones
    3072         to divide by 2 the size in memory. This means that the current
    3073         implementation does not run on a graph of more than 65536 nodes (this
    3074         can be easily changed if necessary, but would require much more
    3075         memory. It may be worth writing two versions). For information, the
    3076         current version of the algorithm on a graph with `65536=2^{16}` nodes
    3077         creates in memory `2` tables on `2^{32}` short elements (2bytes each),
    3078         for a total of `2^{34}` bytes or `16` gigabytes.
    3079 
    3080     EXAMPLE:
    3081 
    3082     On a grid::
    3083 
    3084         sage: g = graphs.Grid2dGraph(10,10)
    3085         sage: from sage.graphs.base.c_graph import all_pairs_shortest_path_BFS
    3086         sage: dist, path = all_pairs_shortest_path_BFS(g)
    3087         sage: all( dist[u][v] == g.distance(u,v) for u in g for v in g )
    3088         True
    3089 
    3090     TESTS:
    3091 
    3092     Too large graphs::
    3093 
    3094         sage: all_pairs_shortest_path_BFS(Graph(65536))
    3095         Traceback (most recent call last):
    3096         ...
    3097         ValueError: The graph backend contains more than 65535 nodes
    3098 
    3099     """
    3100     from sage.rings.infinity import Infinity
    3101 
    3102     cdef CGraph cg = <CGraph> gg._backend._cg
    3103 
    3104     cdef list vertices = cg.verts()
    3105     cdef int n = max(vertices)+1
    3106 
    3107     if n > <unsigned short> -1:
    3108         raise ValueError("The graph backend contains more than "+str(<unsigned short> -1)+" nodes")
    3109 
    3110     # The vertices which have already been visited
    3111     cdef bitset_t seen
    3112     bitset_init(seen, cg.active_vertices.size)
    3113 
    3114     # The list of waiting vertices, the beginning and the end of the list
    3115 
    3116     cdef unsigned short * waiting_list = <unsigned short *> sage_malloc(n*sizeof(short))
    3117     cdef unsigned short waiting_beginning = 0
    3118     cdef unsigned short waiting_end = 0
    3119 
    3120     cdef unsigned short source
    3121     cdef unsigned short v, u
    3122 
    3123     # Dictionaries of dictionaries of distances/predecessors
    3124     cdef dict d_distances = dict()
    3125     cdef dict d_prec      = dict()
    3126 
    3127     # Temporary dictionaries of distances/predecessors
    3128     cdef dict tmp_distances
    3129     cdef dict tmp_prec
    3130 
    3131     # Temporary vectors of distances/predecessors
    3132     cdef unsigned short * v_distances = <unsigned short *> sage_malloc(n*sizeof(short))
    3133     cdef unsigned short * v_prec      = <unsigned short *> sage_malloc(n*sizeof(short))
    3134 
    3135     cdef dict ggbvi = gg._backend.vertex_ints
    3136     cdef dict ggbvl = gg._backend.vertex_labels
    3137 
    3138     cdef int *outneighbors
    3139     cdef int o_n_size
    3140     cdef int i
    3141 
    3142     # Copying the whole graph to obtain the list of neighbors quicker than by
    3143     # calling out_neighbors
    3144 
    3145     # The edges are stored in the vector p_edges. This vector contains, from
    3146     # left to right The list of the first vertex's outneighbors, then the
    3147     # second's, then the third's, ...
    3148     #
    3149     # The outneighbors of vertex i are enumerated from
    3150     #
    3151     # p_vertices[i] to p_vertices[i+1] - 1
    3152     # (if p_vertices[i] is equal to p_vertices[i+1], then i has no outneighbours)
    3153 
    3154     cdef short ** p_vertices = <short **> sage_malloc(n*sizeof(short *))
    3155     cdef short * p_edges = <short *> sage_malloc(cg.num_arcs*sizeof(short))
    3156     cdef short * p_next = p_edges
    3157 
    3158     for v in vertices:
    3159         outneighbors = <int *>sage_malloc(cg.out_degrees[v] * sizeof(int))
    3160         cg.out_neighbors_unsafe(v, outneighbors, cg.out_degrees[v])
    3161         p_vertices[v] = p_next
    3162        
    3163         for 0<= i < cg.out_degrees[v]:
    3164             p_next[0] = <short> outneighbors[i]
    3165             p_next = p_next + 1
    3166 
    3167         sage_free(outneighbors)
    3168 
    3169     # We run n different BFS taking each vertex as a source
    3170     for source in vertices:
    3171 
    3172         # The source is seen
    3173         bitset_set_first_n(seen, 0)
    3174         bitset_add(seen, source)
    3175 
    3176         # Its parameters can already be set
    3177         v_distances[source] = 0
    3178         v_prec[source] = source
    3179 
    3180         # and added to the queue
    3181         waiting_list[0] = source
    3182         waiting_beginning = 0
    3183         waiting_end = 0
    3184 
    3185         # For as long as there are vertices left to explore
    3186         while waiting_beginning <= waiting_end:
    3187 
    3188             # We pick the first one
    3189             v = waiting_list[waiting_beginning]
    3190 
    3191             # Iterating over all the outneighbors u of v
    3192             for 0 <= i < cg.out_degrees[v]:
    3193                 u = p_vertices[v][i]
    3194 
    3195                 # If we notice one of these neighbors is not seen yet, we set
    3196                 # its parameters and add it to the queue to be explored later.
    3197                 if not bitset_in(seen, u):
    3198                     v_distances[u] = v_distances[v]+1
    3199                     v_prec[u] = v
    3200                     bitset_add(seen, u)
    3201                     waiting_end += 1
    3202                     waiting_list[waiting_end] = u
    3203 
    3204             waiting_beginning += 1
    3205 
    3206         # Copying the information to dictionaries (very costly, comparatively)
    3207         # to be returned later
    3208         tmp_distances = dict()
    3209         tmp_prec = dict()
    3210         for v in vertices:
    3211             vv = vertex_label(v, ggbvi, ggbvl, cg)
    3212 
    3213             if bitset_in(seen, v):
    3214                 tmp_prec[vv] = vertex_label(v_prec[v], ggbvi, ggbvl, cg)
    3215                 tmp_distances[vv] = v_distances[v]
    3216             else:
    3217                 tmp_prec[vv] = None
    3218                 tmp_distances[vv] = Infinity
    3219 
    3220         vv = vertex_label(source, ggbvi, ggbvl, cg)
    3221         tmp_prec[vv] = None
    3222         d_prec[vv] = tmp_prec
    3223         d_distances[vv] = tmp_distances
    3224 
    3225     bitset_free(seen)
    3226     sage_free(waiting_list)
    3227     sage_free(v_distances)
    3228     sage_free(v_prec)
    3229     sage_free(p_vertices)
    3230     sage_free(p_edges)
    3231 
    3232     return d_distances, d_prec
    3233        
    3234 
    32352837cdef class Search_iterator:
    32362838    r"""
    32372839    An iterator for traversing a (di)graph.
  • new file sage/graphs/distances_all_pairs.pyx

    diff --git a/sage/graphs/distances_all_pairs.pyx b/sage/graphs/distances_all_pairs.pyx
    new file mode 100644
    - +  
     1r"""
     2Distances/shortest paths between all pairs of vertices
     3
     4This module implements a few functions that deal with the computation of
     5distances or shortest paths between all pairs of vertices.
     6
     7Because these functions involve listing many times the (out)-neighborhoods of
     8(di)-graphs, it is useful in terms of efficiency to build a temporary copy of
     9the graph in a data structure that makes it easy to compute quickly. These
     10functions also work on large volume of data, typically dense matrices of size
     11`n^2`, and are expected to return corresponding dictionaries of size `n^2`,
     12where the integers corresponding to the vertices have first been converted to
     13the vertices' labels. Sadly, this last translating operation turns out to be the
     14most time-consuming, and for this reason it is also nice to have a Cython
     15module, and version of these functions that return C arrays, in order to avoid
     16these operations when they are not necessary.
     17
     18The module's main function
     19--------------------------
     20
     21The C function ``all_pairs_shortest_path_BFS`` actually does all the
     22computations, and all the others (except for ``Floyd_Warshall``) are just
     23wrapping it. This function begins with copying the graph in a data structure
     24that makes it fast to query the out-neighbors of a vertex, then starts one
     25Breadth First Search per vertex of the (di)graph.
     26
     27**What can this function compute ?**
     28
     29    - The matrix of predecessors.
     30
     31      This matrix `P` has size `n^2`, and is such that vertex `P[u,v]` is a
     32      predecessor of `v` on a shortest `uv`-path. Hence, this matrix efficiently
     33      encodes the information of a shortest `uv`-path for any `u,v\in G` :
     34      indeed, to go from `u` to `v` you should first find a shortest
     35      `uP[u,v]`-path, then jump from `P[u,v]` to `v` as it is one of its
     36      outneighbors. Apply recursively and find out what the whole path is !.
     37
     38    - The matrix of distances.
     39
     40      This matrix has size `n^2` and associates to any `uv` the distance from
     41      `u` to `v`.
     42
     43    - The vector of eccentricities.
     44
     45      This vector of size `n` encodes for each vertex `v` the distance to vertex
     46      which is furthest from `v` in the graph. In particular, the diameter of
     47      the graph is the maximum of these values.
     48
     49**What does it take as input ?**
     50
     51    - ``gg`` a (Di)Graph.
     52
     53    - ``unsigned short * predecessors`` -- a pointer toward an array of size
     54      `n^2\text{sizeof(unsigned short)}`. Set to ``NULL`` if you do not want to
     55      compute the predecessors.
     56
     57    - ``unsigned short * distances`` -- a pointer toward an array of size
     58      `n^2\text{sizeof(unsigned short)}`. The computation of the distances is
     59      necessary for the algorithm, so this value can **not** be set to ``NULL``.
     60
     61    - ``unsigned short * eccentricity`` -- a pointer toward an array of size
     62      `n\text{sizeof(unsigned short)}`. Set to ``NULL`` if you do not want to
     63      compute the eccentricity.
     64
     65**Technical details**
     66
     67
     68    - The vertices are encoded as `1, ..., n` as they appear in the ordering of
     69      ``G.vertices()``.
     70
     71    - Because this function works on matrices whose size is quadratic compared
     72      to the number of vertices, it uses short variables to store the vertices'
     73      names instead of long ones to divide by 2 the size in memory. This means
     74      that the current implementation does not run on a graph of more than 65536
     75      nodes (this can be easily changed if necessary, but would require much
     76      more memory. It may be worth writing two versions). For information, the
     77      current version of the algorithm on a graph with `65536=2^{16}` nodes
     78      creates in memory `2` tables on `2^{32}` short elements (2bytes each), for
     79      a total of `2^{34}` bytes or `16` gigabytes.
     80
     81    - In the C version of these functions, a value of ``<unsigned short> -1 =
     82      65535`` for a distance or a predecessor means respectively ``+Infinity``
     83      and ``None``. These case happens when the input is a disconnected graph,
     84      or a non-strongly-connected digraph.
     85
     86    .. WARNING::
     87
     88        The function ``all_pairs_shortest_path_BFS`` has **no reason** to be
     89        called by the user, even though he would be writing his code in Cython
     90        and look for efficiency. This module contains wrappers for this function
     91        that feed it with the good parameters. As the function is inlined, using
     92        those wrappers actually saves time as it should avoid testing the
     93        parameters again and again in the main function's body.
     94
     95AUTHOR:
     96
     97- Nathann Cohen (2011)
     98
     99
     100Functions
     101---------
     102"""
     103
     104##############################################################################
     105#       Copyright (C) 2011 Nathann Cohen <nathann.cohen@gmail.com>
     106#  Distributed under the terms of the GNU General Public License (GPL)
     107#  The full text of the GPL is available at:
     108#                  http://www.gnu.org/licenses/
     109##############################################################################
     110
     111include "../misc/bitset_pxd.pxi"
     112include "../misc/bitset.pxi"
     113from sage.graphs.base.c_graph cimport CGraph
     114from sage.graphs.base.c_graph cimport vertex_label
     115from sage.graphs.base.c_graph cimport get_vertex
     116
     117
     118cdef inline all_pairs_shortest_path_BFS(gg,
     119                                        unsigned short * predecessors,
     120                                        unsigned short * distances,
     121                                        unsigned short * eccentricity):
     122    """
     123    See the module's documentation.
     124    """
     125    from sage.rings.infinity import Infinity
     126
     127    cdef CGraph cg = <CGraph> gg._backend._cg
     128
     129    cdef list int_to_vertex = gg.vertices()
     130    cdef dict vertex_to_int = {}
     131    cdef int i
     132
     133    for i, l in enumerate(int_to_vertex):
     134        int_to_vertex[i] = get_vertex(l, gg._backend.vertex_ints, gg._backend.vertex_labels, cg)
     135        vertex_to_int[int_to_vertex[i]] = i
     136       
     137   
     138    cdef int n = len(int_to_vertex)
     139
     140    if n > <unsigned short> -1:
     141        raise ValueError("The graph backend contains more than "+str(<unsigned short> -1)+" nodes")
     142
     143    # The vertices which have already been visited
     144    cdef bitset_t seen
     145    bitset_init(seen, n)
     146
     147    # The list of waiting vertices, the beginning and the end of the list
     148
     149    cdef unsigned short * waiting_list = <unsigned short *> sage_malloc(n*sizeof(short))
     150    cdef unsigned short waiting_beginning = 0
     151    cdef unsigned short waiting_end = 0
     152
     153    cdef int * degree = <int *> sage_malloc(n*sizeof(int))
     154
     155    cdef unsigned short source
     156    cdef unsigned short v, u
     157
     158    cdef unsigned short * c_predecessors = predecessors
     159    cdef unsigned short * c_eccentricity = eccentricity
     160    cdef unsigned short * c_distances
     161
     162    # If the user does not want to compute the distances, the distances variable
     163    # is set to NULL. We *need* to compute the distances, though, if we want to
     164    # compute anything with this function. In this case we allocate a vector of
     165    # size n (and not a vector of size n^2, that is the size of the distance
     166    # variable when it is not null)
     167
     168    if distances != NULL:
     169        c_distances = distances
     170    else:
     171        c_distances = <unsigned short *> sage_malloc( n * sizeof(unsigned short *))
     172
     173    cdef int * outneighbors
     174    cdef int o_n_size
     175
     176    # Copying the whole graph to obtain the list of neighbors quicker than by
     177    # calling out_neighbors
     178
     179    # The edges are stored in the vector p_edges. This vector contains, from
     180    # left to right The list of the first vertex's outneighbors, then the
     181    # second's, then the third's, ...
     182    #
     183    # The outneighbors of vertex i are enumerated from
     184    #
     185    # p_vertices[i] to p_vertices[i+1] - 1
     186    # (if p_vertices[i] is equal to p_vertices[i+1], then i has no outneighbours)
     187
     188    cdef unsigned short ** p_vertices = <unsigned short **> sage_malloc(n*sizeof(short *))
     189    cdef unsigned short * p_edges = <unsigned short *> sage_malloc(cg.num_arcs*sizeof(short))
     190    cdef unsigned short * p_next = p_edges
     191
     192    for v,i in enumerate(int_to_vertex):
     193        outneighbors = <int *>sage_malloc(cg.out_degrees[i] * sizeof(int))
     194        cg.out_neighbors_unsafe(i, outneighbors, cg.out_degrees[i])
     195        p_vertices[v] = p_next
     196       
     197        for 0<= j < cg.out_degrees[i]:
     198            p_next[0] = <unsigned short> vertex_to_int[outneighbors[j]]
     199            p_next = p_next + 1
     200
     201        degree[v] = p_next - p_vertices[v]
     202
     203        sage_free(outneighbors)
     204
     205
     206    # We run n different BFS taking each vertex as a source
     207    for source in range(n):
     208
     209        # The source is seen
     210        bitset_set_first_n(seen, 0)
     211        bitset_add(seen, source)
     212
     213        # Its parameters can already be set
     214        c_distances[source] = 0
     215
     216        if predecessors != NULL:
     217            c_predecessors[source] = source
     218
     219        # and added to the queue
     220        waiting_list[0] = source
     221        waiting_beginning = 0
     222        waiting_end = 0
     223
     224        # For as long as there are vertices left to explore
     225        while waiting_beginning <= waiting_end:
     226
     227            # We pick the first one
     228            v = waiting_list[waiting_beginning]
     229
     230            # Iterating over all the outneighbors u of v
     231            for 0 <= i < degree[v]:
     232                u = p_vertices[v][i]
     233
     234                # If we notice one of these neighbors is not seen yet, we set
     235                # its parameters and add it to the queue to be explored later.
     236                if not bitset_in(seen, u):
     237                    c_distances[u] = c_distances[v]+1
     238                    if predecessors != NULL:
     239                        c_predecessors[u] = v
     240                    bitset_add(seen, u)
     241                    waiting_end += 1
     242                    waiting_list[waiting_end] = u
     243
     244            waiting_beginning += 1
     245
     246        # If not all the vertices have been met
     247        for v in range(n):
     248            if not bitset_in(seen, v):
     249                c_distances[v] = -1
     250                if predecessors != NULL:
     251                    c_predecessors[v] = -1
     252
     253        if predecessors != NULL:
     254            c_predecessors += n
     255
     256        if eccentricity != NULL:
     257            c_eccentricity[0] = 0
     258            for i in range(n):
     259                c_eccentricity[0] = c_eccentricity[0] if c_eccentricity[0] >= c_distances[i] else c_distances[i]
     260
     261            c_eccentricity += 1
     262
     263        if distances != NULL:
     264            c_distances += n
     265
     266    bitset_free(seen)
     267    sage_free(waiting_list)
     268    sage_free(p_vertices)
     269    sage_free(p_edges)
     270    sage_free(degree)
     271    if distances == NULL:
     272        sage_free(c_distances)
     273
     274
     275################
     276# Predecessors #
     277################
     278
     279cdef unsigned short * c_shortest_path_all_pairs(G):
     280    r"""
     281    Returns the matrix of predecessors in G.
     282
     283    The matrix `P` returned has size `n^2`, and is such that vertex `P[u,v]` is
     284    a predecessor of `v` on a shortest `uv`-path. Hence, this matrix efficiently
     285    encodes the information of a shortest `uv`-path for any `u,v\in G` : indeed,
     286    to go from `u` to `v` you should first find a shortest `uP[u,v]`-path, then
     287    jump from `P[u,v]` to `v` as it is one of its outneighbors.
     288    """
     289
     290    cdef int n = G.order()
     291    cdef unsigned short * distances = <unsigned short *> sage_malloc(n*n*sizeof(unsigned short *))
     292    cdef unsigned short * predecessors = <unsigned short *> sage_malloc(n*n*sizeof(unsigned short *))
     293    all_pairs_shortest_path_BFS(G, predecessors, distances, NULL)
     294
     295    sage_free(distances)
     296
     297    return predecessors
     298
     299def shortest_path_all_pairs(G):
     300    r"""
     301    Returns the matrix of predecessors in G.
     302
     303    The matrix `P` returned has size `n^2`, and is such that vertex `P[u,v]` is
     304    a predecessor of `v` on a shortest `uv`-path. Hence, this matrix efficiently
     305    encodes the information of a shortest `uv`-path for any `u,v\in G` : indeed,
     306    to go from `u` to `v` you should first find a shortest `uP[u,v]`-path, then
     307    jump from `P[u,v]` to `v` as it is one of its outneighbors.
     308
     309    The integer corresponding to a vertex is its index in the list
     310    ``G.vertices()``.
     311
     312    EXAMPLE::
     313
     314        sage: from sage.graphs.distances_all_pairs import shortest_path_all_pairs
     315        sage: g = graphs.PetersenGraph()
     316        sage: shortest_path_all_pairs(g)
     317        {0: {0: None, 1: 0, 2: 1, 3: 4, 4: 0, 5: 0, 6: 1, 7: 5, 8: 5, 9: 4},
     318         1: {0: 1, 1: None, 2: 1, 3: 2, 4: 0, 5: 0, 6: 1, 7: 2, 8: 6, 9: 6},
     319         2: {0: 1, 1: 2, 2: None, 3: 2, 4: 3, 5: 7, 6: 1, 7: 2, 8: 3, 9: 7},
     320         3: {0: 4, 1: 2, 2: 3, 3: None, 4: 3, 5: 8, 6: 8, 7: 2, 8: 3, 9: 4},
     321         4: {0: 4, 1: 0, 2: 3, 3: 4, 4: None, 5: 0, 6: 9, 7: 9, 8: 3, 9: 4},
     322         5: {0: 5, 1: 0, 2: 7, 3: 8, 4: 0, 5: None, 6: 8, 7: 5, 8: 5, 9: 7},
     323         6: {0: 1, 1: 6, 2: 1, 3: 8, 4: 9, 5: 8, 6: None, 7: 9, 8: 6, 9: 6},
     324         7: {0: 5, 1: 2, 2: 7, 3: 2, 4: 9, 5: 7, 6: 9, 7: None, 8: 5, 9: 7},
     325         8: {0: 5, 1: 6, 2: 3, 3: 8, 4: 3, 5: 8, 6: 8, 7: 5, 8: None, 9: 6},
     326         9: {0: 4, 1: 6, 2: 7, 3: 4, 4: 9, 5: 7, 6: 9, 7: 9, 8: 6, 9: None}}
     327    """
     328
     329    cdef int n = G.order()
     330
     331    cdef unsigned short * predecessors = c_shortest_path_all_pairs(G)
     332    cdef unsigned short * c_predecessors = predecessors
     333
     334    cdef dict d = {}
     335    cdef dict d_tmp
     336
     337    cdef CGraph cg = <CGraph> G._backend._cg
     338
     339    cdef list int_to_vertex = G.vertices()
     340    cdef int i, j
     341
     342    for i, l in enumerate(int_to_vertex):
     343        int_to_vertex[i] = get_vertex(l, G._backend.vertex_ints, G._backend.vertex_labels, cg)
     344
     345    for j in range(n):
     346        d_tmp = {}
     347        for i in range(n):
     348            if c_predecessors[i] == <unsigned short> -1:
     349                d_tmp[int_to_vertex[i]] = None
     350            else:
     351                d_tmp[int_to_vertex[i]] = int_to_vertex[c_predecessors[i]]
     352
     353        d_tmp[int_to_vertex[j]] = None
     354        d[int_to_vertex[j]] = d_tmp
     355       
     356        c_predecessors += n
     357       
     358    sage_free(predecessors)
     359    return d
     360
     361#############
     362# Distances #
     363#############
     364
     365cdef unsigned short * c_distances_all_pairs(G):
     366    r"""
     367    Returns the matrix of distances in G.
     368
     369    The matrix `M` returned is of length `n^2`, and the distance between
     370    vertices `u` and `v` is `M[u,v]`. The integer corresponding to a vertex is
     371    its index in the list ``G.vertices()``.
     372    """
     373
     374    cdef int n = G.order()
     375    cdef unsigned short * distances = <unsigned short *> sage_malloc(n*n*sizeof(unsigned short *))
     376    all_pairs_shortest_path_BFS(G, NULL, distances, NULL)
     377
     378    return distances
     379
     380def distances_all_pairs(G):
     381    r"""
     382    Returns the matrix of distances in G.
     383
     384    This function returns a double dictionary ``D`` of vertices, in which the
     385    distance between vertices ``u`` and ``v`` is ``D[u][v]``.
     386
     387    EXAMPLE::
     388
     389        sage: from sage.graphs.distances_all_pairs import distances_all_pairs
     390        sage: g = graphs.PetersenGraph()
     391        sage: distances_all_pairs(g)
     392        {0: {0: 0, 1: 1, 2: 2, 3: 2, 4: 1, 5: 1, 6: 2, 7: 2, 8: 2, 9: 2},
     393        1: {0: 1, 1: 0, 2: 1, 3: 2, 4: 2, 5: 2, 6: 1, 7: 2, 8: 2, 9: 2},
     394        2: {0: 2, 1: 1, 2: 0, 3: 1, 4: 2, 5: 2, 6: 2, 7: 1, 8: 2, 9: 2},
     395        3: {0: 2, 1: 2, 2: 1, 3: 0, 4: 1, 5: 2, 6: 2, 7: 2, 8: 1, 9: 2},
     396        4: {0: 1, 1: 2, 2: 2, 3: 1, 4: 0, 5: 2, 6: 2, 7: 2, 8: 2, 9: 1},
     397        5: {0: 1, 1: 2, 2: 2, 3: 2, 4: 2, 5: 0, 6: 2, 7: 1, 8: 1, 9: 2},
     398        6: {0: 2, 1: 1, 2: 2, 3: 2, 4: 2, 5: 2, 6: 0, 7: 2, 8: 1, 9: 1},
     399        7: {0: 2, 1: 2, 2: 1, 3: 2, 4: 2, 5: 1, 6: 2, 7: 0, 8: 2, 9: 1},
     400        8: {0: 2, 1: 2, 2: 2, 3: 1, 4: 2, 5: 1, 6: 1, 7: 2, 8: 0, 9: 2},
     401        9: {0: 2, 1: 2, 2: 2, 3: 2, 4: 1, 5: 2, 6: 1, 7: 1, 8: 2, 9: 0}}
     402    """
     403
     404    from sage.rings.infinity import Infinity
     405
     406    cdef int n = G.order()
     407    cdef unsigned short * distances = c_distances_all_pairs(G)
     408    cdef unsigned short * c_distances = distances
     409
     410    cdef dict d = {}
     411    cdef dict d_tmp
     412
     413    cdef list int_to_vertex = G.vertices()
     414    cdef int i, j
     415
     416    for j in range(n):
     417        d_tmp = {}
     418        for i in range(n):
     419            if c_distances[i] == <unsigned short> -1:
     420                d_tmp[int_to_vertex[i]] = Infinity
     421            else:
     422                d_tmp[int_to_vertex[i]] = c_distances[i]
     423
     424
     425        d[int_to_vertex[j]] = d_tmp
     426       
     427        c_distances += n
     428       
     429    sage_free(distances)
     430    return d
     431
     432###################################
     433# Both distances and predecessors #
     434###################################
     435
     436def distances_and_predecessors_all_pairs(G):
     437    r"""
     438    Returns the matrix of distances in G and the matrix of predecessors.
     439
     440    Distances : the matrix `M` returned is of length `n^2`, and the distance
     441    between vertices `u` and `v` is `M[u,v]`. The integer corresponding to a
     442    vertex is its index in the list ``G.vertices()``.
     443
     444    Predecessors : the matrix `P` returned has size `n^2`, and is such that
     445    vertex `P[u,v]` is a predecessor of `v` on a shortest `uv`-path. Hence, this
     446    matrix efficiently encodes the information of a shortest `uv`-path for any
     447    `u,v\in G` : indeed, to go from `u` to `v` you should first find a shortest
     448    `uP[u,v]`-path, then jump from `P[u,v]` to `v` as it is one of its
     449    outneighbors.
     450
     451    The integer corresponding to a vertex is its index in the list
     452    ``G.vertices()``.
     453
     454
     455    EXAMPLE::
     456
     457        sage: from sage.graphs.distances_all_pairs import distances_and_predecessors_all_pairs
     458        sage: g = graphs.PetersenGraph()
     459        sage: distances_and_predecessors_all_pairs(g)
     460        ({0: {0: 0, 1: 1, 2: 2, 3: 2, 4: 1, 5: 1, 6: 2, 7: 2, 8: 2, 9: 2},
     461          1: {0: 1, 1: 0, 2: 1, 3: 2, 4: 2, 5: 2, 6: 1, 7: 2, 8: 2, 9: 2},
     462          2: {0: 2, 1: 1, 2: 0, 3: 1, 4: 2, 5: 2, 6: 2, 7: 1, 8: 2, 9: 2},
     463          3: {0: 2, 1: 2, 2: 1, 3: 0, 4: 1, 5: 2, 6: 2, 7: 2, 8: 1, 9: 2},
     464          4: {0: 1, 1: 2, 2: 2, 3: 1, 4: 0, 5: 2, 6: 2, 7: 2, 8: 2, 9: 1},
     465          5: {0: 1, 1: 2, 2: 2, 3: 2, 4: 2, 5: 0, 6: 2, 7: 1, 8: 1, 9: 2},
     466          6: {0: 2, 1: 1, 2: 2, 3: 2, 4: 2, 5: 2, 6: 0, 7: 2, 8: 1, 9: 1},
     467          7: {0: 2, 1: 2, 2: 1, 3: 2, 4: 2, 5: 1, 6: 2, 7: 0, 8: 2, 9: 1},
     468          8: {0: 2, 1: 2, 2: 2, 3: 1, 4: 2, 5: 1, 6: 1, 7: 2, 8: 0, 9: 2},
     469          9: {0: 2, 1: 2, 2: 2, 3: 2, 4: 1, 5: 2, 6: 1, 7: 1, 8: 2, 9: 0}},
     470         {0: {0: None, 1: 0, 2: 1, 3: 4, 4: 0, 5: 0, 6: 1, 7: 5, 8: 5, 9: 4},
     471          1: {0: 1, 1: None, 2: 1, 3: 2, 4: 0, 5: 0, 6: 1, 7: 2, 8: 6, 9: 6},
     472          2: {0: 1, 1: 2, 2: None, 3: 2, 4: 3, 5: 7, 6: 1, 7: 2, 8: 3, 9: 7},
     473          3: {0: 4, 1: 2, 2: 3, 3: None, 4: 3, 5: 8, 6: 8, 7: 2, 8: 3, 9: 4},
     474          4: {0: 4, 1: 0, 2: 3, 3: 4, 4: None, 5: 0, 6: 9, 7: 9, 8: 3, 9: 4},
     475          5: {0: 5, 1: 0, 2: 7, 3: 8, 4: 0, 5: None, 6: 8, 7: 5, 8: 5, 9: 7},
     476          6: {0: 1, 1: 6, 2: 1, 3: 8, 4: 9, 5: 8, 6: None, 7: 9, 8: 6, 9: 6},
     477          7: {0: 5, 1: 2, 2: 7, 3: 2, 4: 9, 5: 7, 6: 9, 7: None, 8: 5, 9: 7},
     478          8: {0: 5, 1: 6, 2: 3, 3: 8, 4: 3, 5: 8, 6: 8, 7: 5, 8: None, 9: 6},
     479          9: {0: 4, 1: 6, 2: 7, 3: 4, 4: 9, 5: 7, 6: 9, 7: 9, 8: 6, 9: None}})
     480    """
     481
     482    from sage.rings.infinity import Infinity
     483    cdef int n = G.order()
     484    cdef unsigned short * distances = <unsigned short *> sage_malloc(n*n*sizeof(unsigned short *))
     485    cdef unsigned short * c_distances = distances
     486    cdef unsigned short * predecessor = <unsigned short *> sage_malloc(n*n*sizeof(unsigned short *))
     487    cdef unsigned short * c_predecessor = predecessor
     488
     489    all_pairs_shortest_path_BFS(G, predecessor, distances, NULL)
     490
     491
     492    cdef dict d_distance = {}
     493    cdef dict d_predecessor = {}
     494    cdef dict t_distance = {}
     495    cdef dict t_predecessor = {}
     496
     497    cdef list int_to_vertex = G.vertices()
     498    cdef int i, j
     499
     500    for j in range(n):
     501        t_distance = {}
     502        t_predecessor = {}
     503
     504        for i in range(n):
     505
     506            if c_distances[i] == <unsigned short> -1:
     507                t_distance[int_to_vertex[i]] = Infinity
     508                t_predecessor[int_to_vertex[i]] = None
     509            else:
     510                t_distance[int_to_vertex[i]] = c_distances[i]
     511                t_predecessor[int_to_vertex[i]] = int_to_vertex[c_predecessor[i]]
     512
     513        t_predecessor[int_to_vertex[j]] = None
     514
     515        d_distance[int_to_vertex[j]] = t_distance
     516        d_predecessor[int_to_vertex[j]] = t_predecessor
     517       
     518        c_distances += n
     519        c_predecessor += n
     520
     521    sage_free(distances)
     522    sage_free(predecessor)
     523
     524    return d_distance, d_predecessor
     525
     526################
     527# Eccentricity #
     528################
     529
     530cdef unsigned short * c_eccentricity(G):
     531    r"""
     532    Returns the vector of eccentricities in G.
     533
     534    The array returned is of length n, and its ith component is the eccentricity
     535    of the ith vertex in ``G.vertices()``.
     536    """
     537    cdef int n = G.order()
     538
     539    cdef unsigned short * ecc = <unsigned short *> sage_malloc(n*sizeof(unsigned short *))
     540    cdef unsigned short * distances = <unsigned short *> sage_malloc(n*n*sizeof(unsigned short *))
     541    all_pairs_shortest_path_BFS(G, NULL, distances, ecc)
     542    sage_free(distances)
     543
     544    return ecc
     545
     546def eccentricity(G):
     547    r"""
     548    Returns the vector of eccentricities in G.
     549
     550    The array returned is of length n, and its ith component is the eccentricity
     551    of the ith vertex in ``G.vertices()``.
     552
     553    EXAMPLE::
     554
     555        sage: from sage.graphs.distances_all_pairs import eccentricity
     556        sage: g = graphs.PetersenGraph()
     557        sage: eccentricity(g)
     558        [2, 2, 2, 2, 2, 2, 2, 2, 2, 2]
     559    """
     560    from sage.rings.infinity import Infinity
     561    cdef int n = G.order()
     562    cdef unsigned short * ecc = c_eccentricity(G)
     563
     564    cdef list l_ecc = []
     565    cdef int i
     566    for i in range(n):
     567        if ecc[i] == <unsigned short> -1:
     568            l_ecc.append(Infinity)
     569        else:
     570            l_ecc.append(ecc[i])
     571
     572    sage_free(ecc)
     573
     574    return l_ecc
     575
     576def diameter(G):
     577    r"""
     578    Returns the diameter of `G`.
     579
     580    EXAMPLE::
     581
     582        sage: from sage.graphs.distances_all_pairs import diameter
     583        sage: g = graphs.PetersenGraph()
     584        sage: diameter(g)
     585        2
     586    """
     587    return max(eccentricity(G))
     588
     589##################
     590# Floyd-Warshall #
     591##################
     592
     593
     594def floyd_warshall(gg, paths = True, distances = False):
     595    r"""
     596    Computes the shortest path/distances between all pairs of vertices.
     597
     598    For more information on the Floyd-Warshall algorithm, see the `Wikipedia
     599    article on Floyd-Warshall
     600    <http://en.wikipedia.org/wiki/Floyd%E2%80%93Warshall_algorithm>`_.
     601
     602    INPUT:
     603
     604        - ``gg`` -- the graph on which to work.
     605
     606        - ``paths`` (boolean) -- whether to return the dictionary of shortest
     607          paths. Set to ``True`` by default.
     608
     609        - ``distances`` (boolean) -- whether to return the dictionary of
     610          distances. Set to ``False`` by default.
     611
     612    OUTPUT:
     613
     614        Depending on the input, this function return the dictionary of paths,
     615        the dictionary of distances, or a pair of dictionaries
     616        ``(distances, paths)`` where ``distance[u][v]`` denotes the distance of a
     617        shortest path from `u` to `v` and ``paths[u][v]`` denotes an inneighbor
     618        `w` of `v` such that `dist(u,v)= 1 + dist(u,w)`.
     619
     620    .. WARNING::
     621
     622        Because this function works on matrices whose size is quadratic compared
     623        to the number of vertices, it uses short variables instead of long ones
     624        to divide by 2 the size in memory. This means that the current
     625        implementation does not run on a graph of more than 65536 nodes (this
     626        can be easily changed if necessary, but would require much more
     627        memory. It may be worth writing two versions). For information, the
     628        current version of the algorithm on a graph with `65536=2^{16}` nodes
     629        creates in memory `2` tables on `2^{32}` short elements (2bytes each),
     630        for a total of `2^{34}` bytes or `16` gigabytes. Let us also remember
     631        that if the memory size is quadratic, the algorithm runs in cubic time.
     632
     633    .. NOTE::
     634
     635        When ``paths = False`` the algorithm saves roughly half of the memory as
     636        it does not have to maintain the matrix of predecessors. However,
     637        setting ``distances=False`` produces no such effect as the algorithm can
     638        not run without computing them. They will not be returned, but they will
     639        be stored while the method is running.
     640
     641    EXAMPLES:
     642
     643    Shortest paths in a small grid ::
     644
     645        sage: g = graphs.Grid2dGraph(2,2)
     646        sage: from sage.graphs.distances_all_pairs import floyd_warshall
     647        sage: print floyd_warshall(g)
     648        {(0, 1): {(0, 1): None, (1, 0): (0, 0), (0, 0): (0, 1), (1, 1): (0, 1)},
     649        (1, 0): {(0, 1): (0, 0), (1, 0): None, (0, 0): (1, 0), (1, 1): (1, 0)},
     650        (0, 0): {(0, 1): (0, 0), (1, 0): (0, 0), (0, 0): None, (1, 1): (0, 1)},
     651        (1, 1): {(0, 1): (1, 1), (1, 0): (1, 1), (0, 0): (0, 1), (1, 1): None}}
     652
     653    Checking the distances are correct ::
     654
     655        sage: g = graphs.Grid2dGraph(5,5)
     656        sage: dist,path = floyd_warshall(g, distances = True)
     657        sage: all( dist[u][v] == g.distance(u,v) for u in g for v in g )
     658        True
     659
     660    Checking a random path is valid ::
     661
     662        sage: u,v = g.random_vertex(), g.random_vertex()
     663        sage: p = [v]
     664        sage: while p[0] != None:
     665        ...     p.insert(0,path[u][p[0]])
     666        sage: len(p) == dist[u][v] + 2
     667        True
     668
     669    Distances for all pairs of vertices in a diamond::
     670
     671        sage: g = graphs.DiamondGraph()
     672        sage: floyd_warshall(g, paths = False, distances = True)
     673        {0: {0: 0, 1: 1, 2: 1, 3: 2},
     674         1: {0: 1, 1: 0, 2: 1, 3: 1},
     675         2: {0: 1, 1: 1, 2: 0, 3: 1},
     676         3: {0: 2, 1: 1, 2: 1, 3: 0}}
     677
     678    TESTS:
     679
     680    Too large graphs::
     681
     682        sage: from sage.graphs.distances_all_pairs import floyd_warshall
     683        sage: floyd_warshall(Graph(65536))
     684        Traceback (most recent call last):
     685        ...
     686        ValueError: The graph backend contains more than 65535 nodes
     687    """
     688
     689    from sage.rings.infinity import Infinity
     690    cdef CGraph g = <CGraph> gg._backend._cg
     691
     692    cdef list gverts = g.verts()
     693
     694    cdef int n = max(gverts) + 1
     695
     696    if n >= <unsigned short> -1:
     697        raise ValueError("The graph backend contains more than "+str(<unsigned short> -1)+" nodes")
     698
     699    # All this just creates two tables prec[n][n] and dist[n][n]
     700    cdef unsigned short * t_prec
     701    cdef unsigned short * t_dist
     702    cdef unsigned short ** prec
     703    cdef unsigned short ** dist
     704
     705    cdef int i
     706    cdef int v_int
     707    cdef int u_int
     708    cdef int w_int
     709
     710    # init dist
     711    t_dist = <unsigned short *> sage_malloc(n*n*sizeof(short))
     712    dist = <unsigned short **> sage_malloc(n*sizeof(short *))
     713    dist[0] = t_dist
     714    for 1 <= i< n:
     715        dist[i] = dist[i-1] + n
     716    memset(t_dist, -1, n*n*sizeof(short))
     717    # Copying the adjacency matrix (vertices at distance 1)
     718    for v_int in gverts:
     719        dist[v_int][v_int] =  0
     720        for u_int in g.out_neighbors(v_int):
     721            dist[v_int][u_int] = 1
     722
     723    if paths:
     724        # init prec
     725        t_prec = <unsigned short *> sage_malloc(n*n*sizeof(short))
     726        prec = <unsigned short **> sage_malloc(n*sizeof(short *))
     727        prec[0] = t_prec
     728        for 1 <= i< n:
     729            prec[i] = prec[i-1] + n
     730        memset(t_prec, 0, n*n*sizeof(short))
     731        # Copying the adjacency matrix (vertices at distance 1)
     732        for v_int in gverts:
     733            prec[v_int][v_int] = v_int
     734            for u_int in g.out_neighbors(v_int):
     735                prec[v_int][u_int] = v_int
     736
     737    # The algorithm itself.
     738    cdef unsigned short *dv, *dw
     739    cdef int dvw
     740    cdef int val
     741
     742    for w_int in gverts:
     743        dw = dist[w_int]
     744        for v_int in gverts:
     745            dv = dist[v_int]
     746            dvw = dv[w_int]
     747            for u_int in gverts:
     748                val = dvw + dw[u_int]
     749                # If it is shorter to go from u to v through w, do it
     750                if dv[u_int] > val:
     751                    dv[u_int] = val
     752                    if paths:
     753                        prec[v_int][u_int] = prec[w_int][u_int]
     754
     755    # Dictionaries of distance, precedent element, and integers
     756    cdef dict d_prec
     757    cdef dict d_dist
     758    cdef dict tmp_prec
     759    cdef dict tmp_dist
     760
     761    cdef dict ggbvi = gg._backend.vertex_ints
     762    cdef dict ggbvl = gg._backend.vertex_labels
     763
     764    if paths: d_prec = {}
     765    if distances: d_dist = {}
     766    for v_int in gverts:
     767        if paths: tmp_prec = {}
     768        if distances: tmp_dist = {}
     769        v = vertex_label(v_int, ggbvi, ggbvl, g)
     770        dv = dist[v_int]
     771        for u_int in gverts:
     772            u = vertex_label(u_int, ggbvi, ggbvl, g)
     773            if paths:
     774                tmp_prec[u] = (None if v == u
     775                               else vertex_label(prec[v_int][u_int], ggbvi, ggbvl, g))
     776            if distances:
     777                tmp_dist[u] = (dv[u_int] if (dv[u_int] != <unsigned short> -1)
     778                               else Infinity)
     779        if paths: d_prec[v] = tmp_prec
     780        if distances: d_dist[v] = tmp_dist
     781
     782    if paths:
     783        sage_free(t_prec)
     784        sage_free(prec)
     785
     786    sage_free(t_dist)
     787    sage_free(dist)
     788
     789    if distances and paths:
     790        return d_dist, d_prec
     791    if paths:
     792        return d_prec
     793    if distances:
     794        return d_dist
  • sage/graphs/generic_graph.py

    diff --git a/sage/graphs/generic_graph.py b/sage/graphs/generic_graph.py
    a b  
    1037610376
    1037710377        A doubly indexed dictionary
    1037810378
     10379        .. NOTE::
     10380
     10381           There is a Cython version of this method that is usually
     10382           much faster for large graphs, as most of the time is
     10383           actually spent building the final double
     10384           dictionary. Everything on the subject is to be found in the
     10385           :mod:`~sage.graphs.distances_all_pairs` module.
     10386
    1037910387        EXAMPLE:
    1038010388
    1038110389        The Petersen Graph::
     
    1039810406                algorithm = "BFS"
    1039910407
    1040010408        if algorithm == "BFS":
    10401             from sage.graphs.base.c_graph import all_pairs_shortest_path_BFS
    10402             return all_pairs_shortest_path_BFS(self)[0]
     10409            from sage.graphs.distances_all_pairs import distances_all_pairs
     10410            return distances_all_pairs(self)
    1040310411
    1040410412        elif algorithm == "Floyd-Warshall":
    10405             from sage.graphs.base.c_graph import floyd_warshall
     10413            from sage.graphs.distances_all_pairs import floyd_warshall
    1040610414            return floyd_warshall(self,paths = False, distances = True)
    1040710415
    1040810416        else:
     
    1045310461            {0: +Infinity, 1: +Infinity}
    1045410462        """
    1045510463        if v is None:
     10464            if dist_dict is None:
     10465                from sage.graphs.distances_all_pairs import eccentricity
     10466
     10467                if with_labels:
     10468                    return dict(zip(self.vertices(), eccentricity(self)))
     10469                else:
     10470                    return eccentricity(self)
     10471
    1045610472            v = self.vertices()
    1045710473        elif not isinstance(v, list):
    1045810474            v = [v]
     
    1054310559        """
    1054410560        Returns the largest distance between any two vertices. Returns
    1054510561        Infinity if the (di)graph is not connected.
    10546        
     10562
    1054710563        EXAMPLES::
    1054810564       
    1054910565            sage: G = graphs.PetersenGraph()
     
    1055910575            sage: G = graphs.EmptyGraph()
    1056010576            sage: G.diameter()
    1056110577            0
    10562         """
    10563         e = self.eccentricity()
    10564         if not isinstance(e, list):
    10565             e = [e]
    10566         if len(e) == 0:
     10578
     10579        """
     10580
     10581        if self.order() > 0:
     10582            return max(self.eccentricity())
     10583        else:
    1056710584            return 0
    10568         return max(e)
    1056910585
    1057010586    def distance_graph(self, dist):
    1057110587        r"""
     
    1110611122
    1110711123    def shortest_paths(self, u, by_weight=False, cutoff=None):
    1110811124        """
    11109         Returns a dictionary d of shortest paths d[v] from u to v, for each
    11110         vertex v connected by a path from u.
    11111        
    11112         INPUT:
    11113        
     11125        Returns a dictionary associating to each vertex v a shortest path from u
     11126        to v, if it exists.
     11127       
     11128        INPUT:
    1111411129       
    1111511130        -  ``by_weight`` - if False, uses a breadth first
    1111611131           search. If True, uses Dijkstra's algorithm to find the shortest
    1111711132           paths by weight.
    1111811133       
    11119         -  ``cutoff`` - integer depth to stop search. Ignored
    11120            if by_weight is True.
    11121        
     11134        -  ``cutoff`` - integer depth to stop search.
     11135
     11136           (ignored if ``by_weight == True``)
    1112211137       
    1112311138        EXAMPLES::
    1112411139       
     
    1114011155            sage: G.shortest_paths(0, by_weight=True)
    1114111156            {0: [0], 1: [0, 1], 2: [0, 1, 2], 3: [0, 1, 2, 3], 4: [0, 4]}
    1114211157        """
    11143         import networkx
     11158
    1114411159        if by_weight:
     11160            import networkx
    1114511161            return networkx.single_source_dijkstra_path(self.networkx_graph(copy=False), u)
    1114611162        else:
    1114711163            try:
     
    1116111177           search. If True, takes edge weightings into account, using
    1116211178           Dijkstra's algorithm.
    1116311179       
    11164         -  ``weight_sums`` - if False, returns the number of
    11165            edges in each path. If True, returns the sum of the weights of
    11166            these edges. Default behavior is to have the same value as
    11167            by_weight.
    11168        
    11169        
    1117011180        EXAMPLES::
    1117111181       
    1117211182            sage: D = graphs.DodecahedralGraph()
     
    1117711187            sage: G.shortest_path_lengths(0, by_weight=True)
    1117811188            {0: 0, 1: 1, 2: 2, 3: 3, 4: 2}
    1117911189        """
    11180         if weight_sums is None:
    11181             weight_sums = by_weight
    1118211190        paths = self.shortest_paths(u, by_weight)
    11183         if weight_sums:
     11191        if by_weight:
    1118411192            weights = {}
    1118511193            for v in paths:
    1118611194                wt = 0
     
    1125111259            1, or equivalently the length of a path is its number of
    1125211260            edges). Besides, they do not deal with graphs larger than 65536
    1125311261            vertices (which already represents 16GB of ram).
     11262
     11263        .. NOTE::
     11264
     11265           There is a Cython version of this method that is usually
     11266           much faster for large graphs, as most of the time is
     11267           actually spent building the final double
     11268           dictionary. Everything on the subject is to be found in the
     11269           :mod:`~sage.graphs.distances_all_pairs` module.
    1125411270       
    1125511271        EXAMPLES::
    1125611272       
     
    1134411360                algorithm = "Floyd-Warshall-Python"
    1134511361
    1134611362        if algorithm == "BFS":
    11347             from sage.graphs.base.c_graph import all_pairs_shortest_path_BFS
    11348             return all_pairs_shortest_path_BFS(self)
     11363            from sage.graphs.distances_all_pairs import distances_and_predecessors_all_pairs
     11364            return distances_and_predecessors_all_pairs(self)
    1134911365
    1135011366        elif algorithm == "Floyd-Warshall-Cython":
    11351             from sage.graphs.base.c_graph import floyd_warshall
     11367            from sage.graphs.distances_all_pairs import floyd_warshall
    1135211368            return floyd_warshall(self, distances = True)
    1135311369
    1135411370        elif algorithm != "Floyd-Warshall-Python":
  • sage/graphs/graph_decompositions/vertex_separation.pyx

    diff --git a/sage/graphs/graph_decompositions/vertex_separation.pyx b/sage/graphs/graph_decompositions/vertex_separation.pyx
    a b  
    100100.. NOTE::
    101101
    102102    Because of its current implementation, this algorithm only works on graphs
    103     on less than 32 vertices. This can be changed to 54 if necessary, but 32
     103    on less than 32 vertices. This can be changed to 64 if necessary, but 32
    104104    vertices already require 4GB of memory.
    105105
    106106**Lower bound on the vertex separation**
     
    141141
    142142def lower_bound(G):
    143143    r"""
    144     Returns a lower bound on the vertex separation of `D`
     144    Returns a lower bound on the vertex separation of `G`
    145145
    146146    INPUT:
    147147
     
    181181    if not isinstance(G, DiGraph):
    182182        raise ValueError("The parameter must be a DiGraph.")
    183183
     184    if G.order() >= 32:
     185        raise ValueError("The graph can have at most 31 vertices.")
     186
    184187    cdef FastDigraph FD = FastDigraph(G)
    185188    cdef int * g = FD.graph
    186189    cdef int n = FD.n
    187190
    188191    # minimums[i] is means to store the value of c'_{i+1}
    189192    minimums = <char *> sage_malloc(sizeof(char)* n)
    190     cdef int i
     193    cdef unsigned int i
    191194
    192195    # They are initialized to n
    193196    for 0<= i< n:
     
    196199    cdef char tmp, tmp_count
    197200
    198201    # We go through all sets
    199     for 1<= i< (1<<n):
     202    for 1<= i< <unsigned int> (1<<n):
    200203        tmp_count = <char> popcount(i)
    201204        tmp = <char> compute_out_neighborhood_cardinality(FD, i)
    202205
     
    320323    if not isinstance(G, DiGraph):
    321324        raise ValueError("The parameter must be a DiGraph.")
    322325
     326    if G.order() >= 32:
     327        raise ValueError("The graph should have at most 31 vertices !")
     328
    323329    cdef FastDigraph g = FastDigraph(G)
    324330
    325331    if verbose:
     
    328334
    329335    _sig_on
    330336
    331     cdef char * neighborhoods = <char *> sage_malloc((1<<g.n))
     337    cdef unsigned int mem = 1 << g.n
     338    cdef char * neighborhoods = <char *> sage_malloc(mem)
    332339
    333     memset(neighborhoods, <char> -1, (1<<g.n))
     340    if neighborhoods == NULL:
     341        raise MemoryError("Error allocating memory. I just tried to allocate "+str(mem>>10)+"MB, could that be too much ?")
     342
     343    memset(neighborhoods, <char> -1, mem)
    334344
    335345    cdef int i,j , k
    336346    for 0 <= k <g.n: