费用流去掉负权边

消去网络中负权边的方法。

首先不能给边势能函数，因为最短路的路径不一定一致。于是在点上做文章，给每个点一个势能函数 \(h(x)\)，满足 \(\forall (x,y)\in E,s.t.h(x)-h(y)+w(x,y)\geqslant0\)，这样跑出来的最短路多的权就是 \(h(s)-h(t)\)。

至于构造方法，每次增广完以后，如果 \((x,y)\in E_r\)，则 \(y\) 一定是 \(x\) 的 successor（不然就不在增广路上），也就是说 \(d_x+w(x,y)+h(x)-h(y)=d_y\)，其中 \(d\) 是 \(s\) 到该点的最短路，可得 \((d_y+h(y))-(d_x+h(x))=w(x,y)\)。

那么令 \(h(x):=h(x)+d(x)\) 即可。

因为不想写指针所以用了 std::vector 管理边的内存……

template::max(),
    const type2 inf2 = std::numeric_limits::max()>
struct mcf_graph {
    const int n;
    struct edge {
        int to, nxt;
        type1 c; type2 f;
        edge(int to, int nxt, type1 c, type2 f) : to(to), nxt(nxt), c(c), f(f) {}
    };
    int *head, *pre, *id, *vis;
    std::vector e;
    type2 *h, *dis;
    mcf_graph(const int n) : n(n),
        head(new int[n+1]), pre(new int[n+1]), id(new int[n+1]), vis(new int[n+1]),
        h(new type2[n+1]), dis(new type2[n+1]) { std::memset(head, -1, (n+1)<<2); }
    void add_edge(const int x, const int y, type1 c, type2 f) {
        assert(1 <= x && x <= n && 1 <= y && y <= n);
        e.emplace_back(y, head[x], c, f),head[x] = int(e.size())-1;
        e.emplace_back(x, head[y], 0, -f),head[y] = int(e.size())-1;
    }
    void initialize_potentials(const int s) {
        std::queue q;
        std::fill(h+1, h+n+1, inf2);
        std::memset(vis+1, 0, n<<2);
        for(q.push(s), vis[s] = 1, h[s] = 0; q.size(); vis[q.front()] = 0, q.pop()) {
            for(int now = q.front(), i = head[now], y; ~i; i = e[i].nxt) {
                if(e[i].c && (h[y = e[i].to] == inf2 || h[y] > h[now]+e[i].f)) {
                    h[y] = h[now]+e[i].f;
                    if(!vis[y]) q.push(y),vis[y] = 1;
                }
            }
        }
    }
    bool augment(const int s, const int t) {
        std::priority_queue, std::vector>, std::greater>> q;
        std::fill(dis+1, dis+n+1, inf2);
        std::memset(vis+1, 0, n<<2);
        for(q.emplace(dis[s] = 0, s); q.size();) {
            int now = q.top().second;
            q.pop();
            if(vis[now]) continue;
            vis[now] = 1;
            for(int i = head[now], y; ~i; i = e[i].nxt) {
                if(e[i].c && (dis[y = e[i].to] == inf2 || dis[y] > dis[now]+e[i].f+h[now]-h[y])) {
                    dis[y] = dis[now]+e[i].f+h[now]-h[y];
                    pre[y] = now,id[y] = i;
                    if(!vis[y]) q.emplace(dis[y], y);
                }
            }
        }
        return dis[t] < inf2;
    }
    std::pair get(const int s, const int t) {
        assert(1 <= s && s <= n && 1 <= t && t <= n);
        type1 res1 = 0; type2 res2 = 0;
        initialize_potentials(s);
        while(augment(s, t)) {
            type1 cur = inf1;
            for(int i = 1; i <= n; ++i) h[i] += dis[i];
            for(int i = t; i != s; i = pre[i]) cmin(cur, e[id[i]].c);
            for(int i = t; i != s; i = pre[i]) e[id[i]].c -= cur,e[id[i]^1].c += cur;
            res1 += cur,res2 += cur*h[t];
        }
        return std::make_pair(res1, res2);
    }
};