图论

图是图论的主要研究对象。图是由若干给定的顶点及连接两顶点的边所构成的图形。

存图方式

用 n 代指图的点数，用 m 代指图的边数。

邻接表

使用一个支持动态增加元素的数据结构构成的数组，如 vector<int> adj[n + 1] 来存边，其中 adj[u] 存储的是点 u 的所有出边的相关信息（终点、边权等）。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
struct edge {
    int v, w;
};

int main {
    int n, m;
    cin >> n >> m;
    vector<vector<edge>> adj(n + 1);
    for (int i = 1; i <= m; ++i) {
        int u, v, w;
        cin >> u >> v >> w;
        adj[u].emplace_back(v, w);
    }
}

拓扑排序

拓扑排序要解决的问题是如何给一个有向无环图的所有节点排序

Kahn 算法

若 tfn 元素数小于节点数则原 DAG 存在环，若任意时刻不存在两个入度为 0 的节点则拓扑序列唯一。

Kahn 算法同样可以为无向图判环，在判断无向图中是否存在环时，是将所有 度 <= 1 的结点入队。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
vector<int> tfn; // 拓扑序列  
auto topu = [&]() {  
    queue<int> q;  
    for (int u = 1; u <= 3 * n; ++u) {  
        if (deg[u]) continue;  
        q.push(u), tfn.push_back(u);  
    }  
    while (!q.empty()) {  
        int u = q.front();  
        q.pop();  
        for (auto v: adj[u]) {  
            if (--deg[v]) continue;  
            q.push(v), tfn.push_back(v);  
        }  
    }  
};

最短路

Floyd 算法

Floyd 算法是一个基于动态规划的全源最短路算法，能确定某张图是否存在负环，并确定那些不受负环影响的顶点对之间的最短路。时间复杂度 $\Theta(n^3)$ ，空间复杂度 $\Theta(n^2)$ 。

若存在任意节点 $i$ 满足 $f[i][i]<0$，则图中存在负环，可达节点 $i$ 的节点均不存在最短路。

1
2
3
4
5
6
7
8
vector<vector<int> > f = adj; // adj 为存储邻居节点最短边权的邻接矩阵
for (int k = 1; k <= n; ++k) {
    for (int i = 1; i <= n; ++i) {
        for (int j = 1; j <= n; ++j) {
            f[i][j] = min(f[i][j], f[i][k] + f[k][j]);
        }
    }
}

Floyd 传递闭包算法

Floyd 传递闭包算法通过动态规划逐步枚举中间节点，在 $\Theta(n^3)$ 时间内计算出有向图中所有点对的可达性。

使用 bitset 优化后复杂度为 $\Theta\left( \frac{n^3}{w} \right)$ ，其中 $w=64$ 为计算机位宽。

初始化时注意自身节点永远可达。

1
2
3
4
5
6
7
constexpr int maxn = 1000;  
vector<bitset<maxn> > reach = adj;
for (int k = 1; k <= n; ++k) {  
    for (int u = 1; u <= n; ++u) {  
        if (reach[u][k]) reach[u] |= reach[k];  
    }  
}

SPFA 算法

SPFA 算法即队列优化的 Bellman–Ford 算法，是一种基于松弛操作的最短路算法，可以求出有负权的图的最短路，并可以对最短路不存在的情况进行判断，最坏情况下时间复杂度为 $O(nm)$ 。

通过记录最短路经过了多少条边，SPFA 可以用于判断从 $s$ 点出发是否能抵达一个负环，当经过了至少 $n$ 条边时，函数返回 false，$s$ 点可以抵达一个负环。

如果需要判断图是否存在负环，需要建立超级源点，即向图上每个节点连一条权值为 $0$ 的边，然后以超级源点为起点执行 SPFA 算法，注意此时需要判断的边数为 $n+1$ 。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
vector<int> dis(n + 1, inf), cnt(n + 1);  
auto spfa = [&](int s) {  
    queue<int> q;  
    dis[s] = 0, q.push(s);  
    while (!q.empty()) {  
        auto u = q.front();  
        q.pop();  
        for (auto [v, w]: adj[u]) {  
            if (dis[u] + w < dis[v]) {  
                dis[v] = dis[u] + w;  
                cnt[v] = cnt[u] + 1;  
                if (cnt[v] >= n) return false;  
                q.push(v);  
            }  
        }  
    }  
    return true;  
};

