ISAP,即改进的SAP算法。
算法基于:如果每次在残量网络中从汇点开始反向BFS, s-t 距离是单调不减的。
流程如下:
首先用BFS从汇点出发计算每个点到汇点的距离标号。
Advance:
每次走满足 $d_v = d_u - 1$ 的路径,即走最短路。维护每个点的当前弧。
(实践证明,只有ISAP需要当前弧优化……原始对偶加了之后反而更慢……Dinic加了之后有时候快很多,有时候慢一倍……)
Retreat:
当沿着允许弧无法增广的时候,说明需要更新距离标号,产生新的允许弧。
于是令 $d_u = \min \{ d_v \} + 1$ . 再从 s 开始走允许弧。特别地,当残量网络中 $u(u \neq t)$ 没有后继时,说明任何允许弧最终都不应该到达 $u$ 点,因为从 $u$ 点继续走无法到 $t$ ,更无法增广。于是此时令 $d_u = n$ ,即可排除 $u$ 点。( $n$ 个点的无权图中最大距离标号小于等于 $n-1$ )
GAP优化:设每个距离标号 $d$ 对应的点有 $num_d$ 个。如果某次 Retreat 操作后,$num_d = 0$ ,说明残量网络中 s-t 不连通,算法立刻结束。原因很简单:在这次修改距离标号之后,$d_s$ 只会增大。在从 $s$ 到 $t$ 的过程中,距离标号逐渐由 $d_s$ 一步一步地减小到0. 如果在这个过程中距离标号出现了断层(“gap”),s-t 最短路长度一定是无穷大,即残量网络中 s-t 不连通。于是当前流即为最大流。
代码
rk1,偷税
为啥第二名打表的都比我慢2333333
/*
* ISAP.cpp - LOJ最大流模板
*
* Copyright (C) 2018 - panda_2134
*
* This is free software; you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation; either version 2 of the License, or
* (at your option) any later version.
*
* This is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with This. If not, see <http://www.gnu.org/licenses/>.
*/
#include <bits/stdc++.h>
using namespace std;
typedef long long int64;
const int MAXN = 100, MAXM = 5000;
const int64 LL_INF = 0x3f3f3f3f3f3f3f3fLL;
namespace Maxflow {
struct Edge {
int u, v;
int64 flow, cap;
int next;
};
int n, s, t, e_ptr = 1, head[MAXN+10]; Edge E[(MAXM+10)<<1];
int p[MAXN+10], cur[MAXN+10], num[MAXN+10], d[MAXN+10], vis[MAXN+10];
void AddEdge(int u, int v, int64 cap) {
E[++e_ptr] = (Edge) { u, v, 0, cap, head[u] }; head[u] = e_ptr;
E[++e_ptr] = (Edge) { v, u, 0, 0, head[v] }; head[v] = e_ptr;
}
int Hd, Tl, Q[MAXN+10];
void BFS() {
Hd = 1, Tl = 0;
memset(d, 0, sizeof(d));
memset(vis, 0, sizeof(vis));
Q[++Tl] = t; vis[t] = true; d[t] = 0;
while(Hd <= Tl) {
int u = Q[Hd++];
for(int j=head[u]; j; j=E[j].next) {
int v = E[j].v; int64 f = E[j^1].flow, c = E[j^1].cap;
if(f < c && !vis[v]) {
vis[v] = true; Q[++Tl] = v;
}
}
}
}
int64 Augment() {
int64 aug = LL_INF;
for(int u = t; u != s; u = E[p[u]].u)
aug = min(aug, E[p[u]].cap - E[p[u]].flow);
for(int u = t; u != s; u = E[p[u]].u) {
E[p[u]].flow += aug;
E[p[u]^1].flow -= aug;
}
return aug;
}
int64 ISAP() {
int u = s; int64 flow = 0;
BFS(); memcpy(cur, head, sizeof(head));
for(int i = 1; i <= n; i++) ++num[d[i]];
while(d[s] < n) {
bool ok = false;
if(u == t) {
flow += Augment();
u = s;
continue;
}
for(int& j=cur[u]; j; j=E[j].next) { // Advance
int v = E[j].v; int64 f = E[j].flow, c = E[j].cap;
if(f < c && d[v] == d[u] - 1) {
p[v] = j; u = v;
ok = true;
break;
}
}
if(!ok) { // Retreat
int slk = n-1;
for(int j=head[u]; j; j=E[j].next) {
int v = E[j].v; int64 f = E[j].flow, c = E[j].cap;
if(f < c) slk = min(slk, d[v]);
}
if(--num[d[u]] == 0) break; // gap
++num[d[u] = slk + 1];
cur[u] = head[u]; // 重置当前弧
if(u != s) u = E[p[u]].u;
}
}
return flow;
}
}
template<typename T>
inline void readint(T& x) {
T f=1, r=0; char c=getchar();
while(!isdigit(c)) { if(c=='-')f=-1; c=getchar(); }
while(isdigit(c)) { r=r*10+c-'0'; c=getchar(); }
x = f*r;
}
int main() {
using namespace Maxflow;
int m, u, v; int64 c;
readint(n); readint(m); readint(s); readint(t);
for(int i = 1; i <= m; i++) {
readint(u); readint(v); readint(c);
AddEdge(u, v, c);
}
cout << ISAP();
}