算法学习笔记(8.1): 网络最大流算法 EK, Dinic, ISAP

算法,学习,笔记,网络,最大,ek,dinic,isap · 浏览次数 : 196

小编点评

```C++ template inline void read(T &t, Args&... args) { read(t), read(args...); } typedef long long Data;using namespace std;const int N = 207, M = 5007;struct Edge { int to; size_t rev; // 反边的位置,用int也没问题 Data flow; Edge(int to, size_t rev, Data f) : to(to), rev(rev), flow(f) {}};class ISAP {public: int n, m, s, t; vector> dep; int q[N * 2], gap[N * 2]; // vector< vector<Edge> > v; vector<Edge>> v[N * 2]; ISAP(int n, int m, int s, int t) : n(n), m(m), s(s), t(t) { input(); } inline void input() { // v.resize(n + 1); for (int x, y, f, i(0); i ^ m; ++i) { read(x, y, f); v[x].push_back(Edge(y, v[y].size() - 1, 0)); } } inline void init() { dep.assign(n + 1, -1); dep[t] = 0, gap[0] = 1; // 如果要用手写队列,要开大一点……避免玄学RE,虽然理论上N就够了 register int qt(0), qf(0); q[qt++] = t; int x, y; while (qf ^ qt) { x = q[qf++]; for (auto &e : v[x]) { if (dep[(y = e.to)] + 1 == dep[x] && e.flow) { f = sap(y, min(rest, e.flow)); if (f) { rest -= f, e.flow -= f; G[e.to][e.rev].flow += f; } if (!rest) return flow; // flow all used } } // bfs优化部分while (dep[s] <= n) { now.assign(n, 0); maxflow += sap(s, INF);}然后……就搞定了QwQ作者有话说一般来说,如果图非常稠密(边数远远大于点数),当前弧优化的力度就非常大了如:Zoj3229 Shoot the Bullet|东方文花帖|【模板】有源汇上下界最大流 - 洛谷这个专题我会放在网络流的其他部分详解,敬请期待……写了当前弧优化的Dinic能轻松过……没写全TLE虽然没写当前弧优化的ISAP能更快的过前三个点,但最后一个点过不了……QwQ我没有试过当前弧优化的ISAP更新:有当前弧优化的ISAP可以过但是如果边数不多,当前弧优化可能就成了负优化了……所以需要根据题目数据合理使用。归纳总结以上内容,生成内容时需要带简单的排版```

正文

网络最大流

前置知识以及更多芝士参考下述链接
网络流合集链接:网络流


最大流,值得是在不超过管道(边)容量的情况下从源点到汇点最多能到达的流量

抽象一点:使 \(\sum_{(S, v) \in E} f(S, v)\) 最大的流函数被称为网络的最大流,此时的流量被称为网络的最大流量


有了最大流量,就可以通过奇奇怪怪的建模解决很多令人摸不着头脑的题

例如二分图

对于一张二分图,经过建模之后我们可以这样画

其中左部点集 \(A = \{1, 2, 3, 4\}\),右部点集 \(B = \{5, 6, 7, 8\}\), 其中源点为 \(0\),汇点为 \(9\)

建模过程:

新增一个源点 \(S\) 和一个汇点 \(T\), 从 \(S\) 到每一个左部点连有向边,从每一个右部点到 \(T\) 连有向边,把原二分图的每条边看作从左部点到右部点的有向边,形成了一张 \(n + 2\) 个点 \(n + m\) 条边的网络。其中每一条边的容量都为 \(1\)

不难发现,二分图的最大匹配数就等于网络的最大流量。求出最大流后,所有有有”流“经过的点,边就是匹配点,匹配边。

进一步的:如果要求二分图多重匹配,依据题目信息改变连接汇点和源点的边的容量即可

计算最大流的算法很多,这里主要讲解 \(EK (Edmonds-Karp)\)\(Dinic\)\(ISAP\) 算法。


EK 增广路算法

这里的增广路与二分图里面的增广路有不一样了 Q^Q

增广路:若一条源点 \(S\)\(T\) 的路径上各边的剩余容量都严格大于 \(0\),则称这条路径为增广路。显然,可以让一个流沿着增广路流过,使得网络的流量增大。

\(EK\)的思路就是不断进行广度优先搜索寻找增广路,直到不存在增广路为止。

而在搜索的时候,我们只考虑图中所有 \(f(x, y) < c(x,y)\) 的边(或者说是有剩余容量的边)。在寻找路径的同时,还要记录其前驱结点以及路径上的最小容量\(minf\), 在找到一条增广路后,则网络的流量可以增加 \(minf\)