Dijkstra 算法

Dijkstra 算法是一种求解非负权图上单源最短路径的算法。

以下代码使用优先队列实现，不使用 vis 数组，时间复杂度 $O((n+m)\log{m})$ 。

注意保证入队判断和更新判断一致且正确，否则会出现时间复杂度错误。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
struct node {  
    int u, d;  
};  
  
vector<int> dis(n + 1, inf);  
auto dijkstra = [&](int s) {  
    constexpr auto cmp = [](auto x, auto y) { return x.d > y.d; };  
    priority_queue<node, vector<node>, decltype(cmp)> q(cmp);  
    dis[s] = 0;  
    q.emplace(s, 0);  
    while (!q.empty()) {  
        auto [u, d] = q.top();  
        q.pop();  
        if (dis[u] < d) continue;  
        for (auto [v, w]: adj[u]) {  
            if (dis[u] + w < dis[v]) {  
                dis[v] = dis[u] + w;  
                q.emplace(v, dis[v]);  
            }  
        }  
    }  
};

Johnson 算法

新建一个虚拟节点（在这里我们就设它的编号为 $0$ ）。从这个点向其他所有点连一条边权为 $0$ 的边，然后用 SPFA 求出 $0$ 号点到其它所有点的最短路作为势能函数，然后根据势能函数重新标注边权使边权非负，然后就可以跑 $n$ 轮 Dijkstra 计算最短路，时间复杂度 $O(nm\log{m})$。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
vector<int> h(n + 1, inf);
vector dis(n + 1, vector(n + 1, inf));
// SPFA and Dijkstra
auto johnson = [&]() {  
    if (!spfa(0)) return false;  
    for (int u = 1; u <= n; u++) {  
        for (auto &[v, w]: adj[u]) {  
            w = w + h[u] - h[v];  
        }  
    }  
    for (int u = 1; u <= n; u++) dijkstra(u);  
    for (int u = 1; u <= n; u++) {  
        for (int v = 1; v <= n; v++) {  
            if (dis[u][v] == inf) continue;  
            dis[u][v] = dis[u][v] + h[v] - h[u];  
        }  
    }  
    return true;  
};

最小生成树

最小生成树（MST）是一个无向连通图中边的子集，它连接所有顶点且边的权值之和最小，同时不形成任何环。

Kruskal 算法

Kruskal 算法是一种贪心算法，通过按边权从小到大排序并依次选择不构成环的边来构造最小生成树。

以下代码使用并查集判环，Kruskal 部分时间复杂度为 $O(m\alpha(N))$，总时间复杂度 $O(m\log{m})$ 。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
vector<tuple<int, int, int>> edges;
auto kruskal = [&]() {  
    ranges::sort(edges, [](auto &a, auto &b) {  
        return get<2>(a) < get<2>(b);  
    });  
  
    DSU dsu(n + 1);  
    int sum = 0, cnt = 0;  
    for (auto &[u, v, w]: edges) {  
        if (dsu.merge(u, v)) {  
            sum += w;  
            if (++cnt == n - 1) break;  
        }  
    }  
    return cnt == n - 1 ? sum : -1;  
};

Prim 算法

Prim 算法是一种贪心算法，其核心思想是从任意节点开始，逐步选择与当前生成树相连且权值最小的边，直到覆盖所有节点。

以下代码使用优先队列实现，时间复杂度为 $O((n+m)\log{m})$ 。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
struct node {  
    int u, d;  
};

auto prim = [&]() {  
    vector<char> vis(n + 1);  
    constexpr auto cmp = [](auto x, auto y) { return x.d > y.d; };  
    priority_queue<node, vector<node>, decltype(cmp)> pq(cmp);  
    pq.emplace(1, 0);  
    int sum = 0, cnt = 0;  
    while (!pq.empty() && cnt < n) {  
        auto [u, d] = pq.top();  
        pq.pop();  
        if (vis[u]) continue;  
        vis[u] = true;  
        sum += d;  
        cnt++;  
        for (auto &[v, w]: adj[u]) {  
            if (vis[v]) continue;  
            pq.emplace(v, w);  
        }  
    }  
    return cnt == n ? sum : -1;  
};

倍增算法 LCA

