ISAP算法

为什么叫ISAP

　　ISAP(Improved Shortest Augment Path)：改进的最短增广路，属于增广路算法

算法

　　Dinic算法中，我们每次都需要BFS出层次图，而在ISAP中，我们只需要初始化时BFS出层次图（从(T)向(S)进行），然后在増广的过程中维护标号(d)（就是到汇点(T)的距离）。

　　増广的过程和Dinic类似，搜索时沿着(d_v=d_u-1)的边走，走出来的肯定是最短路。不同的是，増广完后，并不会马上对标号(d)进行更新，而是増广到没办法再増广时再更新。具体来说，就是不断BFS找增广路进行増广。如果増广过程中发现没有从(S)到(u)再到(T)的增广路，就要对(d_u)进行更新：令(d_u=min_{(u,v)in E'}d_v+1)，其中(E')为残留网络的边集。这样可以保证下一次找到的増广路还是最短路。如果(u)没有出边，就令(d_u=n)，这样以后増广时就不会经过这个点。更新标号后要往回退一步。

　　那什么时候停止呢？当然是不存在(S)到(T)的最短路的时候。因为最短路是无环的，所以如果(S)到(T)有路径那么一定有(d_Sleq n-1)。

　　时间复杂度是(O(n^2m))，然而一般不会到这个上界。这就是为什么増广路算法的时间复杂度被称为玄学。

步骤

• 初始化，BFS出层次图
• 每次沿着(d_v=d_u-1)的边进行搜索，找到一条増广路就増广
• 如果没有从(u)出发到(T)的路径，就更新(d_u)，往回退一步
• 当(d_Sgeq n)时停止

一些优化

• 多次増广。每次找到一条増广路后可以往回退，尝试着找第二条増广路。

• DFS优化。这样每次找到増广路后就不用把(S)到(T)路径的信息全部重新更新一遍了。

• GAP优化。如果不存在标号为(x)的点，说明残留网络出现了断层。修改标号时如果发现当前点是最后一个(h_u=x)的点，就直接退出。可以维护一个(num)数组，(num_i)表示标号为(i)的点有多少个。

• 当前弧优化。如果一个点的标号没有被修改，那么这个点已经遍历过的边就不用再遍历了，因为不可能通过这些边再找到増广路。修改一个点的标号时顺便把遍历过的边设为空。实现时维护没遍历过得第一条边。

• 标号(d)的修改有连续性，即我们不需要找到最小的(d_v)，而是直接把(d_u+1)。

• 如果流量用完就退出。

代码

#include<cstdio>
#include<cstring>
#include<algorithm>
#include<cstdlib>
#include<ctime>
#include<utility>
using namespace std;
typedef long long ll;
typedef unsigned long long ull;
typedef pair<int,int> pii;
int h[1000010];
int v[8000010];
ll c[8000010];
int t[8000010];
int cnt=0;
void add(int x,int y,ll z)
{
	cnt++;
	v[cnt]=y;
	c[cnt]=z;
	t[cnt]=h[x];
	h[x]=cnt;
	cnt++;
	v[cnt]=x;
	c[cnt]=0;
	t[cnt]=h[y];
	h[y]=cnt;
}
int op(int x)
{
	return ((x-1)^1)+1;
}
int S,T;
int d[1000010];
int to[1000010];
int q[1000010];
int num[1000010];
int cur[1000010];
int n,m;
void bfs()
{
	memset(d,-1,sizeof d);
	d[T]=0;
	int head=1,tail=0;
	q[++tail]=T;
	while(tail>=head)
	{
		int i,x=q[head++];
		num[d[x]]++;
		for(i=h[x];i;i=t[i])
			if(c[op(i)]&&d[v[i]]==-1)
			{
				d[v[i]]=d[x]+1;
				q[++tail]=v[i];
			}
	}
}
int stop;
ll dfs(int x,ll flow)
{
	if(x==T)
		return flow;
	ll s=0,u;
	int &i=cur[x];
	for(;i;i=t[i])
		if(d[v[i]]==d[x]-1&&c[i])
		{
			u=dfs(v[i],min(c[i],flow));
			s+=u;
			flow-=u;
			c[i]-=u;
			c[op(i)]+=u;
			if(stop)
				return s;
			if(!flow)
				return s;
		}
	num[d[x]]--;
	if(!num[d[x]])
	{
		d[S]=n;
		stop=1;
		return s;
	}
	d[x]++;
	num[d[x]]++;
	cur[x]=h[x];
	return s;
}
ll maxflow()
{
	stop=0;
	bfs();
	ll ans=0;
	memcpy(cur,h,sizeof h);
	while(d[S]<=n-1)
		ans+=dfs(S,0x7fffffffffffffffll);
	return ans;
}
int main()
{
#ifdef DEBUG
	freopen("isap.in","r",stdin);
	freopen("isap.out","w",stdout);
#endif
	scanf("%d%d%d%d",&n,&m,&S,&T);
	int i,x,y;
	ll z;
	for(i=1;i<=m;i++)
	{
		scanf("%d%d%lld",&x,&y,&z);
		add(x,y,z);
	}
	ll ans=maxflow();
	printf("%lld
",ans);
	return 0;
}