##网络流##-CSDN博客

本文链接：https://blog.csdn.net/chm517/article/details/37879257

HLPP参考：http://blog.sina.com.cn/s/blog_60707c0f01010oqt.html

ISAP参考：http://blog.sina.com.cn/s/blog_6165e6920100ea08.html

参考：http://blog.sina.com.cn/s/blog_476a25110100mag7.html

http://wenku.baidu.com/link?url=Q3y_0x9BaO1L1W5Kr3AJzgCUVJVKNjG8UGw2nt44d6HDWTT3QLWygZ8IT_lb-V6GNjYoV4Fh74bsB88uTMeJA4xWBnbAyNfUr9VRCo-rCFS

//algo1:
//Ford-Fulkerson
//O(n*m*U)
#define N 115
struct arc{
    int c,f;
}e[N][N];
int flag[N],al[N],pre[N];
int n,m,np,nc;
int ford(int s,int t,int n)
{   while (1)
        {   memset(flag,255,sizeof(flag));
            memset(pre,255,sizeof(pre));
            memset(al,255,sizeof(al));
            queue<int > q;
            flag[s]=0; pre[s]=-1; al[s]=inf;
            q.push(s);
            while (!q.empty()&&flag[t]==-1)
                {
                    int v=q.front(); q.pop();
                    for (int i=0;i<n;i++)
                    if (flag[i]==-1)
                        {
                            if (e[v][i].f<e[v][i].c)
                                {
                                    flag[i]=0; pre[i]=v;
                                    al[i]=min(al[v],e[v][i].c-e[v][i].f);
                                    q.push(i);
                                }
                            else if (e[i][v].f>0)
                                {   flag[i]=0; pre[i]=-v;
                                    al[i]=min(al[v],e[i][v].f);
                                    q.push(i);
                                }
                        }
                    flag[v]=1;
                }

            if (flag[t]==-1||al[t]==0) break;
            int k1=t,k2=abs(pre[k1]),alph=al[t];
            while (1)
                {   if (e[k2][k1].c>0)
                         e[k2][k1].f+=alph;
                    else e[k1][k2].f-=alph;
                    if (k2==s) break;
                    k1=k2; k2=abs(pre[k1]);
                }
        }
    int ans=0;
    for (int i=0;i<n;i++)
        ans+=e[s][i].f;
    return ans;
}
void doit()
{   for (int i=0;i<N;i++)
    for (int j=0;j<N;j++)
        e[i][j].c=e[i][j].f=0;
    //add edges
    printf("%d\n",ford(s,t,N));
}

//algo2:  
//Dinic
//O(n*n*m)
struct dinic{
    int g[N],d[N],G[N];
    int np[M],to[M],cap[M];
    int n,s,t,e;
    void AddEdge(int x,int y,int z,int z2=0)
         {  to[e]=y;cap[e]=z;
            np[e]=g[x];g[x]=e++;

            to[e]=x;cap[e]=z2;
            np[e]=g[y];g[y]=e++;
         }
    void init(int n)
         {  this->n=n;e=0;
            memset(g,255,sizeof(g));
         }
    int maxflow(int s,int t)
         {  this->s=s; this->t=t;
            int ans=0;
            while (bfs())
                {for(int i=0;i<n;i++)G[i]=g[i];
                 ans+=dfs(s,inf);
                }
            return ans;
         }
    bool bfs()//建立层次图
    {
        memset(d,-1,sizeof(d));
        queue<int > q;
        q.push(s);
        d[s]=0;
        while(!q.empty())
        {
            int x=q.front(); q.pop();
            for(int i=g[x];~i;i=np[i])
                 {
                     int y=to[i];
                     if(cap[i]&&d[y]==-1)
                     {  d[y]=d[x]+1;
                        if(y==t)return 1;
                        q.push(y);
                     }
                 }
        }
        return 0;
    }

   int dfs(int x,int sum)
    {
        if (x==t|| sum==0) return sum;
        int ss=sum,ret;
        for(int &i=G[x];~i;i=np[i]) //当前弧优化
        {
            int y=to[i];
            if(d[y]==d[x]+1&&(ret=dfs(y,min(sum,cap[i]))))
                {
                    cap[i]-=ret;
                    cap[i^1]+=ret;
                    sum-=ret;
                    if (!sum) break;
                }
        }
        return ss-sum;
    }
}it;