通过 dfs 预处理出每个节点第 $2^{i}$ 祖先节点，可求两节点间最近公共祖先 (LCA)。

树上两点路径的信息可以在倍增算法求 LCA 的过程中同时处理，如路径长度。

倍增算法的预处理时间复杂度为 $O(n\log n)$ ，单次查询时间复杂度为 $O(\log n)$ 。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
constexpr int maxd = (int) bit_width(100000ull); // 视题目范围确定，也可直接取31

vector<int> dep(n + 1);
vector<array<int, maxd> > fa(n + 1);
vector<array<i64, maxd> > cost(n + 1);
auto dfs = [&](auto &&self, int u) -> void {
    dep[u] = dep[fa[u][0]] + 1;
    for (int i = 1; i < maxd; ++i) {
        fa[u][i] = fa[fa[u][i - 1]][i - 1];
        cost[u][i] = cost[fa[u][i - 1]][i - 1] + cost[u][i - 1];
    }
    for (auto [v, w]: adj[u]) {
        if (v == fa[u][0]) continue;
        fa[v][0] = u, cost[v][0] = w;
        self(self, v);
    }
};
dfs(dfs, 1); // 为根节点调用dfs

auto lca = [&](int x, int y) {
    if (dep[x] > dep[y]) swap(x, y);
	i64 ans = 0;
    int tmp = dep[y] - dep[x];
    for (int j = 0; tmp; ++j, tmp >>= 1)
        if (tmp & 1) ans += cost[y][j], y = fa[y][j];
    if (y == x) return ans; // 此时 x 与 y 均为最近公共祖先
    for (int j = maxd - 1; j >= 0; --j) {
        if (fa[x][j] == fa[y][j]) continue;
        ans += cost[x][j] + cost[y][j];
        x = fa[x][j], y = fa[y][j];
    }
    ans += cost[x][0] + cost[y][0];
    return ans; // 此时 x 或 y 的父节点为最近公共祖先
};

重链剖分

重链剖分可以将树上的任意一条路径划分成不超过 $O(\log n)$ 条连续的链，每条链上的点深度互不相同（即是自底向上的一条链，链上所有点的 LCA 为链的一个端点）。

树上每个节点都属于且仅属于一条重链，所有的重链将整棵树完全剖分，在剖分时重边优先遍历。

最后树的 DFS 序上，重链和一颗子树内内的 DFS 序是连续的。按 DFN 排序后的序列即为剖分后的链。

  • fa(x) 表示节点 x 在树上的父亲。
  • dep(x) 表示节点 x 在树上的深度。
  • sz(x) 表示节点 x 的子树的节点个数。
  • son(x) 表示节点 x 的 重儿子。
  • top(x) 表示节点 x 所在 重链 的顶部节点（深度最小）。
  • btn(x) 表示以节点 x 为根节点的子树 DFS 序最靠后的一个节点。
  • dfn(x) 表示节点 x 的 DFS 序，也是其在线段树中的编号。
  • rnk(x) 表示 DFS 序所对应的节点编号，有 rnk(dfn(x))=x。

树链剖分配合线段树、树状数组可以维护树上两点路径上的信息，子树上的信息。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
vector<i64> fa(n + 1), dep(n + 1), sz(n + 1), son(n + 1);  
auto dfs1 = [&](auto &&self, i64 u) -> i64 {  
    son[u] = 0, sz[u] = 1;  
    for (auto [v, w]: adj[u]) {  
        if (v == fa[u]) continue;  
        fa[v] = u, dep[v] = dep[u] + 1;  
        sz[u] += self(self, v);  
        if (sz[v] > sz[son[u]]) son[u] = v;  
    }  
    return sz[u];  
};  
dfs1(dfs1, r);  
  
i64 id = 0;  
vector<i64> top(n + 1), btn(n + 1), dfn(n + 1), rnk(n + 1);  
auto dfs2 = [&](auto &&self, i64 u, i64 t) -> i64 {  
    top[u] = t, id++;  
    dfn[u] = id, rnk[id] = u;  
    if (!son[u]) return btn[u] = u;  
    btn[u] = self(self, son[u], t);  
    for (auto [v, w]: adj[u]) {  
        if (v == fa[u] || v == son[u]) continue;  
        btn[u] = self(self, v, v);  
    }  
    return btn[u];  
};  
dfs2(dfs2, r, r);  
  
