算法学习：费用流
【前置知识】
【问题描述】
给定一个图，图上的边每条边有流量，流量有费用，
求在流量为 F 时费用的最小值
【分析思路】
求费用最小，我们很容易想到求最短路径，我们只需要将费用看作代价然后求最短路不久能求出来最小费用了嘛
但是，问题来了
我们又如何能够在求最小费用的同时，满足流量最终为F的条件
在最大流中提到，每次增广都是在残余网络中通过询问 s 至 t 的残余网络进行增广
所以，我们每次在残余网络上求最短路就可以了
这里注意，残余网络中的反向边的费用应该时原边费用的相反数
这样在计算的时候，才能保证过程是可逆的
【证明】
在残余网络上求最短路是最小费用
证：
1 设有 f 为以以上方式求出的结果， 设 f ’ 和 f 流量相同，但是费用更少
因为 f ‘ - f 的流量流入流出为 0，所以他是成环的
又因为 f ’ - f 的费用是负的，所以在这些环中，存在负环
结论：f 是最小费用流 <====> 残余网络中没有负环
（因为负环显而易见的能够通过转圈减少费用）
2 对于流量为 0 的流 f0 ，其残余网络为原图，原图没有负环，则f0 就是最小费用流
设流量为i 的流 fi 是最小流，并且下一步我们求得流量为 i + 1的流 f[i+1]，按照方法，f[ i + 1 ] - f [ i ]就是 f[ i ]对应的参余网络中 s 到 t 的最短路
假设 f [ i + 1 ] ' 的费用比 f [ i + 1 ]还小，则 f [ i + 1 ] '- f [ i ] 的费用比f[ i + 1 ] - f [ i ] 还小，所以f [ i + 1 ] '- f [ i ] 中有负环
所以这与假设矛盾，于是可以证明这种方法是正确的
【代码】
（基本上就是在最大流的基础上稍微改进了下）
// luogu-judger-enable-o2 #include<cstdio> #include<cstring> #include<algorithm> #include<iostream> #include<queue> #include<cstring> #define ll long long using namespace std; const int MAXN = 50010; const int MAXM = 100010; const ll INF = (1ll << 62) - 1; typedef pair<int, int> P; struct note { int to; int nt; int rev; ll cost; ll cal; }; struct edge { note arr[MAXM]; int st[MAXN]; ll dis[MAXN]; ll h[MAXN]; int cur[MAXN]; int depth[MAXN]; int pre[MAXN]; int pree[MAXN]; int top; int n, m, s, t; edge() { memset(st, -1, sizeof(st)); memset(depth, -1, sizeof(depth)); memset(dis, -1, sizeof(dis)); top = 0; } void read() { top = 0; scanf("%d%d%d%d", &n, &m, &s, &t); for (int i = 1; i <= m; i++) { int x, y; ll z,c; scanf("%d%d%lld%lld", &x, &y, &z,&c); add(x, y, z,c); } } bool dep() { queue<int> q; q.push(s); memset(depth, -1, sizeof(depth)); depth[s] = 0; while (!q.empty()) { int v = q.front(); q.pop(); for (int i = st[v]; i != -1; i = arr[i].nt) { int to = arr[i].to; if (!arr[i].cal) continue; if (depth[to] != -1) continue; depth[to] = depth[v] + 1; q.push(to); } } return (depth[t] != -1); } void add(int x, int y, ll z,ll c) { top++; arr[top] = { y,st[x],top + 1,c,z }; st[x] = top; top++; arr[top] = { x,st[y],top - 1,-c,0 }; st[y] = top; } ll dfs(int now, ll val) { if (now == t || !val) return val; ll flow = 0; for (int& i = cur[now]; i != -1; i = arr[i].nt) { int to = arr[i].to; if (depth[to] != depth[now] + 1) continue; ll f = dfs(to, min(arr[i].cal, val)); if (!f || !arr[i].cal) continue; flow += f; arr[i].cal -= f; arr[arr[i].rev].cal += f; val -= f; if (!val) return flow; } return flow; } ll dinic() { ll flow = 0; ll f; while (dep()) { for (int i = 1; i <= n; i++) cur[i] = st[i]; while (f = dfs(s, INF)) flow += f; } return flow; } ll min_cost(ll f) { ll res = 0; while (f > 0) { fill(dis, dis + 1 + n, INF); priority_queue<P, vector<P>, greater<P>> q; dis[s] = 0; q.push(P(0, s)); while (!q.empty()) { P p = q.top(); q.pop(); int v = p.second; if (dis[v] < p.first) continue; for (int i = st[v]; i != -1; i = arr[i].nt) { note &e = arr[i]; if (e.cal > 0 && dis[e.to] > dis[v] + e.cost + h[v] - h[e.to]) { dis[e.to] = dis[v] + e.cost + h[v] - h[e.to]; pre[e.to] = v; pree[e.to] = i; q.push(P(dis[e.to],e.to)); } } } if (dis[t] == INF) return -1; for (int i = 0; i <= n; i++) h[i] += dis[i]; ll d = f; //printf("2\n"); for (int i = t; i != s; i = pre[i]) { d = min(d, arr[pree[i]].cal); } f -= d; res += d * h[t]; for (int i = t; i != s; i = pre[i]) { arr[pree[i]].cal -= d;; int rev = arr[pree[i]].rev; arr[rev].cal += d; } //printf("3\n"); } return res; } }; edge road; edge tmp; int main() { int T; road.read(); tmp = road; ll f = tmp.dinic(); printf("%lld ", f); printf("%lld", road.min_cost(f)); return 0; }
注意：因为求最大流之后网络会发生变化，所以需要把最开始的原图记录下来
费用流的求取也要在原图上进行
【2020年2月7日 对代码进行修正】
添加 siz ，保留 n ，作为有的拆点之后网络流中结点增加之后遍历的更改
之前是直接改 n 用的，但是存在很多不必要的弊端，所以进行少量修改
主要是遍历和对 dis[] 的填充
【模板】
1 // luogu-judger-enable-o2 2 #include<cstdio> 3 #include<cstring> 4 #include<algorithm> 5 #include<iostream> 6 #include<queue> 7 #include<cstring> 8 #define ll int 9 using namespace std; 10 const int MAXN = 100010; 11 const int MAXM = 20000010; 12 const int INF = 2e9 + 7; 13 typedef pair<int, int> P; 14 struct note 15 { 16 int to; 17 int nt; 18 int rev; 19 ll cost; 20 ll cal; 21 }; 22 struct node 23 { 24 int x, y; 25 ll cost, z; 26 }; 27 int tn, tm, tk, sum; 28 int p[45], a[45][110]; 29 void read() 30 { 31 scanf("%d%d", &tn, &tm); 32 for (int i = 1; i <= tn; i++) 33 scanf("%d", &p[i]), sum += p[i]; 34 for (int i = 1; i <= tn; i++) 35 for (int j = 1; j <= tm; j++) 36 scanf("%d", &a[i][j]); 37 } 38 struct edge 39 { 40 note arr[MAXM]; 41 int st[MAXN]; 42 ll dis[MAXN]; 43 ll h[MAXN]; 44 int cur[MAXN]; 45 int depth[MAXN]; 46 int pre[MAXN]; 47 int pree[MAXN]; 48 int top; 49 int n, m, s, t, siz, k; 50 edge() 51 { 52 memset(st, -1, sizeof(st)); 53 memset(depth, -1, sizeof(depth)); 54 memset(dis, -1, sizeof(dis)); 55 top = 0; 56 } 57 void build(int ts, int tt) 58 { 59 s = ts, t = tt, siz = tt; 60 //增加了siz 61 n = tn; m = tm; 62 for (int i = 1; i <= sum; i++) 63 add(s, i, 1, 0); 64 for (int i = sum + 1; i < t; i++) 65 add(i, t, 1, 0); 66 int cnt = 0; 67 for (int i = 1; i <= n; i++) 68 { 69 for (int j = 1; j <= p[i]; j++) 70 { 71 int pos = (cnt + j); 72 for (int j1 = 1; j1 <= m; j1++) 73 for (int k = 1; k <= sum; k++) 74 { 75 add(pos, sum + (j1 - 1) * sum + k, 1, a[i][j1] * k); 76 } 77 } 78 cnt += p[i]; 79 } 80 } 81 bool dep() 82 { 83 queue<int> q; 84 q.push(s); 85 memset(depth, -1, sizeof(depth)); 86 depth[s] = 0; 87 while (!q.empty()) 88 { 89 int v = q.front(); q.pop(); 90 for (int i = st[v]; i != -1; i = arr[i].nt) 91 { 92 int to = arr[i].to; 93 if (!arr[i].cal) 94 continue; 95 if (depth[to] != -1) 96 continue; 97 depth[to] = depth[v] + 1; 98 q.push(to); 99 } 100 } 101 return (depth[t] != -1); 102 } 103 void add(int x, int y, ll z, ll c) 104 { 105 top++; arr[top] = { y,st[x],top + 1,c,z }; st[x] = top; 106 top++; arr[top] = { x,st[y],top - 1,-c,0 }; st[y] = top; 107 } 108 ll dfs(int now, ll val) 109 { 110 if (now == t || !val) 111 return val; 112 ll flow = 0; 113 for (int& i = cur[now]; i != -1; i = arr[i].nt) 114 { 115 int to = arr[i].to; 116 if (depth[to] != depth[now] + 1) 117 continue; 118 ll f = dfs(to, min(arr[i].cal, val)); 119 if (!f || !arr[i].cal) 120 continue; 121 flow += f; 122 arr[i].cal -= f; 123 arr[arr[i].rev].cal += f; 124 val -= f; 125 if (!val) 126 return flow; 127 } 128 return flow; 129 } 130 ll dinic() 131 { 132 ll flow = 0; 133 ll f; 134 while (dep()) 135 { 136 for (int i = 0; i <= siz; i++) 137 //这里 0代表的是起点 138 cur[i] = st[i]; 139 while (f = dfs(s, INF)) 140 flow += f; 141 } 142 return flow; 143 } 144 ll min_cost(ll f) 145 { 146 ll res = 0; 147 while (f > 0) 148 { 149 fill(dis, dis + 1 + siz, INF); 150 priority_queue<P, vector<P>, greater<P>> q; 151 dis[s] = 0; 152 q.push(P(0, s)); 153 while (!q.empty()) 154 { 155 P p = q.top(); q.pop(); 156 int v = p.second; 157 if (dis[v] < p.first) continue; 158 for (int i = st[v]; i != -1; i = arr[i].nt) 159 { 160 note& e = arr[i]; 161 if (e.cal > 0 && dis[e.to] > dis[v] + e.cost + h[v] - h[e.to]) 162 { 163 dis[e.to] = dis[v] + e.cost + h[v] - h[e.to]; 164 pre[e.to] = v; 165 pree[e.to] = i; 166 q.push(P(dis[e.to], e.to)); 167 } 168 } 169 } 170 if (dis[t] == INF) return -1; 171 for (int i = 0; i <= siz; i++) 172 //这里改了 0指的是s 173 h[i] += dis[i]; 174 ll d = f; 175 //printf("2\n"); 176 for (int i = t; i != s; i = pre[i]) 177 { 178 d = min(d, arr[pree[i]].cal); 179 } 180 f -= d; 181 res += d * h[t]; 182 for (int i = t; i != s; i = pre[i]) 183 { 184 arr[pree[i]].cal -= d;; 185 int rev = arr[pree[i]].rev; 186 arr[rev].cal += d; 187 } 188 } 189 return res; 190 } 191 }; 192 edge road, new_road; 193 edge tmp; 194 int main() 195 { 196 read(); 197 //如果需要求最大流之前记得保存原图 198 //min_cost参数是最大流 199 road.build(0, sum + sum * tm + 1); 200 printf("%d", road.min_cost(sum)); 201 return 0; 202 }
【常见题型】
0.比较憨的模板题
【luogu P4014】分配问题
1.动态枚举，灵活处理
【luogu P2604】 网络扩容
2.暴力拆点
动态加边优化
3.隐式图
【luogu P2153】晨跑
