From 3d04dd3335a2c4c6007ff4e2a58a2855c7a9c52a Mon Sep 17 00:00:00 2001 From: Simon Tatham Date: Sat, 21 Apr 2018 17:03:10 +0100 Subject: [PATCH] Remove maxflow completely. Its ability to solve general network flow problems was never actually used in this code base; it was _always_ used for the restricted problem of finding a matching in an unweighted bipartite graph. So now I've switched all its clients over to the new matching.c, maxflow is no longer needed. --- maxflow.c | 461 ------------------------------------------------------ maxflow.h | 95 ----------- 2 files changed, 556 deletions(-) delete mode 100644 maxflow.c delete mode 100644 maxflow.h diff --git a/maxflow.c b/maxflow.c deleted file mode 100644 index 97ae8c4..0000000 --- a/maxflow.c +++ /dev/null @@ -1,461 +0,0 @@ -/* - * Edmonds-Karp algorithm for finding a maximum flow and minimum - * cut in a network. Almost identical to the Ford-Fulkerson - * algorithm, but apparently using breadth-first search to find the - * _shortest_ augmenting path is a good way to guarantee - * termination and ensure the time complexity is not dependent on - * the actual value of the maximum flow. I don't understand why - * that should be, but it's claimed on the Internet that it's been - * proved, and that's good enough for me. I prefer BFS to DFS - * anyway :-) - */ - -#include -#include -#include - -#include "maxflow.h" - -#include "puzzles.h" /* for snewn/sfree */ - -int maxflow_with_scratch(void *scratch, int nv, int source, int sink, - int ne, const int *edges, const int *backedges, - const int *capacity, int *flow, int *cut) -{ - int *todo = (int *)scratch; - int *prev = todo + nv; - int *firstedge = todo + 2*nv; - int *firstbackedge = todo + 3*nv; - int i, j, head, tail, from, to; - int totalflow; - - /* - * Scan the edges array to find the index of the first edge - * from each node. - */ - j = 0; - for (i = 0; i < ne; i++) - while (j <= edges[2*i]) - firstedge[j++] = i; - while (j < nv) - firstedge[j++] = ne; - assert(j == nv); - - /* - * Scan the backedges array to find the index of the first edge - * _to_ each node. - */ - j = 0; - for (i = 0; i < ne; i++) - while (j <= edges[2*backedges[i]+1]) - firstbackedge[j++] = i; - while (j < nv) - firstbackedge[j++] = ne; - assert(j == nv); - - /* - * Start the flow off at zero on every edge. - */ - for (i = 0; i < ne; i++) - flow[i] = 0; - totalflow = 0; - - /* - * Repeatedly look for an augmenting path, and follow it. - */ - while (1) { - - /* - * Set up the prev array. - */ - for (i = 0; i < nv; i++) - prev[i] = -1; - - /* - * Initialise the to-do list for BFS. - */ - head = tail = 0; - todo[tail++] = source; - - /* - * Now do the BFS loop. - */ - while (head < tail && prev[sink] < 0) { - from = todo[head++]; - - /* - * Try all the forward edges out of node `from'. For a - * forward edge to be valid, it must have flow - * currently less than its capacity. - */ - for (i = firstedge[from]; i < ne && edges[2*i] == from; i++) { - to = edges[2*i+1]; - if (to == source || prev[to] >= 0) - continue; - if (capacity[i] >= 0 && flow[i] >= capacity[i]) - continue; - /* - * This is a valid augmenting edge. Visit node `to'. - */ - prev[to] = 2*i; - todo[tail++] = to; - } - - /* - * Try all the backward edges into node `from'. For a - * backward edge to be valid, it must have flow - * currently greater than zero. - */ - for (i = firstbackedge[from]; - j = backedges[i], i < ne && edges[2*j+1]==from; i++) { - to = edges[2*j]; - if (to == source || prev[to] >= 0) - continue; - if (flow[j] <= 0) - continue; - /* - * This is a valid augmenting edge. Visit node `to'. - */ - prev[to] = 2*j+1; - todo[tail++] = to; - } - } - - /* - * If prev[sink] is non-null, we have found an augmenting - * path. - */ - if (prev[sink] >= 0) { - int max; - - /* - * Work backwards along the path figuring out the - * maximum flow we can add. - */ - to = sink; - max = -1; - while (to != source) { - int spare; - - /* - * Find the edge we're currently moving along. - */ - i = prev[to]; - from = edges[i]; - assert(from != to); - - /* - * Determine the spare capacity of this edge. - */ - if (i & 1) - spare = flow[i / 2]; /* backward edge */ - else if (capacity[i / 2] >= 0) - spare = capacity[i / 2] - flow[i / 2]; /* forward edge */ - else - spare = -1; /* unlimited forward edge */ - - assert(spare != 0); - - if (max < 0 || (spare >= 0 && spare < max)) - max = spare; - - to = from; - } - /* - * Fail an assertion if max is still < 0, i.e. there is - * an entirely unlimited path from source to sink. Also - * max should not _be_ zero, because by construction - * this _should_ be an augmenting path. - */ - assert(max > 0); - - /* - * Now work backwards along the path again, this time - * actually adjusting the flow. - */ - to = sink; - while (to != source) { - /* - * Find the edge we're currently moving along. - */ - i = prev[to]; - from = edges[i]; - assert(from != to); - - /* - * Adjust the edge. - */ - if (i & 1) - flow[i / 2] -= max; /* backward edge */ - else - flow[i / 2] += max; /* forward edge */ - - to = from; - } - - /* - * And adjust the overall flow counter. - */ - totalflow += max; - - continue; - } - - /* - * If we reach here, we have failed to find an augmenting - * path, which means we're done. Output the `cut' array if - * required, and leave. - */ - if (cut) { - for (i = 0; i < nv; i++) { - if (i == source || prev[i] >= 0) - cut[i] = 0; - else - cut[i] = 1; - } - } - return totalflow; - } -} - -int maxflow_scratch_size(int nv) -{ - return (nv * 4) * sizeof(int); -} - -void maxflow_setup_backedges(int ne, const int *edges, int *backedges) -{ - int i, n; - - for (i = 0; i < ne; i++) - backedges[i] = i; - - /* - * We actually can't use the C qsort() function, because we'd - * need to pass `edges' as a context parameter to its - * comparator function. So instead I'm forced to implement my - * own sorting algorithm internally, which is a pest. I'll use - * heapsort, because I like it. - */ - -#define LESS(i,j) ( (edges[2*(i)+1] < edges[2*(j)+1]) || \ - (edges[2*(i)+1] == edges[2*(j)+1] && \ - edges[2*(i)] < edges[2*(j)]) ) -#define PARENT(n) ( ((n)-1)/2 ) -#define LCHILD(n) ( 2*(n)+1 ) -#define RCHILD(n) ( 2*(n)+2 ) -#define SWAP(i,j) do { int swaptmp = (i); (i) = (j); (j) = swaptmp; } while (0) - - /* - * Phase 1: build the heap. We want the _largest_ element at - * the top. - */ - n = 0; - while (n < ne) { - n++; - - /* - * Swap element n with its parent repeatedly to preserve - * the heap property. - */ - i = n-1; - - while (i > 0) { - int p = PARENT(i); - - if (LESS(backedges[p], backedges[i])) { - SWAP(backedges[p], backedges[i]); - i = p; - } else - break; - } - } - - /* - * Phase 2: repeatedly remove the largest element and stick it - * at the top of the array. - */ - while (n > 0) { - /* - * The largest element is at position 0. Put it at the top, - * and swap the arbitrary element from that position into - * position 0. - */ - n--; - SWAP(backedges[0], backedges[n]); - - /* - * Now repeatedly move that arbitrary element down the heap - * by swapping it with the more suitable of its children. - */ - i = 0; - while (1) { - int lc, rc; - - lc = LCHILD(i); - rc = RCHILD(i); - - if (lc >= n) - break; /* we've hit bottom */ - - if (rc >= n) { - /* - * Special case: there is only one child to check. - */ - if (LESS(backedges[i], backedges[lc])) - SWAP(backedges[i], backedges[lc]); - - /* _Now_ we've hit bottom. */ - break; - } else { - /* - * The common case: there are two children and we - * must check them both. - */ - if (LESS(backedges[i], backedges[lc]) || - LESS(backedges[i], backedges[rc])) { - /* - * Pick the more appropriate child to swap with - * (i.e. the one which would want to be the - * parent if one were above the other - as one - * is about to be). - */ - if (LESS(backedges[lc], backedges[rc])) { - SWAP(backedges[i], backedges[rc]); - i = rc; - } else { - SWAP(backedges[i], backedges[lc]); - i = lc; - } - } else { - /* This element is in the right place; we're done. */ - break; - } - } - } - } - -#undef LESS -#undef PARENT -#undef LCHILD -#undef RCHILD -#undef SWAP - -} - -int maxflow(int nv, int source, int sink, - int ne, const int *edges, const int *capacity, - int *flow, int *cut) -{ - void *scratch; - int *backedges; - int size; - int ret; - - /* - * Allocate the space. - */ - size = ne * sizeof(int) + maxflow_scratch_size(nv); - backedges = smalloc(size); - if (!backedges) - return -1; - scratch = backedges + ne; - - /* - * Set up the backedges array. - */ - maxflow_setup_backedges(ne, edges, backedges); - - /* - * Call the main function. - */ - ret = maxflow_with_scratch(scratch, nv, source, sink, ne, edges, - backedges, capacity, flow, cut); - - /* - * Free the scratch space. - */ - sfree(backedges); - - /* - * And we're done. - */ - return ret; -} - -#ifdef TESTMODE - -#define MAXEDGES 256 -#define MAXVERTICES 128 -#define ADDEDGE(i,j) do{edges[ne*2] = (i); edges[ne*2+1] = (j); ne++;}while(0) - -int compare_edge(const void *av, const void *bv) -{ - const int *a = (const int *)av; - const int *b = (const int *)bv; - - if (a[0] < b[0]) - return -1; - else if (a[0] > b[0]) - return +1; - else if (a[1] < b[1]) - return -1; - else if (a[1] > b[1]) - return +1; - else - return 0; -} - -int main(void) -{ - int edges[MAXEDGES*2], ne, nv; - int capacity[MAXEDGES], flow[MAXEDGES], cut[MAXVERTICES]; - int source, sink, p, q, i, j, ret; - - /* - * Use this algorithm to find a maximal complete matching in a - * bipartite graph. - */ - ne = 0; - nv = 0; - source = nv++; - p = nv; - nv += 5; - q = nv; - nv += 5; - sink = nv++; - for (i = 0; i < 5; i++) { - capacity[ne] = 1; - ADDEDGE(source, p+i); - } - for (i = 0; i < 5; i++) { - capacity[ne] = 1; - ADDEDGE(q+i, sink); - } - j = ne; - capacity[ne] = 1; ADDEDGE(p+0,q+0); - capacity[ne] = 1; ADDEDGE(p+1,q+0); - capacity[ne] = 1; ADDEDGE(p+1,q+1); - capacity[ne] = 1; ADDEDGE(p+2,q+1); - capacity[ne] = 1; ADDEDGE(p+2,q+2); - capacity[ne] = 1; ADDEDGE(p+3,q+2); - capacity[ne] = 1; ADDEDGE(p+3,q+3); - capacity[ne] = 1; ADDEDGE(p+4,q+3); - /* capacity[ne] = 1; ADDEDGE(p+2,q+4); */ - qsort(edges, ne, 2*sizeof(int), compare_edge); - - ret = maxflow(nv, source, sink, ne, edges, capacity, flow, cut); - - printf("ret = %d\n", ret); - - for (i = 0; i < ne; i++) - printf("flow %d: %d -> %d\n", flow[i], edges[2*i], edges[2*i+1]); - - for (i = 0; i < nv; i++) - if (cut[i] == 0) - printf("difficult set includes %d\n", i); - - return 0; -} - -#endif diff --git a/maxflow.h b/maxflow.h deleted file mode 100644 index d490f45..0000000 --- a/maxflow.h +++ /dev/null @@ -1,95 +0,0 @@ -/* - * Edmonds-Karp algorithm for finding a maximum flow and minimum - * cut in a network. Almost identical to the Ford-Fulkerson - * algorithm, but apparently using breadth-first search to find the - * _shortest_ augmenting path is a good way to guarantee - * termination and ensure the time complexity is not dependent on - * the actual value of the maximum flow. I don't understand why - * that should be, but it's claimed on the Internet that it's been - * proved, and that's good enough for me. I prefer BFS to DFS - * anyway :-) - */ - -#ifndef MAXFLOW_MAXFLOW_H -#define MAXFLOW_MAXFLOW_H - -/* - * The actual algorithm. - * - * Inputs: - * - * - `scratch' is previously allocated scratch space of a size - * previously determined by calling `maxflow_scratch_size'. - * - * - `nv' is the number of vertices. Vertices are assumed to be - * numbered from 0 to nv-1. - * - * - `source' and `sink' are the distinguished source and sink - * vertices. - * - * - `ne' is the number of edges in the graph. - * - * - `edges' is an array of 2*ne integers, giving a (source, dest) - * pair for each network edge. Edge pairs are expected to be - * sorted in lexicographic order. - * - * - `backedges' is an array of `ne' integers, each a distinct - * index into `edges'. The edges in `edges', if permuted as - * specified by this array, should end up sorted in the _other_ - * lexicographic order, i.e. dest taking priority over source. - * - * - `capacity' is an array of `ne' integers, giving a maximum - * flow capacity for each edge. A negative value is taken to - * indicate unlimited capacity on that edge, but note that there - * may not be any unlimited-capacity _path_ from source to sink - * or an assertion will be failed. - * - * Output: - * - * - `flow' must be non-NULL. It is an array of `ne' integers, - * each giving the final flow along each edge. - * - * - `cut' may be NULL. If non-NULL, it is an array of `nv' - * integers, which will be set to zero or one on output, in such - * a way that: - * + the set of zero vertices includes the source - * + the set of one vertices includes the sink - * + the maximum flow capacity between the zero and one vertex - * sets is achieved (i.e. all edges from a zero vertex to a - * one vertex are at full capacity, while all edges from a - * one vertex to a zero vertex have no flow at all). - * - * - the returned value from the function is the total flow - * achieved. - */ -int maxflow_with_scratch(void *scratch, int nv, int source, int sink, - int ne, const int *edges, const int *backedges, - const int *capacity, int *flow, int *cut); - -/* - * The above function expects its `scratch' and `backedges' - * parameters to have already been set up. This allows you to set - * them up once and use them in multiple invocates of the - * algorithm. Now I provide functions to actually do the setting - * up. - */ -int maxflow_scratch_size(int nv); -void maxflow_setup_backedges(int ne, const int *edges, int *backedges); - -/* - * Simplified version of the above function. All parameters are the - * same, except that `scratch' and `backedges' are constructed - * internally. This is the simplest way to call the algorithm as a - * one-off; however, if you need to call it multiple times on the - * same network, it is probably better to call the above version - * directly so that you only construct `scratch' and `backedges' - * once. - * - * Additional return value is now -1, meaning that scratch space - * could not be allocated. - */ -int maxflow(int nv, int source, int sink, - int ne, const int *edges, const int *capacity, - int *flow, int *cut); - -#endif /* MAXFLOW_MAXFLOW_H */