///|
/// Immutable snapshot of an original edge.
pub struct MfEdge[Cap] {
  from : Int
  to : Int
  cap : Cap
  flow : Cap
} derive(Debug, Eq)

///|
priv struct Residual {
  to : Int
  mut cap : Int64
} derive(Debug)

///|
/// Dinic's maximum-flow graph with Int or Int64 capacities.
pub struct MfGraph[_] {
  priv adjacency : Array[Array[Int]]
  priv residual : Array[Residual]
  priv sources : Array[Int]
} derive(Debug)

///|
pub fn[C] MfGraph::new(n : Int) -> MfGraph[C] {
  guard 0 <= n && n <= 100000000 else { panic() }
  { adjacency: Array::makei(n, _ => []), residual: [], sources: [], }
}

///|
pub fn[C : @algebra.FlowInt] MfGraph::add_edge(
  self : MfGraph[C],
  from : Int,
  to : Int,
  cap : C,
) -> Int {
  let n = self.adjacency.length()
  let cap = @algebra.FlowInt::to_i64(cap)
  guard 0 <= from && from < n && 0 <= to && to < n && cap >= 0L else { panic() }
  let index = self.residual.length()
  self.adjacency[from].push(index)
  self.adjacency[to].push(index + 1)
  self.residual.push({ to, cap, })
  self.residual.push({ to: from, cap: 0L, })
  self.sources.push(from)
  index / 2
}

///|
pub fn[C : @algebra.FlowInt] MfGraph::get_edge(
  self : MfGraph[C],
  i : Int,
) -> MfEdge[C] {
  guard 0 <= i && i < self.sources.length() else { panic() }
  let forward = self.residual[2 * i]
  let reverse = self.residual[2 * i + 1]
  {
    from: self.sources[i],
    to: forward.to,
    cap: @algebra.FlowInt::from_i64(forward.cap + reverse.cap),
    flow: @algebra.FlowInt::from_i64(reverse.cap),
  }
}

///|
pub fn[C : @algebra.FlowInt] MfGraph::edges(
  self : MfGraph[C],
) -> Array[MfEdge[C]] {
  Array::makei(self.sources.length(), i => self.get_edge(i))
}

///|
pub fn[C : @algebra.FlowInt] MfGraph::change_edge(
  self : MfGraph[C],
  i : Int,
  new_cap : C,
  new_flow : C,
) -> Unit {
  guard 0 <= i && i < self.sources.length() else { panic() }
  let cap = @algebra.FlowInt::to_i64(new_cap)
  let flow = @algebra.FlowInt::to_i64(new_flow)
  guard 0L <= flow && flow <= cap else { panic() }
  self.residual[2 * i].cap = cap - flow
  self.residual[2 * i + 1].cap = flow
}

///|
/// Adds up to flow_limit units of flow. May be called repeatedly on the residual graph.
pub fn[C : @algebra.FlowInt] MfGraph::flow(
  self : MfGraph[C],
  s : Int,
  t : Int,
  flow_limit? : C,
) -> C {
  let n = self.adjacency.length()
  let limit = match flow_limit {
    Some(x) => @algebra.FlowInt::to_i64(x)
    None => C::maximum()
  }
  guard 0 <= s && s < n && 0 <= t && t < n && s != t && limit >= 0L else {
    panic()
  }
  let mut total = 0L
  while total < limit {
    let level = Array::make(n, -1)
    let queue = [s]
    level[s] = 0
    let mut head = 0
    while head < queue.length() {
      let v = queue[head]
      head += 1
      for index in self.adjacency[v] {
        let edge = self.residual[index]
        if edge.cap > 0L && level[edge.to] == -1 {
          level[edge.to] = level[v] + 1
          queue.push(edge.to)
        }
      }
    }
    if level[t] == -1 {
      break
    }
    let cursor = Array::make(n, 0)
    while total < limit {
      let vertices = [s]
      let path : Array[Int] = []
      let bottleneck = [limit - total]
      let mut reached = false
      while !vertices.is_empty() {
        let v = vertices[vertices.length() - 1]
        if v == t {
          reached = true
          break
        }
        while cursor[v] < self.adjacency[v].length() {
          let edge = self.residual[self.adjacency[v][cursor[v]]]
          if edge.cap > 0L && level[edge.to] == level[v] + 1 {
            break
          }
          cursor[v] += 1
        }
        if cursor[v] == self.adjacency[v].length() {
          level[v] = -1
          ignore(vertices.pop())
          ignore(bottleneck.pop())
          if !path.is_empty() {
            ignore(path.pop())
            cursor[vertices[vertices.length() - 1]] += 1
          }
        } else {
          let index = self.adjacency[v][cursor[v]]
          let edge = self.residual[index]
          path.push(index)
          vertices.push(edge.to)
          bottleneck.push(
            Int64::min(bottleneck[bottleneck.length() - 1], edge.cap),
          )
        }
      }
      if !reached {
        break
      }
      let amount = bottleneck[bottleneck.length() - 1]
      total += amount
      for index in path {
        self.residual[index].cap -= amount
        self.residual[index ^ 1].cap += amount
      }
    }
  }
  @algebra.FlowInt::from_i64(total)
}

///|
/// Vertices reachable from s in the current residual graph. O(n+m).
pub fn[C] MfGraph::min_cut(self : MfGraph[C], s : Int) -> Array[Bool] {
  guard 0 <= s && s < self.adjacency.length() else { panic() }
  let reached = Array::make(self.adjacency.length(), false)
  reached[s] = true
  let queue = [s]
  let mut head = 0
  while head < queue.length() {
    let v = queue[head]
    head += 1
    for index in self.adjacency[v] {
      let edge = self.residual[index]
      if edge.cap > 0L && !reached[edge.to] {
        reached[edge.to] = true
        queue.push(edge.to)
      }
    }
  }
  reached
}