但是,考虑到斜对称的性质,由于我们需要把增广路上的所有边的剩余容量 \(e(u, v)\) 减去\(minf\),所以要在对其反边容量 \(e(v, u)\) 加上 \(minf\)

初始化的时候,如果是有向边 \((u, v)\),则 \(e(u, v) = c(u, v), e(v, u) = 0\),如果是无向边,则 \(e(u, v) = e(v, u) = c(u, v)\)

?: 为什么会出现无向边,网络流不是有向图吗?

考虑双向道路马路,既可以顺流,又可以逆着。

例如:[ICPC-Beijing 2006] 狼抓兔子 - 洛谷

这道题需要用到最小割,顺便说一下,最小割 = 最大流

?: 为什么使用BFS,而不是DFS?

因为DFS可能会绕圈圈……在讲述DInic的时候我会再提及

复杂度:复杂度上界为 \(O(nm^2)\),然而实际上远远达不到这个上界,效率还行,可以处理 \(10^3 \sim 10^4\) 规模的网络

--《算法竞赛进阶指南》

我不会证明,下面的两个算法也不会 Q^Q

这里给出一种参考代码

【模板】网络最大流 - 洛谷

提交记录:记录详情

#include <iostream>
#include <cstring>
#include <algorithm>
#include <deque>

using std::deque;

const int N = 2e3 + 7, M = 5e5 + 7, INF = 0x7F7F7F7F;

int n, m, s, t;
int to[M], nex[M], wi[M] = {INF};
int head[N], tot = 1;

void add(int u, int v, int w) {
    to[++tot] = v, nex[tot] = head[u], wi[tot] = w, head[u] = tot;
}

void read() {
    scanf("%d %d %d %d", &n, &m, &s, &t);

    int u, v, w;
    for (int i = 0; i < m; ++i) {
        scanf("%d %d %d", &u, &v, &w);
        add(u, v, w);
        add(v, u, 0);
    }
}

#define min(x, y) ((x)<(y)?(x):(y))
#define pop() que.pop_front()
#define top() que.front()
#define push(x) que.push_back(x);
#define empty() que.empty()

static int inq[N], it = 0;
// px记录上一个点,pe记录遍历过来的边
static int px[N], pe[N];
inline int bfs() {
    deque<int> que;

    push(s); inq[s] = ++it;

    int x, y;
    while (!empty()) {
        x = top(); pop();
        for (int i = head[x]; i; i = nex[i]) {
            if ((inq[(y = to[i])] ^ it) && wi[i]) {
                px[y] = x, pe[y] = i;
                if (y == t) return 1; // 找到增广路了,
                inq[y] = it; push(y);
            }
        }
    }
    return 0; // 到不到了,没有增广路了
}

void work(long long & res) {
    while (bfs()) {
        int val = INF;
        for (int x = t; x ^ s; x = px[x]) {
            val = min(val, wi[pe[x]]);
        }

        for (int x = t; x ^ s; x = px[x]) {
            wi[pe[x]] -= val;
            // 处理反边的时候利用了成对变换的方法!
            wi[pe[x] ^ 1] += val;
        }

        res += val;
    }
}

int main() {
    read();

    long long res = 0;
    work(res);
    printf("%lld\n", res);
    return 0;
}

Dinic

考虑到 \(EK\) 算法每一次在残量网络上只找出来的一条增广路,太慢了,所以有了更优化的东西 Dinic?歌姬吧

先引入一点点概念:

深度:在搜索树上的深度(BFS搜索时的层数)

残量网络:网络中所有节点以及剩余容量大于 \(0\) 的边构成的子图

分层图:依据深度分层的一段段图……或者说在残量网络上,所有满足 \(dep[u] + 1 = dep[v]\) 的边 \((u, v)\) 构成的子图。

分层图显然是一张有向无环图

Dinic 算法不断重复下述过程,直到在残量网络中,\(S\) 不能到达 \(T\)

  • 利用BFS求出分层图

  • 在分层图上DFS寻找增广路,在回溯的时候实时更新剩余容量。另外,每个点可以同时流出到多个结点,每个点也可以接收多个点的流。

?: 这里为什么可以使用DFS

由于我们分了层,意味着DFS只会向更深的地方搜索,而不会在同一层乱跳,甚至搜索到前面。这也是为什么EK用BFS更优秀

