eta_graph/algorithms/
dinic.rs

1use eta_algorithms::data_structs::array::Array;
2use eta_algorithms::data_structs::queue::Queue;
3use eta_algorithms::data_structs::stack::Stack;
4use crate::handles::types::{Edge, VHandle, Weight};
5use crate::handles::{pack, set_wgt, vh, vhu, wgt};
6use crate::traits::{StoreVertex, WeightedEdgeManipulate};
7pub struct DinicGraph<'a, VertexType, VertexStorageType, EdgeStorageType>
8where
9    VertexStorageType: StoreVertex<VertexType=VertexType>,
10    EdgeStorageType: WeightedEdgeManipulate,
11{
12    pub vertices: &'a VertexStorageType,
13    pub edge_storage: EdgeStorageType,
14    pub layer_data: Array<Weight>,
15}
16
17impl<'a, VertexType, VertexStorageType, EdgeStorageType> DinicGraph<'a, VertexType, VertexStorageType, EdgeStorageType>
18where
19    VertexStorageType: StoreVertex<VertexType=VertexType>,
20    EdgeStorageType: WeightedEdgeManipulate,
21    VertexType: std::fmt::Debug + std::fmt::Display,
22{
23    pub fn from(vertices: &'a VertexStorageType, edge_storage: &EdgeStorageType, src_handle: VHandle, sink_handle: VHandle) -> Self {
24        let vertices_len = vertices.len();
25        let mut dinic_graph = DinicGraph {
26            vertices,
27            edge_storage: edge_storage.clone(),
28            layer_data: Array::new(vertices_len),
29        };
30
31        dinic_graph.perform_search(src_handle, sink_handle);
32        dinic_graph.finalize_flow_calc(edge_storage);
33        dinic_graph
34    }
35
36    fn finalize_flow_calc(&mut self, original_edges: &EdgeStorageType)
37    where
38        VertexStorageType: StoreVertex<VertexType=VertexType>,
39        EdgeStorageType: WeightedEdgeManipulate,
40        VertexType: std::fmt::Debug + std::fmt::Display,
41    {
42        let zipped_iters = original_edges.iter().zip(self.edge_storage.iter_mut());
43        for edges in zipped_iters {
44            let (original_edge, dinic_edge) = edges;
45            let original_wgt = wgt(*original_edge);
46            let current_wgt = wgt(*dinic_edge);
47            *dinic_edge = set_wgt(*dinic_edge, original_wgt - current_wgt);
48        }
49    }
50
51    pub fn perform_search(&mut self, src_handle: VHandle, sink_handle: VHandle) {
52        let mut stack = Stack::new(self.vertices.len());
53        let mut queue = Queue::<VHandle>::new_pow2_sized(self.vertices.len()); // Direct pointer access is faster than offsets
54        self.layer_data.fill(Weight::MAX);
55
56        while mark_levels(src_handle, sink_handle, &mut self.edge_storage, &mut queue, &mut self.layer_data).is_ok() {
57            loop {
58                let start_edges = self.edge_storage.edges_as_mut_ptr(src_handle);
59                let mut root_edge = pack(src_handle, Weight::MAX);
60                stack.push((start_edges, (&mut root_edge) as *mut Edge));
61
62                let mut augmented_path_found = false;
63                let mut bottleneck_value = Weight::MAX;
64                let mut current_layer;
65
66                while !stack.is_empty() {
67                    let (next_edges_ptr, edge_ptr) = stack.top_mut().unwrap();
68                    let outgoing_edge_val = unsafe { **edge_ptr };
69                    current_layer = self.layer_data[vh(outgoing_edge_val) as usize];
70
71                    if vh(outgoing_edge_val) == sink_handle {
72                        augmented_path_found = true;
73                    }
74
75                    // In case of augmented path found, we need to backtrack and ignore everything else
76                    if augmented_path_found {
77                        unsafe {
78                            let modified_edge = set_wgt(**edge_ptr, wgt(**edge_ptr) - bottleneck_value);
79                            (*edge_ptr).write(modified_edge);
80                        };
81                        stack.pop();
82                        continue;
83                    }
84
85                    if wgt(outgoing_edge_val) < bottleneck_value {
86                        bottleneck_value = wgt(outgoing_edge_val);
87                    }
88
89                    let next = next_edges_ptr.next();
90                    // Backtracking
91                    if next.is_none() {
92                        stack.pop();
93                        continue;
94                    }
95                    let next = next.unwrap();
96
97                    let next_edges = self.edge_storage.edges_as_mut_ptr(vh(*next));
98                    let next_edge_layer = self.layer_data[vh(*next) as usize];
99
100                    // Exploring deeper
101                    if wgt(*next) != 0 && next_edge_layer > current_layer {
102                        stack.push((next_edges, next));
103                    }
104                }
105
106                if !augmented_path_found {
107                    break;
108                }
109            }
110        }
111    }
112}
113
114pub(in crate) fn mark_levels<EdgeStorageType>(
115    src_handle: VHandle,
116    sink_handle: VHandle,
117    edge_storage: &mut EdgeStorageType,
118    queue: &mut Queue<VHandle>,
119    layer_data: &mut Array<Weight>,
120) -> Result<(), &'static str>
121where
122    EdgeStorageType: WeightedEdgeManipulate,
123{
124    let mut found_sink = false;
125    queue.push(src_handle);
126    let mut layer = 0;
127
128    let mut sibling_counter = 0;
129    let mut last_sibling_in_layer = 1;
130    let mut next_last_sibling_in_layer = 1;
131    layer_data[src_handle as usize] = 0;
132
133    while !queue.is_empty() {
134        let v_handle = queue.dequeue().unwrap();
135        if v_handle == sink_handle {
136            found_sink = true;
137        }
138
139        for next_edge in edge_storage.edges_iter_mut(v_handle){
140            let next_edge_layer = unsafe{*layer_data.index_unchecked(vhu(*next_edge))};
141            if next_edge_layer != Weight::MAX {
142                continue;
143            }
144
145            if wgt(*next_edge) == 0 {
146                continue;
147            }
148            layer_data[vh(*next_edge) as usize] = layer + 1;
149
150            queue.push(vh(*next_edge));
151            next_last_sibling_in_layer += 1;
152        }
153        sibling_counter += 1;
154        if sibling_counter == last_sibling_in_layer {
155            last_sibling_in_layer = next_last_sibling_in_layer;
156            layer += 1;
157        }
158    }
159
160    if !found_sink {
161        return Err("Sink not found");
162    }
163    Ok(())
164}
165