食用前请先了解 SPFA + Dinic/EK 求解 MCMF。
Sol.
总所周知,SPFA 牺牲了。于是我们寻求一些更稳定的算法求解 MCMF。
网络流算法的时间属于玄学,暂且判定为混乱中的稳定。那么我们就只能考虑在最短路算法上寻求优化。于是就想到了 Dijkstra。
但 Dijkstra 有一个致命的弱点:无法处理负权边。而我们应用的场景显然含有负权。
开动脑筋想一想可以想到一个“给所有边权加上巨大多权值进而规避负权边”的方法。
但这样在实现中,还需要记录一条最短路目前经过了哪些边之类的奇怪信息。不是说不可做,但确实复杂。
那我们考虑把这样累加的权值放在点上?就有了正解想法:在点上叠加势能。
这当然不是什么基于经典力学与计算机科学的逻辑思想统一能广泛应用的以简化实现提高效率推动人类社会发展促进世界和平为目的最新学科交叉成果。
它只是和“势能”同名,仅此而已。
具体的讲,我们考虑将每一个点附上一个权值 \(h\),并将 \(u, v\) 两点间的边 \((u, v)\) 的权值替换为 \(w_{u, v} + h_u - h_v\)。
在这张图上,对于一条起点为 \(s\) 终点为 \(t\) 的路径边集 \(D\),边权和为 \(\sum \limits _{(u, v) \in D} (w_{u, v} + h_u - h_v) = h_s - h_t + \sum \limits _{(u, v) \in D} w_{u, v}\)。
也就是说我们可以通过这张图上的最短路长度还原原图的最短路长度,且两图的路径一定一一对应。
容易发现,如果所有的 \(w_{u, v} + h_u - h_v\) 均大于等于 \(0\),即 \(w_{u, v} + h_u \geq h_v\),那么我们就可以用 Dijkstra 来解决问题。
好像构造 \(h\) 成为了麻烦事,不过我们貌似可以用一些现成的量来“充数”。
观察上面的不等式,发现是一个类三角不等式,而在最短路问题中也存在这样的三角不等式,即:记 \(f_x\) 表示起点 \(s\) 到点 \(x\) 的最短路长度,则满足 \(f_u + w_{u, v} \geq f_v\)。这和上面的式子形式一模一样!
而 Dinic/EK 是需要跑很多次最短路的,所以我们可以自然想到将上一次的最短路答案当作这一次的 \(h\),注意这里指的最短路是原图中的最短路而不是被势能改过的最短路,要注意区分。
那么此问题就结了。
最后一点就是说,势能的初值通常设 \(0\) ,而这并不一定能满足第一次 Dijkstra 直接就能跑。所以我们需要先跑一个 SPFA 去求出初始的势能(也还是等于最短路)。
只有一只 SPFA 牺牲了不会让整个代码都牺牲。
Code.
「Luogu P3381」这里是 EK 实现。
#include <queue>
#include <cstdio>
using namespace std;
int Abs(int x) { return x < 0 ? -x : x; }
int Max(int x, int y) { return x > y ? x : y; }
int Min(int x, int y) { return x < y ? x : y; }
int read() {
int x = 0, k = 1;
char s = getchar();
while(s < '0' || s > '9') {
if(s == '-')
k = -1;
s = getchar();
}
while('0' <= s && s <= '9') {
x = (x << 3) + (x << 1) + (s ^ 48);
s = getchar();
}
return x * k;
}
void write(int x) {
if(x < 0) {
x = -x;
putchar('-');
}
if(x > 9)
write(x / 10);
putchar(x % 10 + '0');
}
void print(int x, char s) {
write(x);
putchar(s);
}
const int MAXN = 5e3 + 5;
const int MAXM = 5e4 + 5;
const int INF = 2147483647;
struct edge {
int v, nxt, Wei, Cap, Flow;
edge() {}
edge(int V, int Nxt, int C, int W, int F) {
v = V, nxt = Nxt, Cap = C, Wei = W, Flow = F;
}
} e[MAXM << 1];
int head[MAXM << 1], cnt = 0;
void Add_Edge(int u, int v, int c, int w) {
e[cnt] = edge(v, head[u], c, w, 0);
head[u] = cnt++;
e[cnt] = edge(u, head[v], 0, -w, 0);
head[v] = cnt++;
}
queue<int> q;
bool vis[MAXN];
int Dist[MAXN], Aug[MAXN], h[MAXN], n, m;
struct Back {
int Pre, id;
Back() {}
Back(int P, int Id) {
Pre = P, id = Id;
}
} Last[MAXN];
struct node {
int x, dis;
node() {}
node(int X, int Dis) {
x = X, dis = Dis;
}
friend bool operator < (node One, node TheOther) {
return One.dis > TheOther.dis;
}
};
void spfa(int s, int t) {
for(int i = 1; i <= n; i++)
h[i] = INF, vis[i] = false;
h[s] = 0, vis[s] = true;
queue<int> q;
q.push(s);
while(!q.empty()) {
int u = q.front(); q.pop();
vis[u] = false;
for(int i = head[u], v; ~i; i = e[i].nxt) {
v = e[i].v;
if(e[i].Cap - e[i].Flow > 0 && h[v] > h[u] + e[i].Wei) {
h[v] = h[u] + e[i].Wei;
if(!vis[v])
vis[v] = true, q.push(v);
}
}
}
}
bool Dijkstra(int s, int t) {
for(int i = 1; i <= n; i++)
Dist[i] = INF, Last[i] = Back(-1, -1), vis[i] = false, Aug[i] = INF;
priority_queue<node> q;
Dist[s] = 0;
q.push(node(s, Dist[s]));
while(!q.empty()) {
int u = q.top().x; q.pop();
if(vis[u])
continue;
vis[u] = true;
for(int i = head[u], v; ~i; i = e[i].nxt) {
v = e[i].v;
if(e[i].Cap - e[i].Flow > 0 && Dist[v] > Dist[u] + e[i].Wei + h[u] - h[v]) {
Last[v] = Back(u, i);
Dist[v] = Dist[u] + e[i].Wei + h[u] - h[v];
Aug[v] = Min(Aug[u], e[i].Cap - e[i].Flow);
q.push(node(v, Dist[v]));
}
}
}
return Dist[t] != INF;
}
int Flow, Cost;
void EK(int s, int t) {
Flow = 0, Cost = 0;
spfa(s, t);
while(Dijkstra(s, t)) {
Flow += Aug[t];
Cost += Aug[t] * (Dist[t] + h[t]);
int pos = t;
while(pos != s) {
e[Last[pos].id].Flow += Aug[t];
e[Last[pos].id ^ 1].Flow -= Aug[t];
pos = Last[pos].Pre;
}
for(int i = 1; i <= n; i++)
if(h[i] < INF)
h[i] += Dist[i];
}
}
int main() {
n = read(), m = read();
int s = read(), t = read();
for(int i = 1; i <= n; i++)
head[i] = -1;
for(int i = 1, u, v, c, w; i <= m; i++) {
u = read(), v = read(), c = read(), w = read();
Add_Edge(u, v, c, w);
}
EK(s, t);
print(Flow, ' '), print(Cost, '\n');
return 0;
}