graph/mcmf.hpp
| Line | Branch | Exec | Source |
|---|---|---|---|
| 1 | #pragma once | ||
| 2 | |||
| 3 | #include <bits/stdc++.h> | ||
| 4 | // #include<bits/extc++.h> | ||
| 5 | #include <ext/pb_ds/priority_queue.hpp> | ||
| 6 | |||
| 7 | namespace wala { | ||
| 8 | |||
| 9 | // NOTE: This doesn't support negative-cost edges; you can adjust edge weights | ||
| 10 | // (e.g. by precomputing a potential function) to make them positive. | ||
| 11 | |||
| 12 | template <typename flow_t = int, typename cost_t = int64_t> | ||
| 13 | struct MCMF_SSPA { | ||
| 14 | int N; | ||
| 15 | std::vector<std::vector<int>> adj; | ||
| 16 | struct edge_t { | ||
| 17 | int dest; | ||
| 18 | flow_t cap; | ||
| 19 | cost_t cost; | ||
| 20 | }; | ||
| 21 | std::vector<edge_t> edges; | ||
| 22 | |||
| 23 | std::vector<char> seen; | ||
| 24 | std::vector<cost_t> pi; | ||
| 25 | std::vector<int> prv; | ||
| 26 | |||
| 27 | ✗ | explicit MCMF_SSPA(int N_) : N(N_), adj(N), pi(N, 0), prv(N) {} | |
| 28 | |||
| 29 | ✗ | void add_edge(int from, int to, flow_t cap, cost_t cost) { | |
| 30 | ✗ | assert(cap >= 0); | |
| 31 | ✗ | assert(cost + pi[from] - pi[to] >= 0); // TODO: Remove this restriction | |
| 32 | ✗ | int e = int(edges.size()); | |
| 33 | ✗ | edges.emplace_back(edge_t{to, cap, cost}); | |
| 34 | ✗ | edges.emplace_back(edge_t{from, 0, -cost}); | |
| 35 | ✗ | adj[from].push_back(e); | |
| 36 | ✗ | adj[to].push_back(e+1); | |
| 37 | } | ||
| 38 | |||
| 39 | static constexpr cost_t INF_COST = std::numeric_limits<cost_t>::max() / 4; | ||
| 40 | static constexpr flow_t INF_FLOW = std::numeric_limits<flow_t>::max() / 4; | ||
| 41 | std::vector<cost_t> dist; | ||
| 42 | __gnu_pbds::priority_queue<std::pair<cost_t, int>> q; | ||
| 43 | std::vector<typename decltype(q)::point_iterator> its; | ||
| 44 | ✗ | cost_t dijkstra(int s, int t) { | |
| 45 | ✗ | dist.assign(N, INF_COST); | |
| 46 | ✗ | dist[s] = 0; | |
| 47 | |||
| 48 | ✗ | its.assign(N, q.end()); | |
| 49 | ✗ | its[s] = q.push({-(dist[s] - pi[s]), s}); | |
| 50 | |||
| 51 | ✗ | while (!q.empty()) { | |
| 52 | ✗ | int i = q.top().second; q.pop(); | |
| 53 | ✗ | cost_t d = dist[i]; | |
| 54 | ✗ | for (int e : adj[i]) { | |
| 55 | ✗ | if (edges[e].cap) { | |
| 56 | ✗ | int j = edges[e].dest; | |
| 57 | ✗ | cost_t nd = d + edges[e].cost; | |
| 58 | ✗ | if (nd < dist[j]) { | |
| 59 | ✗ | dist[j] = nd; | |
| 60 | ✗ | prv[j] = e; | |
| 61 | ✗ | if (its[j] == q.end()) { | |
| 62 | ✗ | its[j] = q.push({-(dist[j] - pi[j]), j}); | |
| 63 | } else { | ||
| 64 | ✗ | q.modify(its[j], {-(dist[j] - pi[j]), j}); | |
| 65 | } | ||
| 66 | } | ||
| 67 | } | ||
| 68 | } | ||
| 69 | } | ||
| 70 | |||
| 71 | ✗ | swap(pi, dist); | |
| 72 | ✗ | return pi[t]; | |
| 73 | } | ||
| 74 | |||
| 75 | ✗ | flow_t path(int s, int t) { | |
| 76 | ✗ | flow_t cur_flow = std::numeric_limits<flow_t>::max(); | |
| 77 | ✗ | for (int cur = t; cur != s; ) { | |
| 78 | ✗ | int e = prv[cur]; | |
| 79 | ✗ | int nxt = edges[e^1].dest; | |
| 80 | ✗ | cur_flow = std::min(cur_flow, edges[e].cap); | |
| 81 | ✗ | cur = nxt; | |
| 82 | } | ||
| 83 | ✗ | for (int cur = t; cur != s; ) { | |
| 84 | ✗ | int e = prv[cur]; | |
| 85 | ✗ | int nxt = edges[e^1].dest; | |
| 86 | ✗ | edges[e].cap -= cur_flow; | |
| 87 | ✗ | edges[e^1].cap += cur_flow; | |
| 88 | ✗ | cur = nxt; | |
| 89 | } | ||
| 90 | ✗ | return cur_flow; | |
| 91 | } | ||
| 92 | |||
| 93 | ✗ | std::vector<std::pair<flow_t, cost_t>> all_flows(int s, int t, cost_t max_cost = INF_COST - 1) { | |
| 94 | ✗ | assert(s != t); | |
| 95 | ✗ | std::vector<std::pair<flow_t, cost_t>> res; | |
| 96 | ✗ | while (dijkstra(s, t) <= max_cost) { | |
| 97 | ✗ | assert(res.empty() || pi[t] >= res.back().second); | |
| 98 | ✗ | flow_t f = path(s, t); | |
| 99 | ✗ | res.push_back({f, pi[t]}); | |
| 100 | } | ||
| 101 | ✗ | return res; | |
| 102 | } | ||
| 103 | |||
| 104 | ✗ | std::pair<flow_t, cost_t> max_flow(int s, int t, cost_t max_cost = INF_COST - 1) { | |
| 105 | ✗ | assert(s != t); | |
| 106 | ✗ | flow_t tot_flow = 0; cost_t tot_cost = 0; | |
| 107 | ✗ | while (dijkstra(s, t) <= max_cost) { | |
| 108 | ✗ | flow_t cur_flow = path(s, t); | |
| 109 | ✗ | tot_flow += cur_flow; | |
| 110 | ✗ | tot_cost += cur_flow * pi[t]; | |
| 111 | } | ||
| 112 | ✗ | return {tot_flow, tot_cost}; | |
| 113 | } | ||
| 114 | }; | ||
| 115 | |||
| 116 | template <typename flow_t = int, typename cost_t = int64_t> | ||
| 117 | struct MCMF_Dinic { | ||
| 118 | int N; | ||
| 119 | std::vector<std::vector<int>> adj; | ||
| 120 | struct edge_t { | ||
| 121 | int dest; | ||
| 122 | flow_t cap; | ||
| 123 | cost_t cost; | ||
| 124 | }; | ||
| 125 | std::vector<edge_t> edges; | ||
| 126 | |||
| 127 | std::vector<char> seen; | ||
| 128 | std::vector<cost_t> pi; | ||
| 129 | |||
| 130 | ✗ | explicit MCMF_Dinic(int N_) : N(N_), adj(N), pi(N, 0) {} | |
| 131 | |||
| 132 | ✗ | void add_edge(int from, int to, flow_t cap, cost_t cost) { | |
| 133 | ✗ | assert(cap >= 0); | |
| 134 | ✗ | assert(cost + pi[from] - pi[to] >= 0); // TODO: Remove this restriction | |
| 135 | ✗ | int e = int(edges.size()); | |
| 136 | ✗ | edges.emplace_back(edge_t{to, cap, cost}); | |
| 137 | ✗ | edges.emplace_back(edge_t{from, 0, -cost}); | |
| 138 | ✗ | adj[from].push_back(e); | |
| 139 | ✗ | adj[to].push_back(e+1); | |
| 140 | } | ||
| 141 | |||
| 142 | static constexpr cost_t INF_COST = std::numeric_limits<cost_t>::max() / 4; | ||
| 143 | static constexpr flow_t INF_FLOW = std::numeric_limits<flow_t>::max() / 4; | ||
| 144 | std::vector<cost_t> dist; | ||
| 145 | __gnu_pbds::priority_queue<std::pair<cost_t, int>> q; | ||
| 146 | std::vector<typename decltype(q)::point_iterator> its; | ||
| 147 | ✗ | cost_t dijkstra(int s, int t) { | |
| 148 | ✗ | dist.assign(N, INF_COST); | |
| 149 | ✗ | dist[s] = 0; | |
| 150 | |||
| 151 | ✗ | its.assign(N, q.end()); | |
| 152 | ✗ | its[s] = q.push({-(dist[s] - pi[s]), s}); | |
| 153 | |||
| 154 | ✗ | while (!q.empty()) { | |
| 155 | ✗ | int i = q.top().second; q.pop(); | |
| 156 | ✗ | cost_t d = dist[i]; | |
| 157 | ✗ | for (int e : adj[i]) { | |
| 158 | ✗ | if (edges[e].cap) { | |
| 159 | ✗ | int j = edges[e].dest; | |
| 160 | ✗ | cost_t nd = d + edges[e].cost; | |
| 161 | ✗ | if (nd < dist[j]) { | |
| 162 | ✗ | dist[j] = nd; | |
| 163 | ✗ | if (its[j] == q.end()) { | |
| 164 | ✗ | its[j] = q.push({-(dist[j] - pi[j]), j}); | |
| 165 | } else { | ||
| 166 | ✗ | q.modify(its[j], {-(dist[j] - pi[j]), j}); | |
| 167 | } | ||
| 168 | } | ||
| 169 | } | ||
| 170 | } | ||
| 171 | } | ||
| 172 | |||
| 173 | ✗ | std::swap(pi, dist); | |
| 174 | ✗ | return pi[t]; | |
| 175 | } | ||
| 176 | |||
| 177 | std::vector<int> buf; | ||
| 178 | std::vector<int> level; | ||
| 179 | ✗ | flow_t dinic_dfs(int cur, int t, flow_t f) { | |
| 180 | ✗ | if (cur == t) return f; | |
| 181 | ✗ | flow_t cur_f = 0; | |
| 182 | ✗ | assert(f > 0); | |
| 183 | ✗ | for (; buf[cur] < int(adj[cur].size()); buf[cur]++) { | |
| 184 | ✗ | int e = adj[cur][buf[cur]]; | |
| 185 | ✗ | int nxt = edges[e].dest; | |
| 186 | ✗ | if (level[nxt] == level[cur] + 1 && edges[e].cap > 0 && edges[e].cost == pi[nxt] - pi[cur]) { | |
| 187 | ✗ | flow_t v = dinic_dfs(nxt, t, std::min(f, edges[e].cap)); | |
| 188 | ✗ | edges[e].cap -= v; | |
| 189 | ✗ | edges[e^1].cap += v; | |
| 190 | ✗ | f -= v; | |
| 191 | ✗ | cur_f += v; | |
| 192 | ✗ | if (f == 0) break; | |
| 193 | } | ||
| 194 | } | ||
| 195 | ✗ | return cur_f; | |
| 196 | } | ||
| 197 | ✗ | flow_t dinic(int s, int t) { | |
| 198 | ✗ | flow_t tot_flow = 0; | |
| 199 | ✗ | while (true) { | |
| 200 | ✗ | buf.clear(); | |
| 201 | ✗ | buf.reserve(N); | |
| 202 | ✗ | level.assign(N, -1); | |
| 203 | ✗ | buf.push_back(s); | |
| 204 | ✗ | level[s] = 0; | |
| 205 | ✗ | for (int z = 0; z < int(buf.size()); z++) { | |
| 206 | ✗ | int cur = buf[z]; | |
| 207 | ✗ | for (int e : adj[cur]) { | |
| 208 | ✗ | int nxt = edges[e].dest; | |
| 209 | ✗ | if (edges[e].cap > 0 && edges[e].cost == pi[nxt] - pi[cur] && level[nxt] == -1) { | |
| 210 | ✗ | level[nxt] = level[cur] + 1; | |
| 211 | ✗ | buf.push_back(nxt); | |
| 212 | } | ||
| 213 | } | ||
| 214 | } | ||
| 215 | ✗ | if (level[t] == -1) break; | |
| 216 | ✗ | buf.assign(N, 0); | |
| 217 | ✗ | tot_flow += dinic_dfs(s, t, INF_FLOW); | |
| 218 | } | ||
| 219 | ✗ | return tot_flow; | |
| 220 | } | ||
| 221 | |||
| 222 | ✗ | std::vector<std::pair<flow_t, cost_t>> all_flows(int s, int t, cost_t max_cost = INF_COST - 1) { | |
| 223 | ✗ | assert(s != t); | |
| 224 | ✗ | std::vector<std::pair<flow_t, cost_t>> res; | |
| 225 | ✗ | while (dijkstra(s, t) <= max_cost) { | |
| 226 | ✗ | assert(res.empty() || pi[t] > res.back().second); | |
| 227 | ✗ | flow_t f = dinic(s, t); | |
| 228 | ✗ | res.push_back({f, pi[t]}); | |
| 229 | } | ||
| 230 | ✗ | return res; | |
| 231 | } | ||
| 232 | |||
| 233 | ✗ | std::pair<flow_t, cost_t> max_flow(int s, int t, cost_t max_cost = INF_COST - 1) { | |
| 234 | ✗ | assert(s != t); | |
| 235 | ✗ | flow_t tot_flow = 0; cost_t tot_cost = 0; | |
| 236 | ✗ | while (dijkstra(s, t) <= max_cost) { | |
| 237 | ✗ | flow_t cur_flow = dinic(s, t); | |
| 238 | ✗ | tot_flow += cur_flow; | |
| 239 | ✗ | tot_cost += cur_flow * pi[t]; | |
| 240 | } | ||
| 241 | ✗ | return {tot_flow, tot_cost}; | |
| 242 | } | ||
| 243 | }; | ||
| 244 | |||
| 245 | template <typename flow_t = int, typename tot_flow_t = flow_t> | ||
| 246 | struct Dinic { | ||
| 247 | int N; | ||
| 248 | std::vector<std::vector<int>> adj; | ||
| 249 | struct edge_t { | ||
| 250 | int dest; | ||
| 251 | flow_t cap; | ||
| 252 | }; | ||
| 253 | std::vector<edge_t> edges; | ||
| 254 | |||
| 255 | std::vector<char> seen; | ||
| 256 | |||
| 257 | ✗ | explicit Dinic(int N_) : N(N_), adj(N) {} | |
| 258 | |||
| 259 | ✗ | void add_edge(int from, int to, flow_t cap) { | |
| 260 | ✗ | return add_bi_edge(from, to, cap, 0); | |
| 261 | } | ||
| 262 | |||
| 263 | ✗ | void add_bi_edge(int from, int to, flow_t cap, flow_t rev_cap) { | |
| 264 | ✗ | assert(cap >= 0); | |
| 265 | ✗ | assert(rev_cap >= 0); | |
| 266 | ✗ | int e = int(edges.size()); | |
| 267 | ✗ | edges.emplace_back(edge_t{to, cap}); | |
| 268 | ✗ | edges.emplace_back(edge_t{from, rev_cap}); | |
| 269 | ✗ | adj[from].push_back(e); | |
| 270 | ✗ | adj[to].push_back(e+1); | |
| 271 | } | ||
| 272 | |||
| 273 | static constexpr tot_flow_t INF_FLOW = std::numeric_limits<tot_flow_t>::max() / 4; | ||
| 274 | std::vector<int> buf; | ||
| 275 | std::vector<int> level; | ||
| 276 | ✗ | tot_flow_t dinic_dfs(int cur, int t, tot_flow_t f) { | |
| 277 | ✗ | if (cur == t) return f; | |
| 278 | ✗ | tot_flow_t cur_f = 0; | |
| 279 | ✗ | assert(f > 0); | |
| 280 | ✗ | for (; buf[cur] < int(adj[cur].size()); buf[cur]++) { | |
| 281 | ✗ | int e = adj[cur][buf[cur]]; | |
| 282 | ✗ | int nxt = edges[e].dest; | |
| 283 | ✗ | if (level[nxt] == level[cur] + 1 && edges[e].cap > 0) { | |
| 284 | ✗ | flow_t v = flow_t(dinic_dfs(nxt, t, std::min<tot_flow_t>(f, edges[e].cap))); | |
| 285 | ✗ | edges[e].cap -= v; | |
| 286 | ✗ | edges[e^1].cap += v; | |
| 287 | ✗ | f -= v; | |
| 288 | ✗ | cur_f += v; | |
| 289 | ✗ | if (f == 0) break; | |
| 290 | } | ||
| 291 | } | ||
| 292 | ✗ | return cur_f; | |
| 293 | } | ||
| 294 | ✗ | tot_flow_t dinic(int s, int t) { | |
| 295 | ✗ | tot_flow_t tot_flow = 0; | |
| 296 | ✗ | while (true) { | |
| 297 | ✗ | buf.clear(); | |
| 298 | ✗ | buf.reserve(N); | |
| 299 | ✗ | level.assign(N, -1); | |
| 300 | ✗ | buf.push_back(s); | |
| 301 | ✗ | level[s] = 0; | |
| 302 | ✗ | for (int z = 0; z < int(buf.size()); z++) { | |
| 303 | ✗ | int cur = buf[z]; | |
| 304 | ✗ | for (int e : adj[cur]) { | |
| 305 | ✗ | int nxt = edges[e].dest; | |
| 306 | ✗ | if (edges[e].cap > 0 && level[nxt] == -1) { | |
| 307 | ✗ | level[nxt] = level[cur] + 1; | |
| 308 | ✗ | buf.push_back(nxt); | |
| 309 | } | ||
| 310 | } | ||
| 311 | } | ||
| 312 | ✗ | if (level[t] == -1) break; | |
| 313 | ✗ | buf.assign(N, 0); | |
| 314 | ✗ | tot_flow += dinic_dfs(s, t, INF_FLOW); | |
| 315 | } | ||
| 316 | ✗ | return tot_flow; | |
| 317 | } | ||
| 318 | ✗ | tot_flow_t max_flow(int s, int t) { return dinic(s, t); } | |
| 319 | }; | ||
| 320 | |||
| 321 | } // namespace wala | ||
| 322 |