vector<Info> val(n + 1);  
for (int i = 1; i <= n; ++i) val[dfn[i]].val = nd[i];  
SegmentTree<Info, Tag> seg(val);  
  
auto path_apply = [&](i64 u, i64 v, i64 w) {  
    while (top[u] != top[v]) {  
        if (dep[top[u]] < dep[top[v]]) swap(u, v);  
        seg.apply(1, 0, n + 1, dfn[top[u]], dfn[u] + 1, {w});  
        u = fa[top[u]];  
    }  
    if (dep[u] < dep[v]) swap(u, v);  
    seg.apply(1, 0, n + 1, dfn[v], dfn[u] + 1, {w});  
};  
  
auto path_query = [&](i64 u, i64 v) {  
    i64 sum = 0;  
    while (top[u] != top[v]) {  
        if (dep[top[u]] < dep[top[v]]) swap(u, v);  
        sum += seg.query(1, 0, n + 1, dfn[top[u]], dfn[u] + 1).val;  
        u = fa[top[u]];  
    }  
    if (dep[u] < dep[v]) swap(u, v);  
    sum += seg.query(1, 0, n + 1, dfn[v], dfn[u] + 1).val;  
    return sum % p;  
};

图的匹配

图的匹配是指在图中选择一组边，使得没有两个边共享同一个节点。

二分图最大匹配

使用匈牙利算法来进行增广路径搜索，从而寻找最大匹配，时间复杂度为 $O(n_{l}n_{r})$ 。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
vector<int> ml(n + 1), mr(m + 1), vis(m + 1);
auto dfs = [&](auto &&self, int u, int t) -> bool {
    for (int v = 1; v <= m; ++v) {
        if (!adj[u][v]) continue;
        if (vis[v] == t) continue;
        vis[v] = t;
        if (!mr[v] || self(self, mr[v], t)) {
            ml[u] = v, mr[v] = u;
            return true;
        }
    }
    return false;
};

int ans = 0;
for (int u = 1; u <= n; ++u) {
    if (ml[u]) continue;
    ans += dfs(dfs, u, u);
}

网络最大流

令 $G=(V,E)$ 是一个有源汇点的网络，我们希望在 $G$ 上指定合适的流 $f$，以最大化整个网络的流量 $|f|$ （即 $\sum_{x \in V} f(s, x) - \sum_{x \in V} f(x, s)$），这一问题被称作最大流问题。

以下算法建图方法均为：

1
2
3
4
struct edge {  
    int v, r;  
    i64 w, f;  
};

Dinic 算法

Dinic 算法是一种用于解决最大流问题的高效算法，通过分层图和分层深度优先搜索增广流量，优化了增广路径的寻找过程，最坏复杂度 $O(n^{2}m)$ ，但在在一些性质良好的图上，Dinic 算法有更好的时间复杂度。

稠密图 稀疏图 平面图 二分图 有向无环图 单位容量网络
$O(n^{2}m)$ $O(n^{\frac{3}{2}})$ $O(n^{\frac{3}{2}})$ $O(\sqrt{n}m)$ $O(nm)$ $O(m \min(m^\frac{1}{2}, n^{\frac{2}{3}}))$ 
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
constexpr i64 inf = 1e18;  
  
struct Flow {  
    vector<vector<edge> > adj;  
    vector<int> dep, cur;  
  
    explicit Flow(int n) : adj(n), dep(n, -1), cur(n) {  
    }  
  
    void add_edge(int u, int v, i64 w) {  
        adj[u].emplace_back(v, (int) adj[v].size(), w, 0);  
        adj[v].emplace_back(u, (int) adj[u].size() - 1, 0, 0);  
    }  
  
    bool bfs(int s, int t) {  
        ranges::fill(dep, -1);  
        queue<int> q;  
        q.push(s);  
        dep[s] = 0;  
        while (!q.empty()) {  
            int u = q.front();  
            q.pop();  
            for (auto &[v, r, w, f]: adj[u]) {  
                if (~dep[v] || w <= f) continue;  
                dep[v] = dep[u] + 1;  
                if (v == t) return true;  
                q.push(v);  
            }  
        }  
        return false;  
    }  
  
