Add {all_pairs,single_source}_bellman_ford_path_length

graphblas_algorithms/algorithms/exceptions.py

@@ -12,3 +12,7 @@ class EmptyGraphError(GraphBlasAlgorithmException):
 
 class PointlessConcept(GraphBlasAlgorithmException):
     pass
+
+
+class Unbounded(GraphBlasAlgorithmException):
+    pass

graphblas_algorithms/algorithms/shortest_paths/__init__.py

@@ -1,2 +1,3 @@
 from .dense import *
 from .generic import *
+from .weighted import *

graphblas_algorithms/algorithms/shortest_paths/weighted.py

@@ -0,0 +1,45 @@
+from graphblas import Vector, binary, monoid, replace
+from graphblas.semiring import min_plus
+
+from ..exceptions import Unbounded
+
+__all__ = ["single_source_bellman_ford_path_length"]
+
+
+def single_source_bellman_ford_path_length(G, source):
+    # No need for `is_weighted=` keyword, b/c this is assumed to be weighted (I think)
+    index = G._key_to_id[source]
+    A = G._A
+    if A.dtype == bool:
+        # Should we upcast e.g. INT8 to INT64 as well?
+        dtype = int
+    else:
+        dtype = A.dtype
+    n = A.nrows
+    d = Vector(dtype, n, name="single_source_bellman_ford_path_length")
+    d[index] = 0
+    cur = d.dup(name="cur")
+    mask = Vector(bool, n, name="mask")
+    for _i in range(n - 1):
+        # This is a slightly modified Bellman-Ford algorithm.
+        # `cur` is the current frontier of values that improved in the previous iteration.
+        # This means that in this iteration we drop values from `cur` that are not better.
+        cur << min_plus(cur @ A)
+
+        # Mask is True where cur not in d or cur < d
+        mask(cur.S, replace) << True  # or: `mask << unary.one[bool](cur)`
+        mask(binary.second) << binary.lt(cur & d)
+
+        # Drop values from `cur` that didn't improve
+        cur(mask.V, replace) << cur
+        if cur.nvals == 0:
+            break
+        # Update `d` with values that improved
+        d(cur.S) << cur
+    else:
+        # Check for negative cycle when for loop completes without breaking
+        cur << min_plus(cur @ A)
+        mask << binary.lt(cur & d)
+        if mask.reduce(monoid.lor):
+            raise Unbounded("Negative cycle detected.")
+    return d

graphblas_algorithms/interface.py

@@ -60,6 +60,12 @@ class Dispatcher:
         nxapi.shortest_paths.dense.floyd_warshall_predecessor_and_distance
     )
     has_path = nxapi.shortest_paths.generic.has_path
+    all_pairs_bellman_ford_path_length = (
+        nxapi.shortest_paths.weighted.all_pairs_bellman_ford_path_length
+    )
+    single_source_bellman_ford_path_length = (
+        nxapi.shortest_paths.weighted.single_source_bellman_ford_path_length
+    )
     # Simple Paths
     is_simple_path = nxapi.simple_paths.is_simple_path
     # S Metric
@@ -103,13 +109,23 @@ def on_start_tests(items):
             import pytest
         except ImportError:  # pragma: no cover (import)
             return
+
+        def key(testpath):
+            filename, path = testpath.split(":")
+            classname, testname = path.split(".")
+            return (testname, frozenset({classname, filename}))
+
+        # Reasons to skip tests
         multi_attributed = "unable to handle multi-attributed graphs"
         multidigraph = "unable to handle MultiDiGraph"
-        freeze = frozenset
+        multigraph = "unable to handle MultiGraph"
+
+        # Which tests to skip
         skip = {
-            ("test_attributes", freeze({"TestBoruvka", "test_mst.py"})): multi_attributed,
-            ("test_weight_attribute", freeze({"TestBoruvka", "test_mst.py"})): multi_attributed,
-            ("test_zero_weight", freeze({"TestFloyd", "test_dense.py"})): multidigraph,
+            key("test_mst.py:TestBoruvka.test_attributes"): multi_attributed,
+            key("test_mst.py:TestBoruvka.test_weight_attribute"): multi_attributed,
+            key("test_dense.py:TestFloyd.test_zero_weight"): multidigraph,
+            key("test_weighted.py:TestBellmanFordAndGoldbergRadzik.test_multigraph"): multigraph,
         }
         for item in items:
             kset = set(item.keywords)

graphblas_algorithms/nxapi/exception.py

@@ -8,6 +8,9 @@ class NetworkXError(Exception):
     class NetworkXPointlessConcept(Exception):
         pass
 
+    class NetworkXUnbounded(Exception):
+        pass
+
     class NodeNotFound(Exception):
         pass
 
@@ -18,6 +21,7 @@ class PowerIterationFailedConvergence(Exception):
     from networkx import (
         NetworkXError,
         NetworkXPointlessConcept,
+        NetworkXUnbounded,
         NodeNotFound,
         PowerIterationFailedConvergence,
     )

graphblas_algorithms/nxapi/shortest_paths/__init__.py

@@ -1,2 +1,3 @@
 from .dense import *
 from .generic import *
+from .weighted import *

graphblas_algorithms/nxapi/shortest_paths/weighted.py

@@ -0,0 +1,36 @@
+from graphblas_algorithms import algorithms
+from graphblas_algorithms.classes.digraph import to_graph
+
+from ..exception import NetworkXUnbounded, NodeNotFound
+
+__all__ = [
+    "all_pairs_bellman_ford_path_length",
+    "single_source_bellman_ford_path_length",
+]
+
+
+def all_pairs_bellman_ford_path_length(G, weight="weight"):
+    # TODO: what if weight is a function?
+    # How should we implement and call `algorithms.all_pairs_bellman_ford_path_length`?
+    # Should we compute in chunks to expose more parallelism?
+    G = to_graph(G, weight=weight)
+    for source in G:
+        try:
+            d = algorithms.single_source_bellman_ford_path_length(G, source)
+        except algorithms.exceptions.Unbounded as e:
+            raise NetworkXUnbounded(*e.args) from e
+        except KeyError as e:
+            raise NodeNotFound(*e.args) from e
+        yield (source, G.vector_to_nodemap(d))
+
+
+def single_source_bellman_ford_path_length(G, source, weight="weight"):
+    # TODO: what if weight is a function?
+    G = to_graph(G, weight=weight)
+    try:
+        d = algorithms.single_source_bellman_ford_path_length(G, source)
+    except algorithms.exceptions.Unbounded as e:
+        raise NetworkXUnbounded(*e.args) from e
+    except KeyError as e:
+        raise NodeNotFound(*e.args) from e
+    return G.vector_to_nodemap(d)