这篇文章上次修改于 224 天前,可能其部分内容已经发生变化,如有疑问可询问作者。

算法学习笔记(20): 网络流费用流

简介

在基础的网络图上,每条边多了一个属性,即费用。这个费用是单位费用,即当这条边有流量 f 的时候,有 f * w 的 费用。

最小费用最大流:在所有的最大流情况下,费用最小的一个情况,称作最小费用最大流。

SSP算法

SSP算法思想为贪心。即在找最大流的时候额外多加一个费用条件。

注意: 我们需要先排除一个情况,就是最小费用流不能存在负价环,不然SSP无法解决。

复习一下原来的步骤:

  • 建立剩余流量图的分层图
  • 找到分层图的一个堵塞流(注意)
  • 在剩余流量图中进行对应的减少以及建立反向边
  • 当分层图找不到堵塞流停止

在原来的基础上,将寻找堵塞流的算法改为优先寻找单位费用最小的路径。

只需要修改原来的寻找堵塞流方法即可.

其中最小费用路径等同于最短路径,以下使用SPFA完成。

const int N = 1e4;
const int M = 1e5;
const int INF = 0x3f3f3f3f;

namespace CostFlow
{
    struct Edge {
        int next, to, capacity, weight;
    };

    Edge edges[M];
    int head[M], edgeCount = 1;

    void addPath(int from, int to, int capacity, int weight) {
        edges[++edgeCount] = (Edge){head[from], to, capacity, weight};
        head[from] = edgeCount;
    }

    void addFlow(int from, int to, int capacity, int weight) {
        addPath(from, to, capacity, weight);
        addPath(to, from, 0, -weight);
    }

    int dist[N], previous[N], flow[N];
    bool visited[N];

    int source = 0, sink = 1;

    bool spfa() {
        memset(dist, 0x3f, sizeof(dist));
        queue<int> q;
        q.push(source), dist[source] = 0, flow[source] = INF, flow[sink] = 0;
        while (!q.empty()) {
            int u = q.front();
            q.pop();
            visited[u] = false;
            for (int i = head[u]; i; i = edges[i].next) {
                const int &v = edges[i].to, &cap = edges[i].capacity, &cost = edges[i].weight;
                if (cap == 0 || dist[v] <= dist[u] + cost) continue;
                dist[v] = dist[u] + cost, flow[v] = min(cap, flow[u]), previous[v] = i;
                if (!visited[v]) q.push(v), visited[v] = true;
            }
        }
        return flow[sink];
    }

    int maxFlow = 0, minCost = 0;

    void update() {
        maxFlow += flow[sink];
        for (int u = sink; u != source; u = edges[previous[u] ^ 1].to) {
            edges[previous[u]].capacity -= flow[sink];
            edges[previous[u] ^ 1].capacity += flow[sink];
            minCost += flow[sink] * edges[previous[u]].weight;
        }
    }

    void start(int S, int T){
        source = S;
        sink = T;
        while(spfa()){
            update();
        }
    }
};

对比理解

对于原来的最大流算法,在为了摆脱 复杂度和最大流挂钩,使用了EK算法思路,将图看做无权图,并去寻找最短路径,实际这个时候,最短路径就是边最少的情况。

而后来Danic为了加速多次的这种寻找,使用了分层图,效果显著。

此时最小费用流的思路是去找最小费用的路径,也就是将图看做其单价费用作为权的 有权图最短路径。SPFA就可以完成这个步骤,同样,在找到最短路后进行增广反向。

缺点分析

复习一下SPFA,实际上是基于队列优化的贝尔曼福特算法,SPFA算法很容易被卡死,一旦遇到了菊花图,队列一下就需要加入一群点。

Primal-Dual 原始对偶算法

既然用到了最短路算法,除了作为代表的SPFA,我们应当想到dijkstra算法,但是它在本题中,本题的图频繁出现负权边,干扰了dijkstra算法的贪心正确性。