复杂度:一般来说,时间复杂度为 \(O(n^2m)\),可以说是不仅简单,而且容易实现的高效算法之一,一般能够处理 \(10^4 \sim 10^5\) 规模的网络。特别的,用此算法求二分图的最大匹配时只需要 \(O(m\sqrt{n})\), 实际上表现会更好。

题目不变

没有当前弧优化:提交详情

有当前弧优化:记录详情

// 重复内容已省略

int dis[N], vis[N], vt = 0;
int now[N]; // 用于当前弧优化
// return true if exists non-0 road to t
bool bfs() {
    memset(dis, 0, sizeof(dis)); dis[s] = 1;

    deque<int> que;
    que.push_back(s);
    while (que.size()) {
        int x = que.front(); que.pop_front();
        now[x] = head[x]; // 更新当前弧
        for (int y, i = head[x]; i; i = nex[i]) {
            if (!dis[y = to[i]] && wi[i]) {
                dis[y] = dis[x] + 1;
                que.push_back(y);
                if (y == t) return true;
            }
        }
    }
    return false;
}

#define min(x, y) ((x) < (y) ? (x) : (y))

long long dinic(int x, long long maxflow) {
    if (x == t) return maxflow;
    long long rest = maxflow, k;
    for (int y, i = now[x]; i && rest; i = nex[i]) {
        now[x] = i; // 更新当前弧
        // 要在更深的一层,以及需要有剩余流量
        if (dis[y = to[i]] == dis[x] + 1 && wi[i]) {
            k = dinic(y, min(rest, wi[i]));
            if (!k) dis[y] = 0;
            wi[i] -= k, wi[i ^ 1] += k;
            rest -= k;
        }
    }
    return maxflow - rest;
}

int main() {
    read();
    long long maxflow = 0, flow;
    while (bfs()) {
        while (flow = dinic(s, INF)) maxflow += flow;
    } 
    printf("%lld\n", maxflow);
}

?: 当前弧优化是个啥玩意

注意到如果我们每一次遍历后,对于当前边 \((u, v)\),不可能再有流量流过这条边,所以我们可以暂时的删除这条边……注意,只是暂时,每一分层的时候是需要考虑这条边的,因为这条边的剩余流量不一定为 0


ISAP

某位大佬的博客上说这是究极最大流算法之一。还有一个HLPP(最高标记预留推进),思路完全与这几个方法不同,不依赖于增广路,我会把它放在另外的文章中单独讲。我可不会告诉你们是我不会优化,太笨了,看不懂大佬的优化

题目链接:【模板】最大流 加强版 / 预流推进 - 洛谷

这是我的:记录详情 4.77s

这是大佬的:记录详情 185ms

由于Dinic需要多次BFS……所以有些不满足的数学家决定优化常数……于是有了ISAP,只需要一次BFS的东西……

可恶,竟然没有找到不用gap优化的写法 T^T

ISAP算法从某种程度上是SAP算法和Dinic的融合

SAP算法就是所谓的EK算法……ISAP也就是Improved SAP……但是主体怎么跟DInic几乎一模一样!

算法流程如下:

  1. \(T\) 开始进行BFS,直到 \(S\) ,标记深度,同时记录当前深度有多少个

  2. 利用DFS不断寻找增广路,思路与Dinic类似

  3. 每次回溯结束后,将所在节点深度加一(远离 \(T\) 一点),同时更新深度记录。如果出现了断层(有一个深度没有点了)那么结束寻找。

?: 为什么需要深度加一

由于我们在遍历过一次过后,这个点不可能再向更靠近 \(T\) 的点送出流量,所以只能退而求其次,给自己同层的结点送流量。

怎么跟Dinic一摸一样啊,关键是也可以用当前弧优化,只是我用写的是vetor存图……需要一点奇技淫巧

参考代码……

提交题目还是【模板】网络最大流 - 洛谷

记录详情

!! 竟然在最优解第二页 O-O

对于下面代码做出一些解释

?: 为什么终止条件是 dep[s] > n

由于我们只有n个点,意味着初始化的时候最大深度为 n

