4.5.5 Maximum Flow (Push-Relabel)

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151

/*

Given a flow network with integer capacities, find the maximum flow from a given source node to a
given sink node. The flow along each edge may not exceed its capacity, and flow is conserved at
every node other than the source and sink.

Rather than augmenting along whole paths, the push-relabel algorithm maintains a preflow and a
height label on each node, repeatedly pushing excess flow to lower-labeled neighbors and relabeling
(raising) a node whose excess cannot yet be pushed, until no node other than the source and sink
holds excess.

- `push_relabel(source, sink)` returns maximum flow for a global capacity matrix `cap` whose nodes
  are numbered [0, `n`), where `n` is `cap.size()`.
- `min_cut(source)` returns the source side of a minimum cut after `push_relabel()` has been called.

Although the push-relabel algorithm is considered one of the most efficient maximum flow algorithms,
it cannot take advantage of the magnitude of the maximum flow being less than $n^3$ (in which case
the Ford-Fulkerson or Edmonds-Karp algorithms may be more efficient).

Time Complexity:
- O(n^3) per call to `push_relabel()`, where $n$ is the number of nodes.
- O(n^2) per call to `min_cut()`.

Space Complexity:
- O(n^2) for storage of the flow network, where $n$ is the number of nodes.
- O(n) auxiliary heap space for `push_relabel()`.

*/

#include <algorithm>
#include <climits>
#include <cstdint>
#include <queue>
#include <vector>

const int INF = INT_MAX / 2;
std::vector<std::vector<int>> cap;
std::vector<std::vector<int>> f;

int64_t push_relabel(int source, int sink) {
  int n = static_cast<int>(cap.size());
  f.assign(n, std::vector<int>(n, 0));
  std::vector<int> e(n, 0), h(n, 0), maxh(n, 0);
  h[source] = n - 1;
  for (int i = 0; i < n; i++) {
    f[source][i] = cap[source][i];
    f[i][source] = -f[source][i];
    e[i] = cap[source][i];
  }
  int size = 0;
  while (true) {
    if (size == 0) {
      for (int i = 0; i < n; i++) {
        if (i != source && i != sink && e[i] > 0) {
          if (size != 0 && h[i] > h[maxh[0]]) {
            size = 0;
          }
          maxh[size++] = i;
        }
      }
    }
    if (size == 0) {
      break;
    }
    while (size != 0) {
      int i = maxh[size - 1];
      bool pushed = false;
      for (int j = 0; j < n && e[i] != 0; j++) {
        if (h[i] == h[j] + 1 && cap[i][j] - f[i][j] > 0) {
          int df = std::min(cap[i][j] - f[i][j], e[i]);
          f[i][j] += df;
          f[j][i] -= df;
          e[i] -= df;
          e[j] += df;
          if (e[i] == 0) {
            size--;
          }
          pushed = true;
        }
      }
      if (pushed) {
        continue;
      }
      h[i] = INF;
      for (int j = 0; j < n; j++) {
        if (h[i] > h[j] + 1 && cap[i][j] - f[i][j] > 0) {
          h[i] = h[j] + 1;
        }
      }
      if (h[i] > h[maxh[0]]) {
        size = 0;
        break;
      }
    }
  }
  int64_t max_flow = 0;
  for (int i = 0; i < n; i++) {
    max_flow += f[source][i];
  }
  return max_flow;
}

std::vector<bool> min_cut(int source) {
  int n = static_cast<int>(cap.size());
  std::vector<bool> reachable(n);
  std::queue<int> q;
  reachable[source] = true;
  q.push(source);
  while (!q.empty()) {
    int u = q.front();
    q.pop();
    for (int v = 0; v < n; v++) {
      if (!reachable[v] && cap[u][v] - f[u][v] > 0) {
        reachable[v] = true;
        q.push(v);
      }
    }
  }
  return reachable;
}

/*** Example Usage ***/

#include <cassert>

int main() {
  int nodes = 6;
  cap.assign(nodes, std::vector<int>(nodes));
  // Example graph after max flow, with each edge labeled flow/capacity:
  //            2/2
  //       1 --------> 3
  //      / \          |
  // 3/4 /   \ 1/1     | 2/2
  //    /     v        v
  //   0       4 ----> 5
  //    \     ^   3/3
  // 2/3 \   / 2/2
  //      v /
  //       2
  cap[0][1] = 4;
  cap[0][2] = 3;
  cap[1][3] = 2;
  cap[1][4] = 1;
  cap[2][4] = 2;
  cap[3][5] = 2;
  cap[4][5] = 3;
  assert(push_relabel(0, 5) == 5);
  std::vector<bool> cut = min_cut(0);
  assert(cut[0] && !cut[5]);
  return 0;
}