# HG changeset patch
# User Nathann Cohen
# Date 1299446050 3600
# Node ID 954991d687ee27a6abc8cbbbb02cd96c7177563e
# Parent 361a4ad7d52c69b64ae2e658ffd0820af0d87e93
Ticket #10885  Cython implementation of FloydWarshall algorithm
diff r 361a4ad7d52c r 954991d687ee sage/graphs/base/c_graph.pyx
 a/sage/graphs/base/c_graph.pyx Fri Feb 25 18:56:01 2011 +0000
+++ b/sage/graphs/base/c_graph.pyx Sun Mar 06 22:14:10 2011 +0100
@@ 2851,6 +2851,165 @@
else:
return True
+def floyd_warshall(gg, paths = True, distances = False):
+ r"""
+ Computes the shortest path/distances between all pairs of vertices.
+
+ For more information on the FloydWarshall algorithm, see the `Wikipedia
+ article on FloydWarshall
+ `_.
+
+ INPUT:
+
+  ``gg``  the graph on which to work.
+
+  ``paths`` (boolean)  whether to return the dictionary of shortest
+ paths. Set to ``True`` by default.
+
+  ``distances`` (boolean)  whether to return the dictionary of
+ distances. Set to ``False`` by default.
+
+ OUTPUT:
+
+ Depending on the input, this function return the dictionary of paths,
+ the dictionary of distances, or a pair of dictionaries
+ ``(distances,paths)`` where ``distance[u][v]`` denotes the distance of a
+ shortest path from `u` to `v` and ``paths[u][v]`` denotes an inneighbor
+ `w` of `v` such that `dist(u,v)= 1 + dist(u,w)`.
+
+ .. WARNING::
+
+ Because this function works on matrices whose size is quadratic compared
+ to the number of vertices, it uses short variables instead of long ones
+ to divide by 2 the size in memory. This means that the current
+ implementation does not run on a graph of more than 65536 nodes (this
+ can be easily changed if necessary, but would require much more
+ memory. It may be worth writing two versions). For information, the
+ current version of the algorithm on a graph with `65536=2^{16}` nodes
+ creates in memory `2` tables on `2^{32}` short elements (2bytes each),
+ for a total of `2^{34}` bytes or `16` gigabytes. Let us also remember
+ that if the memory size is quadratic, the algorithm runs in cubic time.
+
+ EXAMPLES:
+
+ Shortest paths in a small grid ::
+
+ sage: g = graphs.Grid2dGraph(2,2)
+ sage: from sage.graphs.base.c_graph import floyd_warshall
+ sage: print floyd_warshall(g)
+ {(0, 1): {(0, 1): None, (1, 0): (0, 0), (0, 0): (0, 1), (1, 1): (0, 1)},
+ (1, 0): {(0, 1): (0, 0), (1, 0): None, (0, 0): (1, 0), (1, 1): (1, 0)},
+ (0, 0): {(0, 1): (0, 0), (1, 0): (0, 0), (0, 0): None, (1, 1): (0, 1)},
+ (1, 1): {(0, 1): (1, 1), (1, 0): (1, 1), (0, 0): (0, 1), (1, 1): None}}
+
+ Checking the distances are correct ::
+
+ sage: g = graphs.Grid2dGraph(5,5)
+ sage: dist,path = floyd_warshall(g, distances = True)
+ sage: all( dist[u][v] == g.distance(u,v) for u in g for v in g )
+ True
+
+ Checking a random path is valid ::
+
+ sage: u,v = g.random_vertex(), g.random_vertex()
+ sage: p = [v]
+ sage: while p[0] != None:
+ ... p.insert(0,path[u][p[0]])
+ sage: len(p) == dist[u][v] + 2
+ True
+ """
+ from sage.rings.infinity import Infinity
+ cdef CGraph g = gg._backend._cg
+ cdef int n = max(g.verts())+1
+
+ if n > 1:
+ raise ValueError("The graph backend contains more than "+( 1)+" nodes")
+
+ # Dictionaries of distance, precedent element, and integers
+ cdef dict d_prec = dict()
+ cdef dict d_dist = dict()
+ cdef dict tmp
+
+ cdef int v_int
+ cdef int u_int
+ cdef int w_int
+ cdef int i
+
+ # All this just creates two tables prec[n][n] and dist[n][n]
+ cdef unsigned short * t_prec = sage_malloc(n*n*sizeof(short))
+ cdef unsigned short * t_dist = sage_malloc(n*n*sizeof(short))
+ cdef unsigned short ** prec = sage_malloc(n*sizeof(short *))
+ cdef unsigned short ** dist = sage_malloc(n*sizeof(short *))
+ prec[0] = t_prec
+ dist[0] = t_dist
+
+ for 1 <= i< n:
+ prec[i] = prec[i1] + n
+ dist[i] = dist[i1] + n
+
+ # Initializing prec and dist
+ memset(t_prec, 0, n*n*sizeof(short))
+ memset(t_dist, 1, n*n*sizeof(short))
+
+ # Copying the adjacency matrix (vertices at distance 1)
+ for v_int in g.verts():
+ prec[v_int][v_int] = v_int
+ dist[v_int][v_int] = 0
+ for u_int in g.out_neighbors(v_int):
+ dist[v_int][u_int] = 1
+ prec[v_int][u_int] = v_int
+
+ # The algorithm itself.
+
+ for w_int in g.verts():
+ for v_int in g.verts():
+ for u_int in g.verts():
+
+ # If it is shorter to go from u to v through w, do it
+ if dist[v_int][u_int] > dist[v_int][w_int] + dist[w_int][u_int]:
+ dist[v_int][u_int] = dist[v_int][w_int] + dist[w_int][u_int]
+ prec[v_int][u_int] = prec[w_int][u_int]
+
+ # If the paths are to be returned
+ if paths:
+ for v_int in g.verts():
+ tmp = dict()
+ v = vertex_label(v_int, gg._backend.vertex_ints, gg._backend.vertex_labels, gg._backend._cg)
+
+ for u_int in g.verts():
+ u = vertex_label(u_int, gg._backend.vertex_ints, gg._backend.vertex_labels, gg._backend._cg)
+ w = (None if v == u
+ else vertex_label(prec[v_int][u_int], gg._backend.vertex_ints, gg._backend.vertex_labels, gg._backend._cg))
+ tmp[u] = w
+
+ d_prec[v] = tmp
+
+ sage_free(t_prec)
+ sage_free(prec)
+
+ # If the distances are to be returned
+ if distances:
+ for v_int in g.verts():
+ tmp = dict()
+ v = vertex_label(v_int, gg._backend.vertex_ints, gg._backend.vertex_labels, gg._backend._cg)
+
+ for u_int in g.verts():
+ u = vertex_label(u_int, gg._backend.vertex_ints, gg._backend.vertex_labels, gg._backend._cg)
+
+ tmp[u] = dist[v_int][u_int] if (dist[v_int][u_int] != 1) else Infinity
+
+ d_dist[v] = tmp
+
+ sage_free(t_dist)
+ sage_free(dist)
+
+ if distances and paths:
+ return d_dist, d_prec
+ elif paths:
+ return d_prec
+ else:
+ return d_dist
+
cdef class Search_iterator:
r"""
An iterator for traversing a (di)graph.
diff r 361a4ad7d52c r 954991d687ee sage/graphs/generic_graph.py
 a/sage/graphs/generic_graph.py Fri Feb 25 18:56:01 2011 +0000
+++ b/sage/graphs/generic_graph.py Sun Mar 06 22:14:10 2011 +0100
@@ 10092,10 +10092,22 @@
"""
return self.shortest_path_length(u, v)
 def distance_all_pairs(self):
