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