|  | 
|  | 1 | +use std::collections::VecDeque; | 
|  | 2 | +use std::ops::{Add, AddAssign, Neg, Sub, SubAssign}; | 
|  | 3 | + | 
|  | 4 | +// We assume that graph vertices are numbered from 1 to n. | 
|  | 5 | + | 
|  | 6 | +/// Adjacency matrix | 
|  | 7 | +type Graph = Vec<Vec<usize>>; | 
|  | 8 | + | 
|  | 9 | +/// We assume that T::default() gives "zero" flow and T supports negative values | 
|  | 10 | +pub struct FlowEdge<T> { | 
|  | 11 | +    pub sink: usize, | 
|  | 12 | +    pub capacity: T, | 
|  | 13 | +    pub flow: T, | 
|  | 14 | +} | 
|  | 15 | + | 
|  | 16 | +pub struct FlowResultEdge<T> { | 
|  | 17 | +    pub source: usize, | 
|  | 18 | +    pub sink: usize, | 
|  | 19 | +    pub flow: T, | 
|  | 20 | +} | 
|  | 21 | + | 
|  | 22 | +impl<T: Clone + Copy + Add + AddAssign + Sub<Output = T> + SubAssign + Ord + Neg + Default> | 
|  | 23 | +    FlowEdge<T> | 
|  | 24 | +{ | 
|  | 25 | +    pub fn new(sink: usize, capacity: T) -> Self { | 
|  | 26 | +        FlowEdge { | 
|  | 27 | +            sink, | 
|  | 28 | +            capacity, | 
|  | 29 | +            flow: T::default(), | 
|  | 30 | +        } | 
|  | 31 | +    } | 
|  | 32 | +} | 
|  | 33 | + | 
|  | 34 | +pub struct DinicMaxFlow<T> { | 
|  | 35 | +    /// BFS Level of each vertex. starts from 1 | 
|  | 36 | +    level: Vec<usize>, | 
|  | 37 | + | 
|  | 38 | +    /// The index of the last visited edge connected to each vertex | 
|  | 39 | +    pub last_edge: Vec<usize>, | 
|  | 40 | + | 
|  | 41 | +    /// Holds wether the solution has already been calculated | 
|  | 42 | +    network_solved: bool, | 
|  | 43 | + | 
|  | 44 | +    pub source: usize, | 
|  | 45 | +    pub sink: usize, | 
|  | 46 | + | 
|  | 47 | +    /// Number of edges added to the residual network | 
|  | 48 | +    pub num_edges: usize, | 
|  | 49 | +    pub num_vertices: usize, | 
|  | 50 | + | 
|  | 51 | +    pub adj: Graph, | 
|  | 52 | + | 
|  | 53 | +    /// The list of flow edges | 
|  | 54 | +    pub edges: Vec<FlowEdge<T>>, | 
|  | 55 | +} | 
|  | 56 | + | 
|  | 57 | +impl<T: Clone + Copy + Add + AddAssign + Sub<Output = T> + SubAssign + Neg + Ord + Default> | 
|  | 58 | +    DinicMaxFlow<T> | 
|  | 59 | +{ | 
|  | 60 | +    pub fn new(source: usize, sink: usize, num_vertices: usize) -> Self { | 
|  | 61 | +        DinicMaxFlow { | 
|  | 62 | +            level: vec![0; num_vertices + 1], | 
|  | 63 | +            last_edge: vec![0; num_vertices + 1], | 
|  | 64 | +            network_solved: false, | 
|  | 65 | +            source, | 
|  | 66 | +            sink, | 
|  | 67 | +            num_edges: 0, | 
|  | 68 | +            num_vertices, | 
|  | 69 | +            adj: vec![vec![]; num_vertices + 1], | 
|  | 70 | +            edges: vec![], | 
|  | 71 | +        } | 
|  | 72 | +    } | 
|  | 73 | +    #[inline] | 
|  | 74 | +    pub fn add_edge(&mut self, source: usize, sink: usize, capacity: T) { | 
|  | 75 | +        self.edges.push(FlowEdge::new(sink, capacity)); | 
|  | 76 | +        // Add the reverse edge with zero capacity | 
|  | 77 | +        self.edges.push(FlowEdge::new(source, T::default())); | 
|  | 78 | +        // We inserted the m'th edge from source to sink | 
|  | 79 | +        self.adj[source].push(self.num_edges); | 
|  | 80 | +        self.adj[sink].push(self.num_edges + 1); | 
|  | 81 | +        self.num_edges += 2; | 
|  | 82 | +    } | 
|  | 83 | + | 
|  | 84 | +    fn bfs(&mut self) -> bool { | 
|  | 85 | +        let mut q: VecDeque<usize> = VecDeque::new(); | 
|  | 86 | +        q.push_back(self.source); | 
|  | 87 | + | 
|  | 88 | +        while !q.is_empty() { | 
|  | 89 | +            let v = q.pop_front().unwrap(); | 
|  | 90 | +            for &e in self.adj[v].iter() { | 
|  | 91 | +                if self.edges[e].capacity <= self.edges[e].flow { | 
|  | 92 | +                    continue; | 
|  | 93 | +                } | 
|  | 94 | +                let u = self.edges[e].sink; | 
|  | 95 | +                if self.level[u] != 0 { | 
|  | 96 | +                    continue; | 
|  | 97 | +                } | 
|  | 98 | +                self.level[u] = self.level[v] + 1; | 
|  | 99 | +                q.push_back(u); | 
|  | 100 | +            } | 
|  | 101 | +        } | 
|  | 102 | + | 
|  | 103 | +        self.level[self.sink] != 0 | 
|  | 104 | +    } | 
|  | 105 | + | 
|  | 106 | +    fn dfs(&mut self, v: usize, pushed: T) -> T { | 
|  | 107 | +        // We have pushed nothing, or we are at the sink | 
|  | 108 | +        if v == self.sink { | 
|  | 109 | +            return pushed; | 
|  | 110 | +        } | 
|  | 111 | +        for e_pos in self.last_edge[v]..self.adj[v].len() { | 
|  | 112 | +            let e = self.adj[v][e_pos]; | 
|  | 113 | +            let u = self.edges[e].sink; | 
|  | 114 | +            if (self.level[v] + 1) != self.level[u] || self.edges[e].capacity <= self.edges[e].flow | 
|  | 115 | +            { | 
|  | 116 | +                continue; | 
|  | 117 | +            } | 
|  | 118 | +            let down_flow = self.dfs( | 
|  | 119 | +                u, | 
|  | 120 | +                std::cmp::min(pushed, self.edges[e].capacity - self.edges[e].flow), | 
|  | 121 | +            ); | 
|  | 122 | +            if down_flow == T::default() { | 
|  | 123 | +                continue; | 
|  | 124 | +            } | 
|  | 125 | +            self.last_edge[v] = e_pos; | 
|  | 126 | +            self.edges[e].flow += down_flow; | 
|  | 127 | +            self.edges[e ^ 1].flow -= down_flow; | 
|  | 128 | +            return down_flow; | 
|  | 129 | +        } | 
|  | 130 | +        self.last_edge[v] = self.adj[v].len(); | 
|  | 131 | +        T::default() | 
|  | 132 | +    } | 
|  | 133 | + | 
|  | 134 | +    pub fn find_maxflow(&mut self, infinite_flow: T) -> T { | 
|  | 135 | +        self.network_solved = true; | 
|  | 136 | +        let mut total_flow: T = T::default(); | 
|  | 137 | +        loop { | 
|  | 138 | +            self.level.fill(0); | 
|  | 139 | +            self.level[self.source] = 1; | 
|  | 140 | +            // There is no longer a path from source to sink in the residual | 
|  | 141 | +            // network | 
|  | 142 | +            if !self.bfs() { | 
|  | 143 | +                break; | 
|  | 144 | +            } | 
|  | 145 | +            self.last_edge.fill(0); | 
|  | 146 | +            let mut next_flow = self.dfs(self.source, infinite_flow); | 
|  | 147 | +            while next_flow != T::default() { | 
|  | 148 | +                total_flow += next_flow; | 
|  | 149 | +                next_flow = self.dfs(self.source, infinite_flow); | 
|  | 150 | +            } | 
|  | 151 | +        } | 
|  | 152 | +        total_flow | 
|  | 153 | +    } | 
|  | 154 | + | 
|  | 155 | +    pub fn get_flow_edges(&mut self, infinite_flow: T) -> Vec<FlowResultEdge<T>> { | 
|  | 156 | +        if !self.network_solved { | 
|  | 157 | +            self.find_maxflow(infinite_flow); | 
|  | 158 | +        } | 
|  | 159 | +        let mut result = Vec::new(); | 
|  | 160 | +        for v in 1..self.adj.len() { | 
|  | 161 | +            for &e_ind in self.adj[v].iter() { | 
|  | 162 | +                let e = &self.edges[e_ind]; | 
|  | 163 | +                // Make sure that reverse edges from residual network are not | 
|  | 164 | +                // included | 
|  | 165 | +                if e.flow > T::default() { | 
|  | 166 | +                    result.push(FlowResultEdge { | 
|  | 167 | +                        source: v, | 
|  | 168 | +                        sink: e.sink, | 
|  | 169 | +                        flow: e.flow, | 
|  | 170 | +                    }); | 
|  | 171 | +                } | 
|  | 172 | +            } | 
|  | 173 | +        } | 
|  | 174 | +        result | 
|  | 175 | +    } | 
|  | 176 | +} | 
|  | 177 | + | 
|  | 178 | +#[cfg(test)] | 
|  | 179 | +mod tests { | 
|  | 180 | +    use super::*; | 
|  | 181 | +    #[test] | 
|  | 182 | +    fn small_graph() { | 
|  | 183 | +        let mut flow: DinicMaxFlow<i32> = DinicMaxFlow::new(1, 6, 6); | 
|  | 184 | +        flow.add_edge(1, 2, 16); | 
|  | 185 | +        flow.add_edge(1, 4, 13); | 
|  | 186 | +        flow.add_edge(2, 3, 12); | 
|  | 187 | +        flow.add_edge(3, 4, 9); | 
|  | 188 | +        flow.add_edge(3, 6, 20); | 
|  | 189 | +        flow.add_edge(4, 2, 4); | 
|  | 190 | +        flow.add_edge(4, 5, 14); | 
|  | 191 | +        flow.add_edge(5, 3, 7); | 
|  | 192 | +        flow.add_edge(5, 6, 4); | 
|  | 193 | + | 
|  | 194 | +        let max_flow = flow.find_maxflow(i32::MAX); | 
|  | 195 | +        assert_eq!(max_flow, 23); | 
|  | 196 | + | 
|  | 197 | +        let mut sm_out = vec![0; 7]; | 
|  | 198 | +        let mut sm_in = vec![0; 7]; | 
|  | 199 | + | 
|  | 200 | +        let flow_edges = flow.get_flow_edges(i32::MAX); | 
|  | 201 | +        for e in flow_edges { | 
|  | 202 | +            sm_out[e.source] += e.flow; | 
|  | 203 | +            sm_in[e.sink] += e.flow; | 
|  | 204 | +        } | 
|  | 205 | +        for i in 2..=5 { | 
|  | 206 | +            assert_eq!(sm_in[i], sm_out[i]); | 
|  | 207 | +        } | 
|  | 208 | +        assert_eq!(sm_in[1], 0); | 
|  | 209 | +        assert_eq!(sm_out[1], max_flow); | 
|  | 210 | +        assert_eq!(sm_in[6], max_flow); | 
|  | 211 | +        assert_eq!(sm_out[6], 0); | 
|  | 212 | +    } | 
|  | 213 | +} | 
0 commit comments