    i64 dfs(int u, int t, i64 flow) {  
        if (u == t) return flow;  
        i64 res = 0;  
        for (int &i = cur[u]; i < adj[u].size(); i++) {  
            auto &[v, r, w, f] = adj[u][i];  
            if (dep[v] != dep[u] + 1 || w <= f) continue;  
            i64 d = dfs(v, t, min(flow - res, w - f));  
            if (d <= 0) continue;  
            res += d, f += d;  
            adj[v][r].f -= d;  
            if (res == flow) return res;  
        }  
        return res;  
    }  
  
    i64 dinic(int s, int t) {  
        i64 res = 0;  
        while (bfs(s, t)) {  
            ranges::fill(cur, 0);  
            res += dfs(s, t, inf);  
        }  
        return res;  
    }  
};

ISAP 算法

ISAP 算法是 Dinic 算法的优化版本，通过反向BFS初始化分层，并在增广的过程中完成重分层，从而避免多次运行 BFS 。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
struct Flow {  
    int n;  
    vector<vector<edge> > adj;  
    vector<int> dep, cur, gap;  
  
    explicit Flow(int n) : n(n), adj(n), dep(n, -1), cur(n), gap(n + 2) {  
    }  
  
    void add_edge(int u, int v, i64 w) {  
        adj[u].emplace_back(v, (int) adj[v].size(), w, 0);  
        adj[v].emplace_back(u, (int) adj[u].size() - 1, 0, 0);  
    }  
  
    void init(int t) {  
        queue<int> q;  
        q.push(t);  
        dep[t] = 0;  
        while (!q.empty()) {  
            int u = q.front();  
            q.pop();  
            for (auto &[v, r, w, f]: adj[u]) {  
                if (~dep[v] || adj[v][r].w <= adj[v][r].f) continue;  
                dep[v] = dep[u] + 1;  
                q.push(v);  
            }  
        }  
    }  
  
    i64 dfs(int u, int s, int t, i64 flow) {  
        if (u == t) return flow;  
        i64 res = 0;  
        for (int &i = cur[u]; i < adj[u].size(); i++) {  
            auto &[v, r, w, f] = adj[u][i];  
            if (dep[v] + 1 != dep[u] || w <= f) continue;  
            i64 d = dfs(v, s, t, min(flow - res, w - f));  
            if (d <= 0) continue;  
            res += d, f += d;  
            adj[v][r].f -= d;  
            if (res == flow) return res;  
        }  
        if (--gap[dep[u]] == 0) dep[s] = n;  
        ++gap[++dep[u]], cur[u] = 0;  
        return res;  
    }  
  
    i64 isap(int s, int t) {  
        init(t);  
        if (dep[s] == -1) return 0;  
        for (int u = 0; u < n; ++u) {  
            if (~dep[u]) gap[dep[u]]++;  
        }  
        i64 res = 0;  
        while (dep[s] < n) res += dfs(s, s, t, infl);  
        return res;  
    }  
};

容量缩放算法

容量缩放优化的 Dinic 算法，时间复杂度 $O(nm\log{U})$，其中 $U$ 为所有边的最大容量。$\log{U}$ 部分上界较紧，正常情况不建议使用。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
constexpr i64 inf = 1e18;  
  
struct Flow {  
    vector<vector<edge> > adj;  
    vector<int> dep, cur;  
  
    explicit Flow(int n) : adj(n), dep(n, -1), cur(n) {  
    }  
  
    void add_edge(int u, int v, i64 w) {  
        adj[u].emplace_back(v, (int) adj[v].size(), w, 0);  
        adj[v].emplace_back(u, (int) adj[u].size() - 1, 0, 0);  
    }  
  
    bool bfs(int s, int t, i64 dt) {  
        ranges::fill(dep, -1);  
        queue<int> q;  
        dep[s] = 0, q.push(s);  
        while (!q.empty()) {  
            int u = q.front();  
            q.pop();  
            for (auto &[v, r, w, f]: adj[u]) {  
                if (~dep[v] || w - f < dt) continue;  
                dep[v] = dep[u] + 1;  
                if (v == t) return true;  
                q.push(v);  
            }  
        }  
        return false;  
    }  
  
    i64 dfs(int u, int t, i64 flow, i64 dt) {  
        if (u == t) return flow;  
        i64 res = 0;  
        for (int &i = cur[u]; i < adj[u].size(); ++i) {  
            auto &[v, r, w, f] = adj[u][i];  
            if (dep[v] != dep[u] + 1 || w - f < dt) continue;  
            i64 d = dfs(v, t, min(flow - res, w - f), dt);  
            if (d <= 0) continue;  
            res += d, f += d;  
            adj[v][r].f -= d;  
            if (res == flow) return res;  
        }  
        return res;  
    }  
  
