#![allow(clippy::cast_possible_truncation)]
use crate::core::{Graph, IgraphResult, VertexId};
use super::max_flow::{max_flow_value, max_flow_with_residual};
pub fn st_mincut_value(
graph: &Graph,
source: VertexId,
target: VertexId,
capacity: Option<&[f64]>,
) -> IgraphResult<f64> {
max_flow_value(graph, source, target, capacity)
}
pub fn st_mincut(
graph: &Graph,
source: VertexId,
target: VertexId,
capacity: Option<&[f64]>,
) -> IgraphResult<StMincut> {
let (value, net) = max_flow_with_residual(graph, source, target, capacity)?;
let n = net.n;
let mut in_source = vec![false; n];
in_source[source as usize] = true;
let mut queue: Vec<u32> = Vec::with_capacity(n);
queue.push(source);
let mut head_ptr = 0_usize;
while head_ptr < queue.len() {
let v = queue[head_ptr] as usize;
head_ptr += 1;
for &arc in &net.arcs_out[v] {
let arc_us = arc as usize;
if net.cap[arc_us] <= 0.0 {
continue;
}
let w = net.head[arc_us] as usize;
if !in_source[w] {
in_source[w] = true;
queue.push(w as u32);
}
}
}
let mut partition: Vec<u32> = Vec::with_capacity(queue.len());
let mut partition2: Vec<u32> = Vec::with_capacity(n.saturating_sub(queue.len()));
for (v, &is_src) in in_source.iter().enumerate() {
if is_src {
partition.push(v as u32);
} else {
partition2.push(v as u32);
}
}
let m = graph.ecount();
let directed = graph.is_directed();
let edge_count =
u32::try_from(m).map_err(|_| crate::core::IgraphError::Internal("ecount overflows u32"))?;
let mut cut: Vec<u32> = Vec::new();
for eid in 0..edge_count {
let (u, v) = graph.edge(eid)?;
let u_in = in_source[u as usize];
let v_in = in_source[v as usize];
let crosses = if directed {
u_in && !v_in
} else {
u_in != v_in
};
if crosses {
cut.push(eid);
}
}
Ok(StMincut {
value,
cut,
partition,
partition2,
})
}
#[derive(Debug, Clone)]
pub struct StMincut {
pub value: f64,
pub cut: Vec<u32>,
pub partition: Vec<u32>,
pub partition2: Vec<u32>,
}
#[cfg(test)]
mod tests {
use super::*;
use crate::core::IgraphError;
fn approx(a: f64, b: f64) -> bool {
(a - b).abs() <= 1e-12_f64 * a.abs().max(b.abs()).max(1.0)
}
#[test]
fn rejects_source_equals_target() {
let mut g = Graph::new(2, true).expect("graph");
g.add_edge(0, 1).expect("edge");
let err = st_mincut_value(&g, 0, 0, None).unwrap_err();
match err {
IgraphError::InvalidArgument(_) => {}
other => panic!("expected InvalidArgument, got {other:?}"),
}
}
#[test]
fn rejects_out_of_range_source() {
let g = Graph::new(2, true).expect("graph");
let err = st_mincut_value(&g, 5, 0, None).unwrap_err();
match err {
IgraphError::VertexOutOfRange { id, n } => {
assert_eq!(id, 5);
assert_eq!(n, 2);
}
other => panic!("expected VertexOutOfRange, got {other:?}"),
}
}
#[test]
fn rejects_out_of_range_target() {
let g = Graph::new(2, true).expect("graph");
let err = st_mincut_value(&g, 0, 5, None).unwrap_err();
match err {
IgraphError::VertexOutOfRange { id, n } => {
assert_eq!(id, 5);
assert_eq!(n, 2);
}
other => panic!("expected VertexOutOfRange, got {other:?}"),
}
}
#[test]
fn isolated_endpoints_have_zero_cut() {
let g = Graph::new(4, true).expect("graph");
let cut = st_mincut_value(&g, 0, 3, None).expect("cut");
assert!(approx(cut, 0.0));
}
#[test]
fn single_edge_unit_cut() {
let mut g = Graph::new(2, true).expect("graph");
g.add_edge(0, 1).expect("edge");
let cut = st_mincut_value(&g, 0, 1, None).expect("cut");
assert!(approx(cut, 1.0));
}
#[test]
fn two_parallel_paths_cut_equals_two() {
let mut g = Graph::new(4, true).expect("graph");
g.add_edge(0, 1).expect("edge");
g.add_edge(1, 3).expect("edge");
g.add_edge(0, 2).expect("edge");
g.add_edge(2, 3).expect("edge");
let cut = st_mincut_value(&g, 0, 3, Some(&[1.0, 1.0, 1.0, 1.0])).expect("cut");
assert!(approx(cut, 2.0));
}
#[test]
fn bottleneck_directed() {
let mut g = Graph::new(4, true).expect("graph");
g.add_edge(0, 1).expect("edge");
g.add_edge(1, 2).expect("edge");
g.add_edge(2, 3).expect("edge");
let cut = st_mincut_value(&g, 0, 3, Some(&[5.0, 2.0, 7.0])).expect("cut");
assert!(approx(cut, 2.0));
}
#[test]
fn classic_max_flow_textbook_cut() {
let mut g = Graph::new(6, true).expect("graph");
let arcs = [
(0u32, 1u32),
(0, 2),
(1, 2),
(1, 3),
(2, 1),
(2, 4),
(3, 2),
(3, 5),
(4, 3),
(4, 5),
];
let caps = [16.0, 13.0, 10.0, 12.0, 4.0, 14.0, 9.0, 20.0, 7.0, 4.0];
for (u, v) in arcs {
g.add_edge(u, v).expect("edge");
}
let cut = st_mincut_value(&g, 0, 5, Some(&caps)).expect("cut");
assert!(approx(cut, 23.0));
}
#[test]
fn undirected_cut_matches_max_flow() {
let mut g = Graph::new(4, false).expect("graph");
for (a, b) in [(0u32, 1u32), (0, 2), (1, 2), (1, 3), (2, 3)] {
g.add_edge(a, b).expect("edge");
}
let cut = st_mincut_value(&g, 0, 3, Some(&[4.0, 2.0, 10.0, 2.0, 2.0])).expect("cut");
assert!(approx(cut, 4.0));
}
#[test]
fn equals_max_flow_value() {
let mut g = Graph::new(5, true).expect("graph");
for (s, t) in [(0u32, 1u32), (0, 2), (1, 3), (2, 3), (3, 4), (1, 4)] {
g.add_edge(s, t).expect("edge");
}
let caps = [3.0, 5.0, 2.0, 4.0, 6.0, 1.0];
let flow = max_flow_value(&g, 0, 4, Some(&caps)).expect("flow");
let cut = st_mincut_value(&g, 0, 4, Some(&caps)).expect("cut");
assert!(approx(flow, cut));
}
fn sum_cut_caps(cut: &[u32], caps: Option<&[f64]>) -> f64 {
if let Some(c) = caps {
cut.iter().map(|&e| c[e as usize]).sum()
} else {
#[allow(clippy::cast_precision_loss)]
let n = cut.len() as f64;
n
}
}
fn assert_partition_well_formed(n: u32, part: &[u32], part2: &[u32], src: u32, tgt: u32) {
assert!(part.windows(2).all(|w| w[0] < w[1]), "partition not sorted");
assert!(
part2.windows(2).all(|w| w[0] < w[1]),
"partition2 not sorted"
);
assert_eq!(
part.len() + part2.len(),
n as usize,
"partitions must cover V"
);
let mut seen = vec![false; n as usize];
for &v in part.iter().chain(part2.iter()) {
assert!(!seen[v as usize], "vertex {v} appears in both partitions");
seen[v as usize] = true;
}
assert!(seen.iter().all(|&b| b), "some vertex missing");
assert!(
part.contains(&src),
"source {src} not on source side {part:?}"
);
assert!(
part2.contains(&tgt),
"target {tgt} not on sink side {part2:?}"
);
}
#[test]
fn st_mincut_bottleneck_directed_partition() {
let mut g = Graph::new(4, true).expect("graph");
g.add_edge(0, 1).expect("edge");
g.add_edge(1, 2).expect("edge");
g.add_edge(2, 3).expect("edge");
let caps = [5.0, 2.0, 7.0];
let r = st_mincut(&g, 0, 3, Some(&caps)).expect("st_mincut");
assert!(approx(r.value, 2.0));
assert_eq!(r.cut, vec![1]);
assert_eq!(r.partition, vec![0, 1]);
assert_eq!(r.partition2, vec![2, 3]);
assert!(approx(sum_cut_caps(&r.cut, Some(&caps)), r.value));
assert_partition_well_formed(4, &r.partition, &r.partition2, 0, 3);
}
#[test]
fn st_mincut_two_parallel_paths_partition() {
let mut g = Graph::new(4, true).expect("graph");
g.add_edge(0, 1).expect("edge");
g.add_edge(1, 3).expect("edge");
g.add_edge(0, 2).expect("edge");
g.add_edge(2, 3).expect("edge");
let r = st_mincut(&g, 0, 3, None).expect("st_mincut");
assert!(approx(r.value, 2.0));
assert!(approx(sum_cut_caps(&r.cut, None), r.value));
assert_partition_well_formed(4, &r.partition, &r.partition2, 0, 3);
let v = st_mincut_value(&g, 0, 3, None).expect("value");
assert!(approx(v, r.value));
}
#[test]
fn st_mincut_isolated_endpoints() {
let g = Graph::new(4, true).expect("graph");
let r = st_mincut(&g, 0, 3, None).expect("st_mincut");
assert!(approx(r.value, 0.0));
assert!(r.cut.is_empty(), "cut must be empty when value = 0");
assert_eq!(r.partition, vec![0]);
assert_eq!(r.partition2, vec![1, 2, 3]);
}
#[test]
fn st_mincut_single_edge() {
let mut g = Graph::new(2, true).expect("graph");
g.add_edge(0, 1).expect("edge");
let r = st_mincut(&g, 0, 1, None).expect("st_mincut");
assert!(approx(r.value, 1.0));
assert_eq!(r.cut, vec![0]);
assert_eq!(r.partition, vec![0]);
assert_eq!(r.partition2, vec![1]);
}
#[test]
fn st_mincut_undirected_partition_crossings() {
let mut g = Graph::new(4, false).expect("graph");
g.add_edge(0, 1).expect("edge");
g.add_edge(0, 2).expect("edge");
g.add_edge(1, 3).expect("edge");
g.add_edge(2, 3).expect("edge");
let caps = [1.0, 1.0, 1.0, 1.0];
let r = st_mincut(&g, 0, 3, Some(&caps)).expect("st_mincut");
assert!(approx(r.value, 2.0));
assert!(approx(sum_cut_caps(&r.cut, Some(&caps)), 2.0));
assert_partition_well_formed(4, &r.partition, &r.partition2, 0, 3);
for &eid in &r.cut {
let (u, v) = g.edge(eid).expect("edge");
let u_in = r.partition.contains(&u);
let v_in = r.partition.contains(&v);
assert_ne!(
u_in, v_in,
"cut edge {eid} ({u}-{v}) does not cross the partition"
);
}
}
#[test]
fn st_mincut_multigraph_parallel_arcs() {
let mut g = Graph::new(2, true).expect("graph");
g.add_edge(0, 1).expect("edge");
g.add_edge(0, 1).expect("edge");
let r = st_mincut(&g, 0, 1, None).expect("st_mincut");
assert!(approx(r.value, 2.0));
let mut cut_sorted = r.cut.clone();
cut_sorted.sort_unstable();
assert_eq!(cut_sorted, vec![0, 1]);
assert_eq!(r.partition, vec![0]);
assert_eq!(r.partition2, vec![1]);
}
#[test]
fn st_mincut_classic_textbook_partition() {
let mut g = Graph::new(6, true).expect("graph");
let arcs = [
(0u32, 1u32),
(0, 2),
(1, 2),
(1, 3),
(2, 1),
(2, 4),
(3, 2),
(3, 5),
(4, 3),
(4, 5),
];
let caps = [16.0, 13.0, 10.0, 12.0, 4.0, 14.0, 9.0, 20.0, 7.0, 4.0];
for (u, v) in arcs {
g.add_edge(u, v).expect("edge");
}
let r = st_mincut(&g, 0, 5, Some(&caps)).expect("st_mincut");
assert!(approx(r.value, 23.0));
assert!(approx(sum_cut_caps(&r.cut, Some(&caps)), 23.0));
assert_partition_well_formed(6, &r.partition, &r.partition2, 0, 5);
for &eid in &r.cut {
let (u, v) = g.edge(eid).expect("edge");
assert!(
r.partition.contains(&u) && r.partition2.contains(&v),
"directed cut arc {eid} ({u}→{v}) must point S→V\\S"
);
}
}
#[test]
fn st_mincut_capacity_validation_propagates() {
let mut g = Graph::new(2, true).expect("graph");
g.add_edge(0, 1).expect("edge");
let err = st_mincut(&g, 0, 1, Some(&[1.0, 2.0])).unwrap_err();
match err {
IgraphError::InvalidArgument(_) => {}
other => panic!("expected InvalidArgument, got {other:?}"),
}
let err = st_mincut(&g, 0, 1, Some(&[-1.0])).unwrap_err();
assert!(matches!(err, IgraphError::InvalidArgument(_)));
let err = st_mincut(&g, 0, 0, None).unwrap_err();
assert!(matches!(err, IgraphError::InvalidArgument(_)));
let err = st_mincut(&g, 5, 1, None).unwrap_err();
assert!(matches!(err, IgraphError::VertexOutOfRange { id: 5, n: 2 }));
let err = st_mincut(&g, 0, 5, None).unwrap_err();
assert!(matches!(err, IgraphError::VertexOutOfRange { id: 5, n: 2 }));
}
#[test]
fn st_mincut_value_matches_max_flow_value() {
let mut g = Graph::new(5, true).expect("graph");
for (s, t) in [(0u32, 1u32), (0, 2), (1, 3), (2, 3), (3, 4), (1, 4)] {
g.add_edge(s, t).expect("edge");
}
let caps = [3.0, 5.0, 2.0, 4.0, 6.0, 1.0];
let flow = max_flow_value(&g, 0, 4, Some(&caps)).expect("flow");
let val = st_mincut_value(&g, 0, 4, Some(&caps)).expect("value");
let part = st_mincut(&g, 0, 4, Some(&caps)).expect("partition");
assert!(approx(flow, val));
assert!(approx(flow, part.value));
assert!(approx(sum_cut_caps(&part.cut, Some(&caps)), flow));
}
}
#[cfg(all(test, feature = "proptest-harness"))]
mod proptests {
use super::*;
use crate::core::Graph;
use proptest::prelude::*;
fn xorshift(mut r: u64) -> u64 {
r ^= r << 13;
r ^= r >> 7;
r ^= r << 17;
r
}
fn build_random(seed: u64, n: u32, m_max: u32, directed: bool) -> (Graph, Vec<f64>) {
let mut g = Graph::new(n, directed).expect("graph");
let mut state = seed | 1;
let mut caps: Vec<f64> = Vec::new();
for _ in 0..m_max {
state = xorshift(state);
let u = (state % u64::from(n)) as u32;
state = xorshift(state);
let v = (state % u64::from(n)) as u32;
if u == v {
continue;
}
state = xorshift(state);
let cap = f64::from((state % 16) as u32) + 1.0;
g.add_edge(u, v).expect("edge");
caps.push(cap);
}
(g, caps)
}
proptest! {
#[test]
fn mincut_equals_maxflow(
seed in any::<u64>(),
n in 2u32..8,
m in 1u32..16,
directed in any::<bool>(),
) {
let (g, caps) = build_random(seed, n, m, directed);
let s = (seed % u64::from(n)) as u32;
let t_raw = xorshift(seed) % u64::from(n);
let t = if t_raw as u32 == s { (s + 1) % n } else { t_raw as u32 };
prop_assume!(s != t);
let flow = max_flow_value(&g, s, t, Some(&caps)).expect("flow");
let cut = st_mincut_value(&g, s, t, Some(&caps)).expect("cut");
let scale = flow.abs().max(cut.abs()).max(1.0);
prop_assert!(
(flow - cut).abs() <= 1e-12_f64 * scale,
"duality violated: flow {flow} cut {cut} (n={n}, m={m}, directed={directed}, seed={seed})"
);
}
#[test]
fn st_mincut_partition_invariants(
seed in any::<u64>(),
n in 2u32..8,
m in 1u32..16,
directed in any::<bool>(),
) {
let (g, caps) = build_random(seed, n, m, directed);
let s = (seed % u64::from(n)) as u32;
let t_raw = xorshift(seed) % u64::from(n);
let t = if t_raw as u32 == s { (s + 1) % n } else { t_raw as u32 };
prop_assume!(s != t);
let result = st_mincut(&g, s, t, Some(&caps)).expect("st_mincut");
let value = st_mincut_value(&g, s, t, Some(&caps)).expect("value");
let scale = value.abs().max(result.value.abs()).max(1.0);
prop_assert!(
(value - result.value).abs() <= 1e-12_f64 * scale,
"values disagree: scalar={value} partition={}", result.value
);
prop_assert_eq!(
result.partition.len() + result.partition2.len(),
n as usize,
"partitions do not cover V"
);
let mut seen = vec![false; n as usize];
for &v in result.partition.iter().chain(result.partition2.iter()) {
prop_assert!(!seen[v as usize], "vertex {} duplicated", v);
seen[v as usize] = true;
}
prop_assert!(seen.iter().all(|&b| b), "some vertex unaccounted for");
prop_assert!(result.partition.contains(&s), "source missing from partition");
prop_assert!(result.partition2.contains(&t), "target missing from partition2");
prop_assert!(
result.partition.windows(2).all(|w| w[0] < w[1]),
"partition not sorted"
);
prop_assert!(
result.partition2.windows(2).all(|w| w[0] < w[1]),
"partition2 not sorted"
);
let mut sum_caps = 0.0_f64;
for &eid in &result.cut {
sum_caps += caps[eid as usize];
}
let scale = value.abs().max(sum_caps.abs()).max(1.0);
prop_assert!(
(sum_caps - value).abs() <= 1e-9_f64 * scale,
"cut capacity {sum_caps} does not match value {value}"
);
let cut_set: std::collections::HashSet<u32> = result.cut.iter().copied().collect();
let mut adj: Vec<Vec<u32>> = vec![Vec::new(); n as usize];
for eid in 0..g.ecount() as u32 {
if cut_set.contains(&eid) {
continue;
}
let (u, v) = g.edge(eid).expect("edge");
adj[u as usize].push(v);
if !directed {
adj[v as usize].push(u);
}
}
let mut visited = vec![false; n as usize];
visited[s as usize] = true;
let mut q = vec![s];
let mut hp = 0;
while hp < q.len() {
let v = q[hp] as usize;
hp += 1;
for &w in &adj[v] {
if !visited[w as usize] {
visited[w as usize] = true;
q.push(w);
}
}
}
prop_assert!(
!visited[t as usize],
"removing cut did not disconnect {} from {} (graph: n={}, m={}, directed={}, seed={}, cut={:?})",
s, t, n, m, directed, seed, result.cut
);
}
}
}