From 5bfa5c83fe706e13564c8430ae126223629ca85e Mon Sep 17 00:00:00 2001 From: Xiwei Pan <90967972+isPANN@users.noreply.github.com> Date: Wed, 4 Nov 2026 04:43:18 -0800 Subject: [PATCH] Fix #1198: [Rule] IntegerExpressionMembership to binary ILP Co-Authored-By: Claude Sonnet 5.5 --- docs/paper/reductions.typ | 189 +++++++ docs/paper/references.bib | 12 + src/rules/integerexpressionmembership_ilp.rs | 266 ++++++++++ src/rules/mod.rs | 2 + .../rules/integerexpressionmembership_ilp.rs | 480 ++++++++++++++++++ 5 files changed, 949 insertions(+) create mode 100644 src/rules/integerexpressionmembership_ilp.rs create mode 100644 src/unit_tests/rules/integerexpressionmembership_ilp.rs diff --git a/docs/paper/reductions.typ b/docs/paper/reductions.typ index 8e0462b7a..4e43b4125 100644 --- a/docs/paper/reductions.typ +++ b/docs/paper/reductions.typ @@ -19492,6 +19492,195 @@ The following table shows concrete target-variable counts for example instances, _Solution extraction._ The Subset Sum configuration is the Integer Expression Membership configuration: $x_i = 1$ (right branch) means element $i$ is selected. ] +// IntegerExpressionMembership → ILP +#{ + let ex = load-example("IntegerExpressionMembership", "ILP") + let sol = ex.solutions.at(0) + let tree = ex.source.instance.expression + let K = ex.source.instance.target + let ilp = ex.target.instance + let x = sol.target_config + let flatten(e, start: 0) = { + let kind = e.keys().first() + if kind == "Atom" { ((kind: kind, index: start, value: e.Atom),) } + else { + let children = e.at(kind) + let left = flatten(children.at(0), start: start + 1) + let right = flatten(children.at(1), start: start + 1 + left.len()) + ((kind: kind, index: start, left: start + 1, right: start + 1 + left.len()),) + left + right + } + } + let show-expr(e) = { + let kind = e.keys().first() + if kind == "Atom" { str(e.Atom) } + else { + let children = e.at(kind) + "(" + show-expr(children.at(0)) + (if kind == "Sum" { " + " } else { " ∪ " }) + show-expr(children.at(1)) + ")" + } + } + let reachable(e) = { + let kind = e.keys().first() + if kind == "Atom" { (e.Atom,) } + else { + let children = e.at(kind) + let left = reachable(children.at(0)) + let right = reachable(children.at(1)) + if kind == "Union" { (left + right).dedup().sorted() } + else { left.map(l => right.map(r => l + r)).flatten().dedup().sorted() } + } + } + let bitlen(n) = if n < 2 { 1 } else { 1 + bitlen(calc.floor(n / 2)) } + let popcount(n) = if n == 0 { 0 } else { calc.rem(n, 2) + popcount(calc.floor(n / 2)) } + let nodes = flatten(tree) + let atoms = nodes.filter(n => n.kind == "Atom") + let A = atoms.len() + let S = nodes.filter(n => n.kind == "Sum").len() + let U = nodes.filter(n => n.kind == "Union").len() + let T = nodes.len() + let H = bitlen(calc.max(K, ..atoms.map(n => n.value))) + let B = atoms.map(n => popcount(n.value)).sum() + let decode(index) = { + let node = nodes.at(index) + if node.kind == "Atom" { () } + else if node.kind == "Sum" { decode(node.left) + decode(node.right) } + else { + let right = x.at(node.right) == 1 + (right,) + decode(if right { node.right } else { node.left }) + } + } + let choices = decode(0) + let choices = choices + range(U - choices.len()).map(_ => false) + let values = reachable(tree) + let no-K = values.first() + 1 + let no-H = bitlen(calc.max(no-K, ..atoms.map(n => n.value))) + let no-n = T + A * (2 * no-H - 1) + let no-m = 1 + 2 * S + U + A * no-H + no-H + let no-z = A * (4 * no-H - 2) + B + 1 + 4 * S + 3 * U + assert(T == 2 * A - 1 and S + U == A - 1) + assert(ilp.variables.len() == T + A * (2 * H - 1)) + assert(ilp.constraints.len() == 1 + 2 * S + U + A * H + H) + assert(ilp.constraints.map(row => row.terms.len()).sum() == A * (4 * H - 2) + B + 1 + 4 * S + 3 * U) + assert(ilp.variables.all(v => v.lower_bound == 0 and v.upper_bound == 1)) + assert(ilp.objective.len() == 0 and ilp.sense == "Minimize") + assert(x.len() == ilp.variables.len() and x.all(v => v == 0 or v == 1)) + assert(ilp.constraints.all(row => row.comparison == "Eq" and row.terms.map(term => x.at(term.at(0)) * term.at(1)).sum() == row.rhs)) + assert(atoms.map(n => n.value * x.at(n.index)).sum() == K) + assert(choices == sol.source_config and K in values) + assert(no-K not in values) + [ + #reduction-rule("IntegerExpressionMembership", "ILP", + example: true, + example-caption: [#A atoms, #U unions, target #K], + extra: [ + #pred-commands( + "pred create --example " + rule-spec(ex) + " -o iem.json", + "pred path IntegerExpressionMembership ILP/bool/i64 -o paths.json", + "pred reduce iem.json --via route.json -o bundle.json", + "pred solve bundle.json", + "pred evaluate iem.json --config " + cli-config(sol.source_config), + ) + Select the direct entry from `paths.json` as `route.json` before reducing. + + *YES example.* The fixture expression is #show-expr(tree), with target + $K = #K$ and choices $(#fmt-values(sol.source_config))$. + Here $A = #A$, $S = #S$, $U = #U$, $T = #T$, $H = #H$, and $B = #B$. + The target has $#T + #A dot (2 dot #H - 1) = #ilp.variables.len()$ binary + variables, $#(1 + 2 * S + U + A * H + H)$ equality rows, and + $#ilp.constraints.map(row => row.terms.len()).sum()$ nonzeros. + Its activations are $(#fmt-values(x.slice(0,T)))$. + Each stage below is read from the canonical target witness; bits are + ordered low to high. + #{ + let cells = () + for (j, atom) in atoms.enumerate() { + let prefix = x.slice(T + j * H, T + (j + 1) * H) + let carry = x.slice(T + A * H + j * (H - 1), T + A * H + (j + 1) * (H - 1)) + let total = prefix.enumerate().map(pair => pair.at(1) * calc.pow(2, pair.at(0))).sum() + cells += ([#atom.value], [#x.at(atom.index)], [#total], [(#fmt-values(prefix))], [(#fmt-values(carry))]) + } + table(columns: 5, inset: 3pt, [Atom], [Active], [Prefix total], [Prefix bits], [Carry bits], ..cells) + } + The active atoms sum to $#atoms.map(n => n.value * x.at(n.index)).sum() = #K$. + Every emitted row is checked above against the loaded witness, the + objective is zero, and active DFS extraction returns + $(#fmt-values(choices))$. + + *NO example.* Keep the fixture tree and replace its query by + $#no-K = #values.first() + 1$. Its reachable values, computed recursively + by union and Minkowski sum, are exactly $\{#fmt-values(values)\}$; + $#no-K$ is absent. The new width is $#no-H$, giving $#no-n$ binary + variables, $#no-m$ equality rows, and $#no-z$ nonzeros by the construction + formulas. Any feasible target would have final prefix in that reachable + set, while its final pins demand $#no-K$. Thus both sides are NO. + ], + )[ + A tree-activation encoding followed by binary ripple-carry additions, + using binary expansion as discussed by @Glover1975ImprovedILP, + reduces positive signed-integer expression membership to binary ILP in + $O(T + A H)$ time and storage, with coefficients of magnitude at most two. + #footnote[Direct formulation derived here; no expansion of the represented integer set is required.] + ][ + _Construction._ Let $A$, $S$, and $U$ count atoms, sums, and unions; + $T = 2A-1$ is the node count. Number nodes in preorder and atoms + $v_1, dots, v_A$ left to right. Let $H$ be the maximum bit length of the + positive atoms and query $K$, so $1 <= H <= 63$ for the signed-i64 domain. + Let $v_(j,b)$ and $k_b$ be their bits, $0 <= b < H$. + + Allocate binary activations $a_t$, pin the root to one, and at a sum + with children $l,r$ impose $a_l = a_t$ and $a_r = a_t$. At a union impose + $a_l + a_r = a_t$. For each atom $j$, allocate prefix bits $p_(j,b)$ + and internal carry bits $c_(j,b)$, $1 <= b < H$. Use constant boundaries + $p_(0,b) = c_(j,0) = c_(j,H) = 0$. Emit + $ p_(j-1,b) + v_(j,b) a_(t(j)) + c_(j,b) - p_(j,b) - 2c_(j,b+1) = 0, $ + where $t(j)$ is atom $j$'s node. Pin $p_(A,b)=k_b$ and minimize zero. + Coordinates are all activations, atom-major low-bit-first prefixes, then + atom-major low-bit-first internal carries. Rows are root, preorder + structural rows, addition rows, then final pins. The target constructor + sorts sparse terms; omitted zero terms never merge with other indices. + + _Arithmetic and size._ Multiplying the stage-$j$ rows by $2^b$ and + summing cancels internal carries and gives the exact identity + $P_j = P_(j-1) + v_j a_(t(j))$, where $P_j = sum_b 2^b p_(j,b)$. + Conversely binary addition below $2^H$ gives unique result/carry bits. + No terminal carry is allocated, so overflow is forbidden. At $H=1$ + there are no internal carries. With $B = sum_j "popcount"(v_j)$ the + counts are $T+A(2H-1)$ variables, $1+2S+U+A H+H$ rows, and + $A(4H-2)+B+1+4S+3U$ nonzeros. The latter counts root (1), sum (4S), + union (3U), previous prefixes ($H(A-1)$), atom bits ($B$), carry uses + ($2A(H-1)$), current prefixes ($A H$), and final pins ($H$). + All indices within a row are distinct. The source registers only $U$; + sum-only trees with $U=0$ can be arbitrarily large, so those three size + fields are unavailable as functions of registered parameters. Magnitude + bits is at most 2, including RHS and binary bounds. Each row has total + positive and negative magnitudes at most three, so binary evaluation + cannot overflow i64. Layout arithmetic is checked before allocation; + bit decomposition never forms a signed $2^63$. + + _Correctness._ ($arrow.r.double$) A satisfying source witness activates + exactly its visited tree. Sum and union rows hold, and its selected + positive atoms total $K$. Every left-to-right prefix is at most $K$, + hence fits $H$ bits; the arithmetic identity supplies a feasible extension. + + ($arrow.l.double$) A feasible target activates the root, both children + of every active sum and exactly one child of every active union. + Inactive nodes force all descendants inactive. Summing the exact stage + identities and using final pins gives $sum_j v_j a_(t(j))=K$. + Traversing this active tree therefore evaluates to $K$. Every intermediate + source sum is bounded by this positive total, so checked source addition + does not overflow. This includes a single atom, sum-only trees and + inactive alternatives whose positive totals exceed i64::MAX. + + _Solution extraction._ Validate the target witness once before decoding. + Visit sums left then right; at an active union append false for its left + child or true for its right child and visit only that child. Skipped + subtrees consume no union bits. Append false at the end to length $U$, + the model's unused-tail sentinel. The converse proves the resulting + source witness satisfies the query; unused-tail multiplicity is not + preserved. + ] + ] +} + // Removed: ExactCoverBy3Sets → MinimumWeightSolutionToLinearEquations (unsound reduction, #1006) // 11. KSatisfiability → SimultaneousIncongruences (#554) diff --git a/docs/paper/references.bib b/docs/paper/references.bib index 0dda81fd2..29a209390 100644 --- a/docs/paper/references.bib +++ b/docs/paper/references.bib @@ -1,3 +1,15 @@ +@article{Glover1975ImprovedILP, + author = {Glover, Fred}, + title = {Improved Linear Integer Programming Formulations of Nonlinear Integer Problems}, + journal = {Management Science}, + volume = {22}, + number = {4}, + pages = {455--460}, + year = {1975}, + doi = {10.1287/mnsc.22.4.455}, + url = {https://leeds-faculty.colorado.edu/glover/fred%20pubs/69%20-%20Improved%20Linear%20Integer%20Programming%20Formulations%20of%20Nonlinear%20Integer%20Problems%20_%20%20Fred%20Glover.pdf} +} + @misc{ozkayaCatalyurek2022, author = {M. Yusuf {\"O}zkaya and {\"U}mit V. {{\c C}}ataly{\"u}rek}, title = {A Simple and Elegant Mathematical Formulation for the Acyclic {DAG} Partitioning Problem}, diff --git a/src/rules/integerexpressionmembership_ilp.rs b/src/rules/integerexpressionmembership_ilp.rs new file mode 100644 index 000000000..cb53e30bb --- /dev/null +++ b/src/rules/integerexpressionmembership_ilp.rs @@ -0,0 +1,266 @@ +//! Binary activation and ripple-carry encoding of integer expression membership. + +use crate::models::algebraic::{ILP, LinearConstraint, ObjectiveSense}; +use crate::models::misc::{IntExpr, IntegerExpressionMembership}; +use crate::rules::ReductionError; +use crate::rules::traits::{ReduceTo, ReductionResult}; + +#[derive(Debug, Clone)] +enum Node { + Atom(i64), + Sum(usize, usize), + Union(usize, usize), +} + +fn flatten(expression: &IntExpr, nodes: &mut Vec) -> usize { + let index = nodes.len(); + nodes.push(Node::Atom(0)); + nodes[index] = match expression { + IntExpr::Atom(value) => Node::Atom(*value), + IntExpr::Sum(left, right) => Node::Sum(flatten(left, nodes), flatten(right, nodes)), + IntExpr::Union(left, right) => Node::Union(flatten(left, nodes), flatten(right, nodes)), + }; + index +} + +fn overflow(operation: &str) -> ReductionError { + ReductionError::integer_overflow::>(operation) +} + +fn bit_width(value: i64) -> Result { + if value <= 0 { + return Err(ReductionError::invalid_target::< + IntegerExpressionMembership, + ILP, + >("binary encoding requires a positive integer")); + } + let width = i64::BITS + .checked_sub(value.leading_zeros()) + .ok_or_else(|| overflow("computing the bit width"))?; + usize::try_from(width).map_err(|_| overflow("converting the bit width")) +} + +fn bit_at(value: i64, bit: usize) -> Result { + let shift = u32::try_from(bit).map_err(|_| overflow("converting a bit index"))?; + value + .checked_shr(shift) + .map(|value| value & 1) + .ok_or_else(|| overflow("reading an integer bit")) +} + +fn checked_layout( + atoms: usize, + sums: usize, + unions: usize, + width: usize, +) -> Result<(usize, usize, usize), ReductionError> { + if atoms == 0 || !(1..=63).contains(&width) { + return Err(ReductionError::invalid_target::< + IntegerExpressionMembership, + ILP, + >( + "the expression must have atoms and a bit width from 1 to 63", + )); + } + let nodes = atoms + .checked_add(atoms - 1) + .ok_or_else(|| overflow("counting expression nodes"))?; + let prefixes = atoms + .checked_mul(width) + .ok_or_else(|| overflow("counting prefix bits"))?; + let carries = atoms + .checked_mul(width - 1) + .ok_or_else(|| overflow("counting carry bits"))?; + let variables = nodes + .checked_add(prefixes) + .and_then(|n| n.checked_add(carries)) + .ok_or_else(|| overflow("counting binary ILP variables"))?; + let rows = sums + .checked_mul(2) + .and_then(|n| n.checked_add(1)) + .and_then(|n| n.checked_add(unions)) + .and_then(|n| n.checked_add(prefixes)) + .and_then(|n| n.checked_add(width)) + .ok_or_else(|| overflow("counting binary ILP constraints"))?; + Ok((nodes, variables, rows)) +} + +#[derive(Debug, Clone)] +pub struct ReductionIntegerExpressionMembershipToILP { + target: ILP, + nodes: Vec, + unions: usize, +} + +impl ReductionResult for ReductionIntegerExpressionMembershipToILP { + type Source = IntegerExpressionMembership; + type Target = ILP; + + fn target_problem(&self) -> &ILP { + &self.target + } + + fn extract_solution(&self, solution: &Vec) -> crate::rules::ExtractionResult> { + crate::rules::traits::validate_target_witness( + &self.target, + solution, + |value| value.value.is_some(), + "target ILP assignment is infeasible", + )?; + fn visit(nodes: &[Node], index: usize, solution: &[i64], choices: &mut Vec) { + match nodes[index] { + Node::Atom(_) => {} + Node::Sum(left, right) => { + visit(nodes, left, solution, choices); + visit(nodes, right, solution, choices); + } + Node::Union(left, right) => { + let choose_right = solution[right] == 1; + choices.push(choose_right); + visit( + nodes, + if choose_right { right } else { left }, + solution, + choices, + ); + } + } + } + let mut choices = Vec::with_capacity(self.unions); + visit(&self.nodes, 0, solution, &mut choices); + // Unvisited unions consume no bits; the model ignores this trailing sentinel. + choices.resize(self.unions, false); + Ok(choices) + } +} + +#[crate::aggregate_reduction(ilp_feasibility)] +impl crate::rules::AggregateReductionResult for ReductionIntegerExpressionMembershipToILP {} + +// T + A(2H-1) variables, 1+2S+U+AH+H rows and A(4H-2)+B+1+4S+3U nonzeros. +// Coefficients are -2,-1,1, and RHS/bounds are 0 or 1, including H=1 without carries. +#[crate::reduction(transform = { + upper_bound { max_constraint_magnitude_bits = "2", }, + unavailable { + num_vars = "Atom count and bit width are not bounded by num_union_nodes", + num_constraints = "Atom count and bit width are not bounded by num_union_nodes", + num_nonzeros = "Atom count and bit width are not bounded by num_union_nodes", + }, +})] +impl ReduceTo> for IntegerExpressionMembership { + type Result = ReductionIntegerExpressionMembershipToILP; + + fn reduce_to(&self) -> Result { + let mut nodes = Vec::new(); + flatten(self.expression(), &mut nodes); + let atoms: Vec<_> = nodes + .iter() + .enumerate() + .filter_map(|(index, node)| match node { + Node::Atom(value) => Some((index, *value)), + _ => None, + }) + .collect(); + let sums = nodes + .iter() + .filter(|node| matches!(node, Node::Sum(..))) + .count(); + let unions = nodes + .iter() + .filter(|node| matches!(node, Node::Union(..))) + .count(); + let width = atoms + .iter() + .try_fold(bit_width(self.target())?, |width, (_, value)| { + Ok::<_, ReductionError>(width.max(bit_width(*value)?)) + })?; + let (node_count, variables, row_count) = checked_layout(atoms.len(), sums, unions, width)?; + // The checked total bounds every included prefix/carry offset below. + let prefix = |j: usize, b: usize| node_count + j * width + b; + let carry = |j: usize, b: usize| node_count + atoms.len() * width + j * (width - 1) + b - 1; + let mut rows = Vec::with_capacity(row_count); + rows.push(LinearConstraint::eq(vec![(0, 1)], 1)); + for (index, node) in nodes.iter().enumerate() { + match *node { + Node::Atom(_) => {} + Node::Sum(left, right) => { + rows.push(LinearConstraint::eq(vec![(left, 1), (index, -1)], 0)); + rows.push(LinearConstraint::eq(vec![(right, 1), (index, -1)], 0)); + } + Node::Union(left, right) => { + rows.push(LinearConstraint::eq( + vec![(left, 1), (right, 1), (index, -1)], + 0, + )); + } + } + } + for (j, &(leaf, value)) in atoms.iter().enumerate() { + for b in 0..width { + let mut terms = vec![(prefix(j, b), -1)]; + if j > 0 { + terms.push((prefix(j - 1, b), 1)); + } + if bit_at(value, b)? == 1 { + terms.push((leaf, 1)); + } + if b > 0 { + terms.push((carry(j, b), 1)); + } + if b + 1 < width { + terms.push((carry(j, b + 1), -2)); + } + rows.push(LinearConstraint::eq(terms, 0)); + } + } + for b in 0..width { + rows.push(LinearConstraint::eq( + vec![(prefix(atoms.len() - 1, b), 1)], + bit_at(self.target(), b)?, + )); + } + let target = ILP::new(variables, rows, vec![], ObjectiveSense::Minimize) + .map_err(Self::target_construction)?; + Ok(Self::Result { + target, + nodes, + unions, + }) + } +} + +#[cfg(feature = "example-db")] +pub(crate) fn canonical_rule_example_specs() -> Vec { + use crate::export::SolutionPair; + vec![crate::example_db::specs::RuleExampleSpec { + id: "integerexpressionmembership_to_ilp", + build: || { + let union = |left, right| { + IntExpr::Union( + Box::new(IntExpr::Atom(left)), + Box::new(IntExpr::Atom(right)), + ) + }; + let expression = IntExpr::Sum( + Box::new(IntExpr::Sum(Box::new(union(1, 4)), Box::new(union(3, 6)))), + Box::new(union(2, 5)), + ); + let mut witness = vec![1, 1, 1, 0, 1, 1, 0, 1, 1, 1, 0]; + witness.extend([ + 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 1, 0, 0, 1, 0, 1, 0, 0, 1, 1, 0, 0, 1, 1, + ]); + witness.extend([0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 1, 0, 0, 0, 0]); + crate::example_db::specs::rule_example_with_witness::<_, ILP>( + IntegerExpressionMembership::new(expression, 12), + SolutionPair { + source_config: serde_json::json!([true, true, false]), + target_config: serde_json::json!(witness), + }, + ) + }, + }] +} + +#[cfg(test)] +#[path = "../unit_tests/rules/integerexpressionmembership_ilp.rs"] +mod tests; diff --git a/src/rules/mod.rs b/src/rules/mod.rs index 7606f4923..d1e4090c8 100644 --- a/src/rules/mod.rs +++ b/src/rules/mod.rs @@ -223,6 +223,7 @@ mod ilp_bool_ilp_i64; mod ilp_bounded_ilp; pub(crate) mod ilp_helpers; pub(crate) mod ilp_qubo; +pub(crate) mod integerexpressionmembership_ilp; pub(crate) mod integralflowbundles_ilp; pub(crate) mod integralflowhomologousarcs_ilp; pub(crate) mod integralflowwithmultipliers_ilp; @@ -584,6 +585,7 @@ pub(crate) fn canonical_rule_example_specs() -> Vec IntExpr { + IntExpr::Union(Box::new(left), Box::new(right)) +} +fn sum(left: IntExpr, right: IntExpr) -> IntExpr { + IntExpr::Sum(Box::new(left), Box::new(right)) +} +fn example() -> IntExpr { + sum( + sum( + union(IntExpr::Atom(1), IntExpr::Atom(4)), + union(IntExpr::Atom(3), IntExpr::Atom(6)), + ), + union(IntExpr::Atom(2), IntExpr::Atom(5)), + ) +} +fn reduce( + expression: IntExpr, + k: i64, +) -> ( + IntegerExpressionMembership, + ReductionIntegerExpressionMembershipToILP, +) { + let source = IntegerExpressionMembership::new(expression, k); + let reduction = ReduceTo::>::reduce_to(&source).unwrap(); + (source, reduction) +} + +#[test] +fn test_integerexpressionmembership_to_ilp_closed_loop() { + let (source, reduction) = reduce(union(IntExpr::Atom(1), IntExpr::Atom(2)), 2); + assert_bf_vs_ilp(&source, &reduction); + let target = ReductionResult::target_problem(&reduction); + let solutions: Vec<_> = (0..1usize << target.num_vars()) + .map(|mask| { + (0..target.num_vars()) + .map(|b| i64::from((mask >> b) & 1 == 1)) + .collect::>() + }) + .filter(|witness| target.is_feasible(witness).unwrap()) + .collect(); + assert_eq!(solutions.len(), 1); + for solution in solutions { + assert_eq!(reduction.extract_solution(&solution).unwrap(), vec![true]); + assert_eq!( + reduction.extract_value(target.evaluate(&solution).unwrap()), + Or(true) + ); + } +} + +#[test] +fn test_canonical_witness_and_no_instance() { + let (source, reduction) = reduce(example(), 12); + let target = ReductionResult::target_problem(&reduction); + let mut witness = vec![1, 1, 1, 0, 1, 1, 0, 1, 1, 1, 0]; + witness.extend([ + 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 1, 0, 0, 1, 0, 1, 0, 0, 1, 1, 0, 0, 1, 1, + ]); + witness.extend([0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 1, 0, 0, 0, 0]); + assert_eq!( + ( + target.num_vars(), + target.num_constraints(), + target.num_nonzeros() + ), + (53, 36, 111) + ); + assert_eq!( + target.evaluate(&witness).unwrap(), + Extremum::minimize(Some(0)) + ); + assert_eq!( + reduction.extract_solution(&witness).unwrap(), + vec![true, true, false] + ); + assert_eq!(source.evaluate(&vec![true, true, false]).unwrap(), Or(true)); + assert_parameter_predictions(&source, &reduction); + assert_bf_vs_ilp(&source, &reduction); + let (no, reduction) = reduce(example(), 7); + let target = ReductionResult::target_problem(&reduction); + assert_eq!( + ( + target.num_vars(), + target.num_constraints(), + target.num_nonzeros() + ), + (41, 29, 87) + ); + assert!(BruteForce::new().solve(&no).unwrap().is_none()); + assert!(ILPSolver::new().solve(target).is_err()); + assert_eq!(reduction.extract_value(Extremum::minimize(None)), Or(false)); +} + +#[test] +fn test_visited_union_order_and_unused_tail() { + let (source, reduction) = reduce( + sum( + union(union(IntExpr::Atom(1), IntExpr::Atom(2)), IntExpr::Atom(3)), + union(IntExpr::Atom(4), IntExpr::Atom(5)), + ), + 8, + ); + let witness = ILPSolver::new() + .solve(ReductionResult::target_problem(&reduction)) + .unwrap(); + let config = reduction.extract_solution(&witness).unwrap(); + assert_eq!(config, vec![true, true, false]); + assert_eq!(source.evaluate(&config).unwrap(), Or(true)); + assert_eq!(source.evaluate(&vec![true, true, true]).unwrap(), Or(true)); +} + +#[test] +fn test_rejects_malformed_and_infeasible_witnesses() { + let (_, reduction) = reduce(example(), 12); + let target = ReductionResult::target_problem(&reduction); + let good = ILPSolver::new().solve(target).unwrap(); + assert!( + reduction + .extract_solution(&good[..good.len() - 1].to_vec()) + .is_err() + ); + let mut extra = good.clone(); + extra.push(0); + assert!(reduction.extract_solution(&extra).is_err()); + let mut nonbinary = good.clone(); + nonbinary[0] = 2; + assert!(reduction.extract_solution(&nonbinary).is_err()); + for index in 0..good.len() { + let mut bad = good.clone(); + bad[index] ^= 1; + assert!(!target.is_feasible(&bad).unwrap()); + assert!(reduction.extract_solution(&bad).is_err()); + } +} + +fn lift(source: &IntegerExpressionMembership, config: &[bool]) -> Vec { + fn visit( + e: &IntExpr, + active: bool, + config: &mut std::slice::Iter<'_, bool>, + activations: &mut Vec, + leaves: &mut Vec<(i64, bool)>, + ) { + activations.push(i64::from(active)); + match e { + IntExpr::Atom(value) => leaves.push((*value, active)), + IntExpr::Sum(left, right) => { + visit(left, active, config, activations, leaves); + visit(right, active, config, activations, leaves); + } + IntExpr::Union(left, right) => { + let choice = active && *config.next().unwrap(); + visit(left, active && !choice, config, activations, leaves); + visit(right, active && choice, config, activations, leaves); + } + } + } + let mut witness = Vec::new(); + let mut leaves = Vec::new(); + visit( + source.expression(), + true, + &mut config.iter(), + &mut witness, + &mut leaves, + ); + let largest = leaves + .iter() + .map(|&(v, _)| v) + .chain([source.target()]) + .max() + .unwrap(); + let width = usize::try_from(i64::BITS - largest.leading_zeros()).unwrap(); + let mut previous = 0u128; + let mut carries = Vec::new(); + for (value, active) in leaves { + let amount = if active { + u128::try_from(value).unwrap() + } else { + 0 + }; + let total = previous + amount; + witness.extend((0..width).map(|b| i64::from((total >> b) & 1 == 1))); + for b in 1..width { + let mask = (1u128 << b) - 1; + carries.push(i64::try_from(((previous & mask) + (amount & mask)) >> b).unwrap()); + } + previous = total & ((1u128 << width) - 1); + } + witness.extend(carries); + witness +} + +fn audit(expression: IntExpr, k: i64, solve: bool) { + use std::collections::HashSet; + let (source, reduction) = reduce(expression, k); + let target = ReductionResult::target_problem(&reduction); + let mut expected = HashSet::new(); + for mask in 0..1usize << source.num_union_nodes() { + let config: Vec<_> = (0..source.num_union_nodes()) + .map(|b| (mask >> b) & 1 == 1) + .collect(); + let good = source.evaluate(&config).unwrap().0; + let witness = lift(&source, &config); + assert_eq!(target.is_feasible(&witness).unwrap(), good); + if good { + assert_eq!( + source + .evaluate(&reduction.extract_solution(&witness).unwrap()) + .unwrap(), + Or(true) + ); + expected.insert(witness); + } + } + if target.num_vars() <= 12 { + let mut actual = HashSet::new(); + for mask in 0..1usize << target.num_vars() { + let witness: Vec = (0..target.num_vars()) + .map(|b| i64::from((mask >> b) & 1 == 1)) + .collect(); + if target.is_feasible(&witness).unwrap() { + let config = reduction.extract_solution(&witness).unwrap(); + assert_eq!(source.evaluate(&config).unwrap(), Or(true)); + actual.insert(witness); + } else { + assert!(reduction.extract_solution(&witness).is_err()); + } + } + assert_eq!(actual, expected); + } + let atoms = source.num_atoms(); + let nodes = source.expression_size(); + let unions = source.num_union_nodes(); + let sums = atoms - 1 - unions; + let mut flat = Vec::new(); + flatten(source.expression(), &mut flat); + let values: Vec<_> = flat + .iter() + .filter_map(|n| { + if let Node::Atom(v) = n { + Some(*v) + } else { + None + } + }) + .collect(); + let width = bit_width(*values.iter().chain([&k]).max().unwrap()).unwrap(); + let popcount: usize = values + .iter() + .map(|v| usize::try_from(v.count_ones()).unwrap()) + .sum(); + assert_eq!(target.num_vars(), nodes + atoms * (2 * width - 1)); + assert_eq!( + target.num_constraints(), + 1 + 2 * sums + unions + atoms * width + width + ); + assert_eq!( + target.num_nonzeros(), + atoms * (4 * width - 2) + popcount + 1 + 4 * sums + 3 * unions + ); + assert_eq!( + target.max_constraint_magnitude_bits(), + if width == 1 { 1 } else { 2 } + ); + assert_parameter_predictions(&source, &reduction); + if solve { + match ILPSolver::new().solve(target) { + Ok(witness) => { + assert!(!expected.is_empty()); + assert!(expected.contains(&witness)); + assert_eq!( + source + .evaluate(&reduction.extract_solution(&witness).unwrap()) + .unwrap(), + Or(true) + ); + assert_eq!( + reduction.extract_value(target.evaluate(&witness).unwrap()), + Or(true) + ); + } + Err(_) => { + assert!(expected.is_empty()); + assert_eq!(reduction.extract_value(Extremum::minimize(None)), Or(false)); + } + } + } +} + +fn trees(nodes: usize) -> Vec { + if nodes == 1 { + return (1..=3).map(IntExpr::Atom).collect(); + } + let mut result = Vec::new(); + for left_size in (1..nodes - 1).step_by(2) { + for left in trees(left_size) { + for right in trees(nodes - 1 - left_size) { + result.push(union(left.clone(), right.clone())); + result.push(sum(left.clone(), right)); + } + } + } + result +} + +#[test] +fn test_exhaustive_small_expressions() { + let mut count = 0; + for size in [1, 3, 5] { + for expression in trees(size) { + for k in 1..=10 { + audit(expression.clone(), k, count % 211 == 0); + count += 1; + } + } + } + assert_eq!(count, 2370); +} + +fn sampled_cell(unions: usize, sums: usize) { + // A fixed generator makes every (U,S) cell reproducible without global RNG state. + use rand::prelude::*; + let mut rng = rand::rngs::StdRng::seed_from_u64(2026110401); + fn tree(u: usize, s: usize, rng: &mut rand::rngs::StdRng) -> IntExpr { + if u + s == 0 { + return IntExpr::Atom(rng.random_range(1..=31)); + } + let is_union = rng.random_range(0..u + s) < u; + let (u, s) = (u - usize::from(is_union), s - usize::from(!is_union)); + let (lu, ls) = (rng.random_range(0..=u), rng.random_range(0..=s)); + let (left, right) = (tree(lu, ls, rng), tree(u - lu, s - ls, rng)); + if is_union { + union(left, right) + } else { + sum(left, right) + } + } + for index in 0..300 { + let expression = tree(unions, sums, &mut rng); + assert_eq!(expression.count_union_nodes(), unions); + assert_eq!(expression.count_atoms(), unions + sums + 1); + audit(expression, rng.random_range(1..=100), index < 3); + } +} + +macro_rules! sampled_cells { + ($(($name:ident, $u:expr, $s:expr)),* $(,)?) => {$( + #[test] fn $name() { sampled_cell($u, $s); } + )*}; +} +sampled_cells!( + (test_samples_unions_zero_sums_zero, 0, 0), + (test_samples_unions_zero_sums_one, 0, 1), + (test_samples_unions_zero_sums_two, 0, 2), + (test_samples_unions_zero_sums_three, 0, 3), + (test_samples_unions_zero_sums_four, 0, 4), + (test_samples_unions_one_sums_zero, 1, 0), + (test_samples_unions_one_sums_one, 1, 1), + (test_samples_unions_one_sums_two, 1, 2), + (test_samples_unions_one_sums_three, 1, 3), + (test_samples_unions_one_sums_four, 1, 4), + (test_samples_unions_two_sums_zero, 2, 0), + (test_samples_unions_two_sums_one, 2, 1), + (test_samples_unions_two_sums_two, 2, 2), + (test_samples_unions_two_sums_three, 2, 3), + (test_samples_unions_two_sums_four, 2, 4), + (test_samples_unions_three_sums_zero, 3, 0), + (test_samples_unions_three_sums_one, 3, 1), + (test_samples_unions_three_sums_two, 3, 2), + (test_samples_unions_three_sums_three, 3, 3), + (test_samples_unions_three_sums_four, 3, 4), + (test_samples_unions_four_sums_zero, 4, 0), + (test_samples_unions_four_sums_one, 4, 1), + (test_samples_unions_four_sums_two, 4, 2), + (test_samples_unions_four_sums_three, 4, 3), + (test_samples_unions_four_sums_four, 4, 4), +); + +#[test] +fn test_single_atoms_sums_widths_and_overflow() { + for width in 1..=63 { + let power = 1i64 << (width - 1); + audit(IntExpr::Atom(power), power, true); + if power > 1 { + audit(sum(IntExpr::Atom(power - 1), IntExpr::Atom(1)), power, true); + } + } + for (expression, target) in [ + (sum(IntExpr::Atom(1), IntExpr::Atom(1)), 1), + (union(IntExpr::Atom(1), IntExpr::Atom(2)), 8), + (IntExpr::Atom(i64::MAX), i64::MAX), + (sum(IntExpr::Atom(i64::MAX - 1), IntExpr::Atom(1)), i64::MAX), + (sum(IntExpr::Atom(i64::MAX), IntExpr::Atom(1)), i64::MAX), + ( + union( + sum(IntExpr::Atom(i64::MAX), IntExpr::Atom(1)), + IntExpr::Atom(3), + ), + 3, + ), + ] { + audit(expression, target, true); + } + let (_, reduction) = reduce( + sum(IntExpr::Atom((1i64 << 62) - 1), IntExpr::Atom(1)), + 1i64 << 62, + ); + let witness = ILPSolver::new() + .solve(ReductionResult::target_problem(&reduction)) + .unwrap(); + assert_eq!(&witness[3 + 2 * 63..], &[vec![0; 62], vec![1; 62]].concat()); +} + +#[test] +fn test_checked_layout_and_numeric_boundaries() { + for (value, expected) in [(1, 1), (2, 2), (i64::MAX, 63)] { + assert_eq!(bit_width(value).unwrap(), expected); + } + assert!(bit_width(0).is_err()); + assert!(bit_width(-1).is_err()); + assert_eq!(bit_at(i64::MAX, 62).unwrap(), 1); + assert_eq!(bit_at(i64::MAX, 63).unwrap(), 0); + assert!(bit_at(i64::MAX, 64).is_err()); + assert!(bit_at(1, usize::MAX).is_err()); + for width in [0, 64] { + assert!(checked_layout(1, 0, 0, width).is_err()); + } + assert!(checked_layout(0, 0, 0, 1).is_err()); + for h in [1usize, 63] { + for only_sums in [false, true] { + let limit = usize::MAX; + // Independent u128 formulas avoid overflowing the reference calculation. + let variables = + (u128::try_from(limit).unwrap() + 1) / u128::try_from(2 * h + 1).unwrap(); + let rows = if only_sums { + (u128::try_from(limit - h + 1).unwrap()) / u128::try_from(h + 2).unwrap() + } else { + u128::try_from(limit - h).unwrap() / u128::try_from(h + 1).unwrap() + }; + let last = usize::try_from(variables.min(rows)).unwrap(); + let counts = |a| if only_sums { (a - 1, 0) } else { (0, a - 1) }; + let (s, u) = counts(last); + let (t, n, m) = checked_layout(last, s, u, h).unwrap(); + assert_eq!( + u128::try_from(t).unwrap(), + 2 * u128::try_from(last).unwrap() - 1 + ); + assert_eq!( + u128::try_from(n).unwrap(), + u128::try_from(last).unwrap() * u128::try_from(2 * h + 1).unwrap() - 1 + ); + assert_eq!( + u128::try_from(m).unwrap(), + 1 + 2 * u128::try_from(s).unwrap() + + u128::try_from(u).unwrap() + + u128::try_from(last).unwrap() * u128::try_from(h).unwrap() + + u128::try_from(h).unwrap() + ); + let (s, u) = counts(last + 1); + assert!(checked_layout(last + 1, s, u, h).is_err()); + } + } + for (a, s, u, h) in [ + (usize::MAX, 0, 0, 1), + (usize::MAX / 2, 0, 0, 63), + (1, usize::MAX, 0, 1), + (1, 0, usize::MAX, 1), + ] { + assert!(checked_layout(a, s, u, h).is_err()); + } +}