    i64 scdinic(int s, int t) {  
        i64 res = 0;  
        for (i64 dt = 1ll << 31; dt > 0; dt >>= 1) {  
            while (bfs(s, t, dt)) {  
                ranges::fill(cur, 0);  
                res += dfs(s, t, inf, dt);  
            }  
        }  
        return res;  
    }  
};

MPM 算法

HLPP 算法

在通用的预流推送算法基础上，优先选择高度最高的溢出结点推送，算法复杂度为 $O(n^2\sqrt{m})$。

包含 BFS 优化和 GAP 优化，只处理高度小于 $n$ 的溢出节点，如需获得每条边的流量需处理高度小于 $2\times{n}$ 的溢出节点。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
constexpr int inf = 1e9;  
  
struct Flow {  
    int n, mh;  
    vector<vector<int> > bk;  
    vector<vector<edge> > adj;  
    vector<int> ht, gap;  
    vector<i64> ex;  
  
    explicit Flow(int n) : n(n), mh(0), bk(n), adj(n), ht(n, inf), gap(n), ex(n) {  
    }  
  
    void add_edge(int u, int v, i64 w) {  
        adj[u].emplace_back(v, (int) adj[v].size(), w, 0);  
        adj[v].emplace_back(u, (int) adj[u].size() - 1, 0, 0);  
    }  
  
    int top() {  
        while (~mh && bk[mh].empty()) mh--;  
        return mh == -1 ? 0 : bk[mh].back();  
    }  
  
    bool push(int u, int s, int t) {  
        bool init = u == s;  
        for (auto &[v, r, w, f]: adj[u]) {  
            if (w <= f || ht[v] == inf || !init && ht[u] != ht[v] + 1) continue;  
            if (v != s && v != t && !ex[v]) {  
                bk[ht[v]].push_back(v);  
                mh = max(mh, ht[v]);  
            }  
            i64 d = init ? w : min(ex[u], w - f);  
            f += d, adj[v][r].f -= d, ex[u] -= d, ex[v] += d;  
            if (!ex[u]) return false;  
        }  
        return true;  
    }  
  
    void relabel(int u) {  
        ht[u] = inf;  
        for (auto &[v, r, w, f]: adj[u]) {  
            if (w > f) ht[u] = min(ht[u], ht[v] + 1);  
        }  
        if (ht[u] < n) {  
            bk[ht[u]].push_back(u);  
            mh = max(mh, ht[u]);  
            ++gap[ht[u]];  
        }  
    }  
  
    bool bfs(int s, int t) {  
        queue<int> q;  
        q.push(t);  
        ht[t] = 0;  
        while (!q.empty()) {  
            int u = q.front();  
            q.pop();  
            for (auto &[v, r, w, f]: adj[u]) {  
                if (w || ht[v] <= ht[u] + 1) continue;  
                ht[v] = ht[u] + 1;  
                q.push(v);  
            }  
        }  
        return ht[s] != inf;  
    }  
  
    i64 hlpp(int s, int t) {  
        if (!bfs(s, t)) return 0;  
        for (auto h: ht) {  
            if (h != inf) gap[h]++;  
        }  
        ht[s] = n;  
  
        int u;  
        push(s, s, t);  
        while ((u = top())) {  
            bk[mh].pop_back();  
            if (!push(u, s, t)) continue;  
            if (!--gap[ht[u]]) {  
                for (int i = 0; i < n; i++) {  
                    if (i != s && i != t && ht[i] > ht[u] && ht[i] < n + 1) {  
                        ht[i] = n + 1;  
                    }  
                }  
            }  
            relabel(u);  
        }  
        return ex[t];  
    }  
};

最小费用最大流

最小费用最大流是指在网络流中找到一个流量最大的可行流，同时满足总费用最小的优化问题。

以下模板建图方法均为：

1
2
3
4
struct edge {  
    int v, r;  
    i64 w, f, c;  
};

SSP-EK 算法

