Coverage for /pythoncovmergedfiles/medio/medio/usr/local/lib/python3.9/dist-packages/networkx/algorithms/flow/capacityscaling.py: 7%
151 statements
« prev ^ index » next coverage.py v7.3.2, created at 2023-10-20 07:00 +0000
1"""
2Capacity scaling minimum cost flow algorithm.
3"""
5__all__ = ["capacity_scaling"]
7from itertools import chain
8from math import log
10import networkx as nx
12from ...utils import BinaryHeap, arbitrary_element, not_implemented_for
15def _detect_unboundedness(R):
16 """Detect infinite-capacity negative cycles."""
17 G = nx.DiGraph()
18 G.add_nodes_from(R)
20 # Value simulating infinity.
21 inf = R.graph["inf"]
22 # True infinity.
23 f_inf = float("inf")
24 for u in R:
25 for v, e in R[u].items():
26 # Compute the minimum weight of infinite-capacity (u, v) edges.
27 w = f_inf
28 for k, e in e.items():
29 if e["capacity"] == inf:
30 w = min(w, e["weight"])
31 if w != f_inf:
32 G.add_edge(u, v, weight=w)
34 if nx.negative_edge_cycle(G):
35 raise nx.NetworkXUnbounded(
36 "Negative cost cycle of infinite capacity found. "
37 "Min cost flow may be unbounded below."
38 )
41@not_implemented_for("undirected")
42def _build_residual_network(G, demand, capacity, weight):
43 """Build a residual network and initialize a zero flow."""
44 if sum(G.nodes[u].get(demand, 0) for u in G) != 0:
45 raise nx.NetworkXUnfeasible("Sum of the demands should be 0.")
47 R = nx.MultiDiGraph()
48 R.add_nodes_from(
49 (u, {"excess": -G.nodes[u].get(demand, 0), "potential": 0}) for u in G
50 )
52 inf = float("inf")
53 # Detect selfloops with infinite capacities and negative weights.
54 for u, v, e in nx.selfloop_edges(G, data=True):
55 if e.get(weight, 0) < 0 and e.get(capacity, inf) == inf:
56 raise nx.NetworkXUnbounded(
57 "Negative cost cycle of infinite capacity found. "
58 "Min cost flow may be unbounded below."
59 )
61 # Extract edges with positive capacities. Self loops excluded.
62 if G.is_multigraph():
63 edge_list = [
64 (u, v, k, e)
65 for u, v, k, e in G.edges(data=True, keys=True)
66 if u != v and e.get(capacity, inf) > 0
67 ]
68 else:
69 edge_list = [
70 (u, v, 0, e)
71 for u, v, e in G.edges(data=True)
72 if u != v and e.get(capacity, inf) > 0
73 ]
74 # Simulate infinity with the larger of the sum of absolute node imbalances
75 # the sum of finite edge capacities or any positive value if both sums are
76 # zero. This allows the infinite-capacity edges to be distinguished for
77 # unboundedness detection and directly participate in residual capacity
78 # calculation.
79 inf = (
80 max(
81 sum(abs(R.nodes[u]["excess"]) for u in R),
82 2
83 * sum(
84 e[capacity]
85 for u, v, k, e in edge_list
86 if capacity in e and e[capacity] != inf
87 ),
88 )
89 or 1
90 )
91 for u, v, k, e in edge_list:
92 r = min(e.get(capacity, inf), inf)
93 w = e.get(weight, 0)
94 # Add both (u, v) and (v, u) into the residual network marked with the
95 # original key. (key[1] == True) indicates the (u, v) is in the
96 # original network.
97 R.add_edge(u, v, key=(k, True), capacity=r, weight=w, flow=0)
98 R.add_edge(v, u, key=(k, False), capacity=0, weight=-w, flow=0)
100 # Record the value simulating infinity.
101 R.graph["inf"] = inf
103 _detect_unboundedness(R)
105 return R
108def _build_flow_dict(G, R, capacity, weight):
109 """Build a flow dictionary from a residual network."""
110 inf = float("inf")
111 flow_dict = {}
112 if G.is_multigraph():
113 for u in G:
114 flow_dict[u] = {}
115 for v, es in G[u].items():
116 flow_dict[u][v] = {
117 # Always saturate negative selfloops.
118 k: (
119 0
120 if (
121 u != v or e.get(capacity, inf) <= 0 or e.get(weight, 0) >= 0
122 )
123 else e[capacity]
124 )
125 for k, e in es.items()
126 }
127 for v, es in R[u].items():
128 if v in flow_dict[u]:
129 flow_dict[u][v].update(
130 (k[0], e["flow"]) for k, e in es.items() if e["flow"] > 0
131 )
132 else:
133 for u in G:
134 flow_dict[u] = {
135 # Always saturate negative selfloops.
136 v: (
137 0
138 if (u != v or e.get(capacity, inf) <= 0 or e.get(weight, 0) >= 0)
139 else e[capacity]
140 )
141 for v, e in G[u].items()
142 }
143 flow_dict[u].update(
144 (v, e["flow"])
145 for v, es in R[u].items()
146 for e in es.values()
147 if e["flow"] > 0
148 )
149 return flow_dict
152@nx._dispatch(node_attrs="demand", edge_attrs={"capacity": float("inf"), "weight": 0})
153def capacity_scaling(
154 G, demand="demand", capacity="capacity", weight="weight", heap=BinaryHeap
155):
156 r"""Find a minimum cost flow satisfying all demands in digraph G.
158 This is a capacity scaling successive shortest augmenting path algorithm.
160 G is a digraph with edge costs and capacities and in which nodes
161 have demand, i.e., they want to send or receive some amount of
162 flow. A negative demand means that the node wants to send flow, a
163 positive demand means that the node want to receive flow. A flow on
164 the digraph G satisfies all demand if the net flow into each node
165 is equal to the demand of that node.
167 Parameters
168 ----------
169 G : NetworkX graph
170 DiGraph or MultiDiGraph on which a minimum cost flow satisfying all
171 demands is to be found.
173 demand : string
174 Nodes of the graph G are expected to have an attribute demand
175 that indicates how much flow a node wants to send (negative
176 demand) or receive (positive demand). Note that the sum of the
177 demands should be 0 otherwise the problem in not feasible. If
178 this attribute is not present, a node is considered to have 0
179 demand. Default value: 'demand'.
181 capacity : string
182 Edges of the graph G are expected to have an attribute capacity
183 that indicates how much flow the edge can support. If this
184 attribute is not present, the edge is considered to have
185 infinite capacity. Default value: 'capacity'.
187 weight : string
188 Edges of the graph G are expected to have an attribute weight
189 that indicates the cost incurred by sending one unit of flow on
190 that edge. If not present, the weight is considered to be 0.
191 Default value: 'weight'.
193 heap : class
194 Type of heap to be used in the algorithm. It should be a subclass of
195 :class:`MinHeap` or implement a compatible interface.
197 If a stock heap implementation is to be used, :class:`BinaryHeap` is
198 recommended over :class:`PairingHeap` for Python implementations without
199 optimized attribute accesses (e.g., CPython) despite a slower
200 asymptotic running time. For Python implementations with optimized
201 attribute accesses (e.g., PyPy), :class:`PairingHeap` provides better
202 performance. Default value: :class:`BinaryHeap`.
204 Returns
205 -------
206 flowCost : integer
207 Cost of a minimum cost flow satisfying all demands.
209 flowDict : dictionary
210 If G is a digraph, a dict-of-dicts keyed by nodes such that
211 flowDict[u][v] is the flow on edge (u, v).
212 If G is a MultiDiGraph, a dict-of-dicts-of-dicts keyed by nodes
213 so that flowDict[u][v][key] is the flow on edge (u, v, key).
215 Raises
216 ------
217 NetworkXError
218 This exception is raised if the input graph is not directed,
219 not connected.
221 NetworkXUnfeasible
222 This exception is raised in the following situations:
224 * The sum of the demands is not zero. Then, there is no
225 flow satisfying all demands.
226 * There is no flow satisfying all demand.
228 NetworkXUnbounded
229 This exception is raised if the digraph G has a cycle of
230 negative cost and infinite capacity. Then, the cost of a flow
231 satisfying all demands is unbounded below.
233 Notes
234 -----
235 This algorithm does not work if edge weights are floating-point numbers.
237 See also
238 --------
239 :meth:`network_simplex`
241 Examples
242 --------
243 A simple example of a min cost flow problem.
245 >>> G = nx.DiGraph()
246 >>> G.add_node("a", demand=-5)
247 >>> G.add_node("d", demand=5)
248 >>> G.add_edge("a", "b", weight=3, capacity=4)
249 >>> G.add_edge("a", "c", weight=6, capacity=10)
250 >>> G.add_edge("b", "d", weight=1, capacity=9)
251 >>> G.add_edge("c", "d", weight=2, capacity=5)
252 >>> flowCost, flowDict = nx.capacity_scaling(G)
253 >>> flowCost
254 24
255 >>> flowDict
256 {'a': {'b': 4, 'c': 1}, 'd': {}, 'b': {'d': 4}, 'c': {'d': 1}}
258 It is possible to change the name of the attributes used for the
259 algorithm.
261 >>> G = nx.DiGraph()
262 >>> G.add_node("p", spam=-4)
263 >>> G.add_node("q", spam=2)
264 >>> G.add_node("a", spam=-2)
265 >>> G.add_node("d", spam=-1)
266 >>> G.add_node("t", spam=2)
267 >>> G.add_node("w", spam=3)
268 >>> G.add_edge("p", "q", cost=7, vacancies=5)
269 >>> G.add_edge("p", "a", cost=1, vacancies=4)
270 >>> G.add_edge("q", "d", cost=2, vacancies=3)
271 >>> G.add_edge("t", "q", cost=1, vacancies=2)
272 >>> G.add_edge("a", "t", cost=2, vacancies=4)
273 >>> G.add_edge("d", "w", cost=3, vacancies=4)
274 >>> G.add_edge("t", "w", cost=4, vacancies=1)
275 >>> flowCost, flowDict = nx.capacity_scaling(
276 ... G, demand="spam", capacity="vacancies", weight="cost"
277 ... )
278 >>> flowCost
279 37
280 >>> flowDict
281 {'p': {'q': 2, 'a': 2}, 'q': {'d': 1}, 'a': {'t': 4}, 'd': {'w': 2}, 't': {'q': 1, 'w': 1}, 'w': {}}
282 """
283 R = _build_residual_network(G, demand, capacity, weight)
285 inf = float("inf")
286 # Account cost of negative selfloops.
287 flow_cost = sum(
288 0
289 if e.get(capacity, inf) <= 0 or e.get(weight, 0) >= 0
290 else e[capacity] * e[weight]
291 for u, v, e in nx.selfloop_edges(G, data=True)
292 )
294 # Determine the maximum edge capacity.
295 wmax = max(chain([-inf], (e["capacity"] for u, v, e in R.edges(data=True))))
296 if wmax == -inf:
297 # Residual network has no edges.
298 return flow_cost, _build_flow_dict(G, R, capacity, weight)
300 R_nodes = R.nodes
301 R_succ = R.succ
303 delta = 2 ** int(log(wmax, 2))
304 while delta >= 1:
305 # Saturate Δ-residual edges with negative reduced costs to achieve
306 # Δ-optimality.
307 for u in R:
308 p_u = R_nodes[u]["potential"]
309 for v, es in R_succ[u].items():
310 for k, e in es.items():
311 flow = e["capacity"] - e["flow"]
312 if e["weight"] - p_u + R_nodes[v]["potential"] < 0:
313 flow = e["capacity"] - e["flow"]
314 if flow >= delta:
315 e["flow"] += flow
316 R_succ[v][u][(k[0], not k[1])]["flow"] -= flow
317 R_nodes[u]["excess"] -= flow
318 R_nodes[v]["excess"] += flow
319 # Determine the Δ-active nodes.
320 S = set()
321 T = set()
322 S_add = S.add
323 S_remove = S.remove
324 T_add = T.add
325 T_remove = T.remove
326 for u in R:
327 excess = R_nodes[u]["excess"]
328 if excess >= delta:
329 S_add(u)
330 elif excess <= -delta:
331 T_add(u)
332 # Repeatedly augment flow from S to T along shortest paths until
333 # Δ-feasibility is achieved.
334 while S and T:
335 s = arbitrary_element(S)
336 t = None
337 # Search for a shortest path in terms of reduce costs from s to
338 # any t in T in the Δ-residual network.
339 d = {}
340 pred = {s: None}
341 h = heap()
342 h_insert = h.insert
343 h_get = h.get
344 h_insert(s, 0)
345 while h:
346 u, d_u = h.pop()
347 d[u] = d_u
348 if u in T:
349 # Path found.
350 t = u
351 break
352 p_u = R_nodes[u]["potential"]
353 for v, es in R_succ[u].items():
354 if v in d:
355 continue
356 wmin = inf
357 # Find the minimum-weighted (u, v) Δ-residual edge.
358 for k, e in es.items():
359 if e["capacity"] - e["flow"] >= delta:
360 w = e["weight"]
361 if w < wmin:
362 wmin = w
363 kmin = k
364 emin = e
365 if wmin == inf:
366 continue
367 # Update the distance label of v.
368 d_v = d_u + wmin - p_u + R_nodes[v]["potential"]
369 if h_insert(v, d_v):
370 pred[v] = (u, kmin, emin)
371 if t is not None:
372 # Augment Δ units of flow from s to t.
373 while u != s:
374 v = u
375 u, k, e = pred[v]
376 e["flow"] += delta
377 R_succ[v][u][(k[0], not k[1])]["flow"] -= delta
378 # Account node excess and deficit.
379 R_nodes[s]["excess"] -= delta
380 R_nodes[t]["excess"] += delta
381 if R_nodes[s]["excess"] < delta:
382 S_remove(s)
383 if R_nodes[t]["excess"] > -delta:
384 T_remove(t)
385 # Update node potentials.
386 d_t = d[t]
387 for u, d_u in d.items():
388 R_nodes[u]["potential"] -= d_u - d_t
389 else:
390 # Path not found.
391 S_remove(s)
392 delta //= 2
394 if any(R.nodes[u]["excess"] != 0 for u in R):
395 raise nx.NetworkXUnfeasible("No flow satisfying all demands.")
397 # Calculate the flow cost.
398 for u in R:
399 for v, es in R_succ[u].items():
400 for e in es.values():
401 flow = e["flow"]
402 if flow > 0:
403 flow_cost += flow * e["weight"]
405 return flow_cost, _build_flow_dict(G, R, capacity, weight)