## Mercurial > repos > shellac > sam_consensus_v3

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

Find changesets by keywords (author, files, the commit message), revision
number or hash, or revset expression.

"planemo upload commit 60cee0fc7c0cda8592644e1aad72851dec82c959"

author | shellac |
---|---|

date | Mon, 22 Mar 2021 18:12:50 +0000 |

parents | |

children |

comparison

equal
deleted
inserted
replaced

-1:000000000000 | 0:4f3585e2f14b |
---|---|

1 """ | |

2 Highest-label preflow-push algorithm for maximum flow problems. | |

3 """ | |

4 | |

5 from collections import deque | |

6 from itertools import islice | |

7 import networkx as nx | |

8 from ...utils import arbitrary_element | |

9 from .utils import build_residual_network | |

10 from .utils import CurrentEdge | |

11 from .utils import detect_unboundedness | |

12 from .utils import GlobalRelabelThreshold | |

13 from .utils import Level | |

14 | |

15 __all__ = ["preflow_push"] | |

16 | |

17 | |

18 def preflow_push_impl(G, s, t, capacity, residual, global_relabel_freq, value_only): | |

19 """Implementation of the highest-label preflow-push algorithm. | |

20 """ | |

21 if s not in G: | |

22 raise nx.NetworkXError(f"node {str(s)} not in graph") | |

23 if t not in G: | |

24 raise nx.NetworkXError(f"node {str(t)} not in graph") | |

25 if s == t: | |

26 raise nx.NetworkXError("source and sink are the same node") | |

27 | |

28 if global_relabel_freq is None: | |

29 global_relabel_freq = 0 | |

30 if global_relabel_freq < 0: | |

31 raise nx.NetworkXError("global_relabel_freq must be nonnegative.") | |

32 | |

33 if residual is None: | |

34 R = build_residual_network(G, capacity) | |

35 else: | |

36 R = residual | |

37 | |

38 detect_unboundedness(R, s, t) | |

39 | |

40 R_nodes = R.nodes | |

41 R_pred = R.pred | |

42 R_succ = R.succ | |

43 | |

44 # Initialize/reset the residual network. | |

45 for u in R: | |

46 R_nodes[u]["excess"] = 0 | |

47 for e in R_succ[u].values(): | |

48 e["flow"] = 0 | |

49 | |

50 def reverse_bfs(src): | |

51 """Perform a reverse breadth-first search from src in the residual | |

52 network. | |

53 """ | |

54 heights = {src: 0} | |

55 q = deque([(src, 0)]) | |

56 while q: | |

57 u, height = q.popleft() | |

58 height += 1 | |

59 for v, attr in R_pred[u].items(): | |

60 if v not in heights and attr["flow"] < attr["capacity"]: | |

61 heights[v] = height | |

62 q.append((v, height)) | |

63 return heights | |

64 | |

65 # Initialize heights of the nodes. | |

66 heights = reverse_bfs(t) | |

67 | |

68 if s not in heights: | |

69 # t is not reachable from s in the residual network. The maximum flow | |

70 # must be zero. | |

71 R.graph["flow_value"] = 0 | |

72 return R | |

73 | |

74 n = len(R) | |

75 # max_height represents the height of the highest level below level n with | |

76 # at least one active node. | |

77 max_height = max(heights[u] for u in heights if u != s) | |

78 heights[s] = n | |

79 | |

80 grt = GlobalRelabelThreshold(n, R.size(), global_relabel_freq) | |

81 | |

82 # Initialize heights and 'current edge' data structures of the nodes. | |

83 for u in R: | |

84 R_nodes[u]["height"] = heights[u] if u in heights else n + 1 | |

85 R_nodes[u]["curr_edge"] = CurrentEdge(R_succ[u]) | |

86 | |

87 def push(u, v, flow): | |

88 """Push flow units of flow from u to v. | |

89 """ | |

90 R_succ[u][v]["flow"] += flow | |

91 R_succ[v][u]["flow"] -= flow | |

92 R_nodes[u]["excess"] -= flow | |

93 R_nodes[v]["excess"] += flow | |

94 | |

95 # The maximum flow must be nonzero now. Initialize the preflow by | |

96 # saturating all edges emanating from s. | |

97 for u, attr in R_succ[s].items(): | |

98 flow = attr["capacity"] | |

99 if flow > 0: | |

100 push(s, u, flow) | |

101 | |

102 # Partition nodes into levels. | |

103 levels = [Level() for i in range(2 * n)] | |

104 for u in R: | |

105 if u != s and u != t: | |

