diff --git a/docs/paper/reductions.typ b/docs/paper/reductions.typ index 61b0b422c..07b45a924 100644 --- a/docs/paper/reductions.typ +++ b/docs/paper/reductions.typ @@ -3657,6 +3657,54 @@ where $P$ is a penalty weight large enough that any constraint violation costs m _Solution extraction._ Discard slack variables: return $bold(z)[0..n]$. ] +#let mwc_qubo = load-example("MinimumMultiwayCut", "QUBO") +#let mwc_qubo_sol = mwc_qubo.solutions.at(0) +#let mwc_qubo_edges = mwc_qubo.source.instance.graph.inner.edges.map(e => (e.at(0), e.at(1))) +#let mwc_qubo_weights = mwc_qubo.source.instance.edge_weights +#let mwc_qubo_terminals = mwc_qubo.source.instance.terminals +#let mwc_qubo_n = mwc_qubo.source.instance.graph.inner.nodes.len() +#let mwc_qubo_k = mwc_qubo_terminals.len() +#let mwc_qubo_nq = mwc_qubo_n * mwc_qubo_k +#let mwc_qubo_alpha = mwc_qubo_weights.fold(0, (a, w) => a + w) + 1 +#let mwc_qubo_cut_indices = mwc_qubo_sol.source_config.enumerate().filter(((i, v)) => v == 1).map(((i, _)) => i) +#let mwc_qubo_cut_cost = mwc_qubo_cut_indices.fold(0, (a, i) => a + mwc_qubo_weights.at(i)) +#reduction-rule("MinimumMultiwayCut", "QUBO", + example: true, + example-caption: [$n = #mwc_qubo_n$ vertices, $k = #mwc_qubo_k$ terminals $T = {#mwc_qubo_terminals.map(str).join(", ")}$, $|E| = #mwc_qubo_edges.len()$ edges], + extra: [ + *Step 1 -- Source instance.* The canonical graph has $n = #mwc_qubo_n$ vertices, $m = #mwc_qubo_edges.len()$ edges with weights $(#mwc_qubo_weights.map(str).join(", "))$, and $k = #mwc_qubo_k$ terminals $T = {#mwc_qubo_terminals.map(str).join(", ")}$. + + *Step 2 -- Introduce binary variables.* Assign $k = #mwc_qubo_k$ indicator variables per vertex: $x_(u,t) = 1$ means vertex $u$ belongs to terminal $t$'s component. This gives $n k = #mwc_qubo_n times #mwc_qubo_k = #mwc_qubo_nq$ QUBO variables: + $ underbrace(x_(0,0) x_(0,1) x_(0,2), "vertex 0") #h(4pt) underbrace(x_(1,0) x_(1,1) x_(1,2), "vertex 1") #h(4pt) dots.c #h(4pt) underbrace(x_(4,0) x_(4,1) x_(4,2), "vertex 4") $ + + *Step 3 -- Penalty coefficient.* $alpha = 1 + sum_(e in E) w(e) = 1 + #mwc_qubo_weights.map(str).join(" + ") = #mwc_qubo_alpha$. + + *Step 4 -- Build $H_A$ (constraints).* One-hot: diagonal entries $Q_(u k+t, u k+t) = -#mwc_qubo_alpha$, off-diagonal $Q_(u k+s, u k+t) = #(2 * mwc_qubo_alpha)$ within each vertex's group. Terminal pinning: for each terminal vertex $t_i$, the wrong-position diagonal entries $Q_(t_i k+s, t_i k+s) += #mwc_qubo_alpha$ for $s != i$, effectively canceling the one-hot incentive for those positions.\ + + *Step 5 -- Build $H_B$ (cut cost).* For each edge $(u,v)$ with weight $w$ and each pair $s != t$, add $w$ to $Q_(u k+s, v k+t)$. For example, edge $(0,1)$ with weight $2$ contributes $2$ to positions $(x_(0,0), x_(1,1))$, $(x_(0,0), x_(1,2))$, $(x_(0,1), x_(1,0))$, $(x_(0,1), x_(1,2))$, $(x_(0,2), x_(1,0))$, and $(x_(0,2), x_(1,1))$.\ + + *Step 6 -- Verify a solution.* The QUBO ground state $bold(x) = (#mwc_qubo_sol.target_config.map(str).join(", "))$ decodes to the partition: vertex 0 in component 0, vertices 1--3 in component 1, vertex 4 in component 2. Cut edges: $\{#mwc_qubo_cut_indices.map(i => "(" + str(mwc_qubo_edges.at(i).at(0)) + "," + str(mwc_qubo_edges.at(i).at(1)) + ")").join(", ")\}$ with total weight #mwc_qubo_cut_indices.map(i => str(mwc_qubo_weights.at(i))).join(" + ") $= #mwc_qubo_cut_cost$ #sym.checkmark. + ], +)[ + The multiway cut problem requires a partition of vertices into $k$ components — one per terminal — minimizing the total weight of edges crossing components. The penalty method (@sec:penalty-method) encodes two constraints as QUBO penalties: (1) each vertex belongs to exactly one component (one-hot), and (2) each terminal is pinned to its own component. The cut-cost Hamiltonian counts edge weight across distinct components. Reference: @Heidari2022. +][ + _Construction._ Given $G = (V, E)$ with $n = |V|$, edge weights $w: E -> RR_(>0)$, and $k$ terminals $T = {t_0, ..., t_(k-1)}$. Introduce $n k$ binary variables $x_(u,t) in {0,1}$ (indexed by $u dot k + t$), where $x_(u,t) = 1$ means vertex $u$ is in terminal $t$'s component. Let $alpha = 1 + sum_(e in E) w(e)$. + + The QUBO Hamiltonian is $H = H_A + H_B$ where: + $ H_A = alpha (sum_(u in V) (1 - sum_(t=0)^(k-1) x_(u,t))^2 + sum_(i=0)^(k-1) sum_(s != i) x_(t_i, s)) $ + The first term is a _one-hot constraint_ ensuring each vertex is assigned to exactly one component. The second term _pins_ each terminal $t_i$ to position $i$ by penalizing any other assignment. Expanding the one-hot term using $x^2 = x$: + $ Q_(u k+t, u k+t) = -alpha, quad Q_(u k+s, u k+t) = 2 alpha quad (s < t) $ + Terminal pinning adds $alpha$ to the diagonal $Q_(t_i k+s, t_i k+s)$ for $s != i$, canceling the one-hot incentive. + + The cut-cost Hamiltonian: + $ H_B = sum_((u,v) in E) sum_(s != t) w(u,v) dot x_(u,s) dot x_(v,t) $ + counts the total weight of edges whose endpoints lie in different components. + + _Correctness._ ($arrow.r.double$) A valid multiway cut with cost $C$ maps to a QUBO solution with $H_A = 0$ (valid partition with correct terminal pinning) and $H_B = C$. ($arrow.l.double$) If $H_A > 0$, the penalty $alpha > sum_e w(e)$ exceeds the entire cut-cost range, so any QUBO minimizer has $H_A = 0$, encoding a valid partition. Among valid partitions, $H_B$ equals the cut cost, and the minimizer achieves the minimum multiway cut. + + _Solution extraction._ For each vertex $u$, find terminal position $t$ with $x_(u,t) = 1$. For each edge $(u,v)$, output 1 (cut) if $u$ and $v$ are in different components, 0 otherwise. +] + #let qubo_ilp = load-example("QUBO", "ILP") #let qubo_ilp_sol = qubo_ilp.solutions.at(0) #reduction-rule("QUBO", "ILP", diff --git a/docs/paper/references.bib b/docs/paper/references.bib index f6a25dcb3..f1313a079 100644 --- a/docs/paper/references.bib +++ b/docs/paper/references.bib @@ -798,6 +798,14 @@ @article{papadimitriou1982 doi = {10.1145/322307.322309} } +@techreport{Heidari2022, + author = {Heidari, Shahrokh and Dinneen, Michael J. and Delmas, Patrice}, + title = {An Equivalent {QUBO} Model to the Minimum Multi-Way Cut Problem}, + institution = {Centre for Discrete Mathematics and Theoretical Computer Science, University of Auckland}, + year = {2022}, + number = {CDMTCS-565} +} + @phdthesis{booth1975, author = {Booth, Kellogg S.}, title = {{PQ}-Tree Algorithms}, diff --git a/src/rules/minimummultiwaycut_qubo.rs b/src/rules/minimummultiwaycut_qubo.rs new file mode 100644 index 000000000..610ec7397 --- /dev/null +++ b/src/rules/minimummultiwaycut_qubo.rs @@ -0,0 +1,163 @@ +//! Reduction from MinimumMultiwayCut to QUBO. +//! +//! Variable mapping: k*n binary variables x_{u,t} for each vertex u and +//! terminal position t. x_{u,t} = 1 means vertex u is assigned to terminal t's +//! component. Variable index: u * k + t. +//! +//! QUBO Hamiltonian: H = H_A + H_B +//! +//! H_A enforces valid partition (one-hot per vertex) and terminal pinning. +//! H_B encodes the cut cost objective. +//! +//! Reference: Heidari, Dinneen & Delmas (2022). + +use crate::models::algebraic::QUBO; +use crate::models::graph::MinimumMultiwayCut; +use crate::reduction; +use crate::rules::traits::{ReduceTo, ReductionResult}; +use crate::topology::{Graph, SimpleGraph}; + +/// Result of reducing MinimumMultiwayCut to QUBO. +#[derive(Debug, Clone)] +pub struct ReductionMinimumMultiwayCutToQUBO { + target: QUBO, + num_vertices: usize, + num_terminals: usize, + edges: Vec<(usize, usize)>, +} + +impl ReductionResult for ReductionMinimumMultiwayCutToQUBO { + type Source = MinimumMultiwayCut; + type Target = QUBO; + + fn target_problem(&self) -> &Self::Target { + &self.target + } + + /// Decode one-hot assignment: for each vertex find its terminal, then + /// for each edge check if endpoints are in different terminals. + fn extract_solution(&self, target_solution: &[usize]) -> Vec { + let k = self.num_terminals; + let n = self.num_vertices; + + // For each vertex, find which terminal position it is assigned to + let assignments: Vec = (0..n) + .map(|u| { + (0..k) + .find(|&t| target_solution[u * k + t] == 1) + .unwrap_or(0) + }) + .collect(); + + // For each edge, output 1 (cut) if endpoints differ, 0 (keep) otherwise + self.edges + .iter() + .map(|&(u, v)| { + if assignments[u] != assignments[v] { + 1 + } else { + 0 + } + }) + .collect() + } +} + +#[reduction(overhead = { num_vars = "num_terminals * num_vertices" })] +impl ReduceTo> for MinimumMultiwayCut { + type Result = ReductionMinimumMultiwayCutToQUBO; + + fn reduce_to(&self) -> Self::Result { + let n = self.num_vertices(); + let k = self.num_terminals(); + let edges = self.graph().edges(); + let edge_weights = self.edge_weights(); + let terminals = self.terminals(); + let nq = n * k; + + // Penalty: sum of all edge weights + 1 + let alpha: f64 = edge_weights.iter().map(|&w| (w as f64).abs()).sum::() + 1.0; + + let mut matrix = vec![vec![0.0f64; nq]; nq]; + + // Helper: add value to upper-triangular position + let mut add_upper = |i: usize, j: usize, val: f64| { + let (lo, hi) = if i <= j { (i, j) } else { (j, i) }; + matrix[lo][hi] += val; + }; + + // H_A: one-hot constraint per vertex + // (1 - sum_t x_{u,t})^2 = 1 - sum_t x_{u,t} + 2 * sum_{s Vec { + use crate::export::SolutionPair; + + vec![crate::example_db::specs::RuleExampleSpec { + id: "minimummultiwaycut_to_qubo", + build: || { + use crate::models::algebraic::QUBO; + use crate::models::graph::MinimumMultiwayCut; + use crate::topology::SimpleGraph; + let graph = SimpleGraph::new(5, vec![(0, 1), (1, 2), (2, 3), (3, 4), (0, 4), (1, 3)]); + let source = MinimumMultiwayCut::new(graph, vec![0, 2, 4], vec![2, 3, 1, 2, 4, 5]); + crate::example_db::specs::rule_example_with_witness::<_, QUBO>( + source, + SolutionPair { + source_config: vec![1, 0, 0, 1, 1, 0], + target_config: vec![1, 0, 0, 0, 1, 0, 0, 1, 0, 0, 1, 0, 0, 0, 1], + }, + ) + }, + }] +} + +#[cfg(test)] +#[path = "../unit_tests/rules/minimummultiwaycut_qubo.rs"] +mod tests; diff --git a/src/rules/mod.rs b/src/rules/mod.rs index 135b34c38..4740bd9bd 100644 --- a/src/rules/mod.rs +++ b/src/rules/mod.rs @@ -24,6 +24,7 @@ mod maximumindependentset_triangular; pub(crate) mod maximummatching_maximumsetpacking; mod maximumsetpacking_casts; pub(crate) mod maximumsetpacking_qubo; +pub(crate) mod minimummultiwaycut_qubo; pub(crate) mod minimumvertexcover_maximumindependentset; pub(crate) mod minimumvertexcover_minimumsetcovering; pub(crate) mod sat_circuitsat; @@ -94,6 +95,7 @@ pub(crate) fn canonical_rule_example_specs() -> Vec>::reduce_to(&source); + let qubo = reduction.target_problem(); + + let solver = BruteForce::new(); + let qubo_solutions = solver.find_all_best(qubo); + + assert!(!qubo_solutions.is_empty(), "QUBO solver found no solutions"); + + // All QUBO optimal solutions should extract to valid source solutions with cost 8 + for sol in &qubo_solutions { + let extracted = reduction.extract_solution(sol); + let metric = source.evaluate(&extracted); + assert_eq!(metric, SolutionSize::Valid(8)); + } +} + +#[test] +fn test_minimummultiwaycut_to_qubo_small() { + // 3 vertices, 2 terminals {0,2}, edges [(0,1),(1,2)] with weights [1,1] + let graph = SimpleGraph::new(3, vec![(0, 1), (1, 2)]); + let source = MinimumMultiwayCut::new(graph, vec![0, 2], vec![1, 1]); + + let reduction = ReduceTo::>::reduce_to(&source); + let qubo = reduction.target_problem(); + + let solver = BruteForce::new(); + let qubo_solutions = solver.find_all_best(qubo); + + assert!(!qubo_solutions.is_empty(), "QUBO solver found no solutions"); + + // All solutions should extract to valid cuts + for sol in &qubo_solutions { + let extracted = reduction.extract_solution(sol); + let metric = source.evaluate(&extracted); + // With 2 terminals and path 0-1-2, minimum cut is 1 (cut either edge) + assert_eq!(metric, SolutionSize::Valid(1)); + } +} + +#[test] +fn test_minimummultiwaycut_to_qubo_sizes() { + // 5 vertices, 3 terminals => QUBO has k*n = 3*5 = 15 variables + let graph = SimpleGraph::new(5, vec![(0, 1), (1, 2), (2, 3), (3, 4), (0, 4), (1, 3)]); + let source = MinimumMultiwayCut::new(graph, vec![0, 2, 4], vec![2, 3, 1, 2, 4, 5]); + + let reduction = ReduceTo::>::reduce_to(&source); + assert_eq!(reduction.target_problem().num_variables(), 15); +} + +#[test] +fn test_minimummultiwaycut_to_qubo_terminal_pinning() { + // Verify that in all QUBO optimal solutions, each terminal vertex is + // assigned to its own terminal position. + let graph = SimpleGraph::new(5, vec![(0, 1), (1, 2), (2, 3), (3, 4), (0, 4), (1, 3)]); + let terminals = vec![0, 2, 4]; + let source = MinimumMultiwayCut::new(graph, terminals.clone(), vec![2, 3, 1, 2, 4, 5]); + + let reduction = ReduceTo::>::reduce_to(&source); + let qubo = reduction.target_problem(); + + let solver = BruteForce::new(); + let qubo_solutions = solver.find_all_best(qubo); + + let k = terminals.len(); + for sol in &qubo_solutions { + for (t_pos, &t_vertex) in terminals.iter().enumerate() { + // Terminal vertex should be assigned to its own position + assert_eq!( + sol[t_vertex * k + t_pos], + 1, + "Terminal {} at position {} should be 1", + t_vertex, + t_pos + ); + // And not assigned to any other position + for s in 0..k { + if s != t_pos { + assert_eq!( + sol[t_vertex * k + s], + 0, + "Terminal {} at position {} should be 0", + t_vertex, + s + ); + } + } + } + } +}