rewrite cost flow code

zimpha · zimpha · commit 1e763a49961b · 2016-08-27T20:03:56.000+08:00
diff --git a/graph-utility/CostFlow.cc b/graph-utility/CostFlow.cc
@@ -1,128 +1,173 @@
 #include <bits/stdc++.h>
-// 首先调用MCMF::init()设置点的个数，源点和汇点
-// 调用MCMF::addedge()加边
-// 调用MCMF::solve()计算答案，答案存在MinCost和MaxFlow里面
-namespace MCMF {
-  typedef int flow_t;
-  typedef int cost_t;
-  const int MAXN = 200, MAXM = 100000;
+
+template<typename flow_t, typename cost_t>
+struct MCMF {
+  static const int N = 200, M = 100000;
   const flow_t inf = 1e9;
-  struct Edge {
-    int from, to, nx;
+  struct node {
+    int from, to, nxt;
     flow_t cap, flow;
     cost_t cost;
-    Edge() {}
-    Edge(int _from, int _to, int _nx, flow_t _cap, cost_t _cost):
-      from(_from), to(_to), nx(_nx), cap(_cap), flow(0), cost(_cost) {}
-  } E[MAXM];
-  cost_t dis[MAXN];
-  int G[MAXN], pre[MAXN], vis[MAXN];
-  int sz, N; //S源点, T汇点
+    node() {}
+    node(int from, int to, int nxt, flow_t cap, cost_t cost):
+      from(from), to(to), nxt(nxt), cap(cap), flow(0), cost(cost) {}
+  } E[M];
+  cost_t dis[N];
+  int G[N], pre[N], vis[N], n, m;
   void init(int n) {
-    memset(G, -1, sizeof(G[0]) * n);
-    sz = 0; N = n;
+    this->n = n;
+    this->m = 0;
+    std::fill(G, G + n, -1);
   }
-  void addedge(int u, int v, flow_t f, cost_t c) {//u -> v
-    E[sz] = Edge(u, v, G[u], f, +c); G[u] = sz++;
-    E[sz] = Edge(v, u, G[v], 0, -c); G[v] = sz++;
+  void link(int u, int v, flow_t f, cost_t c) {
+    E[m] = node(u, v, G[u], f, +c); G[u] = m++;
+    E[m] = node(v, u, G[v], 0, -c); G[v] = m++;
   }
   bool extand(int S, int T) {
-    std::queue<int> Q;
-    memset(vis, 0, sizeof(vis[0]) * N);
-    memset(pre, -1, sizeof(pre[0]) * N);
-    for (int i = 0; i < N; ++i) dis[i] = inf;
-    Q.push(S); vis[S] = 1; dis[S] = 0;
-    while (!Q.empty()) {
-      int u = Q.front();
-      Q.pop(); vis[u]=0;
-      for (int it = G[u]; ~it; it = E[it].nx) {
+    std::fill(vis, vis + n, 0);
+    std::fill(dis, dis + n, inf);
+    std::queue<int> queue;
+    dis[S] = 0;
+    queue.push(S);
+    for (; !queue.empty(); queue.pop()) {
+      int u = queue.front();
+      vis[u] = false;
+      for (int it = G[u]; ~it; it = E[it].nxt) {
         int v = E[it].to;
         if (E[it].cap > E[it].flow && dis[v] > dis[u] + E[it].cost) {
-          dis[v] = dis[u] + E[it].cost, pre[v] = it;
-          if (!vis[v]) Q.push(v), vis[v] = true;
+          dis[v] = dis[u] + E[it].cost;
+          pre[v] = it;
+          if (!vis[v]) queue.push(v);
+          vis[v] = true;
         }
       }
     }
     return dis[T] < inf; // 改成dis[T] <= 0 求可行流
   }
-  std::pair<flow_t, cost_t> solve(int S, int T) {
-    flow_t MaxFlow = 0;
-    cost_t MinCost = 0;
+  std::pair<flow_t, cost_t> run(int S, int T) {
+    flow_t max_flow = 0;
+    cost_t min_cost = 0;
     while (extand(S, T)) {
-      flow_t aug = inf;
-      for (int it = pre[T]; ~it; it = pre[E[it].from]) {
-        aug = std::min(aug, E[it].cap - E[it].flow);
+      flow_t delta = inf;
+      for (int u = T; u != S; u = E[pre[u]].from) {
+        delta = std::min(delta, E[pre[u]].cap - E[pre[u]].flow);
       }
-      MinCost += aug * dis[T], MaxFlow += aug;
-      for (int it = pre[T]; ~it; it = pre[E[it].from]) {
-        E[it].flow += aug, E[it ^ 1].flow -= aug;
+      min_cost += delta * dis[T];
+      max_flow += delta;
+      for (int u = T; u != S; u = E[pre[u]].from) {
+        E[pre[u]].flow += delta;
+        E[pre[u] ^ 1].flow -= delta;
       }
     }
-    return std::make_pair(MaxFlow, MinCost);
+    return {max_flow, min_cost};
   }
-}
+};
 
-namespace ZKW {
-  const int MAXN = 50000 + 10, MAXM = 300000, inf = 1e9;
-  typedef int VAL;
-  struct Edge {
-    int v,f;
-    VAL c;
-    int nx;
-    Edge() {}
-    Edge(int v,int f,VAL c,int nx):v(v),f(f),c(c),nx(nx) {}
-  } E[MAXM];
-  int G[MAXN], Q[MAXN];
-  VAL dis[MAXN],mincost;
-  bool vis[MAXN],mark[MAXN];
-  int N,S,T,sz;
-  void init(int _n, int _s, int _t) {
-    N = _n; S = _s; T = _t; sz = 0;
-    memset(G, -1, sizeof(G[0]) * N);
+template<typename flow_t, typename cost_t>
+struct MCMF_ZKW {
+  static const int N = 200, M = 3000, inf = 1e9;
+  struct node {
+    int to, nxt;
+    flow_t cap, flow;
+    cost_t cost;
+    node() {}
+    node(int to, int nxt, flow_t cap, cost_t cost):
+      to(to), nxt(nxt), cap(cap), flow(0), cost(cost) {}
+  } E[M];
+  int G[N], n, m;
+  cost_t min_cost, len;
+  flow_t max_flow;
+  bool done[N];
+  void init(int n) {
+    this->n = n;
+    this->m = 0;
+    std::fill(G, G + n, -1);
   }
-  void link(int u, int v, int f, VAL c) {
-    E[sz]=Edge(v,f,+c,G[u]); G[u]=sz++;
-    E[sz]=Edge(u,0,-c,G[v]); G[v]=sz++;
+  void link(int u, int v, flow_t f, cost_t c) {
+    E[m] = node(v, G[u], f, +c); G[u] = m++;
+    E[m] = node(u, G[v], 0, -c); G[v] = m++;
   }
-  bool spfa() {
-    for (int i=0;i<N;++i) dis[i]=inf,vis[i]=0;
-    dis[T]=0; Q[0]=T; vis[T]=1;
-    for (int h=0,t=1;h!=t; ) {
-      int u=Q[h++],v; vis[u]=false;
-      if (h==MAXN) h=0;
-      for (int it=G[u];~it;it=E[it].nx) {
-        if (dis[u]-E[it].c<dis[v=E[it].v]&&E[it^1].f) {
-          dis[v]=dis[u]-E[it].c;
-          if (!vis[v]) {
-            vis[v]=true; Q[t++]=v;
-            if (t==MAXN) t=0;
-          }
-        }
+  flow_t aug(int now, int T, flow_t max_cap) {
+    if (now == T) {
+      max_flow += max_cap;
+      min_cost += max_cap * len;
+      return max_cap;
+    }
+    done[now] = true;
+    flow_t upp = max_cap;
+    for (int it = G[now]; ~it && upp; it = E[it].nxt) {
+      if (E[it].cap > E[it].flow && !E[it].cost && !done[E[it].to]) {
+        flow_t delta = aug(E[it].to, T, std::min(upp, E[it].cap - E[it].flow));
+        E[it].flow += delta;
+        E[it ^ 1].flow -= delta;
+        upp -= delta;
       }
     }
-    return dis[S]!=inf;
+    return max_cap - upp;
   }
-  int dfs(int u, int low) {
-    mark[u]=true; if (u==T) return low;
-    int ret=0,tmp=0;
-    for (int it=G[u],v;~it&&ret<low;it=E[it].nx) {
-      if (E[it].f&&!mark[v=E[it].v]&&dis[u]-E[it].c==dis[v]) {
-        tmp=dfs(v,std::min(E[it].f,low-ret));
-        E[it].f-=tmp; E[it^1].f+=tmp;
-        mincost+=E[it].c*tmp; ret+=tmp;
+  bool label(int S, int T) {//不能用于负费用
+    cost_t delta = inf;
+    for (int u = 0; u < n; ++u) if (done[u]) {
+      for (int it = G[u]; ~it; it = E[it].nxt) {
+        if (E[it].cap > E[it].flow && !done[E[it].to]) {
+          delta = std::min(delta, E[it].cost);
+        }
+      }
+    }
+    if (delta == inf) return false;
+    for (int u = 0; u < n; ++u) if (done[u]) {
+      for (int it = G[u]; ~it; it = E[it].nxt) {
+        E[it].cost -= delta;
+        E[it ^ 1].cost += delta;
       }
     }
-    return ret;
+    len += delta;
+    return true;
   }
-  VAL solve() {
-    int maxflow = 0; mincost = 0;
-    while (spfa()) {
-      mark[T] = true;
-      while (mark[T]) {
-        memset(mark, 0, sizeof(mark[0]) * N);
-        maxflow += dfs(S, inf);
+  cost_t dis[N];
+  bool label_primal_dual(int S, int T) {
+    for (int i = 0; i < n; ++i) dis[i] = inf;
+    std::fill(done, done + n, 0);
+    dis[T] = 0;
+    std::queue<int> queue;
+    queue.push(T);
+    for (; !queue.empty(); queue.pop()) {
+      int u = queue.front();
+      done[u] = false;
+      for (int it = G[u]; ~it; it = E[it].nxt) {
+        int v = E[it].to;
+        cost_t cost = dis[u] - E[it].cost;
+        if (E[it ^ 1].cap > E[it ^ 1].flow && cost < dis[v]) {
+          dis[v] = cost;
+          if (!done[v]) queue.push(v);
+          done[v] = true;
+        }
+      }
+    }
+    for (int u = 0; u < n; ++u) {
+      for (int it = G[u]; ~it; it = E[it].nxt) {
+        E[it].cost += dis[E[it].to] - dis[u];
       }
     }
-    return mincost;
+    len += dis[S];
+    return dis[S] < inf;
+  }
+  std::pair<flow_t, cost_t> run_primal_dual(int S, int T) {
+    max_flow = min_cost = len = 0;
+    while (label_primal_dual(S, T)) {
+      do {
+        std::fill(done, done + n, 0);
+      } while (aug(S, T, inf));
+    }
+    return {max_flow, min_cost};
+  }
+  std::pair<flow_t, cost_t> run(int S, int T) {
+    max_flow = min_cost = len = 0;
+    do {
+      do {
+        std::fill(done, done + n, 0);
+      } while (aug(S, T, inf));
+    } while (label_primal_dual(S, T));
+    return {max_flow, min_cost};
   }
-}
+};