SAP的两个优化：

GAP优化——如果一次重标号时，出现距离断层，则可以证明ST无可行流，此时则可以直接退出算法。
当前弧优化——为了使每次找增广路的时间变成均摊O(V)，还有一个重要的优化是对于每个点保存“当前弧”：初始时当前弧是邻接表的第一条弧；在邻接表中查找时从当前弧开始查找，找到了一条允许弧，就把这条弧设为当前弧；改变距离标号时，把当前弧重新设为邻接表的第一条弧。

//algo3_1:
//SAP
//O(n*n*m)
struct SAP{
    int g[N],d[N],G[N],pre[N],numd[N];
    int np[M],to[M],cap[M];
    int n,e,s,t;
    void AddEdge(int x,int y,int z,int z2=0)
        {  to[e]=y;cap[e]=z;
           np[e]=g[x];g[x]=e++;

           to[e]=x;cap[e]=z2;
           np[e]=g[y];g[y]=e++;
        }
    void init(int n)
        {   this->n=n;
            e=0;memset(g,255,sizeof(g));
        }
    void bfs()//先用广度优先算出每个点的d值
        {   memset(numd,0,sizeof(numd));
            for(int i=0; i<n; i++)
                numd[d[i]=n]++;
            d[t]=0; numd[n]--; numd[0]++;
            queue<int> Q;Q.push(t);
            while(!Q.empty())
            {   int v=Q.front();Q.pop();
                int i=g[v];
                while(i!=-1)
                    {
                        int u=to[i];
                        if(d[u]<n){i=np[i];continue;}
                        d[u]=d[v]+1;
                        numd[n]--;numd[d[u]]++;
                        Q.push(u);
                        i=np[i];
                    }
            }
        }

    int maxflow(int s,int t)
    {   this->s=s; this->t=t;
        for(int i=0;i<n;i++)
            G[i]=g[i];
        int max_flow=0;
        bfs();
        int u=s ;//从源点搜一条到汇点的增广路
        while(d[s]<n)//就算所有的点连成一条线源点的d值也是最多是n-1
        {
            if(u==t)//如果找到一条增广路径
            {
                    int cur_flow=inf,neck;//找到那条瓶颈边
                    for(int from=s; from!=t; from=to[G[from]])
                        if(cur_flow>cap[G[from]])
                        {   neck=from;cur_flow=cap[G[from]];}


                for(int from=s; from!=t; from=to[G[from]]) //修改增广路上的边的容量
                    {
                    int tmp=G[from];
                    cap[tmp]-=cur_flow;
                    cap[tmp^1]+=cur_flow;
                    }
                max_flow+=cur_flow;//累加计算最大流
                u=neck;//下一次搜索直接从瓶颈边的前一个节点搜起
            }

            int i;
            for(i=G[u];~i; i=np[i]) //从当前点开始找一条允许弧
            if(cap[i]&&d[u]==d[to[i]]+1)//如果找到跳出循环
                break;

            if(i!=-1)//找到一条允许弧
            {
                G[u]=i;//从点u出发的允许弧的地址
                pre[to[i]]=u;//允许弧上下一个点的前驱为u
                u=to[i];//u变成下一个点继续搜直到搜出一条增广路
            }
            else  //如果没有搜到允许弧
            {
                numd[d[u]]--; //d[u]将被修改所以numd[d[u]]减一
                if(!numd[d[u]]) break;  //如果没有点的d值为d[u]则不可能再搜到增广路结束搜索
                G[u]=g[u];  //当前点的允许弧为第一条边
                int tmp=n;
                for(int j=g[u]; ~j; j=np[j]) //搜与u相连的点中d值最小的
                    if(cap[j]&&tmp>d[to[j]])tmp=d[to[j]];

                d[u]=tmp+1; //修改d[u]
                numd[d[u]]++;
                if(u!=s)u=pre[u];//从u的前驱搜，因为从u没有搜到允许弧
            }
        }
    return max_flow;
    }
}it;