106 level = levels[R_nodes[u]["height"]] | |

107 if R_nodes[u]["excess"] > 0: | |

108 level.active.add(u) | |

109 else: | |

110 level.inactive.add(u) | |

111 | |

112 def activate(v): | |

113 """Move a node from the inactive set to the active set of its level. | |

114 """ | |

115 if v != s and v != t: | |

116 level = levels[R_nodes[v]["height"]] | |

117 if v in level.inactive: | |

118 level.inactive.remove(v) | |

119 level.active.add(v) | |

120 | |

121 def relabel(u): | |

122 """Relabel a node to create an admissible edge. | |

123 """ | |

124 grt.add_work(len(R_succ[u])) | |

125 return ( | |

126 min( | |

127 R_nodes[v]["height"] | |

128 for v, attr in R_succ[u].items() | |

129 if attr["flow"] < attr["capacity"] | |

130 ) | |

131 + 1 | |

132 ) | |

133 | |

134 def discharge(u, is_phase1): | |

135 """Discharge a node until it becomes inactive or, during phase 1 (see | |

136 below), its height reaches at least n. The node is known to have the | |

137 largest height among active nodes. | |

138 """ | |

139 height = R_nodes[u]["height"] | |

140 curr_edge = R_nodes[u]["curr_edge"] | |

141 # next_height represents the next height to examine after discharging | |

142 # the current node. During phase 1, it is capped to below n. | |

143 next_height = height | |

144 levels[height].active.remove(u) | |

145 while True: | |

146 v, attr = curr_edge.get() | |

147 if height == R_nodes[v]["height"] + 1 and attr["flow"] < attr["capacity"]: | |

148 flow = min(R_nodes[u]["excess"], attr["capacity"] - attr["flow"]) | |

149 push(u, v, flow) | |

150 activate(v) | |

151 if R_nodes[u]["excess"] == 0: | |

152 # The node has become inactive. | |

153 levels[height].inactive.add(u) | |

154 break | |

155 try: | |

156 curr_edge.move_to_next() | |

157 except StopIteration: | |

158 # We have run off the end of the adjacency list, and there can | |

159 # be no more admissible edges. Relabel the node to create one. | |

160 height = relabel(u) | |

161 if is_phase1 and height >= n - 1: | |

162 # Although the node is still active, with a height at least | |

163 # n - 1, it is now known to be on the s side of the minimum | |

164 # s-t cut. Stop processing it until phase 2. | |

165 levels[height].active.add(u) | |

166 break | |

167 # The first relabel operation after global relabeling may not | |

168 # increase the height of the node since the 'current edge' data | |

169 # structure is not rewound. Use height instead of (height - 1) | |

170 # in case other active nodes at the same level are missed. | |

171 next_height = height | |

172 R_nodes[u]["height"] = height | |

173 return next_height | |

174 | |

175 def gap_heuristic(height): | |

176 """Apply the gap heuristic. | |

177 """ | |

178 # Move all nodes at levels (height + 1) to max_height to level n + 1. | |

179 for level in islice(levels, height + 1, max_height + 1): | |

180 for u in level.active: | |

181 R_nodes[u]["height"] = n + 1 | |

182 for u in level.inactive: | |

183 R_nodes[u]["height"] = n + 1 | |

184 levels[n + 1].active.update(level.active) | |

185 level.active.clear() | |

186 levels[n + 1].inactive.update(level.inactive) | |

187 level.inactive.clear() | |

188 | |

189 def global_relabel(from_sink): | |

190 """Apply the global relabeling heuristic. | |

191 """ | |

192 src = t if from_sink else s | |

193 heights = reverse_bfs(src) | |

194 if not from_sink: | |

195 # s must be reachable from t. Remove t explicitly. | |

196 del heights[t] | |

197 max_height = max(heights.values()) | |

198 if from_sink: | |

199 # Also mark nodes from which t is unreachable for relabeling. This | |

200 # serves the same purpose as the gap heuristic. | |

201 for u in R: | |

202 if u not in heights and R_nodes[u]["height"] < n: | |

203 heights[u] = n + 1 | |

204 else: | |

205 # Shift the computed heights because the height of s is n. | |

206 for u in heights: | |

207 heights[u] += n | |

208 max_height += n | |

209 del heights[src] | |

210 for u, new_height in heights.items(): | |

211 old_height = R_nodes[u]["height"] | |

212 if new_height != old_height: | |

213 if u in levels[old_height].active: | |

214 levels[old_height].active.remove(u) | |

215 levels[new_height].active.add(u) | |

216 else: | |