+ def distance_all_pairs(self, algorithm = "BFS"):
r"""
Returns the distances between all pairs of vertices.
+ INPUT:
+
+  ``"algorithm"`` (string)  two algorithms are available
+
+ * ``algorithm = "BFS"`` in which case the distances are computed
+ through `n1` different breadthfirstsearch.
+
+ * ``algorithm = "FloydWarshall"``, in which case the
+ FloydWarshall algorithm is used.
+
+ The default is ``algorithm = "BFS"``.
+
OUTPUT:
A doubly indexed dictionary
@@ 10115,20 +10127,27 @@
sage: all([g.distance(0,v) == distances[0][v] for v in g])
True
"""

 from sage.rings.infinity import Infinity
 distances = dict([(v, self.shortest_path_lengths(v)) for v in self])

 # setting the necessary +Infinity
 cc = self.connected_components()
 for cc1 in cc:
 for cc2 in cc:
 if cc1 != cc2:
 for u in cc1:
 for v in cc2:
 distances[u][v] = Infinity

 return distances
+ if algorithm == "BFS":
+ from sage.rings.infinity import Infinity
+ distances = dict([(v, self.shortest_path_lengths(v)) for v in self])
+
+ # setting the necessary +Infinity
+ cc = self.connected_components()
+ for cc1 in cc:
+ for cc2 in cc:
+ if cc1 != cc2:
+ for u in cc1:
+ for v in cc2:
+ distances[u][v] = Infinity
+
+ return distances
+
+ elif algorithm == "FloydWarshall":
+ from sage.graphs.base.c_graph import floyd_warshall
+ return floyd_warshall(self,paths = False, distances = True)
+
+ else:
+ raise ValueError("The algorithm keyword can be equal to either \"BFS\" or \"FloydWarshall\"")
def eccentricity(self, v=None, dist_dict=None, with_labels=False):
"""
@@ 10916,34 +10935,49 @@
lengths[v] = len(paths[v])  1
return lengths
 def shortest_path_all_pairs(self, by_weight=True, default_weight=1):
+ def shortest_path_all_pairs(self, by_weight=False, default_weight=1):
"""
Uses the FloydWarshall algorithm to find a shortest weighted path
for each pair of vertices.
 The weights (labels) on the vertices can be anything that can be
 compared and can be summed.

 INPUT:


  ``by_weight``  If False, figure distances by the
 numbers of edges.

  ``default_weight``  (defaults to 1) The default
 weight to assign edges that don't have a weight (i.e., a label).


 OUTPUT: A tuple (dist, pred). They are both dicts of dicts. The
 first indicates the length dist[u][v] of the shortest weighted path
 from u to v. The second is more complicated it indicates the
 predecessor pred[u][v] of v in the shortest path from u to v.