//algo3_2:
//SAP
//O(n*n*m)
#define N 105
int n,m,np,nc,st,ed;
struct SAP{
    int g[N][N];//直接在原图上操作了
    int d[N],gap[N];
    int n,s,t;
    void AddEdge(int x,int y,int z)
        {
            g[x][y]=z;
        }
    void init(int n)
        {   this->n=n;
            memset(g,0,sizeof(g));
            memset(gap,0,sizeof(gap));
            memset(d,0,sizeof(d));
        }
    int maxflow(int s,int t)
        {  this->s=s;this->t=t;
           int ans=0;
           for(gap[0]=n;d[s]<n;)
                ans+=sap(s,inf);
           return ans;
        }
    int sap(int u,int flow)
        {
            if(u==t) return flow;
            int res=flow;
            for(int i=0;i<n;++i)
                if(g[u][i] && d[u]==d[i]+1)
                {
                    int tmp=sap(i,min(res,g[u][i]));
                    g[u][i]-=tmp,g[i][u]+=tmp;
                    if(!(res-=tmp))return flow;
                }
            if(!--gap[d[u]])d[s]=n;
            ++gap[++d[u]];
            return flow-res;
        }
}it;

//algo4:  
//ISAP  
//O(n*n*m)  
//ISAP 的一个所谓 GAP 优化:  
//由于从 s 到 t 的一条最短路径的顶点距离标号单调递减，且相邻顶点标号差严格等于1，  
//因此可以预见如果在当前网络中距离标号为 k (0 <= k < n) 的顶点数为 0，  
//那么可以知道一定不存在一条从 s 到 t 的增广路径，此时可直接跳出主循环。  
struct ISAP  
{   struct Edge  
    {  
        int from,to,cap,flow;  
        Edge(){}  
        Edge(int a,int b,int c,int d):from(a),to(b),cap(c),flow(d){}  
    };  
    int n,m,s,t;//结点数，边数(含反向弧)，源点，汇点  
    vector<Edge> edges;//边表，edges[e]&edges[e^1]互为反向弧  
    vector<int> G[N];//邻接表，G[i][j]表示结点i的第j条边在e数组中的序号  
    bool vis[N];//BFS使用  
    int d[N];//从起点到i的距离  
    int cur[N];//当前弧下标  
    int p[N];//可增广路上的上一条弧  
    int num[N];//距离标号计数  
  
    void AddEdge(int from,int to,int cap)//重边不影响  
    {  
        edges.push_back(Edge(from,to,cap,0));  
        edges.push_back(Edge(to,from,0,0));//容量为0,表示反向弧  
        m=edges.size();  
        G[from].push_back(m-2);  
        G[to].push_back(m-1);  
    }  
  
    void init(int n)  
    {  
        this->n=n;  
        for(int i=0;i<n;++i) G[i].clear();  
        edges.clear();  
    }  
  
    void BFS()//反向  
    {   for (int i=0;i<n;i++)
            d[i]=n+10; //!!!
        memset(vis,0,sizeof(vis));
        queue<int> Q;
        Q.push(t);
        d[t]=0;
        vis[t]=1;
        while(!Q.empty())
        {
            int x=Q.front();
            Q.pop();
            for(int i=0; i<G[x].size(); ++i)
            {
                Edge& e=edges[G[x][i]^1];
                if(!vis[e.from]&&e.cap>e.flow)
                {
                    vis[e.from]=1;
                    d[e.from]=d[x]+1;
                    Q.push(e.from);
                }
            }
        }
    }  
  
    int Augment()  
    {  
        int x=t,a=INF;  
        while(x!=s)  
        {  
            Edge& e=edges[p[x]];  
            a=min(a,e.cap-e.flow);  
            x=edges[p[x]].from;  
        }  
        x=t;  
        while(x!=s)  
        {  
            edges[p[x]].flow+=a;  
            edges[p[x]^1].flow-=a;  
            x=edges[p[x]].from;  
        }  
        return a;  
    }  
  
