problemreductions/rules/
travelingsalesman_qubo.rs

1//! Reduction from TravelingSalesman to QUBO.
2//!
3//! Uses the standard position-based QUBO encoding for TSP:
4//! - Binary variables x_{v,p} = 1 iff vertex v is at position p in the tour
5//! - H_A: each vertex appears exactly once (row constraint)
6//! - H_B: each position has exactly one vertex (column constraint)
7//! - H_C: objective encoding edge costs between consecutive positions
8
9use crate::models::algebraic::QUBO;
10use crate::models::graph::TravelingSalesman;
11use crate::reduction;
12use crate::rules::traits::{ReduceTo, ReductionResult};
13use crate::topology::{Graph, SimpleGraph};
14use std::collections::HashMap;
15
16/// Result of reducing TravelingSalesman to QUBO.
17#[derive(Debug, Clone)]
18pub struct ReductionTravelingSalesmanToQUBO {
19    target: QUBO<f64>,
20    num_vertices: usize,
21    num_edges: usize,
22    edge_index: HashMap<(usize, usize), usize>,
23}
24
25impl ReductionResult for ReductionTravelingSalesmanToQUBO {
26    type Source = TravelingSalesman<SimpleGraph, i32>;
27    type Target = QUBO<f64>;
28
29    fn target_problem(&self) -> &Self::Target {
30        &self.target
31    }
32
33    /// Decode position encoding back to edge-based configuration.
34    ///
35    /// The QUBO solution uses n^2 binary variables x_{v,p} (vertex v at position p).
36    /// We extract the tour order, then map consecutive pairs to edge indices.
37    fn extract_solution(&self, target_solution: &[usize]) -> Vec<usize> {
38        let n = self.num_vertices;
39
40        // For each position p, find the vertex v where x_{v,p} == 1
41        let mut tour = vec![0usize; n];
42        for p in 0..n {
43            for v in 0..n {
44                if target_solution[v * n + p] == 1 {
45                    tour[p] = v;
46                    break;
47                }
48            }
49        }
50
51        // Build edge-based config: for each consecutive pair in the tour, mark the edge
52        let mut config = vec![0usize; self.num_edges];
53        for p in 0..n {
54            let u = tour[p];
55            let v = tour[(p + 1) % n];
56            let key = (u.min(v), u.max(v));
57            if let Some(&idx) = self.edge_index.get(&key) {
58                config[idx] = 1;
59            }
60        }
61
62        config
63    }
64}
65
66#[reduction(
67    overhead = {
68        num_vars = "num_vertices^2",
69    }
70)]
71impl ReduceTo<QUBO<f64>> for TravelingSalesman<SimpleGraph, i32> {
72    type Result = ReductionTravelingSalesmanToQUBO;
73
74    fn reduce_to(&self) -> Self::Result {
75        let n = self.num_vertices();
76        let edges = self.edges();
77
78        // Build edge weight map (both directions for undirected lookup)
79        let mut edge_weight_map: HashMap<(usize, usize), f64> = HashMap::new();
80        let mut weight_sum: f64 = 0.0;
81        for &(u, v, w) in &edges {
82            let wf = w as f64;
83            edge_weight_map.insert((u, v), wf);
84            edge_weight_map.insert((v, u), wf);
85            weight_sum += wf.abs();
86        }
87
88        // Build edge index map: canonical (min, max) → edge index
89        let graph_edges = self.graph().edges();
90        let num_edges = graph_edges.len();
91        let mut edge_index: HashMap<(usize, usize), usize> = HashMap::new();
92        for (idx, &(u, v)) in graph_edges.iter().enumerate() {
93            edge_index.insert((u.min(v), u.max(v)), idx);
94        }
95
96        // Penalty weight: must exceed any possible tour cost
97        let a = 1.0 + weight_sum;
98
99        // Build n^2 x n^2 upper-triangular QUBO matrix
100        let dim = n * n;
101        let mut matrix = vec![vec![0.0f64; dim]; dim];
102
103        // Helper: add value to upper-triangular position
104        let mut add_upper = |i: usize, j: usize, val: f64| {
105            let (lo, hi) = if i <= j { (i, j) } else { (j, i) };
106            matrix[lo][hi] += val;
107        };
108
109        // H_A: each vertex visited exactly once (row constraint)
110        // For each vertex v: (sum_p x_{v,p} - 1)^2
111        // = sum_p x_{v,p}^2 - 2*sum_p x_{v,p} + 1
112        // = -sum_p x_{v,p} + 2*sum_{p1<p2} x_{v,p1}*x_{v,p2} + const
113        for v in 0..n {
114            for p in 0..n {
115                // Diagonal: -A (from expanding (sum - 1)^2, the -2*x + x^2 = -x for binary)
116                add_upper(v * n + p, v * n + p, -a);
117            }
118            for p1 in 0..n {
119                for p2 in (p1 + 1)..n {
120                    // Cross terms: 2*A * x_{v,p1} * x_{v,p2}
121                    add_upper(v * n + p1, v * n + p2, 2.0 * a);
122                }
123            }
124        }
125
126        // H_B: each position has exactly one vertex (column constraint)
127        // For each position p: (sum_v x_{v,p} - 1)^2
128        for p in 0..n {
129            for v in 0..n {
130                add_upper(v * n + p, v * n + p, -a);
131            }
132            for v1 in 0..n {
133                for v2 in (v1 + 1)..n {
134                    add_upper(v1 * n + p, v2 * n + p, 2.0 * a);
135                }
136            }
137        }
138
139        // H_C: distance objective
140        // For each pair (u, v), add cost for x_{u,p} * x_{v,p_next} and x_{v,p} * x_{u,p_next}
141        for u in 0..n {
142            for v in (u + 1)..n {
143                let cost = edge_weight_map.get(&(u, v)).copied().unwrap_or(a);
144                for p in 0..n {
145                    let p_next = (p + 1) % n;
146                    // x_{u,p} * x_{v,p_next}
147                    add_upper(u * n + p, v * n + p_next, cost);
148                    // x_{v,p} * x_{u,p_next}
149                    add_upper(v * n + p, u * n + p_next, cost);
150                }
151            }
152        }
153
154        let target = QUBO::from_matrix(matrix);
155
156        ReductionTravelingSalesmanToQUBO {
157            target,
158            num_vertices: n,
159            num_edges,
160            edge_index,
161        }
162    }
163}
164
165#[cfg(feature = "example-db")]
166pub(crate) fn canonical_rule_example_specs() -> Vec<crate::example_db::specs::RuleExampleSpec> {
167    use crate::export::SolutionPair;
168    use crate::models::algebraic::QUBO;
169
170    vec![crate::example_db::specs::RuleExampleSpec {
171        id: "travelingsalesman_to_qubo",
172        build: || {
173            let source = TravelingSalesman::new(
174                SimpleGraph::new(3, vec![(0, 1), (0, 2), (1, 2)]),
175                vec![1, 2, 3],
176            );
177            crate::example_db::specs::rule_example_with_witness::<_, QUBO<f64>>(
178                source,
179                SolutionPair {
180                    source_config: vec![1, 1, 1],
181                    target_config: vec![0, 0, 1, 1, 0, 0, 0, 1, 0],
182                },
183            )
184        },
185    }]
186}
187
188#[cfg(test)]
189#[path = "../unit_tests/rules/travelingsalesman_qubo.rs"]
190mod tests;
problemreductions/rules/travelingsalesman_qubo.rs

problemreductions/rules/
travelingsalesman_qubo.rs
