rs-graph 0.20.0

// Copyright (c) 2015-2021 Frank Fischer <frank-fischer@shadow-soft.de>
//
// This program 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.
//
// This program 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 this program.  If not, see  <http://www.gnu.org/licenses/>
//

//! This module implements Dinic' max flow algorithm
//!
//! # Example
//!
//! ```
//! use rs_graph::traits::*;
//! use rs_graph::maxflow::dinic;
//! use rs_graph::{EdgeVec, Net};
//! use rs_graph::string::{Data, from_ascii};
//!
//! let Data { graph: g, weights: upper, nodes } = from_ascii::<Net>(r"
//!      a---2-->b
//!     @|\      ^\
//!    / | \     | 4
//!   5  |  \    |  \
//!  /   |   |   |   @
//! s    1   1   2    t
//!  \   |   |   |   @
//!   5  |    \  |  /
//!    \ |     \ | 5
//!     @v      @|/
//!      c---2-->d
//!     ").unwrap();
//!
//! let s = nodes[&'s'];
//! let t = nodes[&'t'];
//! let v1 = nodes[&'a'];
//! let v2 = nodes[&'b'];
//! let v3 = nodes[&'c'];
//! let v4 = nodes[&'d'];
//!
//! let (value, flow, mut mincut) = dinic(&g, s, t, |e| upper[e.index()]);
//!
//! assert_eq!(value, 5);
//! assert!(g.edges().all(|e| flow[e] >= 0 && flow[e] <= upper[e.index()]));
//! assert!(g.nodes().filter(|&u| u != s && u != t).all(|u| {
//!     g.outedges(u).map(|(e,_)| flow[e]).sum::<usize>() ==
//!     g.inedges(u).map(|(e,_)| flow[e]).sum::<usize>()
//! }));
//!
//! mincut.sort_by_key(|u| u.index());
//! assert_eq!(mincut, vec![v1, s, v3]);
//! ```
//!
//! ```
//! use rs_graph::traits::*;
//! use rs_graph::maxflow::dinic;
//! use rs_graph::{EdgeVec, Net};
//! use rs_graph::string::{Data, from_ascii};
//!
//! let Data { graph: g, weights: upper, nodes } = from_ascii::<Net>(r"
//!                ---8-->a---10---
//!               /       |        \
//!              /        1  --3--  |
//!             /         | /     \ |
//!            /          v@       \v
//!      ---->b-----9---->c----8--->d----
//!     /      \         @         @^    \
//!   18        ---6--  /         / |     33
//!   /               \/         /  |      \
//!  /                /\    -----   |       @
//! s           --5--- |   /        |        t
//!  \         /       |  /         |       @
//!   27      |  ----2-|--         /       /
//!    \      | /      |  /----8---       6
//!     \     |/       @  |              /
//!      ---->e----9-->f------6---->g----
//!            \          |        @
//!             \         |       /
//!              --5----->h---4---
//!     ").unwrap();
//!
//! let s = nodes[&'s'];
//! let t = nodes[&'t'];
//!
//! assert_eq!(g.num_edges(), 18);
//!
//! let (value, flow, mut mincut) = dinic(&g, s, t, |e| upper[e.index()]);
//! assert_eq!(value, 29);
//!
//! mincut.sort_by_key(|u| u.index());
//! assert_eq!(mincut, "bcsef".chars().map(|v| nodes[&v]).collect::<Vec<_>>());
//! ```

use crate::maxflow::MaxFlow;

use crate::adapters::{Network, NetworkEdge};
use crate::traits::{Directed, GraphSize, IndexDigraph, IndexGraph};
use crate::vec::{EdgeVec, NodeVec};

use std::cmp::min;
use std::collections::VecDeque;

use crate::num::traits::NumAssign;

/// The dinic max-flow algorithm.
pub struct Dinic<'a, G, F>
where
    G: 'a + IndexDigraph<'a>,
{
    g: Network<'a, G>,
    nodes: NodeVec<'a, &'a G, NodeInfo<G::Edge>>,
    edges: EdgeVec<'a, Network<'a, G>, EdgeInfo<G::Edge, F>>,
    queue: VecDeque<G::Node>,
    value: F,
}

#[derive(Clone)]
struct NodeInfo<E> {
    dist: usize,
    first_lvl: Option<NetworkEdge<E>>,
}