217 levels[old_height].inactive.remove(u) | |

218 levels[new_height].inactive.add(u) | |

219 R_nodes[u]["height"] = new_height | |

220 return max_height | |

221 | |

222 # Phase 1: Find the maximum preflow by pushing as much flow as possible to | |

223 # t. | |

224 | |

225 height = max_height | |

226 while height > 0: | |

227 # Discharge active nodes in the current level. | |

228 while True: | |

229 level = levels[height] | |

230 if not level.active: | |

231 # All active nodes in the current level have been discharged. | |

232 # Move to the next lower level. | |

233 height -= 1 | |

234 break | |

235 # Record the old height and level for the gap heuristic. | |

236 old_height = height | |

237 old_level = level | |

238 u = arbitrary_element(level.active) | |

239 height = discharge(u, True) | |

240 if grt.is_reached(): | |

241 # Global relabeling heuristic: Recompute the exact heights of | |

242 # all nodes. | |

243 height = global_relabel(True) | |

244 max_height = height | |

245 grt.clear_work() | |

246 elif not old_level.active and not old_level.inactive: | |

247 # Gap heuristic: If the level at old_height is empty (a 'gap'), | |

248 # a minimum cut has been identified. All nodes with heights | |

249 # above old_height can have their heights set to n + 1 and not | |

250 # be further processed before a maximum preflow is found. | |

251 gap_heuristic(old_height) | |

252 height = old_height - 1 | |

253 max_height = height | |

254 else: | |

255 # Update the height of the highest level with at least one | |

256 # active node. | |

257 max_height = max(max_height, height) | |

258 | |

259 # A maximum preflow has been found. The excess at t is the maximum flow | |

260 # value. | |

261 if value_only: | |

262 R.graph["flow_value"] = R_nodes[t]["excess"] | |

263 return R | |

264 | |

265 # Phase 2: Convert the maximum preflow into a maximum flow by returning the | |

266 # excess to s. | |

267 | |

268 # Relabel all nodes so that they have accurate heights. | |

269 height = global_relabel(False) | |

270 grt.clear_work() | |

271 | |

272 # Continue to discharge the active nodes. | |

273 while height > n: | |

274 # Discharge active nodes in the current level. | |

275 while True: | |

276 level = levels[height] | |

277 if not level.active: | |

278 # All active nodes in the current level have been discharged. | |

279 # Move to the next lower level. | |

280 height -= 1 | |

281 break | |

282 u = arbitrary_element(level.active) | |

283 height = discharge(u, False) | |

284 if grt.is_reached(): | |

285 # Global relabeling heuristic. | |

286 height = global_relabel(False) | |

287 grt.clear_work() | |

288 | |

289 R.graph["flow_value"] = R_nodes[t]["excess"] | |

290 return R | |

291 | |

292 | |

293 def preflow_push( | |

294 G, s, t, capacity="capacity", residual=None, global_relabel_freq=1, value_only=False | |

295 ): | |