    int maxflow(int s,int t)//结点数  
    {  
        this->s=s,this->t=t;  
        int flow=0;  
        BFS();  
        memset(num,0,sizeof(num));  
        for(int i=0;i<n;++i) if (d[i]!=n+10)++num[d[i]];//!!! 
        int x=s;  
        memset(cur,0,sizeof(cur));  
        while(d[s]<n)  
        {  
            if(x==t)  
            {  
                flow+=Augment();  
                x=s;  
            }  
            int ok=0;  
            for(int i=cur[x];i<G[x].size();++i)  
            {  
                Edge& e=edges[G[x][i]];  
                if(e.cap>e.flow&&d[x]==d[e.to]+1)//Advance  
                {  
                    ok=1;  
                    p[e.to]=G[x][i];  
                    cur[x]=i;  
                    x=e.to;  
                    break;  
                }  
            }  
            if(!ok)//Retreat  
            {  
                int m=n-1;  
                for(int i=0;i<G[x].size();++i)  
                {  
                    Edge& e=edges[G[x][i]];  
                    if(e.cap>e.flow) m=min(m,d[e.to]);  
                }  
                if(--num[d[x]]==0) break;//gap优化  
                num[d[x]=m+1]++;  
                cur[x]=0;  
                if(x!=s) x=edges[p[x]].from;  
            }  
        }  
        return flow;  
    }  
}it;

algorithm Improved-Shortest-Augmenting-Path
  f <-- 0
  从终点 t 开始进行一遍反向 BFS 求得所有顶点的起始距离标号 d(i)
  i <-- s
  while d(s) < n do
      if i = t then    // 找到增广路
          Augment
          i <-- s      // 从源点 s 开始下次寻找
      if Gf 包含从 i 出发的一条允许弧 (i,j)
          then Advance(i)
         else Retreat(i)    // 没有从 i 出发的允许弧则回退
  return f

procedure Advance(i)
  设 (i,j) 为从 i 出发的一条允许弧
  pi(j) <-- i    // 保存一条反向路径，为回退时准备
  i <-- j        // 前进一步，使 j 成为当前结点

procedure Retreat(i)
  d(i) <-- 1 + min{d(j):(i,j)属于残量网络Gf}    // 称为重标号的操作
  if i != s
     then i <-- pi(i)    // 回退一步

procedure Augment
  pi 中记录为当前找到的增广路 P
  delta <-- min{rij:(i,j)属于P}
  沿路径 P 增广 delta 的流量
  更新残量网络 Gf

~~~~~~~~~~~~~~~~~~~~~~~~~~~增广路与预流推进的分割线~~~~~~~~~~~~~~~~~~~~~~~~~~~~

定义e(u)为节点u上的余流，实际上我们并不需要计算各个h(v)，只需将h(s)标为n，其余标为0，算法也能逐步计算出各个点的h。网络流中的节点(除s, t)仅需要满足流入量不小于流出量。其中流入量大于流出量的节点，我们称之为活动节点。预流推进算法就是不断地将活动结点，变为非活动结点，使得预流成为可行流。具体分为两个方法：
①压入操作，选出一个活动顶点u。并依次判断残余网络G中每条边(u, v)，若h(u) = h(v) + 1，则顺着这里条边推流，流量变化量df(u,v)=min(e[u],cf(u,v))，同时更新相应的e(u)，e(v)；
②重标记操作，选出一个活动结点。如果对于残余网络G中每条边(u, v)，都有h(v)大于或等于h(u)，则需要对u进行重新标号，h(u) = min{h(v) + 1}；
具体实现又可以根据对活跃点的推进过程分为：
①一般的预流推进算法（Generic Preflow-Push Algorithm），在残余网络中，维护一个预流，不断对活跃点执行push操作，或者relable操作来重新调整这个预流，直到不能操作。时间复杂度为O(V*E*E)；
②先进先出预流推进算法，在残余网络中，以先进先出队列维护活跃点，时间复杂度为O(V*V*V)；
③最高标号预流推进算法（Highest Preflow-Push Algorithm），在残余网络图中，每次检查最高标号的活跃点，需要用到优先队列，时间复杂度为O(V*V*E^0.5)；

