一些简单的模拟费用流问题

---

祝贺我又学了一个学不明白的东西 (

综述

模拟费用流, 大概是利用流的一些性质, 从而用数据结构高效模拟费用流.

基础理论和模型可以看看文末的参考资料, 这里只是记录我做过的一些题目.

例题

题目顺序随心情.

---

首先是一些 “模拟” 费用流的例子.

CF280D k-Maximum Subsequence Sum

题目链接

• https://codeforces.com/problemset/problem/280/D

解题思路

首先有一个费用流做法. 设源汇点 $S$, $T$, 同时将每个位置拆为入点 $i$ 和出点 $i’$.

1. 连边 $(S, S’)$, 容量为 $k$, 费用为 $0$.
2. 连边 $(S’, i)$, 容量为 $1$, 费用为 $0$.
3. 连边 $(i, i’)$, 容量为 $1$, 费用为 $a_i$.
4. 连边 $(i’, i + 1)$, 容量为 $1$, 费用为 $0$.
5. 连边 $(i’, T)$, 容量为 $1$, 费用为 $0$.

此时将费用取反, 求最小费用流即可. 可惜时间复杂度不够优秀.

注意到每次增广的过程, 本质上是选择序列一段区间中的最大连续子段和, 直接用线段树维护即可.

同时需要维护反向边, 也就是将该区间翻转并将权值取反. 那么记录子段和时需记录该子段对应区间端点, 以及最小子段和, 用于权值取反后快速计算现有区间最大子段和. 这样区间取反可打标记维护.

根据这个思路, 单点修改就没什么意思了. 注意每组询问直接独立, 在区间取反后需要撤销影响.

时间复杂度 $O(n \log n)$.

代码实现

// CF280D
// DeP
#include <cctype>
#include <cstdio>
#include <algorithm>
using namespace std;

namespace IO {
  const int MAXSIZE = 1 << 18 | 1;
  char buf[MAXSIZE], *p1, *p2;

  inline int Gc() {
    return p1 == p2 &&
      (p2 = (p1 = buf) + fread(buf, 1, MAXSIZE, stdin), p1 == p2)? EOF: *p1++;
  }
  template<typename T> inline void read(T& x) {
    x = 0; int f = 0, ch = Gc();
    while (!isdigit(ch)) f |= ch == '-', ch = Gc();
    while (isdigit(ch)) x = x * 10 + ch - '0', ch = Gc();
    if (f) x = -x;
  }
}
using IO::read;

const int MAXN = 1e5 + 5, MAXM = MAXN << 2;

int n, q, top;
int A[MAXN];

struct Item {
  int L, R, v;
  Item() { L = MAXN, R = -MAXN, v = 0; }
  Item(int _L, int _R, int _v): L(_L), R(_R), v(_v) { }
  Item operator + (const Item& rhs) const {
    return Item(min(L, rhs.L), max(R, rhs.R), v + rhs.v);
  }
  bool operator < (const Item& rhs) const {
    return v < rhs.v;
  }
} stk[MAXN];

#define lc (nd << 1)
#define rc (nd << 1 | 1)
namespace SGT {
  struct Node {
    Item s, MxL, MxR, MxM, MnL, MnR, MnM;
    Node() { }
    Node(Item v) { s = MxL = MxR = MxM = MnL = MnR = MnM = v; }
    Node operator + (const Node& rhs) const {
      Node ret;
      ret.s = s + rhs.s;

      ret.MxL = max(MxL, s + rhs.MxL);
      ret.MxR = max(rhs.MxR, rhs.s + MxR);
      ret.MxM = max(MxR + rhs.MxL, max(MxM, rhs.MxM));

      ret.MnL = min(MnL, s + rhs.MnL);
      ret.MnR = min(rhs.MnR, rhs.s + MnR);
      ret.MnM = min(MnR + rhs.MnL, min(MnM, rhs.MnM));
      return ret;
    }
  } dat[MAXM];
  bool res[MAXM];

  inline void maintain(int nd) { dat[nd] = dat[lc] + dat[rc]; }

  inline void push(int nd) {
    swap(dat[nd].MxM, dat[nd].MnM);
    swap(dat[nd].MxL, dat[nd].MnL), swap(dat[nd].MxR, dat[nd].MnR);
    dat[nd].s.v *= -1;
    dat[nd].MxL.v *= -1, dat[nd].MxR.v *= -1, dat[nd].MxM.v *= -1;
    dat[nd].MnL.v *= -1, dat[nd].MnR.v *= -1, dat[nd].MnM.v *= -1;
    res[nd] ^= 1;
  }
  inline void pushdown(int nd) {
    if (res[nd]) push(lc), push(rc), res[nd] = 0;
  }

  void build(int nd, int L, int R) {
    if (L == R)
      return void( dat[nd] = Item(L, R, A[L]) );
    int Mid = (L + R) / 2;
    build(lc, L, Mid), build(rc, Mid + 1, R);
    maintain(nd);
  }

  void Rev(int nd, int L, int R, const int& opL, const int& opR) {
    if (opL <= L && R <= opR) return push(nd);
    int Mid = (L + R) / 2;
    pushdown(nd);
    if (opL <= Mid) Rev(lc, L, Mid, opL, opR);
    if (opR > Mid) Rev(rc, Mid + 1, R, opL, opR);
    maintain(nd);
  }

  void Mdy(int nd, int L, int R, const int& p, const int& v) {
    if (L == R)
      return void( dat[nd] = Item(L, R, v) );
    pushdown(nd);
    int Mid = (L + R) / 2;
    if (p <= Mid) Mdy(lc, L, Mid, p, v); else Mdy(rc, Mid + 1, R, p, v);
    maintain(nd);
  }

  Node Qry(int nd, int L, int R, const int& opL, const int& opR) {
    if (opL <= L && R <= opR) return dat[nd];
    int Mid = (L + R) / 2;
    pushdown(nd);
    if (opR <= Mid) return Qry(lc, L, Mid, opL, opR);
    if (opL > Mid) return Qry(rc, Mid + 1, R, opL, opR);
    return Qry(lc, L, Mid, opL, opR) + Qry(rc, Mid + 1, R, opL, opR);
  }
}
#undef lc
#undef rc

int main() {
#ifndef ONLINE_JUDGE
  freopen("input.in", "r", stdin);
#endif
  read(n);
  for (int i = 1; i <= n; ++i) read(A[i]);

  SGT::build(1, 1, n);

  read(q);
  for (int opt, L, R, k, ans; q; --q) {
    read(opt), read(L);
    switch (opt) {
      case 0:
        read(k), SGT::Mdy(1, 1, n, L, k);
        break;
      case 1:
        read(R), read(k), ans = 0;
        while (k--) {
          Item v = SGT::Qry(1, 1, n, L, R).MxM;
          if (v.v <= 0) break;
          ans += v.v, SGT::Rev(1, 1, n, v.L, v.R), stk[++top] = v;
        }
        while (top > 0) {
          Item v = stk[top--];
          SGT::Rev(1, 1, n, v.L, v.R);
        }
        printf("%d\n", ans);
        break;
      default: fprintf(stderr, "ERR\n");
    }
  }
  return 0;
}

Luogu P6122「NEERC2016」Mole Tunnels

题目链接

• https://www.luogu.com.cn/problem/P6122

解题思路

首先有一个费用流思路. 设源汇点分别为 $S$, $T$, 将每只鼹鼠视作流量.

1. 对于树上每条边, 连边 $(u, \lfloor \frac{u}{2} \rfloor)$, $(\lfloor \frac{u}{2} \rfloor, u)$, 容量为 $\infty$, 费用为 $1$.
2. 对于树上每个点 $u$, 连边 $(u, T)$, 容量为 $c_u$, 费用为 $0$.
3. 对于每只鼹鼠, 依次向其树上所在位置 $u$ 连边 $(S, u)$, 容量为 $1$, 费用为 $0$.

观察每次增广的过程, 都是选择树上一条路径, 为保证费用尽量小, 显然是在要求路径尽量短.

考虑到题中给定树为二叉树, 其深度为 $O(\log n)$. 枚举新增鼹鼠位置 $s$ 和最终匹配位置的 LCA $t$. 同时对每个节点 $u$ 记录位置 $p(u)$, 表示 $u$ 子树内同 $u$ 匹配的最优位置.

此时增广的路径即 $s \rightarrow t \rightarrow {p}(t)$, 更新流量计算费用即可.

${p}(u)$ 可从下向上递推计算. 另外需记录树边中的流量. 为了便于计算路径费用, 将连向祖先的流量视作正数, 相反则视作负数.

时间复杂度 $O(n \log n)$.

代码实现

// Luogu P6122
// DeP
#include <cstdio>
#include <cstring>
using namespace std;

const int MAXN = 1e5 + 5;
const int INF = 0x3f3f3f3f;

template<typename T> inline bool ckmin(T& a, const T& b) {
  return (a > b)? a = b, true: false;
}

int n, m;
int C[MAXN], flow[MAXN], d[MAXN], p[MAXN];
// flow: u --> pre(u) 流量, d: u --> 匹配位置 距离

inline void augment(int u) {
  int lc = u << 1, rc = u << 1 | 1;
  d[u] = INF, p[u] = 0;
  if (C[u] > 0) d[u] = 0, p[u] = u;
  if (lc <= n && ckmin(d[u], d[lc] + ((flow[lc] > 0)? -1: 1)))
    p[u] = p[lc];
  if (rc <= n && ckmin(d[u], d[rc] + ((flow[rc] > 0)? -1: 1)))
    p[u] = p[rc];
}

int main() {
#ifndef ONLINE_JUDGE
  freopen("input.in", "r", stdin);
#endif
  scanf("%d%d", &n, &m);
  for (int i = 1; i <= n; ++i) scanf("%d", C + i);

  memset(d, 0x3f, sizeof d);
  for (int u = n; u; --u) augment(u);

  int ans = 0;
  for (int k = 1; k <= m; ++k) {
    int s, t = 0, cost = INF, ds = 0;
    scanf("%d", &s);
    for (int u = s; u > 0; u /= 2) {
      if (ckmin(cost, ds + d[u])) t = u;
      ds += ((flow[u] < 0)? -1: 1);
    }
    ans += cost;
    for (int u = s; u != t; u /= 2) ++flow[u];
    for (int u = p[t]; u != t; u /= 2) --flow[u];
    --C[p[t]];
    for (int u = p[t]; u != t; u /= 2) augment(u);
    for (int u = s; u > 0; u /= 2) augment(u);
    printf("%d%c", ans, " \n"[k == m]);
  }
  return 0;
}

HDU 6634「2019 Multi-University Training Contest」Salty fish

题目链接

• http://acm.hdu.edu.cn/showproblem.php?pid=6634

解题思路

考虑最大权闭合子图的建图方式. 设源汇点分别为 $S$, $T$.

1. 对于树上每个节点 $i$, 连边 $(S, i)$, 容量为 $a_i$.
2. 对于每个摄像头 $j$, 连边 $(j, T)$, 容量为 $c_j$.
3. 对于摄像头 $j$, 若节点 $i$ 在 $j$ 的监控范围内, 连边 $(i, j)$, 容量为 $\infty$.

此时答案为 $\sum {a_i} - \text{mincut}$.

根据最大流-最小割定理, 我们只需要计算出网络的最大流即可. 对于节点 $i$, 其流量一定流向能够监控到 $i$ 的, 且深度尽量小的摄像头. 因为深度大的摄像头更可能监控到深度大的节点. 换言之, 每个摄像头一定选择其范围内深度尽量大的节点.

那么我们可以对每个节点 $u$, 记录 $f(u, i)$ 表示 $u$ 子树内, 同根节点距离为 $i$ 节点的 $a$ 之和. 由深度从大到小的顺序选择每个摄像头, 从满足要求的最大深度开始向前推进, 更新流量和答案即可.

注意到 $f(u, i)$ 和深度有关, 可利用长链剖分维护. 此时需要在一个数组中快速修改 / 删除 / 加入一个数, 用 map 维护 $f(u, i)$ 就好了. 当然也可以写线段树合并.

时间复杂度 $O\big((n + m) \log n\big)$.

代码实现

// HDU 6634
// DeP
#include <map>
#include <cstdio>
#include <vector>
#include <cstring>
using namespace std;

typedef long long LL;
typedef map<int, LL>::iterator IT;

const int MAXN = 3e5 + 5;

int n, m;
int pre[MAXN], A[MAXN], K[MAXN], C[MAXN];
vector<int> P[MAXN];

LL ans;
int Idx[MAXN], idx;
map<int, LL> M[MAXN]; // depth, flow

namespace Graph {
  struct Edge { int nxt, to; } edges[MAXN << 1];
  int head[MAXN], eidx;

  inline void init() { memset(head, -1, sizeof head), eidx = 1; }
  inline void AddEdge(int from, int to) {
    edges[++eidx] = (Edge){ head[from], to }, head[from] = eidx;
  }

  int depth[MAXN], d[MAXN], son[MAXN];

  void dfs(int u) {
    d[u] = 1, son[u] = -1;
    depth[u] = depth[pre[u]] + 1;
    for (int v, i = head[u]; ~i; i = edges[i].nxt) {
      dfs(v = edges[i].to);
      if (son[u] == -1 || d[son[u]] < d[v])
        son[u] = v, d[u] = d[v] + 1;
    }
    Idx[u] = (~son[u])? Idx[son[u]]: ++idx;
    map<int, LL>& f = M[Idx[u]];
    f[depth[u]] += A[u];
    for (int v, i = head[u]; ~i; i = edges[i].nxt) {
      if ((v = edges[i].to) == son[u]) continue;
      for (IT ite = M[Idx[v]].begin(); ite != M[Idx[v]].end(); ++ite)
        f[ite->first] += ite->second;
    }
    for (size_t k = 0; k < P[u].size(); ++k) {
      const int& i = P[u][k];
      IT ite = f.upper_bound(depth[u] + K[i]);
      if (ite == f.begin()) continue;
      --ite;
      while (C[i] > 0) {
        LL w = min((LL) C[i], ite->second);
        ite->second -= w, C[i] -= w, ans -= w;
        if (ite->second > 0) continue;
        if (ite == f.begin()) { f.erase(ite); break; }
        int dep = ite->first;
        --ite, f.erase(dep);
      }
    }
  }
}

int main() {
#ifndef ONLINE_JUDGE
  freopen("input.in", "r", stdin);
#endif
  int Ti;
  scanf("%d", &Ti);
  while (Ti--) {
    for (int i = 1; i <= n; ++i) P[i].clear();
    for (int i = 1; i <= idx; ++i) M[i].clear();
    ans = idx = 0, Graph::init();

    scanf("%d%d", &n, &m);
    for (int i = 2; i <= n; ++i)
      scanf("%d", pre + i), Graph::AddEdge(pre[i], i);
    for (int i = 1; i <= n; ++i)
      scanf("%d", A + i), ans += A[i];
    for (int u, i = 1; i <= m; ++i)
      scanf("%d%d%d", &u, K + i, C + i), P[u].push_back(i);

    Graph::dfs(1);
    printf("%lld\n", ans);
  }
  return 0;
}

UOJ #318「NOI2017」蔬菜

题目链接

• https://uoj.ac/problem/318

解题思路

首先有一个朴素的网络流思路. 设源汇点分别为 $S$, $T$, 同时将每天视作一个节点, 记 $p$ 为最大天数.

1. 对于每个蔬菜 $i$, 记其最大售卖天数 $k$ 为 $\min\{ p, \lceil \frac{c_i}{x_i} \rceil \}$.
   • 向第 $1 \sim k - 1$ 天依次连边, 容量为 $x_i$, 费用为 $a_i$.
   • 向第 $k$ 天连边, 容量为 $c_i - k x_i - 1$, 费用为 $a_i$.
   • 向第 $k$ 天连边, 容量为 $1$, 费用为 $a_i + s_i$
2. 对于第 $i$ 天
   • 连边 $(i, i - 1)$, 容量为 $\infty$, 费用为 $0$.
   • 连边 $(i, T)$, 容量为 $m$, 费用为 $0$.

观察每次增广时的情况. 在天数不断增加时, 原有流量改变一定不优, 否则在前一天可改变流量使得前一天答案更优. 也就是说, 第 $i$ 天售出的蔬菜一定是第 $i + 1$ 天售出蔬菜的一部分, 且网络中不存在退流.

同时, 注意到天数对应点之间的连边, 对于每个蔬菜, 尽量选择靠后天数是更优秀的. 从建图解释, 选择靠后天数可通过容量为 $\infty$ 的边增广. 从实际意义解释, 将蔬菜尽量放在快过期的时候卖出更优. 同时要选择权值尽量大的蔬菜卖出.

那么用堆维护当前可选的蔬菜, 每天的销售情况视作增广 $m$ 的流量, 选择权值最大的蔬菜即可. 同时需要知道当前天对应节点在向 $T$ 的边满流之后, 向前增广的位置. 可以利用并查集维护.

时间复杂度 $O(n \log n + pm\log n)$.

代码实现

// UOJ #318
// DeP
#include <queue>
#include <cctype>
#include <cstdio>
#include <algorithm>
using namespace std;

namespace IO {
  const int MAXSIZE = 1 << 18 | 1;
  char buf[MAXSIZE], *p1, *p2;

  inline int Gc() {
    return p1 == p2 &&
      (p2 = (p1 = buf) + fread(buf, 1, MAXSIZE, stdin), p1 == p2)? EOF: *p1++;
  }
  template<typename T> inline void read(T& x) {
    x = 0; int f = 0, ch = Gc();
    while (!isdigit(ch)) f |= ch == '-', ch = Gc();
    while (isdigit(ch)) x = x * 10 + ch - '0', ch = Gc();
    if (f) x = -x;
  }
}
using IO::read;

typedef long long LL;
const int MAXN = 1e5 + 5;

int n, m, q;
int A[MAXN], X[MAXN], S[MAXN], C[MAXN], Q[MAXN];

LL Ans[MAXN];
int cnt[MAXN];

struct Item {
  int idx, v;
  Item(int _i, int _v): idx(_i), v(_v) { }
  bool operator < (const Item& rhs) const {
    return v < rhs.v;
  }
};
priority_queue<Item> PQ;

namespace DSU {
  int fa[MAXN];

  inline void init(int _n) {
    for (int i = 0; i <= _n; ++i) fa[i] = i;
  }
  int findfa(int u) {
    return (fa[u] == u)? u: fa[u] = findfa(fa[u]);
  }
}

int main() {
#ifndef ONLINE_JUDGE
  freopen("data/ex_vegetables3.in", "r", stdin);
#endif
  read(n), read(m), read(q);
  for (int i = 1; i <= n; ++i)
    read(A[i]), read(S[i]), read(C[i]), read(X[i]);
  int Mx = 0;
  for (int i = 1; i <= q; ++i)
    read(Q[i]), Mx = max(Mx, Q[i]);

  DSU::init(Mx);
  for (int i = 1; i <= n; ++i)
    PQ.push(Item(i, A[i] + S[i]));

  LL ans = 0;
  for (int p = 1; p <= Mx; ++p) {
    int flow = m;
    while (flow > 0 && !PQ.empty()) {
      int i = PQ.top().idx, v = PQ.top().v; PQ.pop();
      int k = (X[i] == 0)? Mx: min(Mx, (C[i] + X[i] - 1) / X[i]);
      k = DSU::findfa(k);
      if (k == 0) continue;
      ++cnt[k];
      if (cnt[k] == m) DSU::fa[k] = DSU::findfa(k - 1);
      --flow, --C[i], ans += v;
      if (C[i] > 0) PQ.push(Item(i, A[i]));
    }
    Ans[p] = ans;
  }

  for (int i = 1; i <= q; ++i)
    printf("%lld\n", Ans[Q[i]]);
  return 0;
}

UOJ #480「NOI2019」序列

题目链接

• https://uoj.ac/problem/480

解题思路

首先有一个朴素的费用流思路. 设源汇点分别为 $S$, $T$, 对两序列的每个位置分别建点 $i$, $i’$.

1. 新建节点 $S’$, 连边 $(S, S’)$, 容量为 $K$, 费用为 $0$.
2. 对于序列的每个位置 $i$, 连边
   • $(S, i)$, 容量为 $1$, 费用为 $a_i$.
   • $(S, i’)$, 容量为 $1$, 费用为 $b_i$.
   • $(i, i’)$, 容量为 $1$, 费用为 $0$.
3. 新建节点 $C$, $D$. 用于描述 “至少有 $L$ 个下标在两个序列中都被指定” 的限制.
   • 连边 $(C, D)$, 容量为 $K - L$, 费用为 $0$.
   • 对于每个位置 $i$, 连边 $(i, C)$, $(D, i’)$, 容量为 $1$, 费用为 $0$.

现在考虑如何模拟该费用流.

存在至多 $K - L$ 个位置不受在两个序列中都被指定的限制. 此时一定先选择 $a_i$, $b_i$ 中权值大的位置, 体现在网络中即边 $(C, D)$ 满流.

现在还剩下 $L$ 个位置需要选择, 即剩下 $L$ 的流量需要增广. 记录每个位置 $i$ 被选择的情况 $f_i$, 其二进制第 $0$ 位表示 $a_i$ 是否被选择, 或是边 $(S, i)$ 的流量. 第 $1$ 位则表示 $b_i$ 的情况.

为使总权值最大, 每次需要选择权值尽量大的位置. 那么利用数个大根堆 (r0, r1, r), 分别维护恰在序列 $a_i$ (r0) 或 $b_i$ (r1) 中没有选择的位置, 以及都没有选择的位置 (r). 对于前两种情况, 在选择位置时需要另一序列某一元素配对, 再利用两个堆 (p0, p1) 记录可配对的位置.

首先单独考虑 $f_i = 3$ 的情况, 此时简直是理想状态, 满足权值尽量大的同时也满足了限制.

为使总权值尽量大, 则需要将边 $(C, D)$ 的某个流量改流在形同 $(i, i’)$ 的边上. 那么记录此种情况的流量总数 flow, 在允许此种增广时选择一组流经 $(C, D)$ 的流量即可. 也就是在 r0, r1 中选择两个元素.

剩下的情况无外乎 $f_i$ 为 $0, 1, 2$. 根据其特性选择增广路, 也就是不同的位置匹配, 记录其权值并取权值最大的情况即可.

细节可能有些多, 还是建议想清楚之后再动键盘. 注意要根据 $f_i$ 的取值删除堆中的不合法元素, 以及及时更新 $f_i$ 和 flow.

时间复杂度 $O(T \cdot n \log n)$.

代码实现

// UOJ #480
// DeP
#include <queue>
#include <cctype>
#include <cstdio>
#include <algorithm>
using namespace std;

namespace IO {
  const int MAXSIZE = 1 << 18 | 1;
  char buf[MAXSIZE], *p1, *p2;

  inline int Gc() {
    return p1 == p2 &&
      (p2 = (p1 = buf) + fread(buf, 1, MAXSIZE, stdin), p1 == p2)? EOF: *p1++;
  }
  template<typename T> inline void read(T& x) {
    x = 0; int f = 0, ch = Gc();
    while (!isdigit(ch)) f |= ch == '-', ch = Gc();
    while (isdigit(ch)) x = x * 10 + ch - '0', ch = Gc();
    if (f) x = -x;
  }
}
using IO::read;

typedef long long LL;
const int MAXN = 2e5 + 5;

int n, K, L;
int A[MAXN], B[MAXN];

struct ItemA {
  int idx;
  ItemA(int _i): idx(_i) { }
  bool operator < (const ItemA& rhs) const {
    return A[idx] < A[rhs.idx];
  }
};

struct ItemB {
  int idx;
  ItemB(int _i): idx(_i) { }
  bool operator < (const ItemB& rhs) const {
    return B[idx] < B[rhs.idx];
  }
};

struct ItemAB {
  int idx;
  ItemAB(int _i): idx(_i) { }
  bool operator < (const ItemAB& rhs) const {
    return A[idx] + B[idx] < A[rhs.idx] + B[rhs.idx];
  }
};

int idx[MAXN], f[MAXN];

int main() {
#ifndef ONLINE_JUDGE
  freopen("data/ex_sequence3.in", "r", stdin);
#endif
  int Ti; read(Ti);
  while (Ti--) {
    read(n), read(K), read(L);
    for (int i = 1; i <= n; ++i) read(A[i]);
    for (int i = 1; i <= n; ++i) read(B[i]);

    priority_queue<ItemA> r0, p0;
    priority_queue<ItemB> r1, p1;
    priority_queue<ItemAB> r;

    LL ans = 0;
    for (int i = 1; i <= n; ++i) f[i] = 0, idx[i] = i;
    sort(idx + 1, idx + 1 + n, [](const int& x, const int& y) {
        return A[x] > A[y];
      });
    for (int i = 1; i <= K - L; ++i)
      ans += A[idx[i]], f[idx[i]] |= 1;
    sort(idx + 1, idx + 1 + n, [](const int& x, const int& y) {
        return B[x] > B[y];
      });
    for (int i = 1; i <= K - L; ++i)
      ans += B[idx[i]], f[idx[i]] |= 2;

    int flow = 0;
    for (int i = 1; i <= n; ++i) {
      switch (f[i]) {
        case 0: r0.push(i), r1.push(i), r.push(i); break;
        case 1: r1.push(i), p1.push(i); break;
        case 2: r0.push(i), p0.push(i); break;
        case 3: ++flow; break;
      }
    }

    while (L--) {
      while (!r.empty() && f[r.top().idx] > 0) r.pop();
      while (!r0.empty() && (f[r0.top().idx] & 1)) r0.pop();
      while (!p0.empty() && (f[p0.top().idx] ^ 2)) p0.pop();
      while (!r1.empty() && (f[r1.top().idx] & 2)) r1.pop();
      while (!p1.empty() && (f[p1.top().idx] ^ 1)) p1.pop();

      if (flow > 0) {
        int i = r0.top().idx, j = r1.top().idx;
        ans += A[i] + B[j], --flow, f[i] |= 1, f[j] |= 2;
        if (f[i] != 3) p1.push(i);
        if (f[j] != 3) p0.push(j);
        flow += (i == j)? 1: ((f[i] == 3) + (f[j] == 3));
      } else {
        static int i1, j1, i2, j2, p3;
        LL Mx, v1 = 0, v2 = 0, v3 = 0;
        if (!r0.empty() && !p1.empty())
          i1 = r0.top().idx, j1 = p1.top().idx, v1 = A[i1] + B[j1];
        if (!p0.empty() && !r1.empty())
          i2 = p0.top().idx, j2 = r1.top().idx, v2 = A[i2] + B[j2];
        if (!r.empty())
          p3 = r.top().idx, v3 = A[p3] + B[p3];
        Mx = max(v1, max(v2, v3)), ans += Mx;
        if (v1 == Mx) {
          f[i1] |= 1, f[j1] |= 2;
          if (f[i1] != 3) p1.push(i1); else ++flow;
        } else if (v2 == Mx) {
          f[i2] |= 1, f[j2] |= 2;
          if (f[j2] != 3) p0.push(j2); else ++flow;
        } else if (v3 == Mx)
          r.pop(), f[p3] |= 3;
      }
    }

    printf("%lld\n", ans);
  }
  return 0;
}

---

以下就是 laofu 课件中 “老鼠进洞” 的模型了. 或者说是, 二分图带权匹配的模型.

考虑到二分图带权匹配已经是经典问题, 下文就不再赘述朴素的网络流做法.

laofu 课件中利用较多篇幅, 从 DP 凸性的角度入手解决了这些问题… 作为没有梦想的咸鱼选手, 我暂且还是选择感性理解的方式.

BZOJ 4977「Lydsy1708 月赛」跳伞求生

题目链接

• http://lydsy.com/JudgeOnline/problem.php?id=4977
• https://darkbzoj.tk/problem/4977

解题思路

laofu 课件中的 Problem 3.

首先思考简化版的问题.

数轴上有 $n$ 只老鼠和 $m$ 个老鼠洞. 第 $i$ 只老鼠的坐标为 $x_i$, 第 $j$ 个老鼠洞的坐标为 $y_i$.
现在要让老鼠进洞, 也就是寻找老鼠和洞的匹配.

假定每只老鼠只能往左走, 每个洞容量为 $1$. 求所有老鼠都进洞的最小总代价（即行走的最小总距离）.

也就是 laofu 课件中的 Problem 1.

首先将洞的权值视作 $-y_i$, 老鼠的权值视作 $x_i$, 所求即最小权匹配. 将数轴上的所有点排序, 那么用栈记录所有的未匹配洞的位置, 在遇到老鼠时选择栈顶的洞匹配上即可.

再来考虑这个问题. 此时每个洞带权, 且不要求所有老鼠都进洞, 所求为最大权匹配.

对于一组匹配, 其对答案的贡献为 $x_i - y_j + w_j$. 考虑到 $w_j$ 的影响, 对所有点按坐标排序后洞的权值并不单调. 那么将栈换作堆维护最值即可.

但直接选择的匹配并不一定为最优解. 联系费用流中的退流操作, 在选择一组匹配后向堆中增加权值为 $-x_i$ 的元素, 从而体现匹配更改的情况.

时间复杂度 $O(n \log n)$.

代码实现

// BZOJ #4977
// DeP
#include <queue>
#include <cstdio>
#include <algorithm>
using namespace std;

typedef long long LL;
const int MAXN = 1e5 + 5;

struct Point {
  int x, w;
  Point() = default;
  Point(int _x, int _w): x(_x), w(_w) { }
  bool operator < (const Point& rhs) const {
    return x < rhs.x || (x == rhs.x && w < rhs.w);
  }
} P[MAXN << 1];

int n, m, p;
priority_queue<int> d;

int main() {
#ifndef ONLINE_JUDGE
  freopen("input.in", "r", stdin);
#endif
  scanf("%d%d", &n, &m);
  for (int a, i = 1; i <= n; ++i)
    scanf("%d", &a), P[++p] = Point(a, 0);
  for (int b, c, i = 1; i <= m; ++i)
    scanf("%d%d", &b, &c), P[++p] = Point(b, c);

  sort(P + 1, P + 1 + p);
  LL ans = 0;
  for (int i = 1; i <= p; ++i) {
    if (P[i].w == 0) {
      if (d.empty() || d.top() + P[i].x <= 0)
        continue;
      ans += d.top() + P[i].x;
      d.pop(), d.push(-P[i].x);
    } else
      d.push(-P[i].x + P[i].w);
  }

  printf("%lld\n", ans);
  return 0;
}

UOJ #455【UER #8】雪灾与外卖

题目链接

• https://uoj.ac/problem/455

解题思路

laofu 课件中的 Problem 6.

不妨先去除每个洞容量的限制, 考虑容量为 $1$ 的情况.

首先匹配交叉一定不优. 同样将数轴上的点从小到大排序, 对于每只老鼠或是每个洞, 讨论其匹配和反悔的权值即可.

具体的讨论可参考 UOJ 官方题解, 或者是文末的参考资料 (

每个洞容量的限制体现在网络流中即边的容量, 直接往之前的堆中丢一个 pair 记录增广的权值和容量即可. 可以证明, 此时对堆操作次数仍为线性.

在老鼠存在 “分身” 的情况下, 也是可以这样做的.

时间复杂度 $O(n \log n)$.

代码实现

// UOJ #455
// DeP
#include <queue>
#include <cctype>
#include <cstdio>
#include <algorithm>
using namespace std;

namespace IO {
  const int MAXSIZE = 1 << 18 | 1;
  char buf[MAXSIZE], *p1, *p2;

  inline int Gc() {
    return p1 == p2 &&
      (p2 = (p1 = buf) + fread(buf, 1, MAXSIZE, stdin), p1 == p2)? EOF: *p1++;
  }
  template<typename T> inline void read(T& x) {
    x = 0; int f = 0, ch = Gc();
    while (!isdigit(ch)) f |= ch == '-', ch = Gc();
    while (isdigit(ch)) x = x * 10 + ch - '0', ch = Gc();
    if (f) x = -x;
  }
}
using IO::read;

typedef long long LL;
typedef pair<LL, LL> Pll;

const int MAXN = 1e5 + 5;
const LL INFLL = 1e18;

struct Point {
  int x, w, c;
  Point() = default;
  Point(int _x, int _w, int _c): x(_x), w(_w), c(_c) { }
  bool operator < (const Point& rhs) const {
    return x < rhs.x;
  }
} P[MAXN << 1];

int n, m, p;
priority_queue<Pll, vector<Pll>, greater<Pll> > d0, d1;

signed main() {
#ifndef ONLINE_JUDGE
  freopen("data/ex_hole7.in", "r", stdin);
#endif
  read(n), read(m);
  for (int x, i = 1; i <= n; ++i)
    read(x), P[++p] = Point(x, -1, -1);
  LL sc = 0;
  for (int x, w, c, i = 1; i <= m; ++i) {
    read(x), read(w), read(c);
    P[++p] = Point(x, w, c), sc += c;
  }
  if (sc < n)
    return puts("-1"), 0;

  sort(P + 1, P + 1 + p);
  LL ans = 0;
  d1.push(Pll(INFLL, n));
  for (int i = 1; i <= p; ++i) {
    int x = P[i].x, w = P[i].w;
    if (w == -1) {
      Pll t = d1.top(); d1.pop();
      ans += t.first + x, --t.second;
      if (t.second > 0) d1.push(t);
      d0.push(Pll(-2LL * x - t.first, 1));
    } else {
      LL tot = 0, c = P[i].c;
      while (c > 0 && !d0.empty()) {
        Pll t = d0.top();
        if (t.first + x + w >= 0) break;
        d0.pop();
        LL s = min((LL) c, t.second);
        ans += s * (t.first + w + x), t.second -= s, c -= s;
        if (t.second > 0) d0.push(t);
        tot += s, d1.push(Pll(-2LL * x - t.first, s));
      }
      if (tot > 0)
        d0.push(Pll(-x - w, tot));
      if (c > 0)
        d1.push(Pll(-x + w, c));
    }
  }

  printf("%lld\n", ans);
  return 0;
}

LOJ #6405「ICPC World Finals 2018」征服世界

题目链接

• https://loj.ac/problem/6405

解题思路

laofu 课件中的 Problem 7.

上一题的模型上树.

首先每个点都存在军队和需求的限制是假的, 每个位置的军队解决自己的问题显然更优. 把军队视为洞, 军队需求视作老鼠, 可以得到老鼠进洞的模型. 换句话说, 对于 $x_i > y_i$ 的节点, 视作容量为 $x_i - y_i$ 的洞; 对于满足 $x_i < y_i$ 的节点, 视作 $y_i - x_i$ 只老鼠.

现在需要让所有老鼠进洞并使花费最小. 记 $\mathrm{dist}(u)$ 表示树上节点 $u$ 到根的距离, 那么花费形同 $\mathrm{dist}(u) + \mathrm{dist}(v) - 2\cdot \mathrm{dist}(\mathrm{LCA}(u, v))$.

对于每个点 $u$, 将 $\mathrm{dist}(u)$ 视作 $u$ 的权值. 同时, 每只老鼠一定要进洞, 直接将权值减去 $\infty$, 在计算答案时加上就好了.

剩下的情况和上一题没有太大区别, 另外需要可并堆快速从下向上合并节点的信息.

时间复杂度 $O(n \log n)$.

代码实现

// LOJ #6405
// DeP
#include <cctype>
#include <cstdio>
#include <cstring>
#include <algorithm>
using namespace std;

namespace IO {
  const int MAXSIZE = 1 << 18 | 1;
  char buf[MAXSIZE], *p1, *p2;

  inline int Gc() {
    return p1 == p2 &&
      (p2 = (p1 = buf) + fread(buf, 1, MAXSIZE, stdin), p1 == p2)? EOF: *p1++;
  }
  template<typename T> inline void read(T& x) {
    x = 0; int f = 0, ch = Gc();
    while (!isdigit(ch)) f |= ch == '-', ch = Gc();
    while (isdigit(ch)) x = x * 10 + ch - '0', ch = Gc();
    if (f) x = -x;
  }
}
using IO::read;

typedef long long LL;
typedef pair<LL, int> Pli;

const int MAXN = 2.5e5 + 5, MAXM = MAXN << 5;
const LL INFLL = MAXN * 1e6;

int n;
int X[MAXN], Y[MAXN];

LL d[MAXN], ans;
int rt0[MAXN], rt1[MAXN];

namespace Heap {
  Pli val[MAXM];
  int ch[2][MAXM], dist[MAXM], nidx;

  inline int newnode(const Pli& v) { return val[++nidx] = v, nidx; }

  int Mrg(int x, int y) {
    if (!x || !y) return x + y;
    if (val[x] > val[y]) swap(x, y);
    ch[1][x] = Mrg(ch[1][x], y);
    if (dist[ch[0][x]] < dist[ch[1][x]]) swap(ch[0][x], ch[1][x]);
    dist[x] = dist[ch[1][x]] + 1;
    return x;
  }

  inline void Pop(int& u) { u = Mrg(ch[0][u], ch[1][u]); }
  inline void Psh(int& u, const Pli& v) { u = Mrg(u, newnode(v)); }
}

namespace Graph {
  struct Edge { int nxt, to, w; } edges[MAXN << 1];
  int head[MAXN], eidx;

  inline void init() { memset(head, -1, sizeof head), eidx = 1; }
  inline void AddEdge(int from, int to, int w) {
    edges[++eidx] = (Edge){ head[from], to, w }, head[from] = eidx;
  }

  void Merge(int x, int y, const LL& d) {
    while (rt0[x] > 0 && rt1[y] > 0) {
      Pli &t0 = Heap::val[rt0[x]], &t1 = Heap::val[rt1[y]];
      if (t0.first + t1.first - 2 * d >= 0) break;
      int s = min(t0.second, t1.second);
      ans += s * (t0.first + t1.first - 2 * d);
      t0.second -= s, t1.second -= s;
      if (t0.second == 0) Heap::Pop(rt0[x]);
      if (t1.second == 0) Heap::Pop(rt1[y]);
      Heap::Psh(rt0[x], Pli(-t1.first + 2 * d, s));
      Heap::Psh(rt1[y], Pli(-t0.first + 2 * d, s));
    }
  }

  void dfs(int u, int fa) {
    if (X[u] > 0)
      Heap::Psh(rt0[u], Pli(d[u], X[u]));
    if (Y[u] > 0)
      Heap::Psh(rt1[u], Pli(d[u] - INFLL, Y[u])), ans += INFLL * Y[u];
    for (int v, i = head[u]; ~i; i = edges[i].nxt) {
      if ((v = edges[i].to) == fa) continue;
      d[v] = d[u] + edges[i].w, dfs(v, u);
      Merge(u, v, d[u]), Merge(v, u, d[u]);
      rt0[u] = Heap::Mrg(rt0[u], rt0[v]);
      rt1[u] = Heap::Mrg(rt1[u], rt1[v]);
    }
  }
}

int main() {
#ifndef ONLINE_JUDGE
  freopen("input.in", "r", stdin);
#endif
  Graph::init();

  read(n);
  for (int u, v, w, i = 1; i < n; ++i) {
    read(u), read(v), read(w);
    Graph::AddEdge(u, v, w), Graph::AddEdge(v, u, w);
  }
  for (int v, i = 1; i <= n; ++i) {
    read(X[i]), read(Y[i]);
    v = min(X[i], Y[i]), X[i] -= v, Y[i] -= v;
  }

  Graph::dfs(1, 0);

  printf("%lld\n", ans);
  return 0;
}

参考资料

• 陈江伦, <模拟费用流问题>, 2019.1.27.
• boshi, 模拟费用流的基本模型.
• mangoyang, 模拟流的问题小结.
• litble, 模拟费用流模型总结.
• 200815147, LOJ #2306「NOI2017」蔬菜 模拟费用流

---