Max-flow

[pathikrit]

int maxFlow(int[][] cap, int source, int sink) {
  for (int flow = 0; ; ++flow)
    if (!augmentPath(cap, new boolean[cap.length], source, sink))
      return flow;
}

boolean augmentPath(int[][] cap, boolean[] vis, int i, int sink) {
  if (i == sink) return true;
  vis[i] = true;
  for (int j = 0; j < vis.length; j++)
    if (!vis[j] && cap[i][j] > 0 && augmentPath(cap, vis, j, sink)) {
      --cap[i][j];
      ++cap[j][i];
      return true;
    }
  return false;
}

[spluxx]

mind if I work on this?

[pathikrit]

no problem

[pathikrit]

@spluxx : I looked over your PR. It is very complex. Why can't we do this:

def maxFlow[V](capacity: Map[(V, V), Int], source: V, sink: V): Option[Int] = {
  val vertices = capacity.keySet.flatMap({case (u, v) => Seq(u, v)})
  val cap = mutable.Map(capacity.toSeq : _*).withDefaultValue(0)

  def augmentingPath(u: V, visited: mutable.Set[V] = mutable.Set.empty): Boolean =
    vertices exists {v =>
      val edge = u -> v
      val res = visited.add(v) && cap(edge) > 0 && augmentingPath(v, visited)
      if (res) {
        capacity(edge) = capacity(edge) - 1
        capacity(edge.swap) = capacity(edge.swap) + 1
      }
      res
    }

  Iterator.from(0).find(i => !augmentingPath(source))
}

This works for any vertex type V and works for both adjacency matrix vs adjacency list.

[spluxx]

@pathikrit The simplicity of your code is amazing, but I still think having flow 1 augmenting path is limiting the full power of the algorithm. How is (basically your code with some modifications / testing done)

  def maxFlow[V](capacity: Map[(V, V), Int], source: V, sink: V): Int = {
    val vertices = capacity.keySet.flatMap({case (u, v) => Seq(u, v)})
    val cap = mutable.Map(capacity.toSeq : _*).withDefaultValue(0)

    def augmentingPath(u: V, aug: Int, visited: mutable.Set[V]): Int =
      base(u == sink)(aug) {
        vertices.findWithDefault(0)(_ > 0) { v =>
          val edge = u -> v
          val res =
            if (cap(edge) > 0 && visited.add(v)) {
              augmentingPath(v, Math.min(aug, cap(edge)), visited)
            } else 0
          cap(edge) -= res
          cap(edge.swap) += res
          res
        }
      }

    assume(source != sink)

    Iterator.iterate(1) { _ =>
      augmentingPath(source, Int.MaxValue, mutable.Set(source))
    }.takeWhile(_ > 0).sum - 1
  }

[pathikrit]

@spluxx . I actually did not try it out before posting. Some improvement:

def maxFlow[V](capacity: Map[(V, V), Int], source: V, sink: V): Int = {
  val vertices = capacity.keySet.flatMap({case (u, v) => Seq(u, v)})
  val cap = mutable.Map(capacity.toSeq : _*).withDefaultValue(0)

  def augmentingPath(u: V, aug: Int, visited: mutable.Set[V]): Option[Int] = {
    def compute() = for {
      v <- vertices.view if visited.add(v) && u != v
      c <- cap.get(u -> v) if c > 0
      res <- augmentingPath(v, c min aug, visited)
    } yield {
      cap(u -> v) -= res
      cap(v -> u) += res
      res
    }

    u match {
      case `sink` => if (aug < Int.MaxValue) Some(aug) else None
      case _ => compute().headOption
    }
  }

  Iterator.continually(augmentingPath(u = source, aug = Int.MaxValue, visited = mutable.Set.empty))
    .takeWhile(_.isDefined)
    .flatten
    .sum
}

I have a question - if we maintain another table flow(u,v) which denotes the current flow, we don't need to sum in the end and simply return flow(source -> sink)?

[spluxx]