//algo5:
//HLPP
//O(n*n*m^(1/2))
#define N 105
int d[N];
struct HLPP{
    int e[N];
    int g[N][N];
    struct node{
         int index;
         friend bool operator < (const node& a, const node& b)
            {return d[a.index] > d[b.index];}
    };
    int n,s,t;
    void Relabel(int u)
        {
            int m1=inf;
            for(int i=0; i<n; ++i)
                if(g[u][i] && d[i] < m1)m1 =d[i];
            d[u]=m1+1;
        }
    void AddEdge(int x,int y,int z)
        {
            g[x][y]+=z;
        }
    void init(int n)
        {   this->n=n;
            memset(g,0,sizeof(g));
        }
    int maxflow(int s,int t)
        {  node no;
           priority_queue<node> Q;
           memset(e,0,sizeof(e));
           memset(d,0,sizeof(d));
           d[s]=n;
           for(int i=0;i<n;++i)
           if(g[s][i])
                {
                    int df=g[s][i];
                    e[s] -= df;
                    e[i] += df;
                    g[s][i] -= df;
                    g[i][s] += df;
                    no.index = i;
                    if(i!=s&& i!=t) Q.push(no);
                }
           while(!Q.empty())
                {
                 no=Q.top();Q.pop();
                 int u=no.index;
                 while(e[u])
                     {  int i;
                        for(i = 0; i < n; ++i)
                            if(g[u][i] && d[u] > d[i])break;
                        if(i == n) Relabel(u);
                        for(int i=0; i<n; ++i)
                             if(d[u]==d[i] + 1 && g[u][i])
                             {
                                 int df = min(e[u], g[u][i]);
                                 g[u][i] -= df;
                                 g[i][u] += df;
                                 e[u] -= df;
                                 e[i] += df;;
                                 no.index = i;
                                 if(e[i] > 0&&i!=s &&i!=t) Q.push(no);
                             }
                     }
                }
            return e[t];
        }
}it;

~~~~~~~~~~~~~~~~~~~~~~~~~~~最大流与最小费用最大流间的分割线~~~~~~~~~~~~~~~~~~~~~~~~~~~~

struct MCMF
{   struct Edge
    {
        int from,to,cap,flow,cost;
        Edge(){}
        Edge(int a,int b,int c,int d,int e):from(a),to(b),cap(c),flow(d),cost(e){}
    };

    int n,m,s,t;
    vector<Edge> edges;
    vector<int> G[N];
    int inq[N],d[N],p[N],a[N];
    queue<int> Q;

    void init(int n)
    {
        this->n=n;
        for(int i=0;i<n;++i) G[i].clear();
        edges.clear();
    }

    void AddEdge(int from,int to,int cap,int cost)
    {   //printf("%d %d %d %d\n",from,to,cap,cost);
        edges.push_back(Edge(from,to,cap,0,cost));
        edges.push_back(Edge(to,from,0,0,-cost));
        m=edges.size();
        G[from].push_back(m-2);
        G[to].push_back(m-1);
    }

    bool SPFA(int s,int t,int & flow,int & cost)
    {
        memset(d,0x3f,sizeof(d));
        memset(inq,0,sizeof(inq));
        d[s]=0;inq[s]=1;p[s]=0;a[s]=inf;
        Q.push(s);
        while(!Q.empty())
        {
            int u=Q.front();Q.pop();
            inq[u]=0;
            for(int i=0;i<G[u].size();++i)
            {
                Edge& e=edges[G[u][i]];
                if(e.cap>e.flow&&d[e.to]>d[u]+e.cost)
                {
                    d[e.to]=d[u]+e.cost;
                    p[e.to]=G[u][i];
                    a[e.to]=min(a[u],e.cap-e.flow);
                    if(!inq[e.to]) {Q.push(e.to);inq[e.to]=1;}
                }
            }
        }
        if(d[t]==inf) return false;
        flow+=a[t];//固定流量k，检查，若flow+a>=k，只增加k-flow个单位
        cost+=d[t]*a[t];
        int u=t;
        while(u!=s)
        {
            edges[p[u]].flow+=a[t];
            edges[p[u]^1].flow-=a[t];
            u=edges[p[u]].from;
        }
        return true;
    }
    int mincost(int s,int t)
    {
        int flow=0,cost=0;
        while(SPFA(s,t,flow,cost));
        return cost;
    }
}it;