ISAP学习笔记
ISAP是OI中求最大流的常用方法之一。相对于Dinic,ISAP的速度提升了很多,但编码复杂度也上升了不少。
约定
采用邻接表存储图,对于每条弧,增加一条容量为0的逆向边。
d数组代表每个点到终点的距离。
引入
求解最大流有一种基础方法,每次在残量网络中dfs找到任意路径增广,直到不存在这样的道路为止。这就是最短路增广算法。但是这样的方法效率太低。
所以我们比较容易想到的一种方法是每次沿最短路增广。这就是最短路增广算法。而ISAP则是针对最短路增广的一种改进(Improved Shortest Augmenting Path)。
算法概述
ISAP基于这样一个事实,每次增广之后,残量网络中的最短路不会变短,但如果残量网络中存在另一条最短路,那么可以继续增广,直到遇到死路为止,所以每次增广后再次bfs找最短路不是必须的(而Dinic则需要每次增广后bfs)。
那么遇到死路后怎么办呢?对于节点u,如果它周围的任意节点v都不满足d[v]=d[u]+1,那么记min为d[v]的最小值,d[u]=d[v]+1。并退回上次经过的节点继续寻找增广路。
但是对于一种特殊情况,如果d[u]=k的节点已完全消失,那么可以直接终止程序,因为此时d值>k的点以无法经过d值等于k的点到达t(d[t]=0),也就是说s必定无法到达k,这就是大名鼎鼎的gap优化。
说到这里,算法流程已大致清楚了,但还有一些程序实现的细节需要注意。
实现
#include <bits/stdc++.h>
using namespace std;
const int maxn = 10005, maxm = 100005;
const int Inf = 1 << 30;
struct edge {
int from, to, next;
int cap, flow;
}e[maxm * 2];
int h[maxn], tot;
int n, m, s, t, flow;
int cur[maxn]; //当前弧编号
int p[maxn]; //增广路上一条边编号(第i点的入边)
int num[maxn]; //和t最短路距离为i的点数量
int d[maxn]; //残量网络中到t的最短距离
inline void add(int u, int v, int w)
{
e[tot++] = (edge) {u, v, h[u], w, 0};
h[u] = tot - 1;
}
inline int gi()
{
char c = getchar();
while(c < '0' || c > '9') c = getchar();
int sum = 0;
while('0' <= c && c <= '9') sum = sum * 10 + c - 48, c = getchar();
return sum;
}
queue<int> que;
//预处理,反向bfs构造d数组
bool bfs()
{
memset(d, 63, sizeof(int) * (n + 1));
d[t] = 0;
que.push(t);
while(!que.empty()) {
int u = que.front(); que.pop();
++num[d[u]];
for(int v, i = h[u]; ~i; i = e[i].next) {
v = e[i].to;
if(e[i].cap == 0 && d[v] > d[u] + 1) {
d[v] = d[u] + 1;
que.push(v);
}
}
}
return d[s] < n;
}
//增广
int augment()
{
int u = t, Min = Inf;
//从汇点到源点通过p追踪增广路径,Min为一路上的最小流量
while(u != s) {
edge E = e[p[u]];
Min = min(Min, E.cap - E.flow);
u = e[p[u]].from;
}
u = t;
while(u != s) {
e[p[u]].flow += Min;
e[p[u] ^ 1].flow -= Min;
u = e[p[u]].from;
}
return Min;
}
int main()
{
n = gi(); m = gi(); s = gi(); t = gi();
memset(h, -1, sizeof(int) * (n + 1));
for(int u, v, w, i = 1; i <= m; ++i) {
u = gi(); v = gi(); w = gi();
add(u, v, w); add(v, u, 0);
}
if(!bfs()) printf("0\n"), exit(0);
int u = s;
memcpy(cur, h, sizeof(int) * (n + 1));
while(d[s] < n) {
if(u == t) {
flow += augment();
u = s;
}
bool advanced = false;
for(int i = cur[u]; ~i; i = e[i].next) {
int v = e[i].to;
if(e[i].cap > e[i].flow && d[u] == d[v] + 1) {
advanced = true;
p[v] = i; cur[u] = i; u = v;
break;
}
}
if(!advanced) { //retreat操作
int Min = n - 1;
for(int i = h[u]; ~i; i = e[i].next)
if(e[i].cap > e[i].flow) Min = min(Min, d[e[i].to]);
if(--num[d[u]] == 0) break; //gap优化
num[d[u] = Min + 1]++;
cur[u] = h[u];
if(u != s) u = e[p[u]].from;
}
}
printf("%d\n", flow);
return 0;
}