#[derive(Clone)]
struct EdgeInfo<E, F> {
    flow: F,
    next_lvl: Option<NetworkEdge<E>>,
}

impl<'a, G, F> MaxFlow<'a> for Dinic<'a, G, F>
where
    G: IndexDigraph<'a>,
    F: NumAssign + Ord + Copy,
{
    type Graph = G;

    type Flow = F;

    fn new(g: &'a G) -> Self {
        Dinic {
            g: Network::new(g),
            nodes: NodeVec::new(
                g,
                NodeInfo {
                    dist: 0,
                    first_lvl: None,
                },
            ),
            edges: EdgeVec::new(
                Network::new(g),
                EdgeInfo {
                    flow: F::zero(),
                    next_lvl: None,
                },
            ),
            queue: VecDeque::with_capacity(g.num_nodes()),
            value: F::zero(),
        }
    }

    fn as_graph(&self) -> &'a Self::Graph {
        self.g.as_graph()
    }

    fn value(&self) -> F {
        self.value
    }

    fn flow(&self, e: G::Edge) -> F {
        self.edges[self.g.from_edge(e)].flow
    }

    fn solve<Us>(&mut self, src: G::Node, snk: G::Node, upper: Us)
    where
        Us: Fn(G::Edge) -> Self::Flow,
    {
        assert_ne!(
            self.g.node_id(src),
            self.g.node_id(snk),
            "Source and sink node must not be equal"
        );

        // initialize network flow of reverse edges
        for e in self.g.edges() {
            self.edges[e].flow = if e.is_forward() { F::zero() } else { upper(e.edge()) };
        }

        self.value = F::zero();
        while self.search(src, snk) {
            let v = self.augment(src, snk, None);
            self.value += v;
        }
    }

    fn mincut(&self) -> Vec<G::Node> {
        let n = self.g.num_nodes();
        self.g.nodes().filter(|&u| self.nodes[u].dist < n).collect()
    }
}

impl<'a, G, F> Dinic<'a, G, F>
where
    G: IndexDigraph<'a>,
    F: NumAssign + Ord + Copy,
{
    fn search(&mut self, src: G::Node, snk: G::Node) -> bool {
        let n = self.g.num_nodes();

        for u in self.g.nodes() {
            self.nodes[u].dist = n;
            self.nodes[u].first_lvl = None;
        }
        self.nodes[src].dist = 0;

        self.queue.clear();
        self.queue.push_back(src);

        let mut snk_d = n;
        while let Some(u) = self.queue.pop_front() {
            let d = self.nodes[u].dist;

            if d >= snk_d {
                return true;
            }

            for (e, v) in self.g.outedges(u) {
                if self.edges[e.reverse()].flow > F::zero() {
                    if self.nodes[v].dist == n {
                        self.nodes[v].dist = d + 1;
                        self.queue.push_back(v);
                        if v == snk {
                            snk_d = d + 1
                        }
                    } else if self.nodes[v].dist != d + 1 {
                        continue;
                    }
                    self.edges[e].next_lvl = self.nodes[u].first_lvl;
                    self.nodes[u].first_lvl = Some(e);
                }
            }
        }

        false
    }

    fn augment(&mut self, src: G::Node, snk: G::Node, target_flow: Option<F>) -> F {
        if src == snk {
            return target_flow.expect("Unbounded flow (are source and sink the same?)");
        }

        let mut df = F::zero();

        while let Some(e) = self.nodes[src].first_lvl {
            let f = e.reverse();
            let rem_cap = match target_flow {
                Some(target_flow) => min(self.edges[f].flow, target_flow - df),
                None => self.edges[f].flow,
            };
            if rem_cap > F::zero() {
                debug_assert!(self.g.src(e) == src);
                debug_assert!(self.g.snk(e) != src);
                let cf = self.augment(self.g.snk(e), snk, Some(rem_cap));
                self.edges[e].flow += cf;
                self.edges[f].flow -= cf;
                df += cf;
                if target_flow.map(|t| df == t).unwrap_or(false) {
                    break;
                }
            }

            // edge is saturated or blocked
            self.nodes[src].first_lvl = self.edges[e].next_lvl;
        }

        if df.is_zero() {
            // nothing can be sent from this node, delete the node and
            // all adjacent edges (we just remove the outgoing edges
            // so that they won't be seen)
            self.nodes[src].first_lvl = None;
        }

        df
    }
}