Skip to main content

problemreductions/rules/
minimummatrixcover_ilp.rs

1//! Reduction from MinimumMatrixCover to ILP (Integer Linear Programming).
2//!
3//! Uses McCormick linearization to convert the quadratic sign assignment
4//! objective into a linear program with binary variables.
5//!
6//! Binary variables x_i ∈ {0,1} where f(i) = 2x_i - 1.
7//! For i<j, auxiliary variables y_{ij} linearize x_i·x_j via:
8//!   y_{ij} ≤ x_i, y_{ij} ≤ x_j, y_{ij} ≥ x_i + x_j - 1
9
10use crate::models::algebraic::MinimumMatrixCover;
11use crate::models::algebraic::{ObjectiveSense, ILP};
12use crate::reduction;
13use crate::rules::ilp_helpers::mccormick_product;
14use crate::rules::traits::{ReduceTo, ReductionResult};
15
16/// Result of reducing MinimumMatrixCover to ILP.
17#[derive(Debug, Clone)]
18pub struct ReductionMinimumMatrixCoverToILP {
19    target: ILP<bool>,
20    n: usize,
21}
22
23impl ReductionResult for ReductionMinimumMatrixCoverToILP {
24    type Source = MinimumMatrixCover;
25    type Target = ILP<bool>;
26
27    fn target_problem(&self) -> &ILP<bool> {
28        &self.target
29    }
30
31    fn extract_solution(
32        &self,
33        target_solution: &<Self::Target as crate::traits::Problem>::Solution,
34    ) -> crate::rules::ExtractionResult<<Self::Source as crate::traits::Problem>::Solution> {
35        crate::rules::traits::validate_target_solution(self.target_problem(), target_solution)?;
36
37        Ok({
38            // First n variables are the sign variables x_0,...,x_{n-1}
39            target_solution[..self.n]
40                .iter()
41                .map(|&value| value == 1)
42                .collect()
43        })
44    }
45}
46
47/// Map pair (i,j) with i<j to auxiliary variable index.
48fn y_index(n: usize, i: usize, j: usize) -> usize {
49    debug_assert!(i < j);
50    // Index into upper triangle: sum_{k=0}^{i-1} (n-1-k) + (j - i - 1)
51    let offset: usize = (0..i).map(|k| n - 1 - k).sum();
52    n + offset + (j - i - 1)
53}
54
55#[reduction(
56    transform = exact {
57        max_constraint_magnitude_bits = "1",
58        num_vars = "num_rows + num_rows * (num_rows - 1) / 2",
59        num_constraints = "3 * num_rows * (num_rows - 1) / 2",
60        num_nonzeros = "7 * num_rows * (num_rows - 1) / 2",
61    },
62)]
63impl ReduceTo<ILP<bool>> for MinimumMatrixCover {
64    type Result = ReductionMinimumMatrixCoverToILP;
65
66    fn reduce_to(&self) -> Result<Self::Result, crate::rules::ReductionError> {
67        let n = self.num_rows();
68        let num_pairs = n * (n.saturating_sub(1)) / 2;
69        let num_vars = n + num_pairs;
70
71        // Build constraints: 3 per pair (i,j) with i<j
72        let mut constraints = Vec::with_capacity(3 * num_pairs);
73        for i in 0..n {
74            for j in (i + 1)..n {
75                let y = y_index(n, i, j);
76
77                constraints.extend(mccormick_product(y, i, j));
78            }
79        }
80
81        // Build objective coefficients.
82        // f(i)·f(j) = (2x_i-1)(2x_j-1) = 4x_ix_j - 2x_i - 2x_j + 1
83        //
84        // For i≠j (using y_{min(i,j),max(i,j)} for x_i·x_j):
85        //   a_ij · f(i)·f(j) = a_ij · (4·y_{..} - 2·x_i - 2·x_j + 1)
86        //
87        // For diagonal (i=j): f(i)² = 1, so a_ii contributes a_ii (constant).
88        //
89        // Objective = Σ_{i≠j} a_ij·(4·y - 2·x_i - 2·x_j + 1) + Σ_i a_ii
90        //           = Σ_{i<j} 4·(a_ij+a_ji)·y_{ij}
91        //             + Σ_k [-2·(Σ_{j≠k} (a_kj + a_jk))]·x_k
92        //             + constant
93        //
94        // The constant does not affect the minimizing assignment. The mapped
95        // source solution is evaluated by the source problem, so omit it here.
96
97        let matrix = self.matrix();
98        let mut obj_coeffs = vec![0i64; num_vars];
99
100        // y_{ij} coefficients: 4·(a_ij + a_ji) for each i<j
101        for (i, row_i) in matrix.iter().enumerate() {
102            for j in (i + 1)..n {
103                let y = y_index(n, i, j);
104                let coefficient = row_i[j]
105                    .checked_add(matrix[j][i])
106                    .and_then(|value| value.checked_mul(4))
107                    .ok_or_else(|| {
108                        crate::rules::ReductionError::integer_overflow::<
109                            MinimumMatrixCover,
110                            ILP<bool>,
111                        >(
112                            "computing an off-diagonal matrix-cover coefficient"
113                        )
114                    })?;
115                obj_coeffs[y] = coefficient;
116            }
117        }
118
119        // x_k coefficients: -2·Σ_{j≠k} (a_kj + a_jk)
120        for (k, row_k) in matrix.iter().enumerate() {
121            let sum = (0..n).filter(|&j| j != k).try_fold(0_i64, |total, j| {
122                let pair = row_k[j].checked_add(matrix[j][k]).ok_or_else(|| {
123                    crate::rules::ReductionError::integer_overflow::<MinimumMatrixCover, ILP<bool>>(
124                        "adding symmetric matrix-cover entries",
125                    )
126                })?;
127                total.checked_add(pair).ok_or_else(|| {
128                    crate::rules::ReductionError::integer_overflow::<MinimumMatrixCover, ILP<bool>>(
129                        "summing matrix-cover row coefficients",
130                    )
131                })
132            })?;
133            let coefficient = sum.checked_mul(-2).ok_or_else(|| {
134                crate::rules::ReductionError::integer_overflow::<MinimumMatrixCover, ILP<bool>>(
135                    "scaling a matrix-cover row coefficient",
136                )
137            })?;
138            obj_coeffs[k] = coefficient;
139        }
140
141        let objective: Vec<(usize, i64)> = obj_coeffs
142            .into_iter()
143            .enumerate()
144            .filter(|&(_, c)| c != 0)
145            .collect();
146
147        let target = ILP::new(num_vars, constraints, objective, ObjectiveSense::Minimize)
148            .map_err(Self::target_construction)?;
149
150        Ok(ReductionMinimumMatrixCoverToILP { target, n })
151    }
152}
153
154#[cfg(feature = "example-db")]
155pub(crate) fn canonical_rule_example_specs() -> Vec<crate::example_db::specs::RuleExampleSpec> {
156    use crate::export::SolutionPair;
157
158    vec![crate::example_db::specs::RuleExampleSpec {
159        id: "minimum_matrix_cover_to_ilp",
160        build: || {
161            // Use a small 2×2 instance for the rule example
162            let source = MinimumMatrixCover::new(vec![vec![0, 3], vec![2, 0]]);
163            // Config [0,1] → f=(-1,+1) → value = 0·1 + 3·(-1) + 2·(-1) + 0·1 = -5
164            // Config [1,0] → f=(+1,-1) → value = 0·1 + 3·(-1) + 2·(-1) + 0·1 = -5
165            // Config [0,0] → f=(-1,-1) → value = 0+3+2+0 = 5
166            // Config [1,1] → f=(+1,+1) → value = 0+3+2+0 = 5
167            // Optimal is [0,1] or [1,0] with value -5
168            // Source config [0,1], target config: x_0=0, x_1=1, y_{01}=0
169            crate::example_db::specs::rule_example_with_witness::<_, ILP<bool>>(
170                source,
171                SolutionPair {
172                    source_config: serde_json::json!(vec![false, true]),
173                    target_config: serde_json::json!(vec![0, 1, 0]),
174                },
175            )
176        },
177    }]
178}
179
180#[cfg(test)]
181#[path = "../unit_tests/rules/minimummatrixcover_ilp.rs"]
182mod tests;