problemreductions/rules/
sumofsquarespartition_ilp.rs

1//! Reduction from SumOfSquaresPartition to ILP (Integer Linear Programming).
2//!
3//! The objective Σ_g (Σ_{i ∈ g} s_i)^2 is quadratic, so we linearize using McCormick:
4//!
5//! Variables:
6//! - x_{i,g}: binary, element i in group g (index: i * K + g)
7//! - z_{i,j,g}: binary product for x_{i,g} * x_{j,g} (index: n*K + (i*n + j) * K + g)
8//!
9//! Constraints:
10//! - Σ_g x_{i,g} = 1 for each element i (assignment)
11//! - McCormick for each (i,j,g):
12//!   z_{i,j,g} ≤ x_{i,g}, z_{i,j,g} ≤ x_{j,g}, z_{i,j,g} ≥ x_{i,g} + x_{j,g} - 1
13//!
14//! Note: Σ_g (Σ_i s_i * x_{i,g})^2 = Σ_g Σ_{i,j} s_i * s_j * x_{i,g} * x_{j,g}
15//!       which equals Σ_g Σ_{i,j} s_i * s_j * z_{i,j,g} after linearization.
16//!
17//! Objective: Minimize Σ_g Σ_{i,j} s_i * s_j * z_{i,j,g}
18
19use crate::models::algebraic::{LinearConstraint, ObjectiveSense, ILP};
20use crate::models::misc::SumOfSquaresPartition;
21use crate::reduction;
22use crate::rules::traits::{ReduceTo, ReductionResult};
23
24/// Result of reducing SumOfSquaresPartition to ILP.
25///
26/// Variable layout:
27/// - x_{i,g} at index i * K + g  (i ∈ 0..n, g ∈ 0..K)
28/// - z_{i,j,g} at index n*K + (i*n + j) * K + g  (i,j ∈ 0..n, g ∈ 0..K)
29///
30/// Total: n*K + n^2*K variables.
31#[derive(Debug, Clone)]
32pub struct ReductionSSPToILP {
33    target: ILP<bool>,
34    num_elements: usize,
35    num_groups: usize,
36}
37
38impl ReductionSSPToILP {
39    fn x_var(&self, i: usize, g: usize) -> usize {
40        i * self.num_groups + g
41    }
42
43    fn z_var(&self, i: usize, j: usize, g: usize) -> usize {
44        let n = self.num_elements;
45        let k = self.num_groups;
46        n * k + (i * n + j) * k + g
47    }
48}
49
50impl ReductionResult for ReductionSSPToILP {
51    type Source = SumOfSquaresPartition;
52    type Target = ILP<bool>;
53
54    fn target_problem(&self) -> &ILP<bool> {
55        &self.target
56    }
57
58    /// Extract solution: for each element i, find the unique group g where x_{i,g} = 1.
59    fn extract_solution(&self, target_solution: &[usize]) -> Vec<usize> {
60        let num_groups = self.num_groups;
61        (0..self.num_elements)
62            .map(|i| {
63                (0..num_groups)
64                    .find(|&g| {
65                        let idx = i * num_groups + g;
66                        idx < target_solution.len() && target_solution[idx] == 1
67                    })
68                    .unwrap_or(0)
69            })
70            .collect()
71    }
72}
73
74#[reduction(
75    overhead = {
76        num_vars = "num_elements * num_groups + num_elements^2 * num_groups",
77        num_constraints = "num_elements + 3 * num_elements^2 * num_groups",
78    }
79)]
80impl ReduceTo<ILP<bool>> for SumOfSquaresPartition {
81    type Result = ReductionSSPToILP;
82
83    fn reduce_to(&self) -> Self::Result {
84        let n = self.num_elements();
85        let k = self.num_groups();
86        let num_vars = n * k + n * n * k;
87
88        let result = ReductionSSPToILP {
89            target: ILP::empty(),
90            num_elements: n,
91            num_groups: k,
92        };
93
94        let mut constraints = Vec::new();
95
96        // Assignment constraints: for each element i, Σ_g x_{i,g} = 1
97        for i in 0..n {
98            let terms: Vec<(usize, f64)> = (0..k).map(|g| (result.x_var(i, g), 1.0)).collect();
99            constraints.push(LinearConstraint::eq(terms, 1.0));
100        }
101
102        // McCormick linearization for z_{i,j,g} = x_{i,g} * x_{j,g}
103        for i in 0..n {
104            for j in 0..n {
105                for g in 0..k {
106                    let z = result.z_var(i, j, g);
107                    let xi = result.x_var(i, g);
108                    let xj = result.x_var(j, g);
109
110                    // z ≤ x_{i,g}
111                    constraints.push(LinearConstraint::le(vec![(z, 1.0), (xi, -1.0)], 0.0));
112                    // z ≤ x_{j,g}
113                    constraints.push(LinearConstraint::le(vec![(z, 1.0), (xj, -1.0)], 0.0));
114                    // z ≥ x_{i,g} + x_{j,g} - 1  →  -z + x_{i,g} + x_{j,g} ≤ 1
115                    constraints.push(LinearConstraint::le(
116                        vec![(z, -1.0), (xi, 1.0), (xj, 1.0)],
117                        1.0,
118                    ));
119                }
120            }
121        }
122
123        // Objective: Minimize Σ_g Σ_{i,j} s_i * s_j * z_{i,j,g}
124        let sizes = self.sizes();
125        let mut objective: Vec<(usize, f64)> = Vec::new();
126        for i in 0..n {
127            for j in 0..n {
128                for g in 0..k {
129                    let coeff = sizes[i] as f64 * sizes[j] as f64;
130                    if coeff.abs() > 0.0 {
131                        objective.push((result.z_var(i, j, g), coeff));
132                    }
133                }
134            }
135        }
136
137        let target = ILP::new(num_vars, constraints, objective, ObjectiveSense::Minimize);
138
139        ReductionSSPToILP {
140            target,
141            num_elements: n,
142            num_groups: k,
143        }
144    }
145}
146
147#[cfg(feature = "example-db")]
148pub(crate) fn canonical_rule_example_specs() -> Vec<crate::example_db::specs::RuleExampleSpec> {
149    vec![crate::example_db::specs::RuleExampleSpec {
150        id: "sumofsquarespartition_to_ilp",
151        build: || {
152            // 4 elements [1, 2, 3, 4], K=2 groups
153            // Group {1,4}: sum=5, Group {2,3}: sum=5 → sos = 25+25 = 50
154            let source = SumOfSquaresPartition::new(vec![1, 2, 3, 4], 2);
155            crate::example_db::specs::rule_example_via_ilp::<_, bool>(source)
156        },
157    }]
158}
159
160#[cfg(test)]
161#[path = "../unit_tests/rules/sumofsquarespartition_ilp.rs"]
162mod tests;
problemreductions/rules/sumofsquarespartition_ilp.rs

problemreductions/rules/
sumofsquarespartition_ilp.rs
