problemreductions/rules/
longestcommonsubsequence_ilp.rs

1//! Reduction from LongestCommonSubsequence to ILP (Integer Linear Programming).
2//!
3//! The source problem is the optimization version of LCS. The ILP builds a
4//! binary model that maximizes the number of active (non-padding) positions:
5//! - `x_(p,a)` selects symbol `a` at witness position `p` (including padding)
6//! - `y_(r,p,j)` selects the matching position `j` in source string `r`
7//!
8//! The constraints enforce exactly one symbol per position (including the
9//! padding symbol), contiguity of padding, conditional matching for active
10//! positions, and character consistency. The objective maximizes the number of
11//! non-padding positions.
12
13use crate::models::algebraic::{LinearConstraint, ObjectiveSense, ILP};
14use crate::models::misc::LongestCommonSubsequence;
15use crate::reduction;
16use crate::rules::traits::{ReduceTo, ReductionResult};
17
18/// Result of reducing LongestCommonSubsequence to ILP.
19#[derive(Debug, Clone)]
20pub struct ReductionLCSToILP {
21    target: ILP<bool>,
22    alphabet_size: usize,
23    max_length: usize,
24}
25
26impl ReductionResult for ReductionLCSToILP {
27    type Source = LongestCommonSubsequence;
28    type Target = ILP<bool>;
29
30    fn target_problem(&self) -> &ILP<bool> {
31        &self.target
32    }
33
34    fn extract_solution(&self, target_solution: &[usize]) -> Vec<usize> {
35        let num_symbols = self.alphabet_size + 1;
36        let mut witness = Vec::with_capacity(self.max_length);
37        for position in 0..self.max_length {
38            let selected = (0..num_symbols)
39                .find(|&symbol| target_solution.get(position * num_symbols + symbol) == Some(&1))
40                .unwrap_or(self.alphabet_size);
41            witness.push(selected);
42        }
43        witness
44    }
45}
46
47#[reduction(
48    overhead = {
49        num_vars = "max_length * (alphabet_size + 1) + max_length * total_length",
50        num_constraints = "max_length + num_transitions + max_length * num_strings + max_length * total_length + num_transitions * sum_triangular_lengths",
51    }
52)]
53impl ReduceTo<ILP<bool>> for LongestCommonSubsequence {
54    type Result = ReductionLCSToILP;
55
56    fn reduce_to(&self) -> Self::Result {
57        let alphabet_size = self.alphabet_size();
58        let max_length = self.max_length();
59        let strings = self.strings();
60        let total_length = self.total_length();
61        let padding = alphabet_size; // padding symbol index
62        let num_symbols = alphabet_size + 1; // includes padding
63
64        let symbol_var_count = max_length * num_symbols;
65        let mut string_offsets = Vec::with_capacity(strings.len());
66        let mut running_offset = 0usize;
67        for string in strings {
68            string_offsets.push(running_offset);
69            running_offset += string.len();
70        }
71
72        let match_var = |string_index: usize, position: usize, char_index: usize| -> usize {
73            symbol_var_count + position * total_length + string_offsets[string_index] + char_index
74        };
75
76        let mut constraints = Vec::new();
77
78        // (1) Exactly one symbol (including padding) per witness position.
79        for position in 0..max_length {
80            let terms = (0..num_symbols)
81                .map(|symbol| (position * num_symbols + symbol, 1.0))
82                .collect();
83            constraints.push(LinearConstraint::eq(terms, 1.0));
84        }
85
86        // (2) Contiguity: once padding starts, it stays padding.
87        // x_(p+1, padding) >= x_(p, padding)
88        for position in 0..max_length.saturating_sub(1) {
89            constraints.push(LinearConstraint::ge(
90                vec![
91                    (position * num_symbols + padding, -1.0),
92                    ((position + 1) * num_symbols + padding, 1.0),
93                ],
94                0.0,
95            ));
96        }
97
98        // (3) For every string and witness position, the sum of match variables
99        // equals 1 when active and 0 when padding:
100        //   sum_j y_(r,p,j) + x_(p, padding) = 1
101        for (string_index, string) in strings.iter().enumerate() {
102            for position in 0..max_length {
103                let mut terms: Vec<(usize, f64)> = (0..string.len())
104                    .map(|char_index| (match_var(string_index, position, char_index), 1.0))
105                    .collect();
106                terms.push((position * num_symbols + padding, 1.0));
107                constraints.push(LinearConstraint::eq(terms, 1.0));
108            }
109        }
110
111        // (4) A chosen source position can only realize the selected witness symbol.
112        // y_(r, p, j) <= x_(p, string[j])
113        for (string_index, string) in strings.iter().enumerate() {
114            for position in 0..max_length {
115                for (char_index, &symbol) in string.iter().enumerate() {
116                    constraints.push(LinearConstraint::le(
117                        vec![
118                            (match_var(string_index, position, char_index), 1.0),
119                            (position * num_symbols + symbol, -1.0),
120                        ],
121                        0.0,
122                    ));
123                }
124            }
125        }
126
127        // (5) Consecutive active witness positions must map to strictly increasing
128        // source positions.
129        for (string_index, string) in strings.iter().enumerate() {
130            for position in 0..max_length.saturating_sub(1) {
131                for previous in 0..string.len() {
132                    for next in 0..=previous {
133                        constraints.push(LinearConstraint::le(
134                            vec![
135                                (match_var(string_index, position, previous), 1.0),
136                                (match_var(string_index, position + 1, next), 1.0),
137                            ],
138                            1.0,
139                        ));
140                    }
141                }
142            }
143        }
144
145        let num_vars = symbol_var_count + max_length * total_length;
146
147        // Objective: maximize number of non-padding positions.
148        // maximize sum_p sum_{a != padding} x_(p,a)
149        let objective: Vec<(usize, f64)> = (0..max_length)
150            .flat_map(|p| (0..alphabet_size).map(move |a| (p * num_symbols + a, 1.0)))
151            .collect();
152
153        let target = ILP::<bool>::new(num_vars, constraints, objective, ObjectiveSense::Maximize);
154
155        ReductionLCSToILP {
156            target,
157            alphabet_size,
158            max_length,
159        }
160    }
161}
162
163#[cfg(feature = "example-db")]
164pub(crate) fn canonical_rule_example_specs() -> Vec<crate::example_db::specs::RuleExampleSpec> {
165    vec![crate::example_db::specs::RuleExampleSpec {
166        id: "longestcommonsubsequence_to_ilp",
167        build: || {
168            // Source: alphabet {0,1,2}, strings [0,1,2] and [1,0,2], max_length = 3
169            // Optimal LCS: [0,2] (length 2) or [1,2] (length 2)
170            // Config with padding: e.g. [0, 2, 3] (symbol 3 = padding)
171            let source = LongestCommonSubsequence::new(3, vec![vec![0, 1, 2], vec![1, 0, 2]]);
172            crate::example_db::specs::rule_example_via_ilp::<_, bool>(source)
173        },
174    }]
175}
176
177#[cfg(test)]
178#[path = "../unit_tests/rules/longestcommonsubsequence_ilp.rs"]
179mod tests;
problemreductions/rules/longestcommonsubsequence_ilp.rs

problemreductions/rules/
longestcommonsubsequence_ilp.rs
