基于 容量放缩( capacity scaling ) 的 弱多项式 最小费用最大流算法

· · 算法·理论

前言

若以下内容出现了谬误，烦请私信我进行指正，以及可私信提问。

可以使用 Dijkstra 求解最短路，但需要引入类似 Primal-Dual 原始对偶和 Johnson 全源最短路的势能，复杂度 O(m^2\log F \log m)，这里 F 指容量的值域。

由于引入势能后同时增加常数与码长，且 SPFA 在此算法上很难卡满，因此这里使用 SPFA 求解最短路，复杂度 O(nm^2\log F)。

使用 Dijkstra 的算法可自行查看 原文献( by ouuan )。

以及笔者偶然发现了 O(n^3 \log nC) 的 cost scaling 算法，这里 C 指费用的值域，由此得到 O(nm \log nC \log \log F) 的 double scaling，但由于笔者蒟蒻，并不会实现，原文献( by ShanLunjiaJian )。

正文

定义

$capt(u,v): $ 边 $(u,v)$ 的容量。 $cost(u,v): $ 边 $(u,v)$ 的费用。 未加特别说明的情况下，默认边 $(u,v)$ 的长度为 $cost(u,v)$，类似地定义 最短路/负环 等概念。 ### 将 最小费用最大流 转为 无源汇最小费用流 我们从 $T$ 到 $S$ 连边，满足 $capt(T,S)=\inf, cost(T,S)=-\inf$。 容易发现，新图的最小费用流一定是原图的最小费用最大流，具体可考虑反证法证明。 因此，以下内容未加特别说明的情况下，默认讨论无源汇最小费用流，且默认最小费用流即指无源汇。 ### 最小费用流 与 负环 一个流是最小费用流，当且仅当其残量网络上不存在负环。 > 证明: > > 充分性：假设当前流为最小费用流，且残量网络存在负环，则沿着负环增广可以得到更小费用的流，与假设矛盾。 > > 必要性：假设当前流 $f$ 不为最小费用流，且残量网络不存在负环，任取一个最小费用流 $f'$，由于 $f,f'$ 都是合法流，则流 $f'-f$ 也合法且费用为负，那么一定存在负环，又因为 $f'-f$ 相当于在 $f$ 上增广得到 $f'$，与假设矛盾。 ### 容量缩放 ( capacity scaling ) 容量缩放考虑这样一个算法流程： 从高到低考虑容量在二进制下的所有位，$O(\log V)$ 轮迭代，每次迭代加入当前位的下一位后更新答案。 例如一条容量为 13 的边，在 4 轮迭代中容量依次为 1,3,6,13。 算法基于这样一个事实：将图中每条边的容量 *2 后，新图中最小费用流每条边的流量也 *2，容易通过反证法证明。 那么，我们只需要考虑每次将图中一条边的容量加一，更新最小费用流即可。 上文证明了最小费用流的残量网络不存在负环，所以我们只在产生负环时增广。 记当前加入的边为 $(u,v)$，为了减少常数，若当前边已在残量网络，直接将其剩余容量 +1 即可。 否则，新出现的负环一定包含 $(u,v)$，由于加入 $(u,v)$ 前当前不存在负环，我们求出 $v \rightarrow u$ 的最短路，判断若 $dis_{v,u}+cost(u,v)<0$，则存在负环进行增广。 容易发现，我们需要进行 $O(m\log F)$ 次最短路求解，使用 $O(nm)$ 的 SPFA 即为 $O(nm^2 \log F)

代码 ( SPFA )

值得一提的是，这份代码并不能完全通过洛谷的模板，提交记录，但可以通过 UOJ 的模板。

同时，我在模板题测试了原文献的代码 Dijkstra 与 SPFA ，均不能通过，且 SPFA 表现比 Dijkstra 良好。


#include<bits/stdc++.h>
using namespace std;

typedef long long LL;
constexpr int N = 5005, M = 1e5+5;
constexpr LL INF = 0x3f3f3f3f3f3f3f3f;

int n, m;
int he[N], ne[M], to[M];
LL rca[M], ca[M], co[M];
int pre[N], q[N], inq[N];
LL dis[N];

void add(int a, int b, LL c1, LL c2)
{
    static int idx=2;
    rca[idx] = c1, ca[idx] = 0, co[idx] = c2;
    to[idx] = b, ne[idx] = he[a], he[a] = idx ++;
}

LL spfa(int s, int t)
{
    int hh = 0, tt = 1;
    memset(dis + 1, INF, sizeof(LL[n]));
    q[0] = s, dis[s] = 0, pre[s] = 0, inq[s] = 1;

    while(hh != tt)
    {
        int u = q[hh ++];
        if(hh == N) hh = 0;
        inq[u] = 0;

        for(int i = he[u]; i; i = ne[i])
        {
            int v = to[i];
            if(ca[i] && dis[v] > dis[u] + co[i])
            {
                dis[v] = dis[u] + co[i], pre[v] = i;
                if(inq[v]) continue;
                q[tt ++] = v, inq[v] = 1;
                if(tt == N) tt = 0;
            }
        }
    }

    return dis[t];
}

int main()
{
    int S, T;
    scanf("%d%d%d%d", &n, &m, &S, &T);
    for(int i = 1; i <= m; ++ i)
    {
        int a, b, c1, c2;
        scanf("%d%d%d%d", &a, &b, &c1, &c2);
        add(a, b, c1, c2), add(b, a, 0, -c2);
    }
    add(T, S, INF, -INF), add(S, T, 0, INF), ++ m;

    for(int i = 30; i >= 0; -- i)
    {
        for(int j = 2; j <= (m << 1 | 1); ++ j)ca[j]<<=1;
        for(int j = 2; j <= (m << 1 | 1); j += 2)
        {
            if(~rca[j] >> i & 1) continue;
            if(ca[j]) {++ ca[j]; continue;}

            int u = to[j ^ 1], v = to[j];
            if(spfa(v, u) + co[j] >= 0) ++ ca[j];
            else
            {
                ++ ca[j ^ 1];
                while(pre[u])
                {
                    -- ca[pre[u]];
                    ++ ca[pre[u] ^ 1];
                    u = to[pre[u] ^ 1];
                }
            }
        }
    }

    LL res = 0;
    for(int i = 3; i < m<<1; i += 2) res += (rca[i] - ca[i]) * co[i];
    printf("%lld %lld", ca[m << 1 | 1], res);

    return 0;
}