图论算法----网络流----最大流sap算法

一、相关概念

1、流网络

流网络G=(V, E)是一个有向图，其中每条边(u,v)均有一非负容量c(u, v)≥0。

流网络中有两个特殊的顶点: 源点s和汇点t。

假定每个顶点都处于从源点到汇点的某条路径上，就是说，对每个顶点v，存在一条路径s→v→t

G为连通图，且|E| ≥|V|-1

2、流

边的流是一个实值函数f，满足下列三个性质：

容量限制：对所有u,v∈V， 要求f(u, v) ≤c(u, v)

    理解：流量不会超过边的容量

反对称性：对所有u,v∈V， 要求f(u, v) = -f(v, u)

    理解：一个方向的流是其反方向流的相反数

流守恒性：对所有u∈V-{s, t}，要求：

    理解：进入点u的总流量=离开点u的总流量

3、网络的流

网络的流f定义为：

即从源点出发的总流表示网络的流。

最大流问题：给出一个源点为s，汇点为t的流网络G，希望找出从s到t的最大值流。

4、残留网络

边的残留容量：r(u,v)=c(u,v)-f(u,v);

残留网络：给定一流网络G=(V,E)和流f，由f导出的G的残留网络是 G f (V,E f )

其中

这就是说，在残留网络中，每条边（称为残留边）能够容纳一个严格为正的流网络

当0<f(u,v)<c(u,v)，由r(u,v)=c(u,v)-f(u,v)>0

此时边(u,v)在残留网络中，边(v,u)也在残留网络。

5、增广路径

这个名字听起来挺熟悉的，我们在二分图匹配里面见过，但是它的定义与二分图匹配中增广路径有一些差别。

增广路径P为残留网络中从s到t的一条简单路径。

增广路径的残留容量是：

沿着路径增广是沿着路径的每条边发送δ(P)单位的流。相应地，修改流的值和残留容量。

f=f+δ(P);

r(u,v)=r(u,v)-δ(P);

r(v,u)=r(v,u)+δ(P);  for(u,v)∈P;

6、Ford-Fulkerson最大流算法

基本结构：

int Ford_Fulkerson()

{

f=0;

创建残留网络 G(f);

while (在G(f)中存在从 s 到 t 的有向路径)

{

令 P是在G(f)中从 s 到 t  的一条路径

Δ =δ(P);

沿着 P发送Δ单位的流;

更新P上的边的残留容量;

f =f + Δ;

}

return f;  //f是最大流

}

7、最大流和最小割

（1）、割集的定义

在网络G=(V,E)中，割集(S,T)是把V分成两个不相交的子集S和T=V-S的划分，使得s∈S且t∈T

（2）、割集的容量（值）

即穿越割的从S到T的边的容量之和。（不包括反向边）

（3）最大流最小切割定理

最大流值是一个割集的最小值.

8、Ford-Fulkerson算法的缺点及改进

如果M为极大值，在找增广路径时不幸按红色所示线路寻找，那么......

就要进行2M次增广才可以算出最大流，岂不是白白浪费了许多时间，

如果我们一开始就按照s->1->t和s->2->t，只需进行两次增广就可以算出最大流了

于是就有了以下两种改进的方法

最大增广路径算法（dinic）

    令P是在G(f)中从s到t使得δ(P)最大的路径

最短增广路径算法（sap）

    令P是在G(f)中从s到t拥有最少边数的路径 

我们一般选择第二种方法： 最短增广路径算法（sap）

9、sap算法的实现

先要增加一个d[]数组，d[v]表示v到t距离标号。

距离标号是一个函数：d: V →Z+距离标号被称为是有效，如果它满足以下:

d(t)=0

d(i) ≤d(j)+1对每个 (i,j)∈G(f)

如果d(i)=d(j)+1，则边 (i,j) ∈G(f)是可进入的

在实践不用初始化d[]数组，因为它在运行sap的时候，它的值就会自己改变。

还要用一个vd[]数组，vd[i]表示距离标号为i的顶点个数为vd[i]。

以下是邻接矩阵版的sap

#define INF 1000000000
int c[MAXN][MAXN],d[MAXN],vd[MAXN];//c[i][j]表示i到j的距离
int n,flow;
int aug(int i,int augco)
{
    int j,augc=augco,mind=n-1,delta;
    if(i==n)
        return augco;
    for(j=1;j<=n;j++){
        if(c[i][j]>0){
            if(d[i]==d[j]+1){
                delta=min(augc,c[i][j]);
                delta=aug(j,delta);
                c[i][j]-=delta;
                c[j][i]+=delta;         
                augc-=delta;
                if(d[1]>=n)
                    return augco-augc;
                if(augc==0) break;
            }
            if(mind>d[j]) mind=d[j];
        }
    }
    if(augco==augc){
        vd[d[i]]--;
        if(vd[d[i]]==0)
            d[1]=n;
        d[i]=mind+1;
        vd[d[i]]++;
    }
    return augco-augc;
}
void sap()
{
    memset(d,0,sizeof(d));
    memset(vd,0,sizeof(vd));
    vd[0]=n;
    while(d[1]<n)
        flow+=aug(1,INF);
}

以下是邻接表版的sap

#define INF 1000000000
struct enode{
    int v,c;
    enode* next;
    enode* back;
}edge[450005];
typedef enode* elist;
elist adj[55005],ecnt;
int h[55005],vh[55005],s,t,flow;
void addedge(int x,int y,int cap)
{
    elist p;
    p=++ecnt;
    p->v=y;
    p->c=cap;
    p->next=adj[x];
    adj[x]=p;
    p->back=ecnt+1;
     
    p=++ecnt;
    p->v=x;
    p->c=0;
    p->next=adj[y];
    adj[y]=p;
    p->back=ecnt-1;
}
int sap(int i,int delt)
{
    int tmp,minh=t-1,ret=0;
    if(i==t)
        return delt;
    for(elist p=adj[i];p!=NULL;p=p->next){
        int j=p->v,cap=p->c;
        if(cap>0){
            if(h[j]+1==h[i]){
                int k=min(cap,delt);
                tmp=sap(j,k);
                delt-=tmp;
                p->c-=tmp;
                p->back->c+=tmp;
                ret+=tmp;
                if(h[s]>=t) return ret;
                if(delt==0) break;
            }
            if(h[j]<minh) minh=h[j];
        }
    }
    if(ret==0){
        vh[h[i]]--;
        if(vh[h[i]]==0)
            h[s]=t;
        h[i]=minh+1;
        vh[h[i]]++;
    }
    return ret;
}
void f()
{
    memset(h,0,sizeof(h));
    memset(vh,0,sizeof(vh));
    vh[0]=t;
    while(h[s]<t)
        flow+=sap(s,INF);
}