problemreductions/rules/
minimumsummulticenter_ilp.rs

1//! Reduction from MinimumSumMulticenter to ILP (Integer Linear Programming).
2//!
3//! The p-median problem is formulated as a binary ILP.
4//!
5//! Variable layout (all binary):
6//! - `x_j` for each vertex j (1 if vertex j is selected as a center), indices `0..n`
7//! - `y_{i,j}` for each ordered pair (i, j), index `n + i*n + j`
8//!   (1 if vertex i is assigned to center j)
9//!
10//! Constraints:
11//! - Cardinality: Σ_j x_j = k (exactly k centers)
12//! - Assignment: ∀i: Σ_j y_{i,j} = 1 (each vertex assigned to exactly one center)
13//! - Assignment link: ∀i,j: if j is reachable from i then y_{i,j} ≤ x_j,
14//!   otherwise y_{i,j} = 0
15//!
16//! Objective: Minimize Σ_{i,j} w_i · d(i,j) · y_{i,j}
17//!
18//! Extraction: first n variables (x_j).
19//!
20//! Note: All-pairs shortest-path distances are computed using weighted shortest
21//! paths over `edge_lengths`. Unreachable assignment variables are forced to 0.
22
23use crate::models::algebraic::{LinearConstraint, ObjectiveSense, ILP};
24use crate::models::graph::MinimumSumMulticenter;
25use crate::reduction;
26use crate::rules::traits::{ReduceTo, ReductionResult};
27use crate::topology::{Graph, SimpleGraph};
28
29/// Result of reducing MinimumSumMulticenter to ILP.
30#[derive(Debug, Clone)]
31pub struct ReductionMSMCToILP {
32    target: ILP<bool>,
33    num_vertices: usize,
34}
35
36impl ReductionResult for ReductionMSMCToILP {
37    type Source = MinimumSumMulticenter<SimpleGraph, i64>;
38    type Target = ILP<bool>;
39
40    fn target_problem(&self) -> &ILP<bool> {
41        &self.target
42    }
43
44    fn extract_solution(
45        &self,
46        target_solution: &<Self::Target as crate::traits::Problem>::Solution,
47    ) -> crate::rules::ExtractionResult<<Self::Source as crate::traits::Problem>::Solution> {
48        crate::rules::traits::validate_target_solution(self.target_problem(), target_solution)?;
49
50        Ok(target_solution[..self.num_vertices]
51            .iter()
52            .map(|&value| value == 1)
53            .collect())
54    }
55}
56
57/// Compute weighted shortest-path distances from `source` in `graph`.
58///
59/// Returns a vector of length `n`; unreachable vertices remain `None`.
60fn weighted_distances_msmc(
61    graph: &SimpleGraph,
62    edge_lengths: &[i64],
63    source: usize,
64    n: usize,
65) -> Vec<Option<i64>> {
66    let mut adj: Vec<Vec<(usize, i64)>> = vec![Vec::new(); n];
67    for (idx, &(u, v)) in graph.edges().iter().enumerate() {
68        let len = edge_lengths[idx];
69        adj[u].push((v, len));
70        adj[v].push((u, len));
71    }
72
73    let mut dist = vec![None; n];
74    let mut visited = vec![false; n];
75    dist[source] = Some(0);
76
77    for _ in 0..n {
78        let mut next = None;
79        for vertex in 0..n {
80            if visited[vertex] {
81                continue;
82            }
83            let Some(dv) = dist[vertex] else {
84                continue;
85            };
86            match next {
87                None => next = Some(vertex),
88                Some(prev) => {
89                    if dv < dist[prev].expect("selected vertex must have a distance") {
90                        next = Some(vertex);
91                    }
92                }
93            }
94        }
95
96        let Some(u) = next else {
97            break;
98        };
99        visited[u] = true;
100        let du = dist[u].expect("selected vertex must have a distance");
101
102        for &(v, len) in &adj[u] {
103            if visited[v] {
104                continue;
105            }
106            let candidate = du + len;
107            let should_update = match dist[v] {
108                None => true,
109                Some(current) => candidate < current,
110            };
111            if should_update {
112                dist[v] = Some(candidate);
113            }
114        }
115    }
116
117    dist
118}
119
120#[reduction(transform = upper_bound {
121    max_constraint_magnitude_bits = "num_vertices + 1",
122    num_vars = "num_vertices + num_vertices^2",
123    num_constraints = "num_vertices^2 + 2 * num_vertices + 1",
124    num_nonzeros = "(num_vertices + num_vertices^2) * (num_vertices^2 + 2 * num_vertices + 1)",
125})]
126impl ReduceTo<ILP<bool>> for MinimumSumMulticenter<SimpleGraph, i64> {
127    type Result = ReductionMSMCToILP;
128
129    fn reduce_to(&self) -> Result<Self::Result, crate::rules::ReductionError> {
130        let n = self.num_vertices();
131        let k = Self::exact_i64(self.k(), "encoding the number of centers")?;
132        let vertex_weights = self.vertex_weights();
133        let edge_lengths = self.edge_lengths();
134
135        // Precompute all-pairs weighted shortest-path distances.
136        let all_dist: Vec<Vec<Option<i64>>> = (0..n)
137            .map(|s| weighted_distances_msmc(self.graph(), edge_lengths, s, n))
138            .collect();
139
140        // Index helpers.
141        let x_var = |j: usize| j;
142        let y_var = |i: usize, j: usize| n + i * n + j;
143
144        let num_vars = n + n * n;
145        // Capacity: n^2 + 2*n + 1
146        let mut constraints = Vec::with_capacity(n * n + 2 * n + 1);
147
148        // Cardinality constraint: Σ_j x_j = k
149        let center_terms: Vec<(usize, i64)> = (0..n).map(|j| (x_var(j), 1)).collect();
150        constraints.push(LinearConstraint::eq(center_terms, k));
151
152        // Assignment constraints: ∀i: Σ_j y_{i,j} = 1
153        for i in 0..n {
154            let terms: Vec<(usize, i64)> = (0..n).map(|j| (y_var(i, j), 1)).collect();
155            constraints.push(LinearConstraint::eq(terms, 1));
156        }
157
158        // Assignment link constraints:
159        // reachable pairs use y_{i,j} ≤ x_j, unreachable pairs force y_{i,j} = 0.
160        for (i, distances) in all_dist.iter().enumerate() {
161            for (j, distance) in distances.iter().enumerate() {
162                if distance.is_some() {
163                    constraints.push(LinearConstraint::le(
164                        vec![(y_var(i, j), 1), (x_var(j), -1)],
165                        0,
166                    ));
167                } else {
168                    constraints.push(LinearConstraint::eq(vec![(y_var(i, j), 1)], 0));
169                }
170            }
171        }
172
173        // Objective: Minimize Σ_{i,j} w_i · d(i,j) · y_{i,j}
174        let mut objective: Vec<(usize, i64)> = Vec::new();
175        for (i, &w) in vertex_weights.iter().enumerate() {
176            for (j, distance) in all_dist[i].iter().enumerate() {
177                if let Some(distance) = distance {
178                    let weighted_distance = w.checked_mul(*distance).ok_or_else(|| {
179                        crate::rules::ReductionError::integer_overflow::<
180                            MinimumSumMulticenter<SimpleGraph, i64>,
181                            ILP<bool>,
182                        >(
183                            "multiplying a vertex weight by a shortest-path distance"
184                        )
185                    })?;
186                    let coeff = weighted_distance;
187                    if coeff != 0 {
188                        objective.push((y_var(i, j), coeff));
189                    }
190                }
191            }
192        }
193
194        let target = ILP::new(num_vars, constraints, objective, ObjectiveSense::Minimize)
195            .map_err(Self::target_construction)?;
196        Ok(ReductionMSMCToILP {
197            target,
198            num_vertices: n,
199        })
200    }
201}
202
203#[cfg(feature = "example-db")]
204pub(crate) fn canonical_rule_example_specs() -> Vec<crate::example_db::specs::RuleExampleSpec> {
205    vec![crate::example_db::specs::RuleExampleSpec {
206        id: "minimumsummulticenter_to_ilp",
207        build: || {
208            // 3-vertex path: 0 - 1 - 2, unit weights, K=1
209            // Optimal center is vertex 1 with total distance 1+0+1 = 2.
210            let source = MinimumSumMulticenter::new(
211                SimpleGraph::new(3, vec![(0, 1), (1, 2)]),
212                vec![1i64; 3],
213                vec![1i64; 2],
214                1,
215            );
216            crate::example_db::specs::rule_example_via_ilp::<_, bool>(source)
217        },
218    }]
219}
220
221#[cfg(test)]
222#[path = "../unit_tests/rules/minimumsummulticenter_ilp.rs"]
223mod tests;
problemreductions/rules/minimumsummulticenter_ilp.rs

problemreductions/rules/
minimumsummulticenter_ilp.rs