考虑如果是dep[s] <= n 的情况,那么要么是有连续的层,要么断层了(此时我们在DFS中会将dep[s]设为n+1

如果源点深度 \(\gt n\) 所以一定会有一个深度是没有对应的点的,意味着一定出现了断层,也就是流量无法到达了

所以,更新答案之后就可以结束循环了

?: 为什么新建ISAP实例的时候需要用static

在国外某知名网站上有这么一句回复When objects are created, the members of the object cannot be initialized directly and this problem of not being able to initialize data members is known as the problem of initialization.

这避免了需要手动清空内存的需要 no memset

同时,在编译时提前创建好实例,可以避免过多造成运行时开销,以及减少栈空间的消耗,减少爆栈的可能(QwQ 我不是特别确定,似乎static创建的对象是存储在BSS部分的,据说会比较慢,还不如定义在外面

不要管我……

?: 为什么我用这么多的vector

考虑到vector的assign在某些情况下比memset要快,而且更方便(我也不知道为什么

vector的 assign(n, v) 指的是将前 n 个元素设为 v

但是如果只需要使用一次或者不需要清零(如手写栈)的数组,就别用vector了,一个数组就搞定了

?: 为什么要多此一举搞一个Data出来

方便最后的时候把int改成long long

不开long long见祖宗

复杂度:与Dinic同阶,但是常数相对小一点

// 写这个的时候,借鉴了写HLPP最优解的大佬写快读的方法……
template<typename T>
inline void read(T &x) {
    char c, f(0); x = 0;
    do if ((c = getchar()) == '-') f = true; while (isspace(c));
    do x = (x<<3) + (x<<1) + (c ^ 48), c = getchar(); while (isdigit(c));
    if (f) x = -x;
}
template <typename T, typename ...Args> inline void read(T &t, Args&... args) { read(t), read(args...); }

typedef long long Data;
using namespace std;

const int N = 207, M = 5007;

struct Edge {
    int to;
    size_t rev; // 反边的位置,用int也没问题
    Data flow;
    Edge(int to, size_t rev, Data f) : to(to), rev(rev), flow(f) {}
};

class ISAP {
public:
    int n, m, s, t;
    vector<int> dep;
    int q[N * 2], gap[N * 2];

    // vector< vector<Edge> > v;
    vector<Edge> v[N * 2];

    ISAP(int n, int m, int s, int t) : n(n), m(m), s(s), t(t) {
        input();
    } 

    inline void input() {
        // v.resize(n + 1);
        for (int x, y, f, i(0); i ^ m; ++i) {
            read(x, y, f);
            v[x].push_back(Edge(y, v[y].size(), f));
            v[y].push_back(Edge(x, v[x].size() - 1, 0));
        }
    }

    inline void init() {
        dep.assign(n + 1, -1);
        dep[t] = 0, gap[0] = 1;

        // 如果要用手写队列,要开大一点……避免玄学RE,虽然理论上N就够了
        register int qt(0), qf(0);
        q[qt++] = t;
        int x, y;
        while (qf ^ qt) {
            x = q[qf++];
            for (auto &e : v[x]) {
                if (dep[(y = e.to)] == -1) // if dep[y] != -1
                    ++gap[(dep[y] = dep[x] + 1)], q[qt++] = y;
            }
        } // bfs end
    }

    inline Data sap(int x, Data flow) {
        if (x == t) return flow;

        Data rest = flow;
        int y, f;
        for (auto &e : v[x]) {
            if (dep[(y = e.to)] + 1 == dep[x] && e.flow) {
                f = sap(y, min(e.flow, rest));
                if (f) {
                    e.flow -= f, v[e.to][e.rev].flow += f;
                    rest -= f;
                }
                if (!rest) return flow; // flow all used
            }
        }

        // change dep
        if (--gap[dep[x]] == 0) dep[s] = n + 1; // can not reach to t
        ++gap[++dep[x]]; // ++depth
        return flow - rest;
    }

    inline Data calc() {
        Data maxflow(0);
        static const Data INF(numeric_limits<Data>::max());
        // dep[s]最大为n,为一条链的时候
        while (dep[s] <= n) {
            // 如果要当前弧优化,在这里需要重置当前弧的now!
            maxflow += sap(s, INF);
        }
        return maxflow;
    }
};

int main() {
    int n, m, s, t;
    read(n, m, s, t);

    static ISAP isap(n, m, s, t);
    isap.init();
    printf("%lld\n", isap.calc());

    return 0;
}

?: 如果我想用vector存图实现当前弧优化怎么整

在sap函数的主体部分

for (int & i = now[x]; i < G[x].size(); ++i) {
    Edge & e = G[x][i]; // 这里就是当前弧优化
    if (dep[(y = e.to)] + 1 == dep[x] && e.flow) {
        f = sap(y, min(rest, e.flow));
        if (f) {
            rest -= f, e.flow -= f;
            G[e.to][e.rev].flow += f;
        }
        if (!rest) return flow;
    }
}

在calc部分

while (dep[s] <= n) {
    now.assign(n, 0);
    maxflow += sap(s, INF);
}

然后……就搞定了QwQ

作者有话说

一般来说,如果图非常稠密(边数远远大于点数),当前弧优化的力度就非常大了

如:Zoj3229 Shoot the Bullet|东方文花帖|【模板】有源汇上下界最大流 - 洛谷

这个专题我会放在网络流的其他部分详解,敬请期待……

写了当前弧优化的Dinic能轻松过……没写全TLE

虽然没写当前弧优化的ISAP能更快的过前三个点,但最后一个点过不了……QwQ

我没有试过当前弧优化的ISAP

更新:有当前弧优化的ISAP可以过

但是如果边数不多,当前弧优化可能就成了负优化了……所以需要根据题目数据合理使用

与算法学习笔记(8.1): 网络最大流算法 EK, Dinic, ISAP相似的内容:

算法学习笔记(8.1): 网络最大流算法 EK, Dinic, ISAP

网络最大流算法 EK, Dinic, ISAP 详解

算法学习笔记(8): 网络流

# 网络流 > 网络流是一个博大精深的一个话题…… ## 箕(基)畚(本)知识 文章链接:[基本知识](https://www.cnblogs.com/jeefy/p/17050220.html) ## 网络最大流 文章链接:[网络最大流](https://www.cnblogs.com/jeefy

算法学习笔记(8.0): 网络流前置知识

网络流基础 网络流合集链接:网络流 网络 $G = (V, E)$ 实际上是一张有向图 对于图中每一条有向边 $(x, y) \in E$ 都有一个给定的容量 $c(x, y)$ 特别的,若 $(x,y) \notin E$ , 则 $c(x, y) = 0$ 图中还有两个指定的特殊结点,$S, T

算法学习笔记(8.2): 上下界网络流

# 上下界网络流 [TOC] > 前置知识以及更多芝士参考下述链接 > 网络流合集链接:[网络流](https://www.cnblogs.com/jeefy/p/17050215.html) 上下界网络流是普通网络流的一种变体,对于网络流,我们不仅关注其流量的上界,下届同样有所体现。 题型大致有五

算法学习笔记(3): 倍增与ST算法

倍增 目录倍增查找 洛谷P2249重点变式练习快速幂ST表扩展 - 运算扩展 - 区间变式答案倍增更多的用法优化矩形查询优化建图优化 DP作者有话说 倍增,字面意思即”成倍增长“ 他与二分十分类似,都是基于”2“的划分思想 那么具体是怎么样,我们以一个例子来看 查找 洛谷P2249 依据题面,我们知

算法学习笔记(6): 树链剖分

树链剖分 树链剖分是一个很神奇,但是在树上可以完成一些区间操作问题 简单来说,就是把一棵树分成一条条的链,通过维护链上的信息来维护整棵树的信息 基础知识可以参考我的另外一篇博客:算法学习笔记(5): 最近公共祖先(LCA) 这里假设你已经掌握了上述博客中的所有相关知识,并清晰了其背后的原理 性质?发

算法学习笔记(11): 原根

原根 此文相对困难,请读者酌情食用 在定义原根之前,我们先定义其他的一点东西 阶 通俗一点来说,对于 $a$ 在模 $p$ 意义下的阶就是 $a^x \equiv 1 \pmod p$ 的最小正整数解 $x$ 或者说,$a$ 在模 $p$ 意义下生成子群的阶(群的大小) 再或者说,是 $a$ 在模

算法学习笔记(30):Kruskal 重构树

Kruskal 重构树 这是一种用于处理与最大/最小边权相关的一个数据结构。 其与 kruskal 做最小生成树的过程是类似的,我们考虑其过程: 按边权排序,利用并查集维护连通性,进行合并。 如果我们在合并时,新建一个节点,其权值为当前处理的边的权值,并将合并的两个节点都连向新建的节点,那么就可以得

算法学习笔记(3.1): ST算法

ST表 在RMQ(区间最值)问题中,著名的ST算法就是倍增的产物。ST算法可以在 \(O(n \log n)\) 的时间复杂度能预处理后,以 \(O(1)\) 的复杂度在线回答区间 [l, r] 内的最值。 当然,ST表不支持动态修改,如果需要动态修改,线段树是一种良好的解决方案,是 \(O(n)\

C++算法之旅、09 力扣篇 | 常见面试笔试题(上)算法小白专用

算法学习笔记,记录容易忘记的知识点和难题。详解时空复杂度、50道常见面试笔试题,包括数组、单链表、栈、队列、字符串、哈希表、二叉树、递归、迭代、分治类型题目,均带思路与C++题解