library/flow/NondecreasingMCF.hpp
Depends on
Verified with
Code
#pragma once
// 辺の重みが流量に対して単調増加な関数
// 現在の流量を引数として、そこに新たに 1 流す時にかかるコストを返す関数を渡す
#include "library/graph/WeightedGraph.hpp"
#define REP_(i, n) for (int i = 0; i < (n); i++)
template <typename TC> class NondecreasingMCF {
using F = std::function<TC(int)>;
struct EdgeInfo {
int cap, flow, rev;
bool reverse_edge;
F cost_func;
EdgeInfo() = default;
EdgeInfo(int cap, F cost_func, int rev, bool reverse_edge)
: cap(cap), cost_func(cost_func), rev(rev),
reverse_edge(reverse_edge), flow(0) {}
TC cost() const {
if (!reverse_edge)
return cost_func(flow);
return -cost_func(cap - 1);
}
};
int n;
WeightedGraph<EdgeInfo> G;
std::vector<TC> potential, dist;
static constexpr TC INF = std::is_same_v<TC, __int128>
? TC(1e30)
: std::numeric_limits<TC>::max() / 2;
// std::numeric_limits<__int128 >::max() は AOJ でバグった
std::vector<std::pair<int, int>> pre; // pre[v]=[u,i] : G[u][i] で v に来た
std::vector<int> in_deg, out_deg;
std::priority_queue<std::pair<TC, int>, std::vector<std::pair<TC, int>>,
std::greater<std::pair<TC, int>>>
que;
bool negative = false; // 負辺存在するか
template <typename T> bool chmin(T &a, const T &b) {
return (a > b and (a = b, true));
}
bool SP_update(int from, int edge_id) {
const auto &e = G[from][edge_id];
if ((e.weight).cap == 0)
return false;
if (chmin(dist[e.to], dist[from] + (e.weight).cost() + potential[from] -
potential[e.to])) {
pre[e.to] = {from, edge_id};
return true;
}
return false;
}
void dijkstra(int s) { // dist[i]:sから残余グラフで辺の重みによるiへの最短路
// となるようにdistを作る
std::ranges::fill(dist, INF);
dist[s] = 0;
que.emplace(0, s);
while (que.size()) {
const auto [now, v] = que.top();
que.pop();
if (dist[v] < now)
continue;
REP_(i, G[v].size())
if (SP_update(v, i))
que.emplace(dist[G[v][i].to], G[v][i].to);
}
}
void DAG(int s) {
negative = false;
std::ranges::fill(dist, INF);
dist[s] = 0;
std::queue<int> que;
REP_(i, n) if (!in_deg[i]) que.push(i);
while (que.size()) {
int v = que.front();
que.pop();
REP_(i, G[v].size()) {
SP_update(v, i);
if (!--in_deg[G[v][i].to])
que.push(G[v][i].to);
}
}
}
public:
NondecreasingMCF() {}
NondecreasingMCF(int n_)
: n(n_), G(n_), potential(n_, 0), dist(n_), pre(n_), in_deg(n_, 0),
out_deg(n_, 0), negative(false) {}
void add_arc(int u, int v, int cap, F cost_func) {
G.add_arc(u, v, {cap, cost_func, out_deg[v], false});
G.add_arc(v, u, {0, cost_func, out_deg[u], true});
out_deg[v]++;
out_deg[u]++;
if (cap > 0) {
in_deg[v]++;
negative |= cost_func(0) < 0;
}
}
std::pair<TC, bool> flow(int s, int t, int f) {
if (!G.is_prepared())
G.build();
TC res = 0;
std::ranges::fill(potential, 0);
for (int i = 0; i < f; i++) {
if (negative)
DAG(s);
else
dijkstra(s);
if (dist[t] == INF)
return std::make_pair(res, false);
REP_(v, n) if (dist[v] != INF) potential[v] += dist[v];
res += potential[t];
for (int v = t; v != s; v = pre[v].first) {
auto &w = G[pre[v].first][pre[v].second].weight;
w.cap--;
w.flow++;
auto &r = G[v][w.rev].weight;
r.cap++;
r.flow--;
}
}
return std::make_pair(res, true);
}
};
#undef REP_
#line 2 "library/flow/NondecreasingMCF.hpp"
// 辺の重みが流量に対して単調増加な関数
// 現在の流量を引数として、そこに新たに 1 流す時にかかるコストを返す関数を渡す
#line 2 "library/graph/WeightedGraph.hpp"
template <typename T> struct WeightedEdge {
WeightedEdge() = default;
WeightedEdge(int from, int to, T weight)
: from(from), to(to), weight(weight) {}
int from, to;
T weight;
operator int() const { return to; }
};
template <typename T> struct WeightedGraph {
int n;
using weight_type = T;
using edge_type = WeightedEdge<T>;
std::vector<edge_type> edges;
protected:
std::vector<int> in_deg;
bool prepared;
class OutgoingEdges {
WeightedGraph *g;
int l, r;
public:
OutgoingEdges(WeightedGraph *g, int l, int r) : g(g), l(l), r(r) {}
edge_type *begin() { return &(g->edges[l]); }
edge_type *end() { return &(g->edges[r]); }
edge_type &operator[](int i) { return g->edges[l + i]; }
int size() const { return r - l; }
};
class ConstOutgoingEdges {
const WeightedGraph *g;
int l, r;
public:
ConstOutgoingEdges(const WeightedGraph *g, int l, int r)
: g(g), l(l), r(r) {}
const edge_type *begin() const { return &(g->edges[l]); }
const edge_type *end() const { return &(g->edges[r]); }
const edge_type &operator[](int i) const { return g->edges[l + i]; }
int size() const { return r - l; }
};
public:
OutgoingEdges operator[](int v) {
assert(prepared);
return {this, in_deg[v], in_deg[v + 1]};
}
const ConstOutgoingEdges operator[](int v) const {
assert(prepared);
return {this, in_deg[v], in_deg[v + 1]};
}
bool is_prepared() const { return prepared; }
WeightedGraph() : n(0), in_deg(1, 0), prepared(false) {}
WeightedGraph(int n) : n(n), in_deg(n + 1, 0), prepared(false) {}
WeightedGraph(int n, int m, bool directed = false, int indexed = 1)
: n(n), in_deg(n + 1, 0), prepared(false) {
scan(m, directed, indexed);
}
void resize(int n) { n = n; }
void add_arc(int from, int to, T weight) {
assert(!prepared);
assert(0 <= from and from < n and 0 <= to and to < n);
edges.emplace_back(from, to, weight);
in_deg[from + 1]++;
}
void add_edge(int u, int v, T weight) {
add_arc(u, v, weight);
add_arc(v, u, weight);
}
void add_arc(const edge_type &e) { add_arc(e.from, e.to, e.weight); }
void add_edge(const edge_type &e) { add_edge(e.from, e.to, e.weight); }
void scan(int m, bool directed = false, int indexed = 1) {
edges.reserve(directed ? m : 2 * m);
while (m--) {
int u, v;
std::cin >> u >> v;
u -= indexed;
v -= indexed;
T weight;
std::cin >> weight;
if (directed)
add_arc(u, v, weight);
else
add_edge(u, v, weight);
}
build();
}
void build() {
assert(!prepared);
prepared = true;
for (int v = 0; v < n; v++)
in_deg[v + 1] += in_deg[v];
std::vector<edge_type> new_edges(in_deg.back());
auto counter = in_deg;
for (auto &&e : edges)
new_edges[counter[e.from]++] = e;
edges = new_edges;
}
void graph_debug() const {
#ifndef __DEBUG
return;
#endif
assert(prepared);
for (int from = 0; from < n; from++) {
std::cerr << from << ";";
for (int i = in_deg[from]; i < in_deg[from + 1]; i++)
std::cerr << "(" << edges[i].to << "," << edges[i].weight
<< ")";
std::cerr << "\n";
}
}
};
#line 5 "library/flow/NondecreasingMCF.hpp"
#define REP_(i, n) for (int i = 0; i < (n); i++)
template <typename TC> class NondecreasingMCF {
using F = std::function<TC(int)>;
struct EdgeInfo {
int cap, flow, rev;
bool reverse_edge;
F cost_func;
EdgeInfo() = default;
EdgeInfo(int cap, F cost_func, int rev, bool reverse_edge)
: cap(cap), cost_func(cost_func), rev(rev),
reverse_edge(reverse_edge), flow(0) {}
TC cost() const {
if (!reverse_edge)
return cost_func(flow);
return -cost_func(cap - 1);
}
};
int n;
WeightedGraph<EdgeInfo> G;
std::vector<TC> potential, dist;
static constexpr TC INF = std::is_same_v<TC, __int128>
? TC(1e30)
: std::numeric_limits<TC>::max() / 2;
// std::numeric_limits<__int128 >::max() は AOJ でバグった
std::vector<std::pair<int, int>> pre; // pre[v]=[u,i] : G[u][i] で v に来た
std::vector<int> in_deg, out_deg;
std::priority_queue<std::pair<TC, int>, std::vector<std::pair<TC, int>>,
std::greater<std::pair<TC, int>>>
que;
bool negative = false; // 負辺存在するか
template <typename T> bool chmin(T &a, const T &b) {
return (a > b and (a = b, true));
}
bool SP_update(int from, int edge_id) {
const auto &e = G[from][edge_id];
if ((e.weight).cap == 0)
return false;
if (chmin(dist[e.to], dist[from] + (e.weight).cost() + potential[from] -
potential[e.to])) {
pre[e.to] = {from, edge_id};
return true;
}
return false;
}
void dijkstra(int s) { // dist[i]:sから残余グラフで辺の重みによるiへの最短路
// となるようにdistを作る
std::ranges::fill(dist, INF);
dist[s] = 0;
que.emplace(0, s);
while (que.size()) {
const auto [now, v] = que.top();
que.pop();
if (dist[v] < now)
continue;
REP_(i, G[v].size())
if (SP_update(v, i))
que.emplace(dist[G[v][i].to], G[v][i].to);
}
}
void DAG(int s) {
negative = false;
std::ranges::fill(dist, INF);
dist[s] = 0;
std::queue<int> que;
REP_(i, n) if (!in_deg[i]) que.push(i);
while (que.size()) {
int v = que.front();
que.pop();
REP_(i, G[v].size()) {
SP_update(v, i);
if (!--in_deg[G[v][i].to])
que.push(G[v][i].to);
}
}
}
public:
NondecreasingMCF() {}
NondecreasingMCF(int n_)
: n(n_), G(n_), potential(n_, 0), dist(n_), pre(n_), in_deg(n_, 0),
out_deg(n_, 0), negative(false) {}
void add_arc(int u, int v, int cap, F cost_func) {
G.add_arc(u, v, {cap, cost_func, out_deg[v], false});
G.add_arc(v, u, {0, cost_func, out_deg[u], true});
out_deg[v]++;
out_deg[u]++;
if (cap > 0) {
in_deg[v]++;
negative |= cost_func(0) < 0;
}
}
std::pair<TC, bool> flow(int s, int t, int f) {
if (!G.is_prepared())
G.build();
TC res = 0;
std::ranges::fill(potential, 0);
for (int i = 0; i < f; i++) {
if (negative)
DAG(s);
else
dijkstra(s);
if (dist[t] == INF)
return std::make_pair(res, false);
REP_(v, n) if (dist[v] != INF) potential[v] += dist[v];
res += potential[t];
for (int v = t; v != s; v = pre[v].first) {
auto &w = G[pre[v].first][pre[v].second].weight;
w.cap--;
w.flow++;
auto &r = G[v][w.rev].weight;
r.cap++;
r.flow--;
}
}
return std::make_pair(res, true);
}
};
#undef REP_
Back to top page