+ INPUT:
+
+
+  ``by_weight``  Whether to use the labels defined over the edges as
+ weights. If ``False`` (default), the distance between `u` and `v` is
+ the minimum number of edges of a path from `u` to `v`.
+
+  ``default_weight``  (defaults to 1) The default weight to assign
+ edges that don't have a weight (i.e., a label).
+
+ Implies ``by_weight == True``.
+
+ OUTPUT:
+
+ A tuple ``(dist, pred)``. They are both dicts of dicts. The first
+ indicates the length ``dist[u][v]`` of the shortest weighted path
+ from `u` to `v`. The second is a compact representation of all the
+ paths it indicates the predecessor ``pred[u][v]`` of `v` in the
+ shortest path from `u` to `v`.
+
+ .. NOTE::
+
+ Two version of the FloydWarshall algorithm are implemented. One is
+ pure Python, the other one is Cython. The first is used whenever
+ ``by_weight = True``, and can perform its computations using the
+ labels on the edges for as long as they can be added together. The
+ second one, used when ``by_weight = False``, is much faster but only
+ deals with the topological distance (each edge is of weight 1, or
+ equivalently the length of a path is its number of edges). Besides,
+ the current version of this second implementation does not deal with
+ graphs larger than 65536 vertices (which already represents 16GB of
+ ram when running the FloydWarshall algorithm).
EXAMPLES::
sage: G = Graph( { 0: {1: 1}, 1: {2: 1}, 2: {3: 1}, 3: {4: 2}, 4: {0: 2} }, sparse=True )
sage: G.plot(edge_labels=True).show() # long time
 sage: dist, pred = G.shortest_path_all_pairs()
+ sage: dist, pred = G.shortest_path_all_pairs(by_weight = True)
sage: dist
{0: {0: 0, 1: 1, 2: 2, 3: 3, 4: 2}, 1: {0: 1, 1: 0, 2: 1, 3: 2, 4: 3}, 2: {0: 2, 1: 1, 2: 0, 3: 1, 4: 3}, 3: {0: 3, 1: 2, 2: 1, 3: 0, 4: 2}, 4: {0: 2, 1: 3, 2: 3, 3: 2, 4: 0}}
sage: pred
@@ 10951,15 +10985,15 @@
sage: pred[0]
{0: None, 1: 0, 2: 1, 3: 2, 4: 0}
 So for example the shortest weighted path from 0 to 3 is obtained
 as follows. The predecessor of 3 is pred[0][3] == 2, the
 predecessor of 2 is pred[0][2] == 1, and the predecessor of 1 is
 pred[0][1] == 0.
+ So for example the shortest weighted path from `0` to `3` is obtained as
+ follows. The predecessor of `3` is ``pred[0][3] == 2``, the predecessor
+ of `2` is ``pred[0][2] == 1``, and the predecessor of `1` is
+ ``pred[0][1] == 0``.
::
sage: G = Graph( { 0: {1:None}, 1: {2:None}, 2: {3: 1}, 3: {4: 2}, 4: {0: 2} }, sparse=True )
 sage: G.shortest_path_all_pairs(by_weight=False)
+ sage: G.shortest_path_all_pairs()
({0: {0: 0, 1: 1, 2: 2, 3: 2, 4: 1},
1: {0: 1, 1: 0, 2: 1, 3: 2, 4: 2},
2: {0: 2, 1: 1, 2: 0, 3: 1, 4: 2},
@@ 10970,7 +11004,7 @@
2: {0: 1, 1: 2, 2: None, 3: 2, 4: 3},
3: {0: 4, 1: 2, 2: 3, 3: None, 4: 3},
4: {0: 4, 1: 0, 2: 3, 3: 4, 4: None}})
 sage: G.shortest_path_all_pairs()
+ sage: G.shortest_path_all_pairs(by_weight = True)
({0: {0: 0, 1: 1, 2: 2, 3: 3, 4: 2},
1: {0: 1, 1: 0, 2: 1, 3: 2, 4: 3},
2: {0: 2, 1: 1, 2: 0, 3: 1, 4: 3},
@@ 10993,6 +11027,13 @@
3: {0: 4, 1: 2, 2: 3, 3: None, 4: 3},
4: {0: 4, 1: 0, 2: 3, 3: 4, 4: None}})
"""
+ if default_weight != 1:
+ by_weight = True
+
+ if by_weight is False:
+ from sage.graphs.base.c_graph import floyd_warshall
+ return floyd_warshall(self, distances = True)
+
from sage.rings.infinity import Infinity
dist = {}
pred = {}