@pathikrit hmm keeping track of flow(u, v) for all pairs (u, v) seems unnecessary work. Ideas for improvement?:

1. sum(capacity(u->sink) - cap(u -> sink)) for all vertices u directed to sink.
2. maintain flow(u) to keep track of incoming flows on vertex u something like

   cap(u -> v) -= res // increase flow from u to v by res
   flow(v) += res // increase incoming flow on v by res
   cap(v -> u) += res

[spluxx]

  def maxFlow[V](capacity: Map[(V, V), Int], source: V, sink: V): Int = {
    val vertices = capacity.keySet.flatMap({case (u, v) => Seq(u, v)})
    val cap = mutable.Map(capacity.toSeq : _*).withDefaultValue(0)
    val flow = mutable.Map.empty[V, Int].withDefaultValue(0)

    def augmentingPath(u: V, aug: Int, visited: mutable.Set[V]): Option[Int] = {
      def compute() = {
        for {
          v <- vertices.view if cap(u -> v) > 0 && visited.add(v) && u != v
          res <- augmentingPath(v, cap(u -> v) min aug, visited)
        } yield {
          cap(u -> v) -= res
          flow(v) += res
          cap(v -> u) += res
          res
        }
      }

      u match {
        case `sink` => if (aug < Int.MaxValue) Some(aug) else None
        case _ => compute().headOption
      }
    }

    Iterator.continually(augmentingPath(u = source, aug = Int.MaxValue, visited = mutable.Set(source)))
      .takeWhile(_.isDefined).toList // -> we need to have this lazy iterator actually evaluated

    flow(sink)
  }

But I still think summing up the augments more "elegant" in some sense.

[pathikrit]

@spluxx : Yes you are right - the sum version is more elegant (this library prioritizes elegance and readability over constant-factor performance gains). This is the Ford-Fulkerson algorithm right which is O(flow*E). Here is my attempt at the Push-Relabel algorithm (FIFO version) which is O(V*V*sqrt(E)):

def pushRelabel[V](capacity: Map[(V, V), Double], source: V, sink: V) = {
  type Edge = (V, V)

  val vertices = capacity.keySet.flatMap({case (u, v) => Seq(u, v)})
  val heights = mutable.Map.empty[V, Int].withDefaultValue(0)
  val excessFlow = mutable.Map.empty[V, Double].withDefaultValue(0)
  val preFlow = mutable.Map.empty[Edge, Double].withDefaultValue(0)
  val visited = mutable.Map.empty[V, mutable.Set[V]].withDefaultValue(mutable.Set.empty)
  val queue = mutable.Queue.empty[V]

  def residualCapacity(edge: Edge) = capacity(edge) - preFlow(edge)

  def canPush(u: V)(v: V) = heights(u) > heights(v) && !visited(u)(v)

  def push(u: V)(v: V) = {
    val delta = excessFlow(u) min residualCapacity(u -> v)
    preFlow(u -> v) += delta
    preFlow(v -> u) -= delta
    excessFlow(u) -= delta
    excessFlow(v) += delta
    visited(u) += v
  }

  def relabel(u: V) = {
    visited(u).clear()
    val next = vertices.collect({case v if residualCapacity(u -> v) > 0 => heights(v)})
    if (next.nonEmpty) heights(u) = 1 + next.min
  }

  def discharge(u: V) = {
    while(excessFlow(u) > 0) {
      vertices.find(canPush(u)).map(push(u)).getOrElse(relabel(u))
    }
  }

  /*******[Init]******/
  heights(source) = vertices.size
  excessFlow(source) = Double.PositiveInfinity
  for {v <- vertices if v != source} {
    push(source)(v)
    queue += v 
  }

  /*******[Run]******/
  while(queue.nonEmpty) {
    val u = queue.dequeue()
    val h1 = heights(u)
    discharge(u)
    val h2 = heights(u)
    if (h2 > h1) queue.enqueue(u)
  }

  preFlow.toMap
}

This code returns the final network of flows.