ISAP算法

为什么叫ISAP

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

算法

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

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

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

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

步骤

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

一些优化

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

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

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

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

• 标号d的修改有连续性，即我们不需要找到最小的dv，而是直接把du+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\n",ans);    return 0;}