296 r"""Find a maximum single-commodity flow using the highest-label | |

297 preflow-push algorithm. | |

298 | |

299 This function returns the residual network resulting after computing | |

300 the maximum flow. See below for details about the conventions | |

301 NetworkX uses for defining residual networks. | |

302 | |

303 This algorithm has a running time of $O(n^2 \sqrt{m})$ for $n$ nodes and | |

304 $m$ edges. | |

305 | |

306 | |

307 Parameters | |

308 ---------- | |

309 G : NetworkX graph | |

310 Edges of the graph are expected to have an attribute called | |

311 'capacity'. If this attribute is not present, the edge is | |

312 considered to have infinite capacity. | |

313 | |

314 s : node | |

315 Source node for the flow. | |

316 | |

317 t : node | |

318 Sink node for the flow. | |

319 | |

320 capacity : string | |

321 Edges of the graph G are expected to have an attribute capacity | |

322 that indicates how much flow the edge can support. If this | |

323 attribute is not present, the edge is considered to have | |

324 infinite capacity. Default value: 'capacity'. | |

325 | |

326 residual : NetworkX graph | |

327 Residual network on which the algorithm is to be executed. If None, a | |

328 new residual network is created. Default value: None. | |

329 | |

330 global_relabel_freq : integer, float | |

331 Relative frequency of applying the global relabeling heuristic to speed | |

332 up the algorithm. If it is None, the heuristic is disabled. Default | |

333 value: 1. | |

334 | |

335 value_only : bool | |

336 If False, compute a maximum flow; otherwise, compute a maximum preflow | |

337 which is enough for computing the maximum flow value. Default value: | |

338 False. | |

339 | |

340 Returns | |

341 ------- | |

342 R : NetworkX DiGraph | |

343 Residual network after computing the maximum flow. | |

344 | |

345 Raises | |

346 ------ | |

347 NetworkXError | |

348 The algorithm does not support MultiGraph and MultiDiGraph. If | |

349 the input graph is an instance of one of these two classes, a | |

350 NetworkXError is raised. | |

351 | |

352 NetworkXUnbounded | |

353 If the graph has a path of infinite capacity, the value of a | |

354 feasible flow on the graph is unbounded above and the function | |

355 raises a NetworkXUnbounded. | |

356 | |

357 See also | |

358 -------- | |

359 :meth:`maximum_flow` | |

360 :meth:`minimum_cut` | |

361 :meth:`edmonds_karp` | |

362 :meth:`shortest_augmenting_path` | |

363 | |

364 Notes | |

365 ----- | |

366 The residual network :samp:`R` from an input graph :samp:`G` has the | |

367 same nodes as :samp:`G`. :samp:`R` is a DiGraph that contains a pair | |

368 of edges :samp:`(u, v)` and :samp:`(v, u)` iff :samp:`(u, v)` is not a | |

369 self-loop, and at least one of :samp:`(u, v)` and :samp:`(v, u)` exists | |

370 in :samp:`G`. For each node :samp:`u` in :samp:`R`, | |

371 :samp:`R.nodes[u]['excess']` represents the difference between flow into | |

372 :samp:`u` and flow out of :samp:`u`. | |

373 | |

374 For each edge :samp:`(u, v)` in :samp:`R`, :samp:`R[u][v]['capacity']` | |

375 is equal to the capacity of :samp:`(u, v)` in :samp:`G` if it exists | |

376 in :samp:`G` or zero otherwise. If the capacity is infinite, | |

377 :samp:`R[u][v]['capacity']` will have a high arbitrary finite value | |

378 that does not affect the solution of the problem. This value is stored in | |

379 :samp:`R.graph['inf']`. For each edge :samp:`(u, v)` in :samp:`R`, | |

380 :samp:`R[u][v]['flow']` represents the flow function of :samp:`(u, v)` and | |

381 satisfies :samp:`R[u][v]['flow'] == -R[v][u]['flow']`. | |

382 | |

383 The flow value, defined as the total flow into :samp:`t`, the sink, is | |

384 stored in :samp:`R.graph['flow_value']`. Reachability to :samp:`t` using | |

385 only edges :samp:`(u, v)` such that | |

386 :samp:`R[u][v]['flow'] < R[u][v]['capacity']` induces a minimum | |

387 :samp:`s`-:samp:`t` cut. | |

388 | |

389 Examples | |

390 -------- | |

391 >>> from networkx.algorithms.flow import preflow_push | |

392 | |

393 The functions that implement flow algorithms and output a residual | |

394 network, such as this one, are not imported to the base NetworkX | |

395 namespace, so you have to explicitly import them from the flow package. | |

396 | |

397 >>> G = nx.DiGraph() | |

398 >>> G.add_edge("x", "a", capacity=3.0) | |

399 >>> G.add_edge("x", "b", capacity=1.0) | |

400 >>> G.add_edge("a", "c", capacity=3.0) | |

401 >>> G.add_edge("b", "c", capacity=5.0) | |

402 >>> G.add_edge("b", "d", capacity=4.0) | |

403 >>> G.add_edge("d", "e", capacity=2.0) | |

404 >>> G.add_edge("c", "y", capacity=2.0) | |

405 >>> G.add_edge("e", "y", capacity=3.0) | |

406 >>> R = preflow_push(G, "x", "y") | |

407 >>> flow_value = nx.maximum_flow_value(G, "x", "y") | |

408 >>> flow_value == R.graph["flow_value"] | |

409 True | |

410 >>> # preflow_push also stores the maximum flow value | |

411 >>> # in the excess attribute of the sink node t | |

412 >>> flow_value == R.nodes["y"]["excess"] | |

413 True | |

414 >>> # For some problems, you might only want to compute a | |

415 >>> # maximum preflow. | |

416 >>> R = preflow_push(G, "x", "y", value_only=True) | |

417 >>> flow_value == R.graph["flow_value"] | |

418 True | |

419 >>> flow_value == R.nodes["y"]["excess"] | |

420 True | |

421 | |

422 """ | |

423 R = preflow_push_impl(G, s, t, capacity, residual, global_relabel_freq, value_only) | |

424 R.graph["algorithm"] = "preflow_push" | |

425 return R |