费用流 Dijkstra 原始对偶方法(primal-dual method)
简单叙述用Dijkstra求费用流
Dijkstra不能求有负权边的最短路。
类似于Johnson算法，我们也可以设计一个势函数，以满足在与原图等价的新图中的边权非负。
但是这个算法并不能处理有负圈的情况（可能需要消圈算法）。
对网络\(G\)中的每一个点设置一个势函数\(h(u)\)，在任意残留网络G'的任意边\((u, v)\)都需要满足\(w_{u, v} + h(u) - h(v) \ge 0\)。
令图G的对偶图（不知道能不能这么说）为\(G'\)，其对应的边\((u, v)\)的权值为\(w'_{u, v} = w_{u, v} + h(u) - h(v)\)
对于原图中的任意一条路径\((u_1, u_2, \cdots, u_k)\)，它在\(G\)中的权值为\(w_{u_1, u_2} + w_{u_2, u_3} + \cdots + w_{u_{k - 1}, u_k}\)，在\(G'\)中的权值为
\(w'_{u_1, u_2} + w'_{u_2, u_3} + \cdots + w'_{u_{k - 1}, u_k}\)
\(= w_{u_1, u_2} + w_{u_2, u_3} + \cdots + w_{u_{k - 1}, u_k} + h(u_1) - h(u_2) + h(u_2) - h(u_3) + \cdots + h(u_{k - 1}) - h(u_k)\)
\(= w_{u_1, u_2} + w_{u_2, u_3} + \cdots + w_{u_{k - 1}, u_k} + h(u_1) - h(u_k)\)
所以，我们在\(G'\)求出的路径都可以对应到\(G\)上，令\(dist_{u, v}\)为图\(G\)点\(u\)到点\(v\)的最短路径，\(dist'_{u, v}\)为图\(G'\)点\(u\)到点\(v\)的最短路径，显然有\(dist_{u, v} = dist'_{u, v} - h(u) + h(v)\)
所以我们只需要求\(G'\)的最短路径，就能对应回原图的最短路径。
若网络\(G\)初始边权非负，我们可令\(h(u) = 0\)
否则我们令\(h(u) = dist_{s, u}\)，这个可以用Bellman-Ford算法解决。（这与Johnson算法是一模一样的）
我们考虑怎么维护势函数\(h(u)\)。
令网路\(G\)上\(dist_{s, u} = d_u\)，网路\(G'\)上\(dist'_{s, u} = d'_u\)
令残余网路\(G'\)上新的势函数为\(h'(u)\)。
对于残余网络上的一条边\((u, v)\)，有两种可能：
\((u, v) \in G\)，那么有\(d_u + w_{u, v} + h(u) - h(v) \ge d_v\)
移项得\(w_{u, v} + (h(u) + d_u) - (h(v) + d_v) \ge 0\)
\((u, v) \in G\)，那么\((v, u) \in G的增广路\)，就有\(d_v + w_{v, u} + h(v) - h(u) = d_u\)
移项得\(-w_{v, u} + (h(u) + d_u) - (h(v) - d_v) = 0\)
由\(w_{u, v} = -w_{v, u}\)可得\(w_{u, v} + (h(u) + d_u) - (h(v) - d_v) = 0\)
也就有\(w_{u, v} + (h(u) + d_u) - (h(v) - d_v) \ge 0\)
所以我们不妨令\(h'(u) = h(u) + d_u\)，这就维护好了势函数\(h(u)\)
单路增广
#include <bits/stdc++.h>
using namespace std;
typedef pair<int, int> PII;
#define fi first
#define se second
#define mp make_pair
const int N = 5005, M = 50005;
const int INF = ~0u >> 1;
int n, m, s, t;
int tot = -1, head[N];
struct Edge {
int p, nxt, c, w;
Edge(int p = 0, int nxt = 0, int c = 0, int w = 0) : p(p), nxt(nxt), c(c), w(w) {}
} edge[M * 2];
inline void Add_Edge(int u, int v, int c, int w) {
edge[++tot] = Edge(v, head[u], c, w);
head[u] = tot;
return;
}
int h[N], d[N], dc[N], pr[N];
priority_queue<PII> pq;
bool Dijkstra(int s, int t, int &mf, int &mc) {
fill(d + 1, d + n + 1, INF);
dc[s] = INF;
d[s] = 0;
pr[s] = -1;
pq.push(mp(0, s));
while (!pq.empty()) {
PII cur = pq.top();
pq.pop();
int u = cur.se;
if (-cur.fi > d[u]) continue;
for (int i = head[u]; ~i; i = edge[i].nxt) {
int v = edge[i].p, c = edge[i].c, w = edge[i].w + h[u] - h[v];
if (!c) continue;
if (d[v] > d[u] + w) {
d[v] = d[u] + w;
pr[v] = i;
dc[v] = min(dc[u], c);
pq.push(mp(-d[v], v));
}
}
}
if (d[t] == INF) return 0;
for (int i = 1; i <= n; ++i)
if (d[i] < INF) h[i] += d[i];
int c = dc[t];
mf += c;
mc += c * h[t];
for (int x = t; ~pr[x]; x = edge[pr[x] ^ 1].p) {
edge[pr[x]].c -= c;
edge[pr[x] ^ 1].c += c;
}
return 1;
}
int main() {
scanf("%d%d%d%d", &n, &m, &s, &t);
memset(head, -1, sizeof(head));
for (int i = 1; i <= m; ++i) {
int u, v, c, w;
scanf("%d%d%d%d", &u, &v, &c, &w);
Add_Edge(u, v, c, w);
Add_Edge(v, u, 0, -w);
}
int mf = 0, mc = 0;
while (Dijkstra(s, t, mf, mc));
printf("%d %d\n", mf, mc);
return 0;
}
多路增广（复杂度不变）
#include <bits/stdc++.h>
using namespace std;
typedef pair<int, int> PII;
#define fi first
#define se second
#define mp make_pair
const int N = 5005, M = 50005;
const int INF = ~0u >> 1;
int n, m, s, t;
int tot = -1, head[N], cur[N];
struct Edge {
int p, nxt, c, w;
Edge(int p = 0, int nxt = 0, int c = 0, int w = 0) : p(p), nxt(nxt), c(c), w(w) {}
} edge[M * 2];
inline void Add_Edge(int u, int v, int c, int w) {
edge[++tot] = Edge(v, head[u], c, w);
head[u] = tot;
return;
}
int h[N], d[N];
priority_queue<PII> pq;
bool Dijkstra(int s, int t) {
fill(d + 1, d + n + 1, INF);
pq.push(mp(d[s] = 0, s));
while (!pq.empty()) {
PII cur = pq.top();
pq.pop();
int u = cur.se;
if (-cur.fi > d[u]) continue;
for (int i = head[u]; ~i; i = edge[i].nxt) {
int v = edge[i].p, c = edge[i].c, w = edge[i].w + h[u] - h[v];
if (c && d[v] > d[u] + w) pq.push(mp(-(d[v] = d[u] + w), v));
}
}
return d[t] < INF;
}
int v[N];
int DFS(int u, int c, int t) {
if (u == t) return c;
int r = c;
v[u] = 1;
for (int &i = cur[u]; ~i && r; i = edge[i].nxt) {
int v = edge[i].p, c = edge[i].c, w = edge[i].w + h[u] - h[v];
if (!::v[v] && c && d[u] + w == d[v]) {
int x = DFS(v, min(r, c), t);
r -= x;
edge[i].c -= x;
edge[i ^ 1].c += x;
}
}
v[u] = 0;
return c - r;
}
int main() {
scanf("%d%d%d%d", &n, &m, &s, &t);
memset(head, -1, sizeof(head));
for (int i = 1; i <= m; ++i) {
int u, v, c, w;
scanf("%d%d%d%d", &u, &v, &c, &w);
Add_Edge(u, v, c, w);
Add_Edge(v, u, 0, -w);
}
int mf = 0, mc = 0;
while (Dijkstra(s, t)) {
memcpy(cur, head, sizeof(cur));
int c = DFS(s, INF, t);
for (int i = 1; i <= n; ++i)
if (d[i] < INF) h[i] += d[i];
mf += c;
mc += c * h[t];
}
printf("%d %d\n", mf, mc);
return 0;
}
