diff env/lib/python3.9/site-packages/networkx/algorithms/flow/preflowpush.py @ 0:4f3585e2f14b draft default tip

"planemo upload commit 60cee0fc7c0cda8592644e1aad72851dec82c959"
author shellac
date Mon, 22 Mar 2021 18:12:50 +0000
parents
children
line wrap: on
line diff
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/env/lib/python3.9/site-packages/networkx/algorithms/flow/preflowpush.py	Mon Mar 22 18:12:50 2021 +0000
@@ -0,0 +1,425 @@
+"""
+Highest-label preflow-push algorithm for maximum flow problems.
+"""
+
+from collections import deque
+from itertools import islice
+import networkx as nx
+from ...utils import arbitrary_element
+from .utils import build_residual_network
+from .utils import CurrentEdge
+from .utils import detect_unboundedness
+from .utils import GlobalRelabelThreshold
+from .utils import Level
+
+__all__ = ["preflow_push"]
+
+
+def preflow_push_impl(G, s, t, capacity, residual, global_relabel_freq, value_only):
+    """Implementation of the highest-label preflow-push algorithm.
+    """
+    if s not in G:
+        raise nx.NetworkXError(f"node {str(s)} not in graph")
+    if t not in G:
+        raise nx.NetworkXError(f"node {str(t)} not in graph")
+    if s == t:
+        raise nx.NetworkXError("source and sink are the same node")
+
+    if global_relabel_freq is None:
+        global_relabel_freq = 0
+    if global_relabel_freq < 0:
+        raise nx.NetworkXError("global_relabel_freq must be nonnegative.")
+
+    if residual is None:
+        R = build_residual_network(G, capacity)
+    else:
+        R = residual
+
+    detect_unboundedness(R, s, t)
+
+    R_nodes = R.nodes
+    R_pred = R.pred
+    R_succ = R.succ
+
+    # Initialize/reset the residual network.
+    for u in R:
+        R_nodes[u]["excess"] = 0
+        for e in R_succ[u].values():
+            e["flow"] = 0
+
+    def reverse_bfs(src):
+        """Perform a reverse breadth-first search from src in the residual
+        network.
+        """
+        heights = {src: 0}
+        q = deque([(src, 0)])
+        while q:
+            u, height = q.popleft()
+            height += 1
+            for v, attr in R_pred[u].items():
+                if v not in heights and attr["flow"] < attr["capacity"]:
+                    heights[v] = height
+                    q.append((v, height))
+        return heights
+
+    # Initialize heights of the nodes.
+    heights = reverse_bfs(t)
+
+    if s not in heights:
+        # t is not reachable from s in the residual network. The maximum flow
+        # must be zero.
+        R.graph["flow_value"] = 0
+        return R
+
+    n = len(R)
+    # max_height represents the height of the highest level below level n with
+    # at least one active node.
+    max_height = max(heights[u] for u in heights if u != s)
+    heights[s] = n
+
+    grt = GlobalRelabelThreshold(n, R.size(), global_relabel_freq)
+
+    # Initialize heights and 'current edge' data structures of the nodes.
+    for u in R:
+        R_nodes[u]["height"] = heights[u] if u in heights else n + 1
+        R_nodes[u]["curr_edge"] = CurrentEdge(R_succ[u])
+
+    def push(u, v, flow):
+        """Push flow units of flow from u to v.
+        """
+        R_succ[u][v]["flow"] += flow
+        R_succ[v][u]["flow"] -= flow
+        R_nodes[u]["excess"] -= flow
+        R_nodes[v]["excess"] += flow
+
+    # The maximum flow must be nonzero now. Initialize the preflow by
+    # saturating all edges emanating from s.
+    for u, attr in R_succ[s].items():
+        flow = attr["capacity"]
+        if flow > 0:
+            push(s, u, flow)
+
+    # Partition nodes into levels.
+    levels = [Level() for i in range(2 * n)]
+    for u in R:
+        if u != s and u != t:
+            level = levels[R_nodes[u]["height"]]
+            if R_nodes[u]["excess"] > 0:
+                level.active.add(u)
+            else:
+                level.inactive.add(u)
+
+    def activate(v):
+        """Move a node from the inactive set to the active set of its level.
+        """
+        if v != s and v != t:
+            level = levels[R_nodes[v]["height"]]
+            if v in level.inactive:
+                level.inactive.remove(v)
+                level.active.add(v)
+
+    def relabel(u):
+        """Relabel a node to create an admissible edge.
+        """
+        grt.add_work(len(R_succ[u]))
+        return (
+            min(
+                R_nodes[v]["height"]
+                for v, attr in R_succ[u].items()
+                if attr["flow"] < attr["capacity"]
+            )
+            + 1
+        )
+
+    def discharge(u, is_phase1):
+        """Discharge a node until it becomes inactive or, during phase 1 (see
+        below), its height reaches at least n. The node is known to have the
+        largest height among active nodes.
+        """
+        height = R_nodes[u]["height"]
+        curr_edge = R_nodes[u]["curr_edge"]
+        # next_height represents the next height to examine after discharging
+        # the current node. During phase 1, it is capped to below n.
+        next_height = height
+        levels[height].active.remove(u)
+        while True:
+            v, attr = curr_edge.get()
+            if height == R_nodes[v]["height"] + 1 and attr["flow"] < attr["capacity"]:
+                flow = min(R_nodes[u]["excess"], attr["capacity"] - attr["flow"])
+                push(u, v, flow)
+                activate(v)
+                if R_nodes[u]["excess"] == 0:
+                    # The node has become inactive.
+                    levels[height].inactive.add(u)
+                    break
+            try:
+                curr_edge.move_to_next()
+            except StopIteration:
+                # We have run off the end of the adjacency list, and there can
+                # be no more admissible edges. Relabel the node to create one.
+                height = relabel(u)
+                if is_phase1 and height >= n - 1:
+                    # Although the node is still active, with a height at least
+                    # n - 1, it is now known to be on the s side of the minimum
+                    # s-t cut. Stop processing it until phase 2.
+                    levels[height].active.add(u)
+                    break
+                # The first relabel operation after global relabeling may not
+                # increase the height of the node since the 'current edge' data
+                # structure is not rewound. Use height instead of (height - 1)
+                # in case other active nodes at the same level are missed.
+                next_height = height
+        R_nodes[u]["height"] = height
+        return next_height
+
+    def gap_heuristic(height):
+        """Apply the gap heuristic.
+        """
+        # Move all nodes at levels (height + 1) to max_height to level n + 1.
+        for level in islice(levels, height + 1, max_height + 1):
+            for u in level.active:
+                R_nodes[u]["height"] = n + 1
+            for u in level.inactive:
+                R_nodes[u]["height"] = n + 1
+            levels[n + 1].active.update(level.active)
+            level.active.clear()
+            levels[n + 1].inactive.update(level.inactive)
+            level.inactive.clear()
+
+    def global_relabel(from_sink):
+        """Apply the global relabeling heuristic.
+        """
+        src = t if from_sink else s
+        heights = reverse_bfs(src)
+        if not from_sink:
+            # s must be reachable from t. Remove t explicitly.
+            del heights[t]
+        max_height = max(heights.values())
+        if from_sink:
+            # Also mark nodes from which t is unreachable for relabeling. This
+            # serves the same purpose as the gap heuristic.
+            for u in R:
+                if u not in heights and R_nodes[u]["height"] < n:
+                    heights[u] = n + 1
+        else:
+            # Shift the computed heights because the height of s is n.
+            for u in heights:
+                heights[u] += n
+            max_height += n
+        del heights[src]
+        for u, new_height in heights.items():
+            old_height = R_nodes[u]["height"]
+            if new_height != old_height:
+                if u in levels[old_height].active:
+                    levels[old_height].active.remove(u)
+                    levels[new_height].active.add(u)
+                else:
+                    levels[old_height].inactive.remove(u)
+                    levels[new_height].inactive.add(u)
+                R_nodes[u]["height"] = new_height
+        return max_height
+
+    # Phase 1: Find the maximum preflow by pushing as much flow as possible to
+    # t.
+
+    height = max_height
+    while height > 0:
+        # Discharge active nodes in the current level.
+        while True:
+            level = levels[height]
+            if not level.active:
+                # All active nodes in the current level have been discharged.
+                # Move to the next lower level.
+                height -= 1
+                break
+            # Record the old height and level for the gap heuristic.
+            old_height = height
+            old_level = level
+            u = arbitrary_element(level.active)
+            height = discharge(u, True)
+            if grt.is_reached():
+                # Global relabeling heuristic: Recompute the exact heights of
+                # all nodes.
+                height = global_relabel(True)
+                max_height = height
+                grt.clear_work()
+            elif not old_level.active and not old_level.inactive:
+                # Gap heuristic: If the level at old_height is empty (a 'gap'),
+                # a minimum cut has been identified. All nodes with heights
+                # above old_height can have their heights set to n + 1 and not
+                # be further processed before a maximum preflow is found.
+                gap_heuristic(old_height)
+                height = old_height - 1
+                max_height = height
+            else:
+                # Update the height of the highest level with at least one
+                # active node.
+                max_height = max(max_height, height)
+
+    # A maximum preflow has been found. The excess at t is the maximum flow
+    # value.
+    if value_only:
+        R.graph["flow_value"] = R_nodes[t]["excess"]
+        return R
+
+    # Phase 2: Convert the maximum preflow into a maximum flow by returning the
+    # excess to s.
+
+    # Relabel all nodes so that they have accurate heights.
+    height = global_relabel(False)
+    grt.clear_work()
+
+    # Continue to discharge the active nodes.
+    while height > n:
+        # Discharge active nodes in the current level.
+        while True:
+            level = levels[height]
+            if not level.active:
+                # All active nodes in the current level have been discharged.
+                # Move to the next lower level.
+                height -= 1
+                break
+            u = arbitrary_element(level.active)
+            height = discharge(u, False)
+            if grt.is_reached():
+                # Global relabeling heuristic.
+                height = global_relabel(False)
+                grt.clear_work()
+
+    R.graph["flow_value"] = R_nodes[t]["excess"]
+    return R
+
+
+def preflow_push(
+    G, s, t, capacity="capacity", residual=None, global_relabel_freq=1, value_only=False
+):
+    r"""Find a maximum single-commodity flow using the highest-label
+    preflow-push algorithm.
+
+    This function returns the residual network resulting after computing
+    the maximum flow. See below for details about the conventions
+    NetworkX uses for defining residual networks.
+
+    This algorithm has a running time of $O(n^2 \sqrt{m})$ for $n$ nodes and
+    $m$ edges.
+
+
+    Parameters
+    ----------
+    G : NetworkX graph
+        Edges of the graph are expected to have an attribute called
+        'capacity'. If this attribute is not present, the edge is
+        considered to have infinite capacity.
+
+    s : node
+        Source node for the flow.
+
+    t : node
+        Sink node for the flow.
+
+    capacity : string
+        Edges of the graph G are expected to have an attribute capacity
+        that indicates how much flow the edge can support. If this
+        attribute is not present, the edge is considered to have
+        infinite capacity. Default value: 'capacity'.
+
+    residual : NetworkX graph
+        Residual network on which the algorithm is to be executed. If None, a
+        new residual network is created. Default value: None.
+
+    global_relabel_freq : integer, float
+        Relative frequency of applying the global relabeling heuristic to speed
+        up the algorithm. If it is None, the heuristic is disabled. Default
+        value: 1.
+
+    value_only : bool
+        If False, compute a maximum flow; otherwise, compute a maximum preflow
+        which is enough for computing the maximum flow value. Default value:
+        False.
+
+    Returns
+    -------
+    R : NetworkX DiGraph
+        Residual network after computing the maximum flow.
+
+    Raises
+    ------
+    NetworkXError
+        The algorithm does not support MultiGraph and MultiDiGraph. If
+        the input graph is an instance of one of these two classes, a
+        NetworkXError is raised.
+
+    NetworkXUnbounded
+        If the graph has a path of infinite capacity, the value of a
+        feasible flow on the graph is unbounded above and the function
+        raises a NetworkXUnbounded.
+
+    See also
+    --------
+    :meth:`maximum_flow`
+    :meth:`minimum_cut`
+    :meth:`edmonds_karp`
+    :meth:`shortest_augmenting_path`
+
+    Notes
+    -----
+    The residual network :samp:`R` from an input graph :samp:`G` has the
+    same nodes as :samp:`G`. :samp:`R` is a DiGraph that contains a pair
+    of edges :samp:`(u, v)` and :samp:`(v, u)` iff :samp:`(u, v)` is not a
+    self-loop, and at least one of :samp:`(u, v)` and :samp:`(v, u)` exists
+    in :samp:`G`. For each node :samp:`u` in :samp:`R`,
+    :samp:`R.nodes[u]['excess']` represents the difference between flow into
+    :samp:`u` and flow out of :samp:`u`.
+
+    For each edge :samp:`(u, v)` in :samp:`R`, :samp:`R[u][v]['capacity']`
+    is equal to the capacity of :samp:`(u, v)` in :samp:`G` if it exists
+    in :samp:`G` or zero otherwise. If the capacity is infinite,
+    :samp:`R[u][v]['capacity']` will have a high arbitrary finite value
+    that does not affect the solution of the problem. This value is stored in
+    :samp:`R.graph['inf']`. For each edge :samp:`(u, v)` in :samp:`R`,
+    :samp:`R[u][v]['flow']` represents the flow function of :samp:`(u, v)` and
+    satisfies :samp:`R[u][v]['flow'] == -R[v][u]['flow']`.
+
+    The flow value, defined as the total flow into :samp:`t`, the sink, is
+    stored in :samp:`R.graph['flow_value']`. Reachability to :samp:`t` using
+    only edges :samp:`(u, v)` such that
+    :samp:`R[u][v]['flow'] < R[u][v]['capacity']` induces a minimum
+    :samp:`s`-:samp:`t` cut.
+
+    Examples
+    --------
+    >>> from networkx.algorithms.flow import preflow_push
+
+    The functions that implement flow algorithms and output a residual
+    network, such as this one, are not imported to the base NetworkX
+    namespace, so you have to explicitly import them from the flow package.
+
+    >>> G = nx.DiGraph()
+    >>> G.add_edge("x", "a", capacity=3.0)
+    >>> G.add_edge("x", "b", capacity=1.0)
+    >>> G.add_edge("a", "c", capacity=3.0)
+    >>> G.add_edge("b", "c", capacity=5.0)
+    >>> G.add_edge("b", "d", capacity=4.0)
+    >>> G.add_edge("d", "e", capacity=2.0)
+    >>> G.add_edge("c", "y", capacity=2.0)
+    >>> G.add_edge("e", "y", capacity=3.0)
+    >>> R = preflow_push(G, "x", "y")
+    >>> flow_value = nx.maximum_flow_value(G, "x", "y")
+    >>> flow_value == R.graph["flow_value"]
+    True
+    >>> # preflow_push also stores the maximum flow value
+    >>> # in the excess attribute of the sink node t
+    >>> flow_value == R.nodes["y"]["excess"]
+    True
+    >>> # For some problems, you might only want to compute a
+    >>> # maximum preflow.
+    >>> R = preflow_push(G, "x", "y", value_only=True)
+    >>> flow_value == R.graph["flow_value"]
+    True
+    >>> flow_value == R.nodes["y"]["excess"]
+    True
+
+    """
+    R = preflow_push_impl(G, s, t, capacity, residual, global_relabel_freq, value_only)
+    R.graph["algorithm"] = "preflow_push"
+    return R