problemreductions/rules/
eulerianpath_ilp.rs

1//! Reduction from EulerianPath to ILP (Integer Linear Programming).
2//!
3//! Encodes the directed Eulerian-trail witness structure as an ILP feasibility
4//! instance with integer variables. Given a directed multigraph
5//! `D = (V, A)` with arc occurrences `A = {a_1, ..., a_m}`:
6//!
7//! - For every compatible ordered pair `(a, b) in P` (where
8//!   `P = { (a, b) : a != b and head(a) = tail(b) }`), introduce an integer
9//!   successor variable `y_{a,b}` (intended `0/1`: `1` iff `b` immediately
10//!   follows `a` in the trail).
11//! - For every arc `a in A`, introduce integer variables `s_a`, `e_a`
12//!   (`0/1`: `s_a = 1` iff `a` is first, `e_a = 1` iff `a` is last) and an
13//!   integer position variable `u_a` (intended value `0..m-1`).
14//! - The predecessor and successor equalities, together with the unique start
15//!   and unique end constraints and Miller--Tucker--Zemlin-style ordering
16//!   constraints, force any feasible solution to encode a directed trail that
17//!   uses every arc occurrence exactly once.
18//!
19//! The empty-arc instance (`m = 0`) maps to the empty ILP with no variables
20//! and no constraints, which is vacuously feasible.
21//!
22//! References: Ebert, "Computing Eulerian trails," IPL 28(2):93--97 (1988);
23//! Bang-Jensen and Gutin, *Digraphs: Theory, Algorithms and Applications*,
24//! 2nd ed., Springer (2009).
25
26use crate::models::algebraic::{LinearConstraint, ObjectiveSense, ILP};
27use crate::models::graph::EulerianPath;
28use crate::reduction;
29use crate::rules::traits::{ReduceTo, ReductionResult};
30
31/// Result of reducing EulerianPath to `ILP<i64>`.
32///
33/// Variable layout (all in the non-negative integer domain, with explicit
34/// upper bounds enforcing the intended `0/1` and `0..m-1` ranges):
35/// - `y_{a,b}` at index `k` for the `k`-th compatible pair in `pairs` order
36///   (a single sweep over `(a, b)` with `a, b in 0..m`, `a != b`, in
37///   row-major order),
38/// - `s_a` at index `p + a` for `a in 0..m`,
39/// - `e_a` at index `p + m + a` for `a in 0..m`,
40/// - `u_a` at index `p + 2 * m + a` for `a in 0..m`,
41///
42/// where `p = pairs.len()` is the number of compatible ordered pairs.
43#[derive(Debug, Clone)]
44pub struct ReductionEulerianPathToILP {
45    target: ILP<i64>,
46    /// Compatible ordered pairs `(a, b)` in the order their `y_{a,b}` variables
47    /// appear in the ILP, for `m > 0`. Empty when `m = 0`.
48    pairs: Vec<(usize, usize)>,
49    /// Number of arc occurrences in the source instance.
50    num_arcs: usize,
51}
52
53impl ReductionEulerianPathToILP {
54    fn s_idx(&self, a: usize) -> usize {
55        self.pairs.len() + a
56    }
57}
58
59impl ReductionResult for ReductionEulerianPathToILP {
60    type Source = EulerianPath;
61    type Target = ILP<i64>;
62
63    fn target_problem(&self) -> &ILP<i64> {
64        &self.target
65    }
66
67    /// Decode an ILP assignment into a source arc ordering.
68    ///
69    /// Reads the unique active start arc (`s_a = 1`) and walks the active
70    /// successor relation (`y_{a,b} = 1`) one step at a time, producing an arc
71    /// permutation of length `m`. Malformed assignments return an extraction
72    /// error instead of fabricating an ordering.
73    fn extract_solution(
74        &self,
75        target_solution: &<Self::Target as crate::traits::Problem>::Solution,
76    ) -> crate::rules::ExtractionResult<<Self::Source as crate::traits::Problem>::Solution> {
77        crate::rules::traits::validate_target_solution(self.target_problem(), target_solution)?;
78
79        Ok({
80            let m = self.num_arcs;
81            if m == 0 {
82                return Ok(Vec::new());
83            }
84
85            // Find the unique active start arc.
86            let mut current = match (0..m).find(|&a| target_solution[self.s_idx(a)] == 1) {
87                Some(a) => a,
88                None => {
89                    return Err(crate::rules::ExtractionError::invalid(
90                        "ILP witness has no active Eulerian-path start arc",
91                    ));
92                }
93            };
94
95            // Walk the active successor relation, recording each visited arc.
96            let mut order = Vec::with_capacity(m);
97            let mut visited = vec![false; m];
98            order.push(current);
99            visited[current] = true;
100
101            for _ in 1..m {
102                let next = self
103                    .pairs
104                    .iter()
105                    .enumerate()
106                    .find(|&(k, &(a, _))| a == current && target_solution[k] == 1)
107                    .map(|(_, &(_, b))| b);
108
109                match next {
110                    Some(b) if !visited[b] => {
111                        order.push(b);
112                        visited[b] = true;
113                        current = b;
114                    }
115                    _ => {
116                        return Err(crate::rules::ExtractionError::invalid(format!(
117                            "ILP witness has no unvisited successor for arc {current}",
118                        )));
119                    }
120                }
121            }
122            order
123        })
124    }
125}
126
127/// Enumerate compatible ordered pairs `(a, b)` with `a != b` and
128/// `head(a) = tail(b)`. The order is `a`-major then `b`-major, matching the
129/// nested-loop construction below.
130fn compatible_pairs(arcs: &[(usize, usize)]) -> Vec<(usize, usize)> {
131    let mut pairs = Vec::new();
132    for (a, &(_, head_a)) in arcs.iter().enumerate() {
133        for (b, &(tail_b, _)) in arcs.iter().enumerate() {
134            if a != b && tail_b == head_a {
135                pairs.push((a, b));
136            }
137        }
138    }
139    pairs
140}
141
142#[reduction(
143    transform = upper_bound {
144        num_vars = "3 * num_arcs + num_arcs * num_arcs",
145        num_constraints = "5 * num_arcs + 2 * num_arcs * num_arcs + 2",
146    },
147    unavailable = {
148        num_nonzeros = "the exact target parameter is not represented by this reduction's symbolic transform",
149    }
150)]
151impl ReduceTo<ILP<i64>> for EulerianPath {
152    type Result = ReductionEulerianPathToILP;
153
154    fn reduce_to(&self) -> Result<Self::Result, crate::rules::ReductionError> {
155        let arcs = self.graph().arcs();
156        let m = arcs.len();
157
158        // Empty-arc instance: vacuously feasible empty ILP.
159        if m == 0 {
160            let target = ILP::new(0, Vec::new(), Vec::new(), ObjectiveSense::Minimize)
161                .map_err(Self::target_construction)?;
162            return Ok(ReductionEulerianPathToILP {
163                target,
164                pairs: Vec::new(),
165                num_arcs: 0,
166            });
167        }
168
169        let pairs = compatible_pairs(&arcs);
170        let p = pairs.len();
171        let num_vars = p + 3 * m;
172
173        // Index helpers (mirroring the struct's accessors, but we need them
174        // before the struct exists).
175        let y_idx = |k: usize| -> usize { k };
176        let s_idx = |a: usize| -> usize { p + a };
177        let e_idx = |a: usize| -> usize { p + m + a };
178        let u_idx = |a: usize| -> usize { p + 2 * m + a };
179
180        // Inverse index: for each arc `a`, list `k` indices of pairs ending in
181        // `a` (`pairs[k].1 == a`) and pairs starting at `a` (`pairs[k].0 == a`).
182        let mut incoming: Vec<Vec<usize>> = vec![Vec::new(); m];
183        let mut outgoing: Vec<Vec<usize>> = vec![Vec::new(); m];
184        for (k, &(a, b)) in pairs.iter().enumerate() {
185            outgoing[a].push(k);
186            incoming[b].push(k);
187        }
188
189        let mut constraints: Vec<LinearConstraint> = Vec::new();
190
191        // (1) Predecessor equality: s_a + sum_{(b,a) in P} y_{b,a} = 1.
192        // (2) Successor equality:   e_a + sum_{(a,b) in P} y_{a,b} = 1.
193        let m_i64 = Self::exact_i64(m, "encoding the arc order")?;
194        for a in 0..m {
195            let mut pred_terms: Vec<(usize, i64)> = vec![(s_idx(a), 1)];
196            for &k in &incoming[a] {
197                pred_terms.push((y_idx(k), 1));
198            }
199            constraints.push(LinearConstraint::eq(pred_terms, 1));
200
201            let mut succ_terms: Vec<(usize, i64)> = vec![(e_idx(a), 1)];
202            for &k in &outgoing[a] {
203                succ_terms.push((y_idx(k), 1));
204            }
205            constraints.push(LinearConstraint::eq(succ_terms, 1));
206        }
207
208        // (3) Binary upper bounds on start / end variables, and position
209        //     upper bound on `u_a`.
210        for a in 0..m {
211            constraints.push(LinearConstraint::le(vec![(s_idx(a), 1)], 1));
212            constraints.push(LinearConstraint::le(vec![(e_idx(a), 1)], 1));
213            constraints.push(LinearConstraint::le(vec![(u_idx(a), 1)], m_i64 - 1));
214        }
215
216        // (4) Binary upper bounds on successor variables.
217        // (5) Order consistency (MTZ): u_b >= u_a + 1 - m * (1 - y_{a,b})
218        //     i.e.  u_a - u_b + m * y_{a,b} <= m - 1.
219        for (k, &(a, b)) in pairs.iter().enumerate() {
220            constraints.push(LinearConstraint::le(vec![(y_idx(k), 1)], 1));
221            constraints.push(LinearConstraint::le(
222                vec![(u_idx(a), 1), (u_idx(b), -1), (y_idx(k), m_i64)],
223                m_i64 - 1,
224            ));
225        }
226
227        // (6) Unique start: sum_a s_a = 1.
228        // (7) Unique end:   sum_a e_a = 1.
229        let start_sum: Vec<(usize, i64)> = (0..m).map(|a| (s_idx(a), 1)).collect();
230        let end_sum: Vec<(usize, i64)> = (0..m).map(|a| (e_idx(a), 1)).collect();
231        constraints.push(LinearConstraint::eq(start_sum, 1));
232        constraints.push(LinearConstraint::eq(end_sum, 1));
233
234        let target = ILP::new(num_vars, constraints, Vec::new(), ObjectiveSense::Minimize)
235            .map_err(Self::target_construction)?;
236
237        Ok(ReductionEulerianPathToILP {
238            target,
239            pairs,
240            num_arcs: m,
241        })
242    }
243}
244
245#[cfg(feature = "example-db")]
246pub(crate) fn canonical_rule_example_specs() -> Vec<crate::example_db::specs::RuleExampleSpec> {
247    use crate::topology::DirectedGraph;
248    vec![crate::example_db::specs::RuleExampleSpec {
249        id: "eulerianpath_to_ilp",
250        build: || {
251            // Canonical issue #1025 instance: V = {0,1,2},
252            // A = [(0,1), (0,1), (1,2), (2,0)] (parallel arcs a_0, a_1).
253            // Witness ordering (a_0, a_2, a_3, a_1) traces 0->1->2->0->1.
254            let source =
255                EulerianPath::new(DirectedGraph::new(3, vec![(0, 1), (0, 1), (1, 2), (2, 0)]));
256            crate::example_db::specs::rule_example_via_ilp::<_, i64>(source)
257        },
258    }]
259}
260
261#[cfg(test)]
262#[path = "../unit_tests/rules/eulerianpath_ilp.rs"]
263mod tests;
problemreductions/rules/eulerianpath_ilp.rs

problemreductions/rules/
eulerianpath_ilp.rs