以下代码是一个基于 EK 算法的 SSP 最小费用最大流实现，使用 SPFA 算法寻找增广路径，时间复杂度为 $O(nmf)$，其中 $f$ 是最大流量。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
constexpr i64 inf = 1e18;  
  
struct Flow {  
    vector<vector<edge> > adj;  
    vector<int> pn, pe;  
    vector<i64> dis;  
    i64 cost;  
  
    explicit Flow(int n) : adj(n), pn(n), pe(n), dis(n), cost(0) {  
    }  
  
    void add_edge(int u, int v, i64 w, i64 c) {  
        adj[u].emplace_back(v, (int) adj[v].size(), w, 0, c);  
        adj[v].emplace_back(u, (int) adj[u].size() - 1, 0, 0, -c);  
    }  
  
    bool spfa(int s, int t) {  
        ranges::fill(dis, inf);  
        queue<int> q;  
        dis[s] = 0, pn[t] = -1, q.push(s);  
        while (!q.empty()) {  
            int u = q.front();  
            q.pop();  
            for (int i = 0; i < adj[u].size(); i++) {  
                auto &[v, r, w, f, c] = adj[u][i];  
                if (w > f && dis[v] > dis[u] + c) {  
                    dis[v] = dis[u] + c;  
                    pn[v] = u, pe[v] = i;  
                    q.push(v);  
                }  
            }  
        }  
        return pn[t] != -1;  
    }  
  
    i64 ssp(int s, int t) {  
        i64 res = 0;  
        cost = 0;  
        while (spfa(s, t)) {  
            i64 d = inf;  
            for (int v = t; v != s; v = pn[v]) {  
                int u = pn[v], i = pe[v];  
                d = min(d, adj[u][i].w - adj[u][i].f);  
            }  
            for (int v = t; v != s; v = pn[v]) {  
                int u = pn[v], i = pe[v];  
                adj[u][i].f += d;  
                adj[adj[u][i].v][adj[u][i].r].f -= d;  
            }  
            res += d;  
            cost += d * dis[t];  
        }  
        return res;  
    }  
};

SSP-zkw 算法

反向求最短路，多路增广优化，时间复杂度仍为 $O(nmf)$，但在大多数图尤其是层次小，流量大，费用小的图更具优势。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
constexpr i64 inf = 1e18;  

struct Flow {  
    vector<vector<edge> > adj;  
    vector<char> vis;  
    vector<i64> dis;  
    i64 cost;  
  
    explicit Flow(int n) : adj(n), vis(n), dis(n), cost(0) {  
    }  
  
    void add_edge(int u, int v, i64 w, i64 c) {  
        adj[u].emplace_back(v, (int) adj[v].size(), w, 0, c);  
        adj[v].emplace_back(u, (int) adj[u].size() - 1, 0, 0, -c);  
    }  
  
    bool spfa(int s, int t) {  
        ranges::fill(vis, false);  
        ranges::fill(dis, inf);  
        queue<int> q;  
        dis[t] = 0, vis[t] = true, q.push(t);  
        while (!q.empty()) {  
            int u = q.front();  
            vis[u] = false, q.pop();  
            for (auto &[v, r, w, f, c]: adj[u]) {  
                if (adj[v][r].w <= adj[v][r].f || dis[v] <= dis[u] - c) continue;  
                dis[v] = dis[u] - c;  
                if (!vis[v]) vis[v] = true, q.push(v);  
            }  
        }  
        return dis[s] != inf;  
    }  
  
    i64 dfs(int u, int t, i64 flow) {  
        vis[u] = true;  
        if (u == t) return flow;  
        i64 res = 0;  
        for (auto &[v, r, w, f, c]: adj[u]) {  
            if (vis[v] || dis[v] + c != dis[u] || w <= f) continue;  
            i64 d = dfs(v, t, min(flow - res, w - f));  
            if (d <= 0) continue;  
            res += d, f += d;  
            adj[v][r].f -= d;  
            cost += d * c;  
            if (res == flow) return res;  
        }  
        return res;  
    }  
  
    i64 ssp(int s, int t) {  
        i64 res = 0;  
        cost = 0;  
        while (spfa(s, t)) {  
            vis[t] = true;  
            while (vis[t]) {  
                ranges::fill(vis, false);  
                res += dfs(s, t, inf);  
            }  
        }  
        return res;  
    }  
};

Primal-Dual 原始对偶算法

网络单纯形算法