Ticket #10905: trac_10905.patch

File trac_10905.patch, 12.3 KB (added by ncohen, 10 years ago)
  • sage/graphs/base/c_graph.pyx

    # HG changeset patch
    # User Nathann Cohen <nathann.cohen@gmail.com>
    # Date 1299773000 -28800
    # Node ID d026f22e0d50bfc9ca1b06b0944bc5638a5c735e
    # Parent  b0a5118e9d0af16c871460c37f02d2c1dc28cc28
    trac 10905 -- shortest path all pairs through BFS computations
    
    diff -r b0a5118e9d0a -r d026f22e0d50 sage/graphs/base/c_graph.pyx
    a b  
    30103010    else:
    30113011        return d_dist
    30123012
     3013cpdef all_pairs_shortest_path_BFS(gg):
     3014    """
     3015    Computes the shortest path and the distances between each pair of vertices
     3016    through successive breadth-first-searches
     3017
     3018    OUTPUT:
     3019
     3020        This function return a pair of dictionaries ``(distances,paths)`` where
     3021        ``distance[u][v]`` denotes the distance of a shortest path from `u` to
     3022        `v` and ``paths[u][v]`` denotes an inneighbor `w` of `v` such that
     3023        `dist(u,v)= 1 + dist(u,w)`.
     3024
     3025    .. WARNING::
     3026
     3027        Because this function works on matrices whose size is quadratic compared
     3028        to the number of vertices, it uses short variables instead of long ones
     3029        to divide by 2 the size in memory. This means that the current
     3030        implementation does not run on a graph of more than 65536 nodes (this
     3031        can be easily changed if necessary, but would require much more
     3032        memory. It may be worth writing two versions). For information, the
     3033        current version of the algorithm on a graph with `65536=2^{16}` nodes
     3034        creates in memory `2` tables on `2^{32}` short elements (2bytes each),
     3035        for a total of `2^{34}` bytes or `16` gigabytes.
     3036
     3037    EXAMPLE:
     3038
     3039    On a grid::
     3040
     3041        sage: g = graphs.Grid2dGraph(10,10)
     3042        sage: from sage.graphs.base.c_graph import all_pairs_shortest_path_BFS
     3043        sage: dist, path = all_pairs_shortest_path_BFS(g)
     3044        sage: all( dist[u][v] == g.distance(u,v) for u in g for v in g )
     3045        True
     3046    """
     3047    from sage.rings.infinity import Infinity
     3048
     3049    cdef CGraph cg = <CGraph> gg._backend._cg
     3050
     3051    cdef list vertices = cg.verts()
     3052    cdef int n = max(vertices)+1
     3053
     3054    if n > <unsigned short> -1:
     3055        raise ValueError("The graph backend contains more than "+(<unsigned short> -1)+" nodes")
     3056
     3057    # The vertices which have already been visited
     3058    cdef bitset_t seen
     3059    bitset_init(seen, cg.active_vertices.size)
     3060   
     3061    # The list of waiting vertices, the beginning and the end of the list
     3062
     3063    cdef unsigned short * waiting_list = <unsigned short *> sage_malloc(n*sizeof(short))
     3064    cdef unsigned short waiting_beginning = 0
     3065    cdef unsigned short waiting_end = 0
     3066
     3067    cdef unsigned short source
     3068    cdef unsigned short v, u
     3069
     3070    # Dictionaries of dictionaries of distances/predecessors
     3071    cdef dict d_distances = dict()
     3072    cdef dict d_prec      = dict()
     3073
     3074    # Temporary dictionaries of distances/predecessors
     3075    cdef dict tmp_distances
     3076    cdef dict tmp_prec
     3077
     3078    # Temporary vectors of distances/predecessors
     3079    cdef unsigned short * v_distances = <unsigned short *> sage_malloc(n*sizeof(short))
     3080    cdef unsigned short * v_prec      = <unsigned short *> sage_malloc(n*sizeof(short))
     3081
     3082    for source in vertices:
     3083        bitset_set_first_n(seen, 0)
     3084        bitset_add(seen, source)
     3085
     3086        v_distances[source] = 0
     3087        v_prec[source] = source
     3088        waiting_list[0] = source
     3089        waiting_beginning = 0
     3090        waiting_end = 0
     3091
     3092        while waiting_beginning <= waiting_end:
     3093            v = waiting_list[waiting_beginning]
     3094
     3095            for u in cg.out_neighbors(v):
     3096                if not bitset_in(seen, u):
     3097                    v_distances[u] = v_distances[v]+1
     3098                    v_prec[u] = v
     3099                    bitset_add(seen, u)
     3100                    waiting_end += 1
     3101                    waiting_list[waiting_end] = u
     3102
     3103            waiting_beginning += 1
     3104
     3105        tmp_distances = dict()
     3106        tmp_prec = dict()
     3107        for v in vertices:
     3108            vv = vertex_label(v, gg._backend.vertex_ints, gg._backend.vertex_labels, gg._backend._cg)
     3109
     3110            if bitset_in(seen, v):
     3111                tmp_prec[vv] = vertex_label(v_prec[v], gg._backend.vertex_ints, gg._backend.vertex_labels, gg._backend._cg)
     3112                tmp_distances[vv] = v_distances[v]
     3113            else:
     3114                tmp_prec[vv] = None
     3115                tmp_distances[vv] = Infinity
     3116
     3117        vv = vertex_label(source, gg._backend.vertex_ints, gg._backend.vertex_labels, gg._backend._cg)
     3118        tmp_prec[vv] = None
     3119        d_prec[vv] = tmp_prec
     3120        d_distances[vv] = tmp_distances
     3121
     3122    bitset_free(seen)
     3123    sage_free(waiting_list)
     3124    sage_free(v_distances)
     3125    sage_free(v_prec)
     3126
     3127    return d_distances, d_prec
     3128       
     3129
    30133130cdef class Search_iterator:
    30143131    r"""
    30153132    An iterator for traversing a (di)graph.
  • sage/graphs/generic_graph.py

    diff -r b0a5118e9d0a -r d026f22e0d50 sage/graphs/generic_graph.py
    a b  
    1009210092        """
    1009310093        return self.shortest_path_length(u, v)
    1009410094
    10095     def distance_all_pairs(self, algorithm = "BFS"):
     10095    def distance_all_pairs(self, algorithm = "auto"):
    1009610096        r"""
    1009710097        Returns the distances between all pairs of vertices.
    1009810098
     
    1010110101            - ``"algorithm"`` (string) -- two algorithms are available
    1010210102
    1010310103                * ``algorithm = "BFS"`` in which case the distances are computed
    10104                   through `n-1` different breadth-first-search.
     10104                  through `n` different breadth-first-search.
    1010510105
    1010610106                * ``algorithm = "Floyd-Warshall"``, in which case the
    1010710107                  Floyd-Warshall algorithm is used.
    1010810108
     10109                * ``algorithm = "auto"``, in which case the Floyd-Warshall
     10110                  algorithm is used for graphs on less than 20 vertices, and BFS
     10111                  otherwise.
     10112
    1010910113                The default is ``algorithm = "BFS"``.
    1011010114
    1011110115        OUTPUT:
     
    1012710131            sage: all([g.distance(0,v) == distances[0][v] for v in g])
    1012810132            True
    1012910133        """
     10134        if algorithm == "auto":
     10135            if self.order() <= 20:
     10136                algorithm = "Floyd-Warshall"
     10137            else:
     10138                algorithm = "BFS"
     10139
    1013010140        if algorithm == "BFS":
    10131             from sage.rings.infinity import Infinity
    10132             distances = dict([(v, self.shortest_path_lengths(v)) for v in self])
    10133        
    10134             # setting the necessary +Infinity
    10135             cc = self.connected_components()
    10136             for cc1 in cc:
    10137                 for cc2 in cc:
    10138                     if cc1 != cc2:
    10139                         for u in cc1:
    10140                             for v in cc2:
    10141                                 distances[u][v] = Infinity
    10142 
    10143             return distances
     10141            from sage.graphs.base.c_graph import all_pairs_shortest_path_BFS
     10142            return all_pairs_shortest_path_BFS(self)[0]
    1014410143
    1014510144        elif algorithm == "Floyd-Warshall":
    1014610145            from sage.graphs.base.c_graph import floyd_warshall
    1014710146            return floyd_warshall(self,paths = False, distances = True)
    1014810147
    1014910148        else:
    10150             raise ValueError("The algorithm keyword can be equal to either \"BFS\" or \"Floyd-Warshall\"")
     10149            raise ValueError("The algorithm keyword can be equal to either \"BFS\" or \"Floyd-Warshall\" or \"auto\"")
    1015110150
    1015210151    def eccentricity(self, v=None, dist_dict=None, with_labels=False):
    1015310152        """
     
    1093510934                lengths[v] = len(paths[v]) - 1
    1093610935            return lengths
    1093710936
    10938     def shortest_path_all_pairs(self, by_weight=False, default_weight=1):
    10939         """
    10940         Uses the Floyd-Warshall algorithm to find a shortest weighted path
    10941         for each pair of vertices.
    10942        
     10937    def shortest_path_all_pairs(self, by_weight=False, default_weight=1, algorithm = "auto"):
     10938        """
     10939        Computes a shortest path between each pair of vertices.
     10940
    1094310941        INPUT:
    1094410942       
    1094510943       
     
    1095110949           edges that don't have a weight (i.e., a label).
    1095210950
    1095310951           Implies ``by_weight == True``.
     10952
     10953        - ``algorithm`` -- four options :
     10954       
     10955           * ``"BFS"`` -- the computation is done through a BFS
     10956             centered on each vertex successively. Only implemented
     10957             when ``default_weight = 1`` and ``by_weight = False``.
     10958
     10959           * ``"Floyd-Warshall-Cython"`` -- through the Cython implementation of
     10960             the Floyd-Warshall algorithm.
     10961
     10962           * ``"Floyd-Warshall-Python"`` -- through the Python implementation of
     10963             the Floyd-Warshall algorithm.
     10964
     10965           * ``"auto"`` -- use the fastest algorithm depending on the input
     10966             (``"BFS"`` if possible, and ``"Floyd-Warshall-Python"`` otherwise)
     10967
     10968             This is the default value.
    1095410969       
    1095510970        OUTPUT:
    1095610971
     
    1096210977
    1096310978        .. NOTE::
    1096410979
    10965             Two version of the Floyd-Warshall algorithm are implemented. One is
    10966             pure Python, the other one is Cython. The first is used whenever
    10967             ``by_weight = True``, and can perform its computations using the
    10968             labels on the edges for as long as they can be added together. The
    10969             second one, used when ``by_weight = False``, is much faster but only
    10970             deals with the topological distance (each edge is of weight 1, or
    10971             equivalently the length of a path is its number of edges). Besides,
    10972             the current version of this second implementation does not deal with
    10973             graphs larger than 65536 vertices (which already represents 16GB of
    10974             ram when running the Floyd-Warshall algorithm).
     10980            Three different implementations are actually available through this method :
     10981
     10982                * BFS (Cython)
     10983                * Floyd-Warshall (Cython)
     10984                * Floyd-Warshall (Python)
     10985
     10986            The BFS algorithm is the fastest of the three, then comes the Cython
     10987            implementation of Floyd-Warshall, and last the Python
     10988            implementation. The first two implementations, however, only compute
     10989            distances based on the topological distance (each edge is of weight
     10990            1, or equivalently the length of a path is its number of
     10991            edges). Besides, they do not deal with graphs larger than 65536
     10992            vertices (which already represents 16GB of ram).
    1097510993       
    1097610994        EXAMPLES::
    1097710995       
     
    1102611044            2: {0: 4, 1: 2, 2: None, 3: 2, 4: 3},
    1102711045            3: {0: 4, 1: 2, 2: 3, 3: None, 4: 3},
    1102811046            4: {0: 4, 1: 0, 2: 3, 3: 4, 4: None}})
     11047
     11048        Checking the distances are equal regardless of the algorithm used::
     11049
     11050            sage: g = graphs.Grid2dGraph(5,5)
     11051            sage: d1, _ = g.shortest_path_all_pairs(algorithm="BFS")
     11052            sage: d2, _ = g.shortest_path_all_pairs(algorithm="Floyd-Warshall-Cython")
     11053            sage: d3, _ = g.shortest_path_all_pairs(algorithm="Floyd-Warshall-Python")
     11054            sage: d1 == d2 == d3
     11055            True
     11056
     11057        Checking a random path is valid ::
     11058   
     11059            sage: dist, path = g.shortest_path_all_pairs(algorithm="BFS")
     11060            sage: u,v = g.random_vertex(), g.random_vertex()
     11061            sage: p = [v]
     11062            sage: while p[0] != None:
     11063            ...     p.insert(0,path[u][p[0]])
     11064            sage: len(p) == dist[u][v] + 2
     11065            True
    1102911066        """
    1103011067        if default_weight != 1:
    1103111068            by_weight = True
    1103211069
    11033         if by_weight is False:
     11070        if algorithm == "auto":
     11071            if by_weight is False:
     11072                algorithm = "BFS"
     11073            else:
     11074                algorithm = "Floyd-Warshall-Python"
     11075
     11076        if algorithm == "BFS":
     11077            from sage.graphs.base.c_graph import all_pairs_shortest_path_BFS
     11078            return all_pairs_shortest_path_BFS(self)
     11079
     11080        elif algorithm == "Floyd-Warshall-Cython":
    1103411081            from sage.graphs.base.c_graph import floyd_warshall
    1103511082            return floyd_warshall(self, distances = True)
    1103611083
     11084        elif algorithm != "Floyd-Warshall-Python":
     11085            raise ValueError("The algorithm keyword can only be set to "+
     11086                             "\"auto\","+
     11087                             " \"BFS\", "+
     11088                             "\"Floyd-Warshall-Cython\" or "+
     11089                             "\"Floyd-Warshall-Cython\"")
     11090
    1103711091        from sage.rings.infinity import Infinity
    1103811092        dist = {}
    1103911093        pred = {}