diff --git a/docs/paper/reductions.typ b/docs/paper/reductions.typ index a412a3a29..4cc56197d 100644 --- a/docs/paper/reductions.typ +++ b/docs/paper/reductions.typ @@ -4737,6 +4737,48 @@ where $P$ is a penalty weight large enough that any constraint violation costs m _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 gp_qubo = load-example("GraphPartitioning", "QUBO") +#let gp_qubo_sol = gp_qubo.solutions.at(0) +#let gp_qubo_edges = gp_qubo.source.instance.graph.edges.map(e => (e.at(0), e.at(1))) +#let gp_qubo_n = gp_qubo.source.instance.graph.num_vertices +#let gp_qubo_m = gp_qubo_edges.len() +#let gp_qubo_penalty = gp_qubo_m + 1 +#let gp_qubo_diag = range(0, gp_qubo_n).map(i => gp_qubo.target.instance.matrix.at(i).at(i)) +#let gp_qubo_cut_edges = gp_qubo_edges.filter(e => gp_qubo_sol.source_config.at(e.at(0)) != gp_qubo_sol.source_config.at(e.at(1))) +#let gp_qubo_cut_size = gp_qubo_cut_edges.len() +#reduction-rule("GraphPartitioning", "QUBO", + example: true, + example-caption: [6-vertex balanced partition instance ($n = #gp_qubo_n$, $|E| = #gp_qubo_m$)], + extra: [ + *Step 1 -- Binary partition variables.* Introduce one binary variable per vertex: $x_i = 0$ means vertex $i$ is in the left block, $x_i = 1$ means it is in the right block. For the canonical instance, this gives $n = #gp_qubo_n$ QUBO variables: + $ x_0, x_1, x_2, x_3, x_4, x_5 $ + + *Step 2 -- Choose the balance penalty.* The source graph has $m = #gp_qubo_m$ edges, so the construction uses $P = m + 1 = #gp_qubo_penalty$. Any imbalance contributes at least $P$, which is already larger than the maximum possible cut size of any balanced partition. + + *Step 3 -- Fill the QUBO matrix.* The diagonal entries are $Q_(i i) = deg(i) + P(1 - n)$, which evaluates here to $(#gp_qubo_diag.map(str).join(", "))$. For every pair $i < j$, start from $Q_(i j) = 2P = #(2 * gp_qubo_penalty)$, then subtract $2$ when $(i,j)$ is an edge. Hence edge coefficients become $18$ while non-edge coefficients stay $20$; the exported upper-triangular matrix matches the issue example exactly.\ + + *Step 4 -- Verify a solution.* The exported witness is $bold(x) = (#gp_qubo_sol.target_config.map(str).join(", "))$, which is also the source partition encoding. The cut edges are #gp_qubo_cut_edges.map(e => "(" + str(e.at(0)) + "," + str(e.at(1)) + ")").join(", "), so the cut size is $#gp_qubo_cut_size$ and the balance penalty vanishes because exactly #(gp_qubo_n / 2) vertices are assigned to the right block #sym.checkmark. + ], +)[ + Graph Partitioning (minimum bisection) asks for a balanced bipartition minimizing the number of crossing edges. Lucas's Ising formulation @lucas2014 translates directly to a QUBO by combining a cut-counting quadratic objective with a quadratic equality penalty enforcing $sum_i x_i = n / 2$. The reduction uses one binary variable per source vertex, so the QUBO has exactly $n$ variables. +][ + _Construction._ Given an undirected graph $G = (V, E)$ with even $n = |V|$ and $m = |E|$, introduce binary variables $x_i in {0,1}$ for each vertex $i in V$. Interpret $x_i = 0$ as $i in A$ and $x_i = 1$ as $i in B$. The cut objective is: + $ H_"cut" = sum_((u,v) in E) (x_u + x_v - 2 x_u x_v) $ + because the term equals $1$ exactly when edge $(u,v)$ crosses the partition. To enforce balance, add: + $ H_"bal" = P (sum_i x_i - n/2)^2 $ + with penalty $P = m + 1$. The QUBO objective is $H = H_"cut" + H_"bal"$. + + Expanding $H_"bal"$ with $x_i^2 = x_i$ gives: + $ H_"bal" = P (1 - n) sum_i x_i + 2 P sum_(i < j) x_i x_j + P n^2 / 4 $ + so the upper-triangular QUBO coefficients are: + $ Q_(i i) = deg(i) + P (1 - n) $ + and for $i < j$, $Q_(i j) = 2 P$ for every pair, then subtract $2$ whenever $(i,j) in E$. Equivalently, edge pairs have coefficient $2P - 2$ and non-edge pairs have coefficient $2P$. The additive constant $P n^2 / 4$ does not affect the minimizer. + + _Correctness._ ($arrow.r.double$) If $bold(x)$ encodes a balanced partition, then $sum_i x_i = n/2$ and $H_"bal" = 0$, so the QUBO objective equals the cut size exactly. ($arrow.l.double$) If $bold(x)$ is imbalanced, then $|sum_i x_i - n/2| >= 1$, hence $H_"bal" >= P = m + 1$. Every balanced partition has cut size at most $m$, so any imbalanced assignment has objective strictly larger than at least one balanced assignment. Therefore every QUBO minimizer is balanced, and among balanced assignments minimizing $H$ is identical to minimizing the cut size. + + _Solution extraction._ Return the QUBO bit-vector directly: the same binary assignment already records the source partition. +] + #let qubo_ilp = load-example("QUBO", "ILP") #let qubo_ilp_sol = qubo_ilp.solutions.at(0) #reduction-rule("QUBO", "ILP", diff --git a/src/rules/graphpartitioning_qubo.rs b/src/rules/graphpartitioning_qubo.rs new file mode 100644 index 000000000..8e32846c9 --- /dev/null +++ b/src/rules/graphpartitioning_qubo.rs @@ -0,0 +1,99 @@ +//! Reduction from GraphPartitioning to QUBO. +//! +//! Uses the penalty-method QUBO +//! H = sum_(u,v in E) (x_u + x_v - 2 x_u x_v) + P (sum_i x_i - n/2)^2 +//! with P = |E| + 1 so any imbalanced partition is dominated by a balanced one. + +use crate::models::algebraic::QUBO; +use crate::models::graph::GraphPartitioning; +use crate::reduction; +use crate::rules::traits::{ReduceTo, ReductionResult}; +use crate::topology::{Graph, SimpleGraph}; + +/// Result of reducing GraphPartitioning to QUBO. +#[derive(Debug, Clone)] +pub struct ReductionGraphPartitioningToQUBO { + target: QUBO, +} + +impl ReductionResult for ReductionGraphPartitioningToQUBO { + type Source = GraphPartitioning; + type Target = QUBO; + + fn target_problem(&self) -> &Self::Target { + &self.target + } + + fn extract_solution(&self, target_solution: &[usize]) -> Vec { + target_solution.to_vec() + } +} + +#[reduction(overhead = { num_vars = "num_vertices" })] +impl ReduceTo> for GraphPartitioning { + type Result = ReductionGraphPartitioningToQUBO; + + fn reduce_to(&self) -> Self::Result { + let n = self.num_vertices(); + let penalty = self.num_edges() as f64 + 1.0; + let mut matrix = vec![vec![0.0f64; n]; n]; + let mut degrees = vec![0usize; n]; + let edges = self.graph().edges(); + + for &(u, v) in &edges { + degrees[u] += 1; + degrees[v] += 1; + } + + for (i, row) in matrix.iter_mut().enumerate() { + row[i] = degrees[i] as f64 + penalty * (1.0 - n as f64); + for value in row.iter_mut().skip(i + 1) { + *value = 2.0 * penalty; + } + } + + for (u, v) in edges { + let (lo, hi) = if u < v { (u, v) } else { (v, u) }; + matrix[lo][hi] -= 2.0; + } + + ReductionGraphPartitioningToQUBO { + target: QUBO::from_matrix(matrix), + } + } +} + +#[cfg(feature = "example-db")] +pub(crate) fn canonical_rule_example_specs() -> Vec { + use crate::export::SolutionPair; + + vec![crate::example_db::specs::RuleExampleSpec { + id: "graphpartitioning_to_qubo", + build: || { + crate::example_db::specs::rule_example_with_witness::<_, QUBO>( + GraphPartitioning::new(SimpleGraph::new( + 6, + vec![ + (0, 1), + (0, 2), + (1, 2), + (1, 3), + (2, 3), + (2, 4), + (3, 4), + (3, 5), + (4, 5), + ], + )), + SolutionPair { + source_config: vec![0, 0, 0, 1, 1, 1], + target_config: vec![0, 0, 0, 1, 1, 1], + }, + ) + }, + }] +} + +#[cfg(test)] +#[path = "../unit_tests/rules/graphpartitioning_qubo.rs"] +mod tests; diff --git a/src/rules/mod.rs b/src/rules/mod.rs index 4b55c9aae..439860c97 100644 --- a/src/rules/mod.rs +++ b/src/rules/mod.rs @@ -12,6 +12,7 @@ pub(crate) mod coloring_qubo; pub(crate) mod factoring_circuit; mod graph; pub(crate) mod graphpartitioning_maxcut; +pub(crate) mod graphpartitioning_qubo; mod kcoloring_casts; mod knapsack_qubo; mod ksatisfiability_casts; @@ -100,6 +101,7 @@ pub(crate) fn canonical_rule_example_specs() -> Vec MaxCut -> SpinGlass -> QUBO is better + ( + "GraphPartitioning {graph: \"SimpleGraph\"}", + "QUBO {weight: \"f64\"}", + ), // Knapsack -> ILP -> QUBO is better than the direct penalty reduction ("Knapsack", "QUBO {weight: \"f64\"}"), // MaxMatching → MaxSetPacking → ILP is better than direct MaxMatching → ILP diff --git a/src/unit_tests/rules/graphpartitioning_qubo.rs b/src/unit_tests/rules/graphpartitioning_qubo.rs new file mode 100644 index 000000000..779eea84f --- /dev/null +++ b/src/unit_tests/rules/graphpartitioning_qubo.rs @@ -0,0 +1,82 @@ +use super::*; +use crate::models::algebraic::QUBO; +use crate::rules::test_helpers::assert_optimization_round_trip_from_optimization_target; +use crate::topology::SimpleGraph; + +fn example_problem() -> GraphPartitioning { + GraphPartitioning::new(SimpleGraph::new( + 6, + vec![ + (0, 1), + (0, 2), + (1, 2), + (1, 3), + (2, 3), + (2, 4), + (3, 4), + (3, 5), + (4, 5), + ], + )) +} + +#[test] +fn test_graphpartitioning_to_qubo_closed_loop() { + let source = example_problem(); + let reduction = ReduceTo::>::reduce_to(&source); + + assert_optimization_round_trip_from_optimization_target( + &source, + &reduction, + "GraphPartitioning->QUBO closed loop", + ); +} + +#[test] +fn test_graphpartitioning_to_qubo_matrix_matches_issue_example() { + let source = example_problem(); + let reduction = ReduceTo::>::reduce_to(&source); + let qubo = reduction.target_problem(); + + assert_eq!(qubo.num_vars(), 6); + + let expected_diagonal = [-48.0, -47.0, -46.0, -46.0, -47.0, -48.0]; + for (index, expected) in expected_diagonal.into_iter().enumerate() { + assert_eq!(qubo.get(index, index), Some(&expected)); + } + + let edge_pairs = [ + (0, 1), + (0, 2), + (1, 2), + (1, 3), + (2, 3), + (2, 4), + (3, 4), + (3, 5), + (4, 5), + ]; + for &(u, v) in &edge_pairs { + assert_eq!(qubo.get(u, v), Some(&18.0), "edge ({u}, {v})"); + } + + let non_edge_pairs = [(0, 3), (0, 4), (0, 5), (1, 4), (1, 5), (2, 5)]; + for &(u, v) in &non_edge_pairs { + assert_eq!(qubo.get(u, v), Some(&20.0), "non-edge ({u}, {v})"); + } +} + +#[cfg(feature = "example-db")] +#[test] +fn test_graphpartitioning_to_qubo_canonical_example_spec() { + let spec = canonical_rule_example_specs() + .into_iter() + .find(|spec| spec.id == "graphpartitioning_to_qubo") + .expect("missing canonical GraphPartitioning -> QUBO example spec"); + let example = (spec.build)(); + + assert_eq!(example.source.problem, "GraphPartitioning"); + assert_eq!(example.target.problem, "QUBO"); + assert_eq!(example.target.instance["num_vars"], 6); + assert!(!example.solutions.is_empty()); +}