引入
一个优秀的代码, 时间复杂度一定是很优的, SPFA + EK/dinic 已经满足不了我们的需求了, 所以吃饱了撑着的善于思考的人类不断地探索发现, 一个更加优化的算法就此诞生。
详解
考虑之前的 SPFA + EK/dinic 算法, 我们发现让我们被卡的飞起的地方就是 SPFA 那个**的 \(O(nm)\) 时间复杂度, 所以我们迫切地需要一个更快, 更高效的算法来解决这个问题。 考虑到比 SPFA 快的单元最短路算法 dijkstra 不能跑负边, 所以我们需要对其进行改进。
考虑对于每一个点都附上一个权值 \(pe(i)\) 再将每条边 (假设为 \(u->v\) ) 的权值变为 \(w(u, v) + pe(u) - pe(v)\) 我们可以发现, 在中间点的权值就被抵消掉了, 只剩下了起点和终点, 也就是说最短路选的路径没有变, 但是最短路的值加上了 \(pe(s)\) 减去了 \(pe(t)\) , 当然, 这对网络流里求增广路没有影响。如果我们的 \(w(u, v) + pe(u) - pe(v)\) 都非负, 我们的 dijkstra 就活过来了。
那么, 我们如何让其非负呢?通过大量实验经过严谨的推理后我们发现将起点到某个点 \(v\) , 的 \(pe(v)\) 是一个很好的解决方案。因为对于一个点 \(v\) 来说, 所有指向它的点 \(u\) 都有 \(pe(u) + w(u, v) \le pe(v)\) 否则 \(pe(v)\) 就不是最短路了。
但这以后又有了一个新的问题, 我们的网络流每次增广后都会产生图的改变, 那么我们所求的最开始的最短路也就可能出现问题, 如果这样的话, 我们每次都要用可以跑负边的 SPFA 去跑一次最短路, 那么 dijkstra 也就失去了意义。
这又如何解决呢?
其实也并不难。 我们最开始需要跑一次 SPFA 来找到最短路(相比于网络流, 一次的时间还是不算什么的)。接下来我们考虑第一次增广的 \(pe(i)\) 是已经找出来的, 那么直接进行此次增广后可能会出现有边反向, 但此时我们先不管它, 因为这些反向的边都是在最短路上的, 所以如果存在 \(u->v\) 反向, 那么因为刚才的保证了 \(pe(u) - pe(v) + w (u, v) \ge 0\) 的, 那么 \(pe(v) + (-w(u, v)) \ge pe(u)\) 也是成立的, 也就是说经过这条反边是不会改变 \(pe(u)\) 的, 那么 \(pe(u)\) 可以继续沿用。 也就是说我们每次寻找最短路的 \(pe(i)\) 都可以直接用上一次增广时使用的最短路, 而上一次增广路是直接用 \(dijkstra\) 求出来的, 也就是增广时顺便求出的。
这样我们就完成了整个算法的讲解
补充
在各类资料中, 我们都将这里的 \(pe(i)\) 称作势能, 这个东西其实挺牛逼的, \(Johnson 全源最短路径算法\) 也是通过这个方法将边全变为正, 然后再用 dijkstra 跑 \(n\)次单元最短路, 最后做出了 \(O(n^2mlogm)\) 的多元最短路。
对于这个费用流的算法, 其实实际测出来并不算快, 而且莫名的特别吃 \(O2\) , 不过比较稳定, 很少能被数据直接卡升天。 虽然还是有数据会将其卡飞, 不过毕竟少了。
最后, 这里再引出一个算法, 后面可能也会写, 就是原始对偶算法, 这个算法听说是网络流的最快算法
代码
废话了这么多, 我知道这应该是你们最期待的
开 \(O2\) 跑飞快, 不开 T 1点的 dinic 版本:
#include <cstdio>
#include <algorithm>
#include <queue>
using namespace std;
#define MAXN 50000
#define MAXM 500000
template <typename type, typename costtype, int maxn, int maxm, type inf>
class SSP {
public:
type sv[maxm];
private:
int H, T;
struct node {
int ed;
node *next;
int lo;
}se[maxm];
costtype scost[maxm];
node *sn[maxn];
node *si[maxn];
costtype pl[maxn];
costtype pe[maxn];
int tot;
bool vis[maxn];
private:
void push_edge (int u, int v, type value, costtype cost) {
node *p = &se[tot];
p->ed = v;
scost[tot] = cost;
sv[tot] = value;
p->lo = tot;
p->next = sn[u];
sn[u] = p;
tot ++;
}
type dfs (int now, type cap) {
vis[now] = 1;
if (now == T) {
vis[now] = 0;
return cap;
}
type num = 0;
for (; si[now]; si[now] = si[now]->next) {
if (pl[si[now]->ed] + scost[si[now]->lo] + pe[si[now]->ed] - pe[now] == pl[now] && sv[si[now]->lo] && !vis[si[now]->ed]) {
type sum = dfs (si[now]->ed, min (cap - num, sv[si[now]->lo]));
sv[si[now]->lo] -= sum;
sv[si[now]->lo ^ 1] += sum;
num += sum;
if (num == cap) {
vis[now] = 1;
return cap;
}
}
}
vis[now] = 0;
return num;
}
pair <type, costtype> dinic (int n) {
if (H == T) {
return make_pair(inf, inf);
}
priority_queue <pair <costtype, int>, vector <pair <costtype, int> >, greater<pair <costtype, int> > > q;
for (int i = 1; i <= n; i ++) {
pl[i] = inf;
vis[i] = 0;
}
pl[T] = 0;
q.push(make_pair (0, T));
while (!q.empty()) {
int l = q.top().second;
q.pop();
if (vis[l]) {
continue;
}
si[l] = sn[l];
vis[l] = 1;
for (node *i = sn[l]; i; i = i->next) {
if (sv[i->lo ^ 1] == 0) {
continue;
}
if (scost[i->lo ^ 1] + pl[l] + pe[l] - pe[i->ed] < pl[i->ed]) {
pl[i->ed] = scost[i->lo ^ 1] + pl[l] + pe[l] - pe[i->ed];
q.push(make_pair (pl[i->ed], i->ed));
}
}
}
if (pl[H] < inf) {
for (int i = 1; i <= n; i ++) {
vis[i] = 0;
}
type cap = dfs (H, inf);
costtype cost = cap * (pl[H] + pe[H]);
return make_pair (cap, cost);
}
return make_pair (0, 0);
}
public:
void clear () {
tot = 0;
for (int i = 1; i <= maxn; i ++) {
sn[i] = 0;
}
}
void push (int u, int v, type value, costtype cost) {
push_edge (u, v, value, cost);
push_edge (v, u, 0, -cost);
}
pair<type, costtype> maxflow (int h, int t, int n){
H = h;
T = t;
pair <type, costtype> ans = make_pair (0, 0);
queue <int> q;
for (int i = 1; i <= n; i ++) {
pe[i] = inf;
vis[i] = 0;
}
pe[T] = 0;
q.push (T);
vis[T] = 1;
while (!q.empty ()) {
int now = q.front ();
q.pop ();
vis[now] = 0;
for (node *i = sn[now]; i; i = i->next) {
if (sv[i->lo ^ 1] && pe[i->ed] > pe[now] + scost[i->lo ^ 1]) {
pe[i->ed] = pe[now] + scost[i->lo ^ 1];
if (!vis[i->ed]) {
vis[i->ed] = 1;
q.push (i->ed);
}
}
}
}
while (1) {
pair<type, costtype> num = dinic (n);
ans.first += num.first;
ans.second += num.second;
if (!num.first) {
break;
}
for (int i = 1; i <= n; i ++) {
pe[i] += pl[i];
}
}
return ans;
}
type limited_maxflow (int h, int t, int n, costtype limit) {
H = h;
T = t;
pair <type, costtype> ans = make_pair (0, 0);
queue <int> q;
for (int i = 1; i <= n; i ++) {
pe[i] = inf;
vis[i] = 0;
}
pe[T] = 0;
q.push (T);
vis[T] = 1;
while (!q.empty ()) {
int now = q.front ();
q.pop ();
vis[now] = 0;
for (node *i = sn[now]; i; i = i->next) {
if (sv[i->lo ^ 1] && pe[i->ed] > pe[now] + scost[i->lo ^ 1]) {
pe[i->ed] = pe[now] + scost[i->lo ^ 1];
if (!vis[i->ed]) {
vis[i->ed] = 1;
q.push (i->ed);
}
}
}
}
while (1) {
pair<type, costtype> num = dinic (n);
if (ans.second + num.second > limit) {
ans.first += (limit - ans.second) / (num.second / num.first);
break;
}
ans.first += num.first;
ans.second += num.second;
if (!num.first) {
break;
}
for (int i = 1; i <= n; i ++) {
pe[i] += pl[i];//注意, pl跑出来是最短路减去刚才的pe, 所以直接加
}
}
return ans.first;
}
};
SSP <long long, long long, MAXN + 5, MAXM + 5, 0x3f3f3f3f3f3f3f3f> ssp;
int main () {
int n, m, h, t;
scanf ("%d %d %d %d", &n, &m, &h, &t);
for (int i = 1; i <= m; i ++) {
int u, v;
long long w, c;
scanf ("%d %d %lld %lld", &u, &v, &w, &c);
ssp.push(u, v, w, c);
}
pair <long long, long long> ans = ssp.maxflow(h, t, n);
printf ("%lld %lld", ans.first, ans.second);
}
直接跑就能过, 开 \(O2\) 却没有那么快的 EK 版本:
#include <cstdio>
#include <algorithm>
#include <queue>
using namespace std;
#define MAXN 50000
#define MAXM 500000
template <typename type, typename costtype, int maxn, int maxm, type inf>
class SSP {
public:
type sv[maxm];
private:
int H, T;
struct node {
int ed;
node *next;
int lo;
}se[maxm];
struct Last {
int ed;
int nd;
}last[maxn];
costtype scost[maxm];
node *sn[maxn];
node *si[maxn];
costtype pl[maxn];
costtype pe[maxn];
int tot;
bool vis[maxn];
private:
void push_edge (int u, int v, type value, costtype cost) {
node *p = &se[tot];
p->ed = v;
scost[tot] = cost;
sv[tot] = value;
p->lo = tot;
p->next = sn[u];
sn[u] = p;
tot ++;
}
pair <type, costtype> EK (int n) {
if (H == T) {
return make_pair(inf, inf);
}
priority_queue <pair <costtype, int>, vector <pair <costtype, int> >, greater<pair <costtype, int> > > q;
for (int i = 1; i <= n; i ++) {
pl[i] = inf;
vis[i] = 0;
last[i].ed = 0;
last[i].nd = 0;
}
pl[T] = 0;
q.push(make_pair (0, T));
while (!q.empty()) {
int l = q.top().second;
q.pop();
if (vis[l]) {
continue;
}
si[l] = sn[l];
vis[l] = 1;
for (node *i = sn[l]; i; i = i->next) {
if (sv[i->lo ^ 1] == 0) {
continue;
}
if (scost[i->lo ^ 1] + pl[l] + pe[l] - pe[i->ed] < pl[i->ed]) {
last[i->ed].nd = l;
last[i->ed].ed = (i->lo ^ 1);
pl[i->ed] = scost[i->lo ^ 1] + pl[l] + pe[l] - pe[i->ed];
q.push(make_pair (pl[i->ed], i->ed));
}
}
}
if (pl[H] < inf) {
for (int i = 1; i <= n; i ++) {
vis[i] = 0;
}
type cap = inf;
for (int i = H; i != T; i = last[i].nd) {
cap = min (cap, sv[last[i].ed]);
}
for (int i = H; i != T; i = last[i].nd) {
sv[last[i].ed] -= cap;
sv[last[i].ed ^ 1] += cap;
}
costtype cost = cap * (pl[H] + pe[H]);
return make_pair (cap, cost);
}
return make_pair (0, 0);
}
public:
void clear () {
tot = 0;
for (int i = 1; i <= maxn; i ++) {
sn[i] = 0;
}
}
void push (int u, int v, type value, costtype cost) {
push_edge (u, v, value, cost);
push_edge (v, u, 0, -cost);
}
pair<type, costtype> maxflow (int h, int t, int n){
H = h;
T = t;
pair <type, costtype> ans = make_pair (0, 0);
queue <int> q;
for (int i = 1; i <= n; i ++) {
pe[i] = inf;
vis[i] = 0;
}
pe[T] = 0;
q.push (T);
vis[T] = 1;
while (!q.empty ()) {
int now = q.front ();
q.pop ();
vis[now] = 0;
for (node *i = sn[now]; i; i = i->next) {
if (sv[i->lo ^ 1] && pe[i->ed] > pe[now] + scost[i->lo ^ 1]) {
pe[i->ed] = pe[now] + scost[i->lo ^ 1];
if (!vis[i->ed]) {
vis[i->ed] = 1;
q.push (i->ed);
}
}
}
}
while (1) {
pair<type, costtype> num = EK (n);
ans.first += num.first;
ans.second += num.second;
if (!num.first) {
break;
}
for (int i = 1; i <= n; i ++) {
pe[i] += pl[i];
}
}
return ans;
}
type limited_maxflow (int h, int t, int n, costtype limit) {
H = h;
T = t;
pair <type, costtype> ans = make_pair (0, 0);
queue <int> q;
for (int i = 1; i <= n; i ++) {
pe[i] = inf;
vis[i] = 0;
}
pe[T] = 0;
q.push (T);
vis[T] = 1;
while (!q.empty ()) {
int now = q.front ();
q.pop ();
vis[now] = 0;
for (node *i = sn[now]; i; i = i->next) {
if (sv[i->lo ^ 1] && pe[i->ed] > pe[now] + scost[i->lo ^ 1]) {
pe[i->ed] = pe[now] + scost[i->lo ^ 1];
if (!vis[i->ed]) {
vis[i->ed] = 1;
q.push (i->ed);
}
}
}
}
while (1) {
pair<type, costtype> num = EK (n);
if (ans.second + num.second > limit) {
ans.first += (limit - ans.second) / (num.second / num.first);
break;
}
ans.first += num.first;
ans.second += num.second;
if (!num.first) {
break;
}
for (int i = 1; i <= n; i ++) {
pe[i] += pl[i];
}
}
return ans.first;
}
};
SSP <long long, long long, MAXN + 5, MAXM + 5, 0x3f3f3f3f3f3f3f3f> ssp;
int main () {
int n, m, h, t;
scanf ("%d %d %d %d", &n, &m, &h, &t);
for (int i = 1; i <= m; i ++) {
int u, v;
long long w, c;
scanf ("%d %d %lld %lld", &u, &v, &w, &c);
ssp.push(u, v, w, c);
}
pair <long long, long long> ans = ssp.maxflow(h, t, n);
printf ("%lld %lld", ans.first, ans.second);
}