我们已经知道了费用流的基本步骤,不同于最大流的唯一之处就是寻找堵塞流。而这个寻找也就是模板问题:

在有向带负权图中找到最短路。

这也是为什么我们使用了SPFA,SPFA可以很好的解决这个问题。

但如果结合 Johnson 全源最短路径算法,又可以通过Dijstra来找到算法。

tips:推荐在 本博客的以前文章 中复习一遍再继续看。

原算法中,有一个势能的概念,而在本题中增广路会在残量中增加边权,因此,需要在每次修改势能,设增广后 s 点到 i 点的最短路径是d。

则H[ i ] += d;

const int maxN = 1e4;
const int maxE = 1e5;
const int INF = 0x3fffffff;
namespace CostFlow {
    struct Edge {
        int to, flow, cost, next;
    } edges[maxE];

    struct Parent {
        int vertex, edge;
    } parent[maxE];

    struct MyPair {
        int distance, vertex;

        bool operator<(const MyPair& other) const {
            return distance > other.distance;
        }

        MyPair(int d, int v) { distance = d, vertex = v; }
    };

    int head[maxE], dist[maxE], visited[maxN], potential[maxN];
    int vertexCount, source, sink, edgeIndex = 1, maxFlow = 0, minCost = 0;

    void addEdge(int from, int to, int flow, int cost) {
        edges[++edgeIndex].to = to;
        edges[edgeIndex].flow = flow;
        edges[edgeIndex].cost = cost;
        edges[edgeIndex].next = head[from];
        head[from] = edgeIndex;
    }

    void addFlow(int from, int to, int flow, int cost) {
        addEdge(from, to, flow, cost);
        addEdge(to, from, 0, -cost);
    }

    bool dijkstra() {
        priority_queue<MyPair> pq;
        for (int i = 1; i <= vertexCount; i++) dist[i] = INF;
        memset(visited, 0, sizeof(visited));
        dist[source] = 0;
        pq.push(MyPair(0, source));
        while (!pq.empty()) {
            int u = pq.top().vertex;
            pq.pop();
            if (visited[u]) continue;
            visited[u] = 1;
            for (int i = head[u]; i; i = edges[i].next) {
                int v = edges[i].to,
                    newCost = edges[i].cost + potential[u] - potential[v];
                if (edges[i].flow && dist[v] > dist[u] + newCost) {
                    dist[v] = dist[u] + newCost;
                    parent[v].vertex = u;
                    parent[v].edge = i;
                    if (!visited[v]) pq.push(MyPair(dist[v], v));
                }
            }
        }
        return dist[sink] != INF;
    }

    void spfa() {
        queue<int> q;
        memset(potential, 63, sizeof(potential));
        potential[source] = 0, visited[source] = 1;
        q.push(source);
        while (!q.empty()) {
            int u = q.front();
            q.pop();
            visited[u] = 0;
            for (int i = head[u]; i; i = edges[i].next) {
                int v = edges[i].to;
                if (edges[i].flow &&
                    potential[v] > potential[u] + edges[i].cost) {
                    potential[v] = potential[u] + edges[i].cost;
                    if (!visited[v]) {
                        visited[v] = 1;
                        q.push(v);
                    }
                }
            }
        }
    }

    void start(int s, int t, int n){
        vertexCount = n;
        source = s;
        sink = t;
        spfa();
        while (dijkstra()) {
            int minFlow = INF;
            for (int i = 1; i <= vertexCount; i++) potential[i] += dist[i];
            for (int i = sink; i != source; i = parent[i].vertex)
                minFlow = min(minFlow, edges[parent[i].edge].flow);
            for (int i = sink; i != source; i = parent[i].vertex) {
                edges[parent[i].edge].flow -= minFlow;
                edges[parent[i].edge ^ 1].flow += minFlow;
            }
            maxFlow += minFlow;
            minCost += minFlow * potential[sink];
        }
    }
};