diff options
Diffstat (limited to 'test/monniaux/glpk-4.65/src/misc/ffalg.c')
| -rw-r--r-- | test/monniaux/glpk-4.65/src/misc/ffalg.c | 221 |
1 files changed, 221 insertions, 0 deletions
diff --git a/test/monniaux/glpk-4.65/src/misc/ffalg.c b/test/monniaux/glpk-4.65/src/misc/ffalg.c new file mode 100644 index 00000000..4ea2913d --- /dev/null +++ b/test/monniaux/glpk-4.65/src/misc/ffalg.c @@ -0,0 +1,221 @@ +/* ffalg.c (Ford-Fulkerson algorithm) */ + +/*********************************************************************** +* This code is part of GLPK (GNU Linear Programming Kit). +* +* Copyright (C) 2009-2013 Andrew Makhorin, Department for Applied +* Informatics, Moscow Aviation Institute, Moscow, Russia. All rights +* reserved. E-mail: <mao@gnu.org>. +* +* GLPK is free software: you can redistribute it and/or modify it +* under the terms of the GNU General Public License as published by +* the Free Software Foundation, either version 3 of the License, or +* (at your option) any later version. +* +* GLPK is distributed in the hope that it will be useful, but WITHOUT +* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY +* or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public +* License for more details. +* +* You should have received a copy of the GNU General Public License +* along with GLPK. If not, see <http://www.gnu.org/licenses/>. +***********************************************************************/ + +#include "env.h" +#include "ffalg.h" + +/*********************************************************************** +* NAME +* +* ffalg - Ford-Fulkerson algorithm +* +* SYNOPSIS +* +* #include "ffalg.h" +* void ffalg(int nv, int na, const int tail[], const int head[], +* int s, int t, const int cap[], int x[], char cut[]); +* +* DESCRIPTION +* +* The routine ffalg implements the Ford-Fulkerson algorithm to find a +* maximal flow in the specified flow network. +* +* INPUT PARAMETERS +* +* nv is the number of nodes, nv >= 2. +* +* na is the number of arcs, na >= 0. +* +* tail[a], a = 1,...,na, is the index of tail node of arc a. +* +* head[a], a = 1,...,na, is the index of head node of arc a. +* +* s is the source node index, 1 <= s <= nv. +* +* t is the sink node index, 1 <= t <= nv, t != s. +* +* cap[a], a = 1,...,na, is the capacity of arc a, cap[a] >= 0. +* +* NOTE: Multiple arcs are allowed, but self-loops are not allowed. +* +* OUTPUT PARAMETERS +* +* x[a], a = 1,...,na, is optimal value of the flow through arc a. +* +* cut[i], i = 1,...,nv, is 1 if node i is labelled, and 0 otherwise. +* The set of arcs, whose one endpoint is labelled and other is not, +* defines the minimal cut corresponding to the maximal flow found. +* If the parameter cut is NULL, the cut information are not stored. +* +* REFERENCES +* +* L.R.Ford, Jr., and D.R.Fulkerson, "Flows in Networks," The RAND +* Corp., Report R-375-PR (August 1962), Chap. I "Static Maximal Flow," +* pp.30-33. */ + +void ffalg(int nv, int na, const int tail[], const int head[], + int s, int t, const int cap[], int x[], char cut[]) +{ int a, delta, i, j, k, pos1, pos2, temp, + *ptr, *arc, *link, *list; + /* sanity checks */ + xassert(nv >= 2); + xassert(na >= 0); + xassert(1 <= s && s <= nv); + xassert(1 <= t && t <= nv); + xassert(s != t); + for (a = 1; a <= na; a++) + { i = tail[a], j = head[a]; + xassert(1 <= i && i <= nv); + xassert(1 <= j && j <= nv); + xassert(i != j); + xassert(cap[a] >= 0); + } + /* allocate working arrays */ + ptr = xcalloc(1+nv+1, sizeof(int)); + arc = xcalloc(1+na+na, sizeof(int)); + link = xcalloc(1+nv, sizeof(int)); + list = xcalloc(1+nv, sizeof(int)); + /* ptr[i] := (degree of node i) */ + for (i = 1; i <= nv; i++) + ptr[i] = 0; + for (a = 1; a <= na; a++) + { ptr[tail[a]]++; + ptr[head[a]]++; + } + /* initialize arc pointers */ + ptr[1]++; + for (i = 1; i < nv; i++) + ptr[i+1] += ptr[i]; + ptr[nv+1] = ptr[nv]; + /* build arc lists */ + for (a = 1; a <= na; a++) + { arc[--ptr[tail[a]]] = a; + arc[--ptr[head[a]]] = a; + } + xassert(ptr[1] == 1); + xassert(ptr[nv+1] == na+na+1); + /* now the indices of arcs incident to node i are stored in + * locations arc[ptr[i]], arc[ptr[i]+1], ..., arc[ptr[i+1]-1] */ + /* initialize arc flows */ + for (a = 1; a <= na; a++) + x[a] = 0; +loop: /* main loop starts here */ + /* build augmenting tree rooted at s */ + /* link[i] = 0 means that node i is not labelled yet; + * link[i] = a means that arc a immediately precedes node i */ + /* initially node s is labelled as the root */ + for (i = 1; i <= nv; i++) + link[i] = 0; + link[s] = -1, list[1] = s, pos1 = pos2 = 1; + /* breadth first search */ + while (pos1 <= pos2) + { /* dequeue node i */ + i = list[pos1++]; + /* consider all arcs incident to node i */ + for (k = ptr[i]; k < ptr[i+1]; k++) + { a = arc[k]; + if (tail[a] == i) + { /* a = i->j is a forward arc from s to t */ + j = head[a]; + /* if node j has been labelled, skip the arc */ + if (link[j] != 0) continue; + /* if the arc does not allow increasing the flow through + * it, skip the arc */ + if (x[a] == cap[a]) continue; + } + else if (head[a] == i) + { /* a = i<-j is a backward arc from s to t */ + j = tail[a]; + /* if node j has been labelled, skip the arc */ + if (link[j] != 0) continue; + /* if the arc does not allow decreasing the flow through + * it, skip the arc */ + if (x[a] == 0) continue; + } + else + xassert(a != a); + /* label node j and enqueue it */ + link[j] = a, list[++pos2] = j; + /* check for breakthrough */ + if (j == t) goto brkt; + } + } + /* NONBREAKTHROUGH */ + /* no augmenting path exists; current flow is maximal */ + /* store minimal cut information, if necessary */ + if (cut != NULL) + { for (i = 1; i <= nv; i++) + cut[i] = (char)(link[i] != 0); + } + goto done; +brkt: /* BREAKTHROUGH */ + /* walk through arcs of the augmenting path (s, ..., t) found in + * the reverse order and determine maximal change of the flow */ + delta = 0; + for (j = t; j != s; j = i) + { /* arc a immediately precedes node j in the path */ + a = link[j]; + if (head[a] == j) + { /* a = i->j is a forward arc of the cycle */ + i = tail[a]; + /* x[a] may be increased until its upper bound */ + temp = cap[a] - x[a]; + } + else if (tail[a] == j) + { /* a = i<-j is a backward arc of the cycle */ + i = head[a]; + /* x[a] may be decreased until its lower bound */ + temp = x[a]; + } + else + xassert(a != a); + if (delta == 0 || delta > temp) delta = temp; + } + xassert(delta > 0); + /* increase the flow along the path */ + for (j = t; j != s; j = i) + { /* arc a immediately precedes node j in the path */ + a = link[j]; + if (head[a] == j) + { /* a = i->j is a forward arc of the cycle */ + i = tail[a]; + x[a] += delta; + } + else if (tail[a] == j) + { /* a = i<-j is a backward arc of the cycle */ + i = head[a]; + x[a] -= delta; + } + else + xassert(a != a); + } + goto loop; +done: /* free working arrays */ + xfree(ptr); + xfree(arc); + xfree(link); + xfree(list); + return; +} + +/* eof */ |