diff --git a/.claude/CLAUDE.md b/.claude/CLAUDE.md index 289937f0..72a8feef 100644 --- a/.claude/CLAUDE.md +++ b/.claude/CLAUDE.md @@ -10,6 +10,7 @@ make clippy # Lint make export-graph # Regenerate reduction graph make paper # Build Typst paper make coverage # Generate coverage report (>95% required) +make qubo-testdata # Regenerate QUBO ground truth JSON ``` ## Verify Changes @@ -28,6 +29,9 @@ make test clippy export-graph # Must pass before PR - `src/registry/` - Compile-time reduction metadata collection - `src/unit_tests/` - Unit test files (mirroring `src/` structure, referenced via `#[path]`) - `tests/main.rs` - Integration tests (modules in `tests/suites/`) +- `tests/data/` - Ground truth JSON for integration tests +- `scripts/` - Python test data generation scripts (managed with `uv`) +- `docs/plans/` - Implementation plans ### Trait Hierarchy diff --git a/.claude/rules/adding-reductions.md b/.claude/rules/adding-reductions.md index d9c01cb2..20411ac8 100644 --- a/.claude/rules/adding-reductions.md +++ b/.claude/rules/adding-reductions.md @@ -5,19 +5,59 @@ paths: # Adding a Reduction Rule (A → B) +## 0. Brainstorm & Generate Test Data First + +Before writing any Rust code, follow this workflow: + +1. **Brainstorm the reduction** — use `superpowers:brainstorming` to discuss with the user: + - Research the mathematical formulation (paper, textbook, or derive it) + - Understand the variable mapping and constraint encoding + - Discuss implementation approach: penalty values, matrix construction, solution extraction + - Read reference implementations in the codebase (e.g., `src/rules/spinglass_qubo.rs`) to understand conventions + - Agree on scope (weighted vs unweighted, specific graph types, const generics) +2. **Generate ground truth test data** — use an existing library (e.g., Python with qubogen, qubovert, or networkx) to create small instances, reduce them, brute-force solve both sides, and export as JSON to `tests/data//`. It is recommended to download the relevant package and check the existing tests to understand how to construct tests. To generate the test data, you can use the following command: + ```bash + # Example: generate QUBO test data + cd scripts && uv run python generate_qubo_tests.py + ``` +3. **Create a practical example** — design a small, explainable instance for `examples/` (e.g., "wireless tower placement" for IndependentSet, "map coloring" for Coloring). This example will also appear in the `docs/paper/reductions.typ`. +4. **Write the implementation plan** — save to `docs/plans/` using `superpowers:writing-plans`. The plan must include implementation details from the brainstorming session (formulas, penalty terms, matrix construction, variable indexing). + ## 1. Implementation -Create `src/rules/_.rs`: + +Create `src/rules/_.rs` following the pattern in `src/rules/spinglass_qubo.rs`: ```rust -use problemreductions::reduction; +use crate::reduction; + +#[derive(Debug, Clone)] +pub struct ReductionSourceToTarget { + target: TargetProblem<...>, + source_size: ProblemSize, + // + any metadata needed for extract_solution +} + +impl ReductionResult for ReductionSourceToTarget { + type Source = SourceProblem<...>; + type Target = TargetProblem<...>; + + fn target_problem(&self) -> &Self::Target { &self.target } + fn extract_solution(&self, target_solution: &[usize]) -> Vec { ... } + fn source_size(&self) -> ProblemSize { self.source_size.clone() } + fn target_size(&self) -> ProblemSize { self.target.problem_size() } +} #[reduction( overhead = { ReductionOverhead::new(vec![...]) } )] -impl ReduceTo> for SourceProblem { +impl ReduceTo> for SourceProblem<...> { type Result = ReductionSourceToTarget; fn reduce_to(&self) -> Self::Result { ... } } + +#[cfg(test)] +#[path = "../unit_tests/rules/_.rs"] +mod tests; ``` The `#[reduction]` macro auto-generates the `inventory::submit!` call. Optional attributes: `source_graph`, `target_graph`, `source_weighted`, `target_weighted`. @@ -28,25 +68,39 @@ mod source_target; pub use source_target::ReductionSourceToTarget; ``` -## 2. Closed-Loop Test (Required) +## 2. Tests (Required) + +- **Unit tests** in `src/unit_tests/rules/_.rs` — closed-loop + edge cases. See `rules/testing.md`. +- **Integration tests** in `tests/suites/reductions.rs` — compare against JSON ground truth from step 0. +- Test name: `test__to__closed_loop` -See `rules/testing.md` for the full pattern. Test name: `test__to__closed_loop`. +## 3. Example Program + +Add a round-trip demo to `examples/` showing a practical, explainable instance: +1. Create source problem with a real-world story +2. Reduce to target, solve, extract solution +3. Print human-readable explanation + +## 4. Documentation -## 3. Documentation Update `docs/paper/reductions.typ` (see `rules/documentation.md` for the pattern): - Add theorem + proof sketch -- Add code example +- Add Rust code example from the example program - Add to summary table with overhead and citation +The goal is to 1. prove the correctness of the reduction to human beings. 2. provide a minimal working example to the readers. + Citations must be verifiable. Use `[Folklore]` or `—` for trivial reductions. -## 4. Regenerate Reduction Graph +## 5. Regenerate Reduction Graph ```bash make export-graph ``` ## Anti-patterns +- Don't write Rust code before understanding the math and having test data - Don't create reductions without closed-loop tests - Don't forget `inventory::submit!` registration (reduction graph won't update) - Don't hardcode weights - use generic `W` parameter - Don't skip overhead polynomial specification +- Don't skip the example program — every reduction needs an explainable demo diff --git a/.claude/skills/issue-to-pr.md b/.claude/skills/issue-to-pr.md index f4946055..3e532175 100644 --- a/.claude/skills/issue-to-pr.md +++ b/.claude/skills/issue-to-pr.md @@ -19,13 +19,15 @@ Convert a GitHub issue into an actionable PR with a plan that auto-triggers Clau digraph issue_to_pr { "Receive issue number" [shape=box]; "Fetch issue with gh" [shape=box]; + "Check the rules to follow" [shape=box]; "Brainstorm with user" [shape=box]; "Write plan file" [shape=box]; "Create branch and PR" [shape=box]; "PR triggers [action]" [shape=doublecircle]; "Receive issue number" -> "Fetch issue with gh"; - "Fetch issue with gh" -> "Brainstorm with user"; + "Fetch issue with gh" -> "Check the rules to follow"; + "Check the rules to follow" -> "Brainstorm with user"; "Brainstorm with user" -> "Write plan file"; "Write plan file" -> "Create branch and PR"; "Create branch and PR" -> "PR triggers [action]"; @@ -53,36 +55,18 @@ Present issue summary to user. **REQUIRED:** Invoke `superpowers:brainstorming` skill with the issue context (if superpowers plugin is available). Otherwise, conduct a manual brainstorming discussion with the user. -This ensures: -- User intent is clarified -- Multiple approaches are explored -- Requirements are understood before planning +Brainstorming must cover: +- **User intent** — clarify what the issue is asking for +- **Multiple approaches** — explore 2-3 different implementation strategies +- **Implementation details** — discuss the mathematical formulation, data structures, variable mappings, constraint encodings, and any non-obvious design choices +- **Existing patterns** — read reference implementations in the codebase (e.g., `spinglass_qubo.rs` for reductions) to understand the conventions +- **Scope** — agree on which variants to implement (e.g., unweighted only, specific K values) Do NOT skip brainstorming. Do NOT write a plan without user discussion. ### 4. Write Plan -After brainstorming concludes, write plan to `issue--.md` in the repo root: - -```markdown -# - -Issue: #<number> - -## Context -<Brief problem statement> - -## Approach -<Chosen approach from brainstorming> - -## Tasks -1. <Specific implementation task> -2. <Another task> -... - -## Acceptance Criteria -- <Criteria from issue/brainstorming> -``` +After brainstorming concludes, write plan to `docs/plans/YYYY-MM-DD-<slug>.md` using `superpowers:writing-plans`: ### 5. Create PR @@ -90,8 +74,8 @@ Issue: #<number> # Create branch git checkout -b issue-<number>-<slug> -# Stage only the plan file -git add issue-<number>-<slug>.md +# Stage the plan file +git add docs/plans/<plan-file>.md # Commit git commit -m "Add plan for #<number>: <title>" @@ -117,19 +101,19 @@ User: /issue-to-pr 42 Claude: Let me fetch issue #42... -[Fetches issue: "Add dark mode support"] +[Fetches issue: "Add IndependentSet → QUBO reduction"] -I'll use superpowers:brainstorming to explore this with you. +I'll read the rules to follow in .claude/rules/adding-reductions.md and use superpowers:brainstorming to explore this with you. [Invokes brainstorming - discusses approaches, user preferences, scope] -Based on our discussion, I'll create the plan... +Based on our discussion, I'll create the plan with superpowers:writing-plans... -[Writes docs/plans/issue-42-dark-mode.md] +[Writes docs/plans/2026-02-09-independentset-to-qubo.md] [Creates branch, commits, pushes] [Creates PR with body starting with "[action]"] -Created PR #45: Fix #42: Add dark mode support +Created PR #45: Fix #42: Add IndependentSet → QUBO reduction The [action] trigger will automatically execute the plan. ``` @@ -139,5 +123,5 @@ The [action] trigger will automatically execute the plan. |---------|-----| | Skipping brainstorming | Always use superpowers:brainstorming (or manual discussion) first | | `[action]` not at start | PR body must BEGIN with `[action]` | -| Including code in PR | Only commit the plan file | +| Including implementation code in initial PR | First PR: plan only | | Generic plan | Use specifics from brainstorming | diff --git a/.gitignore b/.gitignore index 64f2d8c3..756e66cc 100644 --- a/.gitignore +++ b/.gitignore @@ -38,11 +38,12 @@ quickcheck-tests.json /dist/ /build/ -# Python (for pre-commit) +# Python __pycache__/ *.py[cod] *$py.class .Python +.venv/ venv/ env/ ENV/ @@ -67,3 +68,6 @@ tarpaulin-report.html # Generated documents *.pdf + +# Reference packages (downloaded for comparison) +pkgref/ diff --git a/Makefile b/Makefile index bce574f7..47fee9fe 100644 --- a/Makefile +++ b/Makefile @@ -1,6 +1,6 @@ # Makefile for problemreductions -.PHONY: help build test fmt clippy doc mdbook paper clean coverage rust-export compare +.PHONY: help build test fmt clippy doc mdbook paper clean coverage rust-export compare qubo-testdata # Default target help: @@ -18,6 +18,7 @@ help: @echo " check - Quick check (fmt + clippy + test)" @echo " rust-export - Generate Rust mapping JSON exports" @echo " compare - Generate and compare Rust mapping exports" + @echo " qubo-testdata - Regenerate QUBO test data (requires uv)" # Build the project build: @@ -65,6 +66,10 @@ clean: check: fmt-check clippy test @echo "✅ All checks passed!" +# Regenerate QUBO test data from Python (requires uv) +qubo-testdata: + cd scripts && uv run python generate_qubo_tests.py + # Generate Rust mapping JSON exports for all graphs and modes GRAPHS := diamond bull house petersen MODES := unweighted weighted triangular diff --git a/README.md b/README.md index 5b803a11..bc285e68 100644 --- a/README.md +++ b/README.md @@ -72,6 +72,21 @@ make check # Quick check before commit (fmt + clippy + test) ``` +## Acknowledgments + +This project draws inspiration from the following packages: + +- **[ProblemReductions.jl](https://github.com/GiggleLiu/ProblemReductions.jl)** — Julia library for computational problem reductions. Our problem trait hierarchy, reduction interface (`ReduceTo`/`ReductionResult`), and graph-based reduction registry are directly inspired by this package. +- **[UnitDiskMapping.jl](https://github.com/QuEraComputing/UnitDiskMapping.jl)** — Julia package for mapping problems to unit disk graphs. Our unit disk graph (King's subgraph / triangular lattice) reductions and the copy-line method are based on this implementation. +- **[qubogen](https://github.com/tamuhey/qubogen)** — Python library for generating QUBO matrices from combinatorial problems. Our QUBO reduction formulas (Vertex Cover, Graph Coloring, Set Packing, Max-2-SAT, binary ILP) reference the implementations in this package. + +## Related Projects + +- **[Karp](https://github.com/REA1/karp)** — A DSL (built on Racket) for writing and testing Karp reductions between NP-complete problems ([PLDI 2022 paper](https://dl.acm.org/doi/abs/10.1145/3519939.3523732)). Focused on education and proof verification rather than a solver pipeline. +- **[Complexity Zoo](https://complexityzoo.net/)** — Comprehensive catalog of 550+ computational complexity classes (Scott Aaronson). +- **[A Compendium of NP Optimization Problems](https://www.csc.kth.se/tcs/compendium/)** — Online catalog of NP optimization problems with approximability results (Crescenzi & Kann). +- **Computers and Intractability** (Garey & Johnson, 1979) — The classic reference cataloging 300+ NP-complete problems with reductions. The most cited book in computer science. + ## License MIT License - see [LICENSE](LICENSE) for details. diff --git a/docs/paper/reduction_graph.json b/docs/paper/reduction_graph.json index 99dce284..399dcbab 100644 --- a/docs/paper/reduction_graph.json +++ b/docs/paper/reduction_graph.json @@ -73,6 +73,14 @@ }, "category": "graph" }, + { + "name": "IndependentSet", + "variant": { + "graph": "SimpleGraph", + "weight": "i32" + }, + "category": "graph" + }, { "name": "KColoring", "variant": {}, @@ -87,6 +95,14 @@ }, "category": "graph" }, + { + "name": "KColoring", + "variant": { + "graph": "SimpleGraph", + "weight": "Unweighted" + }, + "category": "graph" + }, { "name": "KColoring", "variant": { @@ -108,6 +124,14 @@ }, "category": "satisfiability" }, + { + "name": "KSatisfiability", + "variant": { + "graph": "SimpleGraph", + "weight": "i32" + }, + "category": "satisfiability" + }, { "name": "Matching", "variant": {}, @@ -186,6 +210,14 @@ }, "category": "set" }, + { + "name": "SetPacking", + "variant": { + "graph": "SimpleGraph", + "weight": "i32" + }, + "category": "set" + }, { "name": "SpinGlass", "variant": {}, @@ -219,6 +251,14 @@ "weight": "Unweighted" }, "category": "graph" + }, + { + "name": "VertexCovering", + "variant": { + "graph": "SimpleGraph", + "weight": "i32" + }, + "category": "graph" } ], "edges": [ @@ -273,6 +313,40 @@ }, "bidirectional": false }, + { + "source": { + "name": "ILP", + "variant": { + "graph": "SimpleGraph", + "weight": "Unweighted" + } + }, + "target": { + "name": "QUBO", + "variant": { + "graph": "SimpleGraph", + "weight": "f64" + } + }, + "bidirectional": false + }, + { + "source": { + "name": "IndependentSet", + "variant": { + "graph": "SimpleGraph", + "weight": "i32" + } + }, + "target": { + "name": "QUBO", + "variant": { + "graph": "SimpleGraph", + "weight": "f64" + } + }, + "bidirectional": false + }, { "source": { "name": "KColoring", @@ -291,6 +365,40 @@ }, "bidirectional": false }, + { + "source": { + "name": "KColoring", + "variant": { + "graph": "SimpleGraph", + "weight": "Unweighted" + } + }, + "target": { + "name": "QUBO", + "variant": { + "graph": "SimpleGraph", + "weight": "f64" + } + }, + "bidirectional": false + }, + { + "source": { + "name": "KSatisfiability", + "variant": { + "graph": "SimpleGraph", + "weight": "i32" + } + }, + "target": { + "name": "QUBO", + "variant": { + "graph": "SimpleGraph", + "weight": "f64" + } + }, + "bidirectional": false + }, { "source": { "name": "Matching", @@ -393,6 +501,23 @@ }, "bidirectional": true }, + { + "source": { + "name": "SetPacking", + "variant": { + "graph": "SimpleGraph", + "weight": "i32" + } + }, + "target": { + "name": "QUBO", + "variant": { + "graph": "SimpleGraph", + "weight": "f64" + } + }, + "bidirectional": false + }, { "source": { "name": "SpinGlass", @@ -460,6 +585,23 @@ } }, "bidirectional": false + }, + { + "source": { + "name": "VertexCovering", + "variant": { + "graph": "SimpleGraph", + "weight": "i32" + } + }, + "target": { + "name": "QUBO", + "variant": { + "graph": "SimpleGraph", + "weight": "f64" + } + }, + "bidirectional": false } ] } \ No newline at end of file diff --git a/docs/paper/reductions.typ b/docs/paper/reductions.typ index b8e468fa..3e86717a 100644 --- a/docs/paper/reductions.typ +++ b/docs/paper/reductions.typ @@ -64,7 +64,7 @@ In all graph problems below, $G = (V, E)$ denotes an undirected graph with $|V| #definition("Independent Set (IS)")[ Given $G = (V, E)$ with vertex weights $w: V -> RR$, find $S subset.eq V$ maximizing $sum_(v in S) w(v)$ such that no two vertices in $S$ are adjacent: $forall u, v in S: (u, v) in.not E$. - _Reduces to:_ Set Packing (@def:set-packing). + _Reduces to:_ Set Packing (@def:set-packing), QUBO (@def:qubo). _Reduces from:_ Vertex Cover (@def:vertex-cover), SAT (@def:satisfiability), Set Packing (@def:set-packing). @@ -87,7 +87,7 @@ In all graph problems below, $G = (V, E)$ denotes an undirected graph with $|V| #definition("Vertex Cover (VC)")[ Given $G = (V, E)$ with vertex weights $w: V -> RR$, find $S subset.eq V$ minimizing $sum_(v in S) w(v)$ such that every edge has at least one endpoint in $S$: $forall (u, v) in E: u in S or v in S$. - _Reduces to:_ Independent Set (@def:independent-set), Set Covering (@def:set-covering). + _Reduces to:_ Independent Set (@def:independent-set), Set Covering (@def:set-covering), QUBO (@def:qubo). _Reduces from:_ Independent Set (@def:independent-set). @@ -132,7 +132,7 @@ In all graph problems below, $G = (V, E)$ denotes an undirected graph with $|V| #definition("Graph Coloring")[ Given $G = (V, E)$ and $k$ colors, find $c: V -> {1, ..., k}$ minimizing $|{(u, v) in E : c(u) = c(v)}|$. - _Reduces to:_ ILP (@def:ilp). + _Reduces to:_ ILP (@def:ilp), QUBO (@def:qubo). _Reduces from:_ SAT (@def:satisfiability). @@ -204,7 +204,7 @@ In all graph problems below, $G = (V, E)$ denotes an undirected graph with $|V| #definition("Set Packing")[ Given universe $U$, collection $cal(S) = {S_1, ..., S_m}$ with $S_i subset.eq U$, weights $w: cal(S) -> RR$, find $cal(P) subset.eq cal(S)$ maximizing $sum_(S in cal(P)) w(S)$ s.t. $forall S_i, S_j in cal(P): S_i inter S_j = emptyset$. - _Reduces to:_ Independent Set (@def:independent-set). + _Reduces to:_ Independent Set (@def:independent-set), QUBO (@def:qubo). _Reduces from:_ Independent Set (@def:independent-set), Matching (@def:matching). @@ -273,11 +273,11 @@ In all graph problems below, $G = (V, E)$ denotes an undirected graph with $|V| ] <def:spin-glass> #definition("QUBO")[ - Given $n$ binary variables $x_i in {0, 1}$, matrix $Q in RR^(n times n)$, minimize $f(bold(x)) = bold(x)^top Q bold(x)$. + Given $n$ binary variables $x_i in {0, 1}$, upper-triangular matrix $Q in RR^(n times n)$, minimize $f(bold(x)) = sum_(i=1)^n Q_(i i) x_i + sum_(i < j) Q_(i j) x_i x_j$ (using $x_i^2 = x_i$ for binary variables). _Reduces to:_ Spin Glass (@def:spin-glass). - _Reduces from:_ Spin Glass (@def:spin-glass). + _Reduces from:_ Spin Glass (@def:spin-glass), Independent Set (@def:independent-set), Vertex Cover (@def:vertex-cover), Graph Coloring (@def:coloring), Set Packing (@def:set-packing), $k$-SAT (@def:k-sat), ILP (@def:ilp). ```rust pub struct QUBO<W = f64> { @@ -298,6 +298,8 @@ In all graph problems below, $G = (V, E)$ denotes an undirected graph with $|V| #definition("Integer Linear Programming (ILP)")[ Given $n$ integer variables $bold(x) in ZZ^n$, constraint matrix $A in RR^(m times n)$, bounds $bold(b) in RR^m$, and objective $bold(c) in RR^n$, find $bold(x)$ minimizing $bold(c)^top bold(x)$ subject to $A bold(x) <= bold(b)$ and variable bounds. + _Reduces to:_ QUBO (@def:qubo). + _Reduces from:_ Graph Coloring (@def:coloring), Factoring (@def:factoring). ```rust @@ -359,7 +361,7 @@ In all graph problems below, $G = (V, E)$ denotes an undirected graph with $|V| #definition([$k$-SAT])[ SAT with exactly $k$ literals per clause. - _Reduces to:_ SAT (@def:satisfiability). + _Reduces to:_ SAT (@def:satisfiability), QUBO (@def:qubo). _Reduces from:_ SAT (@def:satisfiability). @@ -502,6 +504,184 @@ let sg_solution = result.extract_solution(&qubo_solutions[0]); assert_eq!(sg_solution.len(), 2); ``` +== Penalty-Method QUBO Reductions <sec:penalty-method> + +The _penalty method_ @glover2019 @lucas2014 converts a constrained optimization problem into an unconstrained QUBO by adding quadratic penalty terms. Given an objective $"obj"(bold(x))$ to minimize and constraints $g_k (bold(x)) = 0$, construct: +$ f(bold(x)) = "obj"(bold(x)) + P sum_k g_k (bold(x))^2 $ +where $P$ is a penalty weight large enough that any constraint violation costs more than the entire objective range. Since $g_k (bold(x))^2 >= 0$ with equality iff $g_k (bold(x)) = 0$, minimizers of $f$ are feasible and optimal for the original problem. Because binary variables satisfy $x_i^2 = x_i$, the resulting $f$ is a quadratic in $bold(x)$, i.e.\ a QUBO. + +#theorem[ + *(IS $arrow.r$ QUBO)* Given $G = (V, E)$ with weights $w$, construct upper-triangular $Q in RR^(n times n)$ with $Q_(i i) = -w_i$ and $Q_(i j) = P$ for $(i,j) in E$ ($i < j$), where $P = 1 + sum_i w_i$. Then minimizing $f(bold(x)) = sum_i Q_(i i) x_i + sum_(i<j) Q_(i j) x_i x_j$ is equivalent to maximizing the IS objective. [_Problems:_ @def:independent-set, @def:qubo.] +] + +#proof[ + _Construction._ The IS objective is: maximize $sum_i w_i x_i$ subject to $x_i x_j = 0$ for $(i,j) in E$. Applying the penalty method (@sec:penalty-method): + $ f(bold(x)) = -sum_i w_i x_i + P sum_((i,j) in E) x_i x_j $ + Reading off the QUBO coefficients: diagonal $Q_(i i) = -w_i$ (linear terms), off-diagonal $Q_(i j) = P$ for edges $i < j$ (quadratic penalty). + + _Correctness._ If $bold(x)$ has any adjacent pair $(x_i = 1, x_j = 1)$ with $(i,j) in E$, the penalty $P > sum_i w_i >= -sum_i Q_(i i) x_i$ exceeds the maximum objective gain, so $bold(x)$ is not a minimizer. Among independent sets ($x_i x_j = 0$ for all edges), $f(bold(x)) = -sum_(i in S) w_i$, minimized exactly when $S$ is a maximum-weight IS. +] + +```rust +// Minimal example: IS -> QUBO -> extract solution +let is = IndependentSet::<SimpleGraph, i32>::new(4, vec![(0, 1), (1, 2), (2, 3)]); +let result = ReduceTo::<QUBO>::reduce_to(&is); +let qubo = result.target_problem(); + +let solver = BruteForce::new(); +let solutions = solver.find_best(qubo); +let is_solution = result.extract_solution(&solutions[0]); +assert!(is.solution_size(&is_solution).is_valid); +``` + +#theorem[ + *(VC $arrow.r$ QUBO)* Given $G = (V, E)$ with weights $w$, construct upper-triangular $Q$ with $Q_(i i) = w_i - P dot "deg"(i)$ and $Q_(i j) = P$ for $(i,j) in E$ ($i < j$), where $P = 1 + sum_i w_i$ and $"deg"(i)$ is the degree of vertex $i$. [_Problems:_ @def:vertex-cover, @def:qubo.] +] + +#proof[ + _Construction._ The VC objective is: minimize $sum_i w_i x_i$ subject to $x_i + x_j >= 1$ for $(i,j) in E$. Applying the penalty method (@sec:penalty-method), the constraint $x_i + x_j >= 1$ is violated iff $x_i = x_j = 0$, with penalty $(1 - x_i)(1 - x_j)$: + $ f(bold(x)) = sum_i w_i x_i + P sum_((i,j) in E) (1 - x_i)(1 - x_j) $ + Expanding: $(1 - x_i)(1 - x_j) = 1 - x_i - x_j + x_i x_j$. + Summing over all edges, each vertex $i$ appears in $"deg"(i)$ terms. The QUBO coefficients are: diagonal $Q_(i i) = w_i - P dot "deg"(i)$ (objective plus linear penalty), off-diagonal $Q_(i j) = P$ for edges. The constant $P |E|$ does not affect the minimizer. +] + +```rust +// Minimal example: VC -> QUBO -> extract solution +let vc = VertexCovering::<SimpleGraph, i32>::new(4, vec![(0, 1), (1, 2), (2, 3), (0, 3)]); +let result = ReduceTo::<QUBO>::reduce_to(&vc); +let qubo = result.target_problem(); + +let solver = BruteForce::new(); +let solutions = solver.find_best(qubo); +let vc_solution = result.extract_solution(&solutions[0]); +assert!(vc.solution_size(&vc_solution).is_valid); +``` + +#theorem[ + *(KColoring $arrow.r$ QUBO)* Given $G = (V, E)$ with $k$ colors, construct upper-triangular $Q in RR^(n k times n k)$ using one-hot encoding $x_(v,c) in {0,1}$ ($n k$ variables indexed by $v dot k + c$). [_Problems:_ @def:coloring, @def:qubo.] +] + +#proof[ + _Construction._ Applying the penalty method (@sec:penalty-method), the QUBO objective combines a one-hot constraint penalty and an edge conflict penalty: + $ f(bold(x)) = P_1 sum_(v in V) (1 - sum_(c=1)^k x_(v,c))^2 + P_2 sum_((u,v) in E) sum_(c=1)^k x_(u,c) x_(v,c) $ + + _One-hot expansion._ For each vertex $v$, using $x_(v,c)^2 = x_(v,c)$: + $ (1 - sum_c x_(v,c))^2 = 1 - sum_c x_(v,c) + 2 sum_(c_1 < c_2) x_(v,c_1) x_(v,c_2) $ + This yields diagonal $Q_(v k+c, v k+c) = -P_1$ and intra-vertex off-diagonal $Q_(v k+c_1, v k+c_2) = 2 P_1$ for $c_1 < c_2$. + + _Edge penalty._ For each edge $(u,v)$ and color $c$, the term $P_2 x_(u,c) x_(v,c)$ contributes to $Q_(u k+c, v k+c) += P_2$ (with appropriate index ordering). + + In our implementation, $P_1 = P = 1 + n$ and $P_2 = P\/2$. + + _Solution extraction._ For each vertex $v$, find $c$ with $x_(v,c) = 1$. +] + +```rust +// Minimal example: KColoring -> QUBO -> extract solution +let kc = KColoring::<3, SimpleGraph, i32>::new(3, vec![(0, 1), (1, 2), (0, 2)]); +let result = ReduceTo::<QUBO>::reduce_to(&kc); +let qubo = result.target_problem(); + +let solver = BruteForce::new(); +let solutions = solver.find_best(qubo); +let kc_solution = result.extract_solution(&solutions[0]); +assert_eq!(solutions.len(), 6); // 3! valid 3-colorings of K3 +``` + +#theorem[ + *(SetPacking $arrow.r$ QUBO)* Equivalent to IS on the intersection graph: $Q_(i i) = -w_i$ and $Q_(i j) = P$ for overlapping sets $i, j$ ($i < j$), where $P = 1 + sum_i w_i$. [_Problems:_ @def:set-packing, @def:qubo.] +] + +#proof[ + Two sets conflict iff they share an element. The intersection graph has sets as vertices and edges between conflicting pairs. Applying the penalty method (@sec:penalty-method) yields the same QUBO as IS on this graph: diagonal rewards selection, off-diagonal penalizes overlap. Correctness follows from the IS→QUBO proof. +] + +```rust +// Minimal example: SetPacking -> QUBO -> extract solution +let sp = SetPacking::<i32>::new(vec![vec![0, 1], vec![1, 2], vec![2, 3, 4]]); +let result = ReduceTo::<QUBO>::reduce_to(&sp); +let qubo = result.target_problem(); + +let solver = BruteForce::new(); +let solutions = solver.find_best(qubo); +let sp_solution = result.extract_solution(&solutions[0]); +assert!(sp.solution_size(&sp_solution).is_valid); +``` + +#theorem[ + *(2-SAT $arrow.r$ QUBO)* Given a Max-2-SAT instance with $m$ clauses over $n$ variables, construct upper-triangular $Q in RR^(n times n)$ where each clause $(ell_i or ell_j)$ contributes a penalty gadget encoding its unique falsifying assignment. [_Problems:_ @def:k-sat, @def:qubo.] +] + +#proof[ + _Construction._ Applying the penalty method (@sec:penalty-method), each 2-literal clause has exactly one falsifying assignment (both literals false). The penalty for that assignment is a quadratic function of $x_i, x_j$: + + #table( + columns: (auto, auto, auto, auto), + inset: 4pt, + align: left, + table.header([*Clause*], [*Falsified when*], [*Penalty*], [*QUBO contributions*]), + [$x_i or x_j$], [$x_i=0, x_j=0$], [$(1-x_i)(1-x_j)$], [$Q_(i i) -= 1, Q_(j j) -= 1, Q_(i j) += 1$], + [$overline(x_i) or x_j$], [$x_i=1, x_j=0$], [$x_i(1-x_j)$], [$Q_(i i) += 1, Q_(i j) -= 1$], + [$x_i or overline(x_j)$], [$x_i=0, x_j=1$], [$(1-x_i)x_j$], [$Q_(j j) += 1, Q_(i j) -= 1$], + [$overline(x_i) or overline(x_j)$], [$x_i=1, x_j=1$], [$x_i x_j$], [$Q_(i j) += 1$], + ) + + Summing over all clauses, $f(bold(x)) = sum_j "penalty"_j (bold(x))$ counts falsified clauses. Minimizers of $f$ maximize satisfied clauses. +] + +```rust +// Minimal example: 2-SAT -> QUBO -> extract solution +let ksat = KSatisfiability::<2, i32>::new(3, vec![ + CNFClause::new(vec![1, 2]), // x1 OR x2 + CNFClause::new(vec![-1, 3]), // NOT x1 OR x3 + CNFClause::new(vec![2, -3]), // x2 OR NOT x3 +]); +let result = ReduceTo::<QUBO>::reduce_to(&ksat); +let qubo = result.target_problem(); + +let solver = BruteForce::new(); +let solutions = solver.find_best(qubo); +let sat_solution = result.extract_solution(&solutions[0]); +assert!(ksat.solution_size(&sat_solution).is_valid); +``` + +#theorem[ + *(Binary ILP $arrow.r$ QUBO)* Given binary ILP: maximize $bold(c)^top bold(x)$ subject to $A bold(x) = bold(b)$, $bold(x) in {0,1}^n$, construct upper-triangular $Q = -"diag"(bold(c) + 2P bold(b)^top A) + P A^top A$ where $P = 1 + ||bold(c)||_1 + ||bold(b)||_1$. [_Problems:_ @def:ilp, @def:qubo.] +] + +#proof[ + _Step 1: Normalize constraints._ Convert inequalities to equalities using slack variables: $bold(a)_k^top bold(x) <= b_k$ becomes $bold(a)_k^top bold(x) + sum_(s=0)^(S_k - 1) 2^s y_(k,s) = b_k$ where $S_k = ceil(log_2 (b_k + 1))$ slack bits. For $>=$ constraints, the slack has a negative sign. The extended system is $A' bold(x)' = bold(b)$ with $bold(x)' = (bold(x), bold(y)) in {0,1}^(n')$. For minimization, negate $bold(c)$ to convert to maximization. + + _Step 2: QUBO construction._ Applying the penalty method (@sec:penalty-method), combine objective and penalty: + $ f(bold(x)') = -bold(c')^top bold(x)' + P sum_(k=1)^m (bold(a)'_k^(top) bold(x)' - b_k)^2 $ + where $bold(c)' = (bold(c), bold(0))$. Expanding the quadratic penalty: + $ = bold(x)'^(top) A'^(top) A' bold(x)' - 2 bold(b)^top A' bold(x)' + ||bold(b)||_2^2 $ + Combining with $-bold(c')^top bold(x)'$ and dropping constants: + $ Q = -"diag"(bold(c)' + 2P bold(b)^top A') + P A'^(top) A' $ + The diagonal contains linear terms; the upper triangle of $A'^(top) A'$ gives quadratic terms (doubled for upper-triangular convention). + + _Solution extraction._ Discard slack variables: return $bold(x)' [0..n]$. +] + +```rust +// Minimal example: binary ILP -> QUBO -> extract solution +let ilp = ILP::binary(3, + vec![ + LinearConstraint::le(vec![(0, 1.0), (1, 1.0)], 1.0), + LinearConstraint::le(vec![(1, 1.0), (2, 1.0)], 1.0), + ], + vec![(0, 1.0), (1, 2.0), (2, 3.0)], + ObjectiveSense::Maximize, +); +let result = ReduceTo::<QUBO<f64>>::reduce_to(&ilp); +let qubo = result.target_problem(); + +let solver = BruteForce::new(); +let solutions = solver.find_best(qubo); +let ilp_solution = result.extract_solution(&solutions[0]); +assert_eq!(ilp_solution, vec![1, 0, 1]); // obj = 4 +``` + == Non-Trivial Reductions #theorem[ @@ -859,7 +1039,12 @@ assert_eq!(p * q, 15); // e.g., (3, 5) or (5, 3) table.cell(fill: gray)[IS $arrow.r$ SetPacking], table.cell(fill: gray)[$O(|V| + |E|)$], table.cell(fill: gray)[—], table.cell(fill: gray)[Matching $arrow.r$ SetPacking], table.cell(fill: gray)[$O(|E|)$], table.cell(fill: gray)[—], table.cell(fill: gray)[VC $arrow.r$ SetCovering], table.cell(fill: gray)[$O(|V| + |E|)$], table.cell(fill: gray)[—], - table.cell(fill: gray)[QUBO $arrow.l.r$ SpinGlass], table.cell(fill: gray)[$O(n^2)$], table.cell(fill: gray)[—], + [IS $arrow.r$ QUBO], [$O(n)$], [@lucas2014 @glover2019], + [VC $arrow.r$ QUBO], [$O(n)$], [@lucas2014 @glover2019], + [KColoring $arrow.r$ QUBO], [$O(n dot k)$], [@lucas2014 @glover2019], + [SetPacking $arrow.r$ QUBO], [$O(n)$], [@glover2019], + [2-SAT $arrow.r$ QUBO], [$O(n)$], [@glover2019], + [Binary ILP $arrow.r$ QUBO], [$O(n)$], [@lucas2014 @glover2019], [SAT $arrow.r$ IS], [$O(sum_j |C_j|^2)$], [@karp1972], [SAT $arrow.r$ 3-Coloring], [$O(n + sum_j |C_j|)$], [@garey1979], [SAT $arrow.r$ DominatingSet], [$O(3n + m)$], [@garey1979], @@ -871,7 +1056,7 @@ assert_eq!(p * q, 15); // e.g., (3, 5) or (5, 3) table.cell(fill: gray)[Factoring $arrow.r$ ILP], table.cell(fill: gray)[$O(m n)$], table.cell(fill: gray)[—], [IS $arrow.r$ GridGraph IS], [$O(n^2)$], [@nguyen2023], ), - caption: [Summary of reductions. Gray rows indicate trivial reductions.] + caption: [Summary of reductions. Gray rows indicate trivial (complement/isomorphism) reductions.] ) <tab:summary> #bibliography("references.bib", style: "ieee") diff --git a/docs/paper/references.bib b/docs/paper/references.bib index 0225594c..fe632703 100644 --- a/docs/paper/references.bib +++ b/docs/paper/references.bib @@ -22,6 +22,16 @@ @book{garey1979 year = {1979} } +@article{glover2019, + author = {Fred Glover and Gary Kochenberger and Yu Du}, + title = {Quantum Bridge Analytics {I}: a tutorial on formulating and using {QUBO} models}, + journal = {4OR}, + volume = {17}, + pages = {335--371}, + year = {2019}, + doi = {10.1007/s10288-019-00424-y} +} + @article{lucas2014, author = {Andrew Lucas}, title = {Ising formulations of many NP problems}, diff --git a/docs/plans/2026-02-08-qubo-reductions-plan.md b/docs/plans/2026-02-08-qubo-reductions-plan.md new file mode 100644 index 00000000..7c1a53c0 --- /dev/null +++ b/docs/plans/2026-02-08-qubo-reductions-plan.md @@ -0,0 +1,414 @@ +# Problem-to-QUBO Reductions Implementation Plan + +> **For Claude:** REQUIRED SUB-SKILL: Use superpowers:executing-plans to implement this plan task-by-task. + +**Goal:** Implement 7 reductions from NP-hard problems to QUBO, with tests, examples, and paper documentation (Issue #18). + +**Architecture:** Each reduction creates a `QUBO<f64>` matrix encoding the source problem's objective + constraints as penalty terms. All reductions follow the existing pattern in `src/rules/spinglass_qubo.rs`: a result struct implementing `ReductionResult`, a `ReduceTo` impl with `#[reduction]` macro, unit tests via `#[path]`, and integration tests in `tests/suites/reductions.rs`. + +**Tech Stack:** Rust, `#[reduction]` proc macro, `inventory` for registration, `BruteForce` solver for tests. Ground truth JSON in `tests/data/qubo/` (already generated via PR #29). + +**Branch:** `issue-18-qubo-reductions` (already exists, PR #29) + +--- + +### Task 1: IndependentSet → QUBO + +Maximize weighted IS = minimize `-Σ w_i·x_i + P·Σ_{(i,j)∈E} x_i·x_j` where `P > Σ w_i`. + +**Files:** +- Create: `src/rules/independentset_qubo.rs` +- Create: `src/unit_tests/rules/independentset_qubo.rs` +- Modify: `src/rules/mod.rs` — add `mod independentset_qubo;` + `pub use` + +**Step 1: Write unit test** + +File: `src/unit_tests/rules/independentset_qubo.rs` + +```rust +use super::*; +use crate::solvers::{BruteForce, Solver}; + +#[test] +fn test_independentset_to_qubo_closed_loop() { + // Path graph: 0-1-2-3 (4 vertices, 3 edges) + // Maximum IS = {0, 2} or {1, 3} (size 2) + let is = IndependentSet::<SimpleGraph, Unweighted>::new(4, vec![(0, 1), (1, 2), (2, 3)]); + let reduction = ReduceTo::<QUBO<f64>>::reduce_to(&is); + let qubo = reduction.target_problem(); + + let solver = BruteForce::new(); + let qubo_solutions = solver.find_best(qubo); + + for sol in &qubo_solutions { + let extracted = reduction.extract_solution(sol); + assert!(is.solution_size(&extracted).is_valid); + // IS of size 2 + assert_eq!(extracted.iter().filter(|&&x| x == 1).count(), 2); + } +} + +#[test] +fn test_independentset_to_qubo_triangle() { + // Triangle: 0-1-2 (complete graph K3) + // Maximum IS = any single vertex (size 1) + let is = IndependentSet::<SimpleGraph, Unweighted>::new(3, vec![(0, 1), (1, 2), (0, 2)]); + let reduction = ReduceTo::<QUBO<f64>>::reduce_to(&is); + let qubo = reduction.target_problem(); + + let solver = BruteForce::new(); + let qubo_solutions = solver.find_best(qubo); + + for sol in &qubo_solutions { + let extracted = reduction.extract_solution(sol); + assert!(is.solution_size(&extracted).is_valid); + assert_eq!(extracted.iter().filter(|&&x| x == 1).count(), 1); + } +} +``` + +**Step 2: Run test, verify it fails** + +Run: `cargo test --all-features test_independentset_to_qubo` +Expected: compilation error (module not found) + +**Step 3: Write reduction implementation** + +File: `src/rules/independentset_qubo.rs` + +```rust +//! Reduction from IndependentSet to QUBO. +//! +//! Maximize Σ w_i·x_i s.t. x_i·x_j = 0 for (i,j) ∈ E +//! = Minimize -Σ w_i·x_i + P·Σ_{(i,j)∈E} x_i·x_j +//! +//! Q[i][i] = -w_i, Q[i][j] = P for edges. P = 1 + Σ w_i. + +use crate::models::graph::IndependentSet; +use crate::models::optimization::QUBO; +use crate::poly; +use crate::reduction; +use crate::rules::registry::ReductionOverhead; +use crate::rules::traits::{ReduceTo, ReductionResult}; +use crate::topology::SimpleGraph; +use crate::traits::Problem; +use crate::types::{ProblemSize, Unweighted}; + +#[derive(Debug, Clone)] +pub struct ReductionISToQUBO { + target: QUBO<f64>, + source_size: ProblemSize, +} + +impl ReductionResult for ReductionISToQUBO { + type Source = IndependentSet<SimpleGraph, Unweighted>; + type Target = QUBO<f64>; + + fn target_problem(&self) -> &Self::Target { &self.target } + fn extract_solution(&self, target_solution: &[usize]) -> Vec<usize> { + target_solution.to_vec() + } + fn source_size(&self) -> ProblemSize { self.source_size.clone() } + fn target_size(&self) -> ProblemSize { self.target.problem_size() } +} + +#[reduction( + source_graph = "SimpleGraph", + overhead = { ReductionOverhead::new(vec![("num_vars", poly!(num_vertices))]) } +)] +impl ReduceTo<QUBO<f64>> for IndependentSet<SimpleGraph, Unweighted> { + type Result = ReductionISToQUBO; + + fn reduce_to(&self) -> Self::Result { + let n = self.num_vertices(); + let edges = self.edges(); + let penalty = 1.0 + n as f64; // P > sum of unit weights + + let mut matrix = vec![vec![0.0; n]; n]; + for i in 0..n { + matrix[i][i] = -1.0; // -w_i (unit weight) + } + for (u, v) in &edges { + let (i, j) = if u < v { (*u, *v) } else { (*v, *u) }; + matrix[i][j] += penalty; + } + + ReductionISToQUBO { + target: QUBO::from_matrix(matrix), + source_size: self.problem_size(), + } + } +} + +#[cfg(test)] +#[path = "../unit_tests/rules/independentset_qubo.rs"] +mod tests; +``` + +**Step 4: Register in `src/rules/mod.rs`** + +Add after `mod spinglass_qubo;`: +```rust +mod independentset_qubo; +``` + +Add after `pub use spinglass_qubo::...`: +```rust +pub use independentset_qubo::ReductionISToQUBO; +``` + +**Step 5: Run tests** + +Run: `cargo test --all-features test_independentset_to_qubo` +Expected: PASS + +**Step 6: Run clippy + full test suite** + +Run: `make test clippy` +Expected: all pass, no warnings + +**Step 7: Commit** + +```bash +git add src/rules/independentset_qubo.rs src/unit_tests/rules/independentset_qubo.rs src/rules/mod.rs +git commit -m "feat: add IndependentSet → QUBO reduction" +``` + +--- + +### Task 2: VertexCovering → QUBO + +Minimize `Σ w_i·x_i + P·Σ_{(i,j)∈E} (1-x_i)(1-x_j)`. Expanding: `Q[i][i] = w_i - P·deg(i)`, `Q[i][j] = P`. + +**Files:** +- Create: `src/rules/vertexcovering_qubo.rs` +- Create: `src/unit_tests/rules/vertexcovering_qubo.rs` +- Modify: `src/rules/mod.rs` + +Same pattern as Task 1. Key differences: +- VC minimizes (same as QUBO), so no sign flip on objective +- Penalty enforces: every edge has at least one endpoint selected +- `Q[i][i] = 1.0 - penalty * degree(i)`, `Q[i][j] = penalty` for edges +- Penalty `P = 1 + n` (unit weights) +- Test: cycle graph C4 (4 vertices, 4 edges) → min VC = 2 vertices + +**Step 1: Write test** (same structure as Task 1) +**Step 2: Verify fails** +**Step 3: Implement** — struct `ReductionVCToQUBO`, same boilerplate +**Step 4: Register in mod.rs** +**Step 5-6: Test + clippy** +**Step 7: Commit** `"feat: add VertexCovering → QUBO reduction"` + +--- + +### Task 3: MaxCut → QUBO + +Maximize cut = Σ_{(i,j)∈E} w_ij·(x_i⊕x_j). Minimize negative: `Q[i][i] = -Σ_j w_ij`, `Q[i][j] = 2·w_ij` (upper triangular). + +Note: MaxCut edges carry weights. Use `self.edges()` which returns `Vec<(usize, usize, W)>`. + +**Files:** +- Create: `src/rules/maxcut_qubo.rs` +- Create: `src/unit_tests/rules/maxcut_qubo.rs` +- Modify: `src/rules/mod.rs` + +Key: MaxCut is `MaxCut<SimpleGraph, W>` with edge weights. For unweighted, use `MaxCut::unweighted(n, edges)`. + +- `Q[i][j] = 2·w_ij` for i < j (upper triangular; the `w_ij(x_i + x_j - 2x_ix_j)` formula) +- `Q[i][i] = -Σ_{j:(i,j)∈E} w_ij` +- Test: cycle C4 → max cut = 4 (all edges cut by bipartition) +- No penalty needed — MaxCut is unconstrained + +**Step 1-7:** Same flow. Commit: `"feat: add MaxCut → QUBO reduction"` + +--- + +### Task 4: Coloring (KColoring) → QUBO + +One-hot encoding: `x_{v,c} = 1` iff vertex v gets color c. QUBO index: `v*K + c`. + +- One-hot penalty: `P₁·Σ_v (1 - Σ_c x_{v,c})²` +- Edge penalty: `P₂·Σ_{(u,v)∈E} Σ_c x_{u,c}·x_{v,c}` +- QUBO has `n·K` variables + +**Special:** `KColoring<const K: usize, G, W>` uses const generic. For the reduction, we implement for a specific K (e.g., `K=3`). Or better: implement for generic K using the existing pattern. + +Actually, looking at `coloring_ilp.rs`, there are two reductions: +- `ReductionColoringToILP` for `Coloring<SimpleGraph, W>` (deprecated Coloring type?) +- `ReductionKColoringToILP<const K: usize, W>` for `KColoring<K, SimpleGraph, W>` + +We should implement for `KColoring<K, SimpleGraph, Unweighted>`. The `extract_solution` decodes one-hot: for each vertex, find which color bit is 1. + +The struct needs to store `num_vertices` and `K` for extraction. + +**Files:** +- Create: `src/rules/coloring_qubo.rs` +- Create: `src/unit_tests/rules/coloring_qubo.rs` +- Modify: `src/rules/mod.rs` + +**Test:** Triangle K3, 3 colors → exactly 6 valid colorings (3! permutations). + +**Step 1-7:** Same flow. Commit: `"feat: add KColoring → QUBO reduction"` + +--- + +### Task 5: SetPacking → QUBO + +Same structure as IS on intersection graph: `Q[i][i] = -w_i`, `Q[i][j] = P` for overlapping pairs. + +Use `self.overlapping_pairs()` to get conflicting set pairs. + +**Files:** +- Create: `src/rules/setpacking_qubo.rs` +- Create: `src/unit_tests/rules/setpacking_qubo.rs` +- Modify: `src/rules/mod.rs` + +**Test:** 3 sets with some overlaps → verify max packing found. + +**Step 1-7:** Same flow. Commit: `"feat: add SetPacking → QUBO reduction"` + +--- + +### Task 6: KSatisfiability (K=2) → QUBO + +Max-2-SAT penalty formulation. Each clause contributes to Q based on literal signs. + +For clause `(l₁ ∨ l₂)` where `l = x` or `l = ¬x`: +- `(x_i ∨ x_j)`: penalty `(1-x_i)(1-x_j)` = `1 - x_i - x_j + x_ix_j` +- `(¬x_i ∨ x_j)`: penalty `x_i(1-x_j)` = `x_i - x_ix_j` +- `(x_i ∨ ¬x_j)`: penalty `(1-x_i)x_j` = `x_j - x_ix_j` +- `(¬x_i ∨ ¬x_j)`: penalty `x_ix_j` + +CNFClause uses 1-indexed signed integers: positive = variable, negative = negated. E.g., `[1, -2]` = `(x₁ ∨ ¬x₂)`. + +**Files:** +- Create: `src/rules/ksatisfiability_qubo.rs` +- Create: `src/unit_tests/rules/ksatisfiability_qubo.rs` +- Modify: `src/rules/mod.rs` + +**Test:** 3 vars, 4 clauses → verify all clauses satisfied by extracted solution. + +**Step 1-7:** Same flow. Commit: `"feat: add KSatisfiability(K=2) → QUBO reduction"` + +--- + +### Task 7: ILP (binary) → QUBO + +Binary ILP: `min c^T x s.t. Ax ≤ b`. Feature-gated behind `ilp`. + +Formulation: `Q[i][i] += c_i` (objective) + `P·Σ_k (Σ_j a_{kj}·x_j - b_k)²` (constraint penalties). + +Expanding the quadratic penalty for constraint k: +- `Q[i][j] += P·a_{ki}·a_{kj}` for i ≤ j +- `Q[i][i] += P·a_{ki}·(a_{ki} - 2·b_k)` (diagonal adjustment) + +ILP fields are public: `self.constraints`, `self.objective`, `self.sense`, `self.bounds`, `self.num_vars`. + +Only valid for binary ILP (all bounds = [0,1]). Should assert this. + +For Maximize objectives, negate the objective coefficients (QUBO minimizes). + +**Files:** +- Create: `src/rules/ilp_qubo.rs` (with `#[cfg(feature = "ilp")]`) +- Create: `src/unit_tests/rules/ilp_qubo.rs` +- Modify: `src/rules/mod.rs` + +**Test:** Binary ILP with 3 vars, 2 constraints → verify feasible optimal found. + +**Step 1-7:** Same flow. Commit: `"feat: add ILP (binary) → QUBO reduction"` + +--- + +### Task 8: Integration Tests + +Add integration tests in `tests/suites/reductions.rs` that load JSON ground truth from `tests/data/qubo/` and compare against Rust reductions. + +**Files:** +- Modify: `tests/suites/reductions.rs` + +For each reduction, add a module like: +```rust +mod is_qubo_reductions { + use super::*; + + #[test] + fn test_is_to_qubo_ground_truth() { + // Load JSON, create source problem, reduce, verify QUBO matrix and optimal match + } +} +``` + +**Commit:** `"test: add integration tests for QUBO reductions against ground truth"` + +--- + +### Task 9: Example Program + +Create `examples/qubo_reductions.rs` demonstrating all 7 reductions with practical stories. + +**File:** Create `examples/qubo_reductions.rs` + +Each demo: +1. Create a small practical instance (e.g., "Find the largest non-conflicting set of wireless towers") +2. Reduce to QUBO +3. Solve with BruteForce +4. Extract and explain the solution + +Run: `cargo run --example qubo_reductions --features ilp` + +**Commit:** `"docs: add QUBO reductions example program"` + +--- + +### Task 10: Paper Documentation + +Update `docs/paper/reductions.typ` with 7 new theorems. + +**File:** Modify `docs/paper/reductions.typ` + +For each reduction: +1. Add theorem in Section 3.1 (trivial reductions — these are standard penalty formulations) +2. Add proof with the QUBO formulation +3. Add Rust code example (from `examples/qubo_reductions.rs`) +4. Update summary table with overhead and reference + +Also update `@def:qubo` to list new "Reduces from" links. + +Run: `make export-graph && make paper` + +**Commit:** `"docs: add QUBO reduction theorems and examples to paper"` + +--- + +### Task 11: Final Verification + +```bash +make test # All tests pass +make clippy # No warnings +make export-graph # Reduction graph updated +make paper # Paper compiles +make coverage # >95% for new code +``` + +**Commit:** any final fixups + +--- + +## Key Reference Files + +| Purpose | Path | +|---------|------| +| Model pattern | `src/rules/spinglass_qubo.rs` | +| Test pattern | `src/unit_tests/rules/spinglass_qubo.rs` | +| Module registry | `src/rules/mod.rs` | +| Integration tests | `tests/suites/reductions.rs` | +| ILP feature gate | `src/rules/mod.rs:28-45` (example) | +| Ground truth JSON | `tests/data/qubo/*.json` | +| Paper | `docs/paper/reductions.typ` | +| IS model | `src/models/graph/independent_set.rs` | +| VC model | `src/models/graph/vertex_covering.rs` | +| MaxCut model | `src/models/graph/max_cut.rs` | +| KColoring model | `src/models/graph/kcoloring.rs` | +| SetPacking model | `src/models/set/set_packing.rs` | +| KSat model | `src/models/satisfiability/ksat.rs` | +| ILP model | `src/models/optimization/ilp.rs` | diff --git a/examples/qubo_reductions.rs b/examples/qubo_reductions.rs new file mode 100644 index 00000000..6d0103e9 --- /dev/null +++ b/examples/qubo_reductions.rs @@ -0,0 +1,216 @@ +//! Demonstrates 6 problem-to-QUBO reductions with practical stories. +//! +//! Run with: `cargo run --example qubo_reductions --features ilp` + +use problemreductions::prelude::*; +use problemreductions::topology::SimpleGraph; + +fn main() { + println!("=== Problem-to-QUBO Reductions ===\n"); + + demo_independent_set(); + demo_vertex_covering(); + demo_coloring(); + demo_set_packing(); + demo_ksat(); + demo_ilp(); +} + +/// Wireless tower placement: find the largest set of non-interfering towers. +fn demo_independent_set() { + println!("--- 1. IndependentSet -> QUBO ---"); + println!("Story: Place wireless towers on a 4-site grid. Adjacent towers interfere."); + println!(" Find the maximum set of non-interfering towers.\n"); + + // Path graph: sites 0-1-2-3, adjacent sites interfere + let is = IndependentSet::<SimpleGraph, i32>::new(4, vec![(0, 1), (1, 2), (2, 3)]); + let reduction = ReduceTo::<QUBO>::reduce_to(&is); + let qubo = reduction.target_problem(); + + let solver = BruteForce::new(); + let solutions = solver.find_best(qubo); + + println!(" QUBO variables: {}", qubo.num_variables()); + println!(" Optimal solutions:"); + for sol in &solutions { + let extracted = reduction.extract_solution(sol); + let sites: Vec<usize> = extracted + .iter() + .enumerate() + .filter(|(_, &x)| x == 1) + .map(|(i, _)| i) + .collect(); + println!(" Tower sites: {:?} (size {})", sites, sites.len()); + } + println!(); +} + +/// Security camera placement: cover all corridors with minimum cameras. +fn demo_vertex_covering() { + println!("--- 2. VertexCovering -> QUBO ---"); + println!("Story: Place security cameras at intersections to cover all corridors."); + println!(" Minimize the number of cameras needed.\n"); + + // Cycle C4: 4 intersections, 4 corridors + let vc = VertexCovering::<SimpleGraph, i32>::new(4, vec![(0, 1), (1, 2), (2, 3), (0, 3)]); + let reduction = ReduceTo::<QUBO>::reduce_to(&vc); + let qubo = reduction.target_problem(); + + let solver = BruteForce::new(); + let solutions = solver.find_best(qubo); + + println!(" QUBO variables: {}", qubo.num_variables()); + println!(" Optimal solutions:"); + for sol in &solutions { + let extracted = reduction.extract_solution(sol); + let cameras: Vec<usize> = extracted + .iter() + .enumerate() + .filter(|(_, &x)| x == 1) + .map(|(i, _)| i) + .collect(); + println!( + " Camera positions: {:?} ({} cameras)", + cameras, + cameras.len() + ); + } + println!(); +} + +/// Map coloring: color a triangle map with 3 colors so no neighbors share a color. +fn demo_coloring() { + println!("--- 3. KColoring -> QUBO ---"); + println!("Story: Color 3 countries on a map with 3 colors so no neighbors match.\n"); + + // Triangle K3: 3 countries, all share borders + let kc = KColoring::<3, SimpleGraph, i32>::new(3, vec![(0, 1), (1, 2), (0, 2)]); + let reduction = ReduceTo::<QUBO>::reduce_to(&kc); + let qubo = reduction.target_problem(); + + let solver = BruteForce::new(); + let solutions = solver.find_best(qubo); + + let colors = ["Red", "Green", "Blue"]; + println!( + " QUBO variables: {} (one-hot: 3 countries x 3 colors)", + qubo.num_variables() + ); + println!(" Valid colorings: {}", solutions.len()); + let extracted = reduction.extract_solution(&solutions[0]); + println!( + " Example: Country0={}, Country1={}, Country2={}", + colors[extracted[0]], colors[extracted[1]], colors[extracted[2]] + ); + println!(); +} + +/// Warehouse selection: pick maximum non-overlapping delivery zones. +fn demo_set_packing() { + println!("--- 4. SetPacking -> QUBO ---"); + println!("Story: Select delivery zones that don't overlap. Maximize coverage.\n"); + + // 3 zones covering different areas + let sp = SetPacking::<i32>::new(vec![ + vec![0, 1], // Zone A covers areas 0,1 + vec![1, 2], // Zone B covers areas 1,2 + vec![2, 3, 4], // Zone C covers areas 2,3,4 + ]); + let reduction = ReduceTo::<QUBO>::reduce_to(&sp); + let qubo = reduction.target_problem(); + + let solver = BruteForce::new(); + let solutions = solver.find_best(qubo); + + println!(" QUBO variables: {}", qubo.num_variables()); + println!(" Optimal solutions:"); + for sol in &solutions { + let extracted = reduction.extract_solution(sol); + let zones: Vec<&str> = extracted + .iter() + .enumerate() + .filter(|(_, &x)| x == 1) + .map(|(i, _)| ["Zone-A", "Zone-B", "Zone-C"][i]) + .collect(); + println!(" Selected: {:?}", zones); + } + println!(); +} + +/// Satisfiability: find a boolean assignment satisfying maximum clauses. +fn demo_ksat() { + println!("--- 5. KSatisfiability(K=2) -> QUBO ---"); + println!("Story: Configure 3 switches (on/off) to satisfy maximum rules.\n"); + + // 4 rules over 3 switches + let ksat = KSatisfiability::<2, i32>::new( + 3, + vec![ + CNFClause::new(vec![1, 2]), // switch1 OR switch2 + CNFClause::new(vec![-1, 3]), // NOT switch1 OR switch3 + CNFClause::new(vec![2, -3]), // switch2 OR NOT switch3 + CNFClause::new(vec![-2, -3]), // NOT switch2 OR NOT switch3 + ], + ); + let reduction = ReduceTo::<QUBO>::reduce_to(&ksat); + let qubo = reduction.target_problem(); + + let solver = BruteForce::new(); + let solutions = solver.find_best(qubo); + + println!(" QUBO variables: {}", qubo.num_variables()); + for sol in &solutions { + let extracted = reduction.extract_solution(sol); + let switches: Vec<&str> = extracted.iter().map(|&x| if x == 1 { "ON" } else { "OFF" }).collect(); + let satisfied = ksat.solution_size(&extracted).size; + println!( + " Switches: [{}] -> {}/{} rules satisfied", + switches.join(", "), + satisfied, + ksat.clauses().len() + ); + } + println!(); +} + +/// Resource allocation: maximize value under budget constraints. +fn demo_ilp() { + println!("--- 6. ILP (binary) -> QUBO ---"); + println!("Story: Select projects to maximize profit under resource constraints.\n"); + + // 3 projects: values 1, 2, 3 + // Constraint 1: projects 0 and 1 share a team (at most one) + // Constraint 2: projects 1 and 2 share equipment (at most one) + let ilp = ILP::binary( + 3, + vec![ + LinearConstraint::le(vec![(0, 1.0), (1, 1.0)], 1.0), + LinearConstraint::le(vec![(1, 1.0), (2, 1.0)], 1.0), + ], + vec![(0, 1.0), (1, 2.0), (2, 3.0)], + ObjectiveSense::Maximize, + ); + let reduction = ReduceTo::<QUBO>::reduce_to(&ilp); + let qubo = reduction.target_problem(); + + let solver = BruteForce::new(); + let solutions = solver.find_best(qubo); + + let names = ["Alpha", "Beta", "Gamma"]; + println!(" QUBO variables: {}", qubo.num_variables()); + for sol in &solutions { + let extracted = reduction.extract_solution(sol); + let selected: Vec<&str> = extracted + .iter() + .enumerate() + .filter(|(_, &x)| x == 1) + .map(|(i, _)| names[i]) + .collect(); + let value = ilp.solution_size(&extracted).size; + println!( + " Selected projects: {:?} (total value: {:.0})", + selected, value + ); + } + println!(); +} diff --git a/scripts/.python-version b/scripts/.python-version new file mode 100644 index 00000000..e4fba218 --- /dev/null +++ b/scripts/.python-version @@ -0,0 +1 @@ +3.12 diff --git a/scripts/generate_qubo_tests.py b/scripts/generate_qubo_tests.py new file mode 100644 index 00000000..6f5df890 --- /dev/null +++ b/scripts/generate_qubo_tests.py @@ -0,0 +1,240 @@ +"""Generate QUBO test datasets using qubogen. + +For each supported problem type, creates a small instance, reduces it to QUBO +via qubogen, brute-force solves both sides, and exports JSON ground truth +to tests/data/qubo/. + +Usage: + uv run python scripts/generate_qubo_tests.py +""" + +import json +from itertools import product +from pathlib import Path + +import numpy as np + +# Monkey-patch for qubogen compatibility with numpy >= 1.24 +np.float = np.float64 +np.int = np.int_ +np.bool = np.bool_ + +import qubogen + + +def brute_force_qubo(Q: np.ndarray) -> dict: + """Brute-force solve a QUBO: minimize x^T Q x over binary x.""" + n = Q.shape[0] + best_val = float("inf") + best_configs = [] + for bits in product(range(2), repeat=n): + x = np.array(bits, dtype=float) + val = float(x @ Q @ x) + if val < best_val - 1e-9: + best_val = val + best_configs = [list(bits)] + elif abs(val - best_val) < 1e-9: + best_configs.append(list(bits)) + return {"value": best_val, "configs": best_configs} + + +def save_test(name: str, data: dict, outdir: Path): + """Save test data as compact JSON.""" + path = outdir / f"{name}.json" + with open(path, "w") as f: + json.dump(data, f, separators=(",", ":")) + print(f" wrote {path} ({path.stat().st_size} bytes)") + + +def generate_vertex_covering(outdir: Path): + """Minimum Vertex Cover on a small graph (4 nodes, 5 edges).""" + edges = [(0, 1), (1, 2), (2, 3), (0, 3), (0, 2)] + n_nodes = 4 + penalty = 8.0 + g = qubogen.Graph(edges=np.array(edges), n_nodes=n_nodes) + Q = qubogen.qubo_mvc(g, penalty=penalty) + + qubo_result = brute_force_qubo(Q) + + save_test("vertexcovering_to_qubo", { + "problem": "VertexCovering", + "source": {"num_vertices": n_nodes, "edges": edges, "penalty": penalty}, + "qubo_matrix": Q.tolist(), + "qubo_num_vars": int(Q.shape[0]), + "qubo_optimal": qubo_result, + }, outdir) + + +def generate_independent_set(outdir: Path): + """Independent Set on a small graph. + + IndependentSet is the complement of VertexCover: maximize |S| s.t. no + two adjacent vertices are in S. We formulate as QUBO by negating the + linear terms of MVC (minimize -|S| + penalty * constraint violations). + """ + edges = [(0, 1), (1, 2), (2, 3), (0, 3)] + n_nodes = 4 + penalty = 8.0 + + # Independent set QUBO: maximize sum(x_i) s.t. x_i*x_j = 0 for edges + # = minimize -sum(x_i) + P * sum_{(i,j)} x_i*x_j + Q = np.zeros((n_nodes, n_nodes)) + for i in range(n_nodes): + Q[i][i] = -1.0 + for i, j in edges: + Q[i][j] += penalty + + qubo_result = brute_force_qubo(Q) + + save_test("independentset_to_qubo", { + "problem": "IndependentSet", + "source": {"num_vertices": n_nodes, "edges": edges, "penalty": penalty}, + "qubo_matrix": Q.tolist(), + "qubo_num_vars": int(Q.shape[0]), + "qubo_optimal": qubo_result, + }, outdir) + + +def generate_graph_coloring(outdir: Path): + """Graph Coloring on a small graph (3 nodes triangle, 3 colors).""" + edges = [(0, 1), (1, 2), (0, 2)] + n_nodes = 3 + n_color = 3 + penalty = 10.0 + g = qubogen.Graph(edges=np.array(edges), n_nodes=n_nodes) + Q = qubogen.qubo_graph_coloring(g, n_color=n_color, penalty=penalty) + + qubo_result = brute_force_qubo(Q) + + # QUBO variables: n_nodes * n_color (one-hot encoding) + save_test("coloring_to_qubo", { + "problem": "Coloring", + "source": { + "num_vertices": n_nodes, + "edges": edges, + "num_colors": n_color, + "penalty": penalty, + }, + "qubo_matrix": Q.tolist(), + "qubo_num_vars": int(Q.shape[0]), + "qubo_optimal": qubo_result, + }, outdir) + + +def generate_set_packing(outdir: Path): + """Set Packing: select maximum-weight non-overlapping sets.""" + # 3 sets over 4 elements + # set 0: {0, 2} + # set 1: {1, 2} + # set 2: {0, 3} + sets = [[0, 2], [1, 2], [0, 3]] + n_elements = 4 + n_sets = len(sets) + weights = [1.0, 2.0, 1.5] + penalty = 8.0 + + # Build incidence matrix (elements x sets) + a = np.zeros((n_elements, n_sets)) + for j, s in enumerate(sets): + for i in s: + a[i][j] = 1 + + Q = qubogen.qubo_set_pack(a, np.array(weights), penalty=penalty) + + qubo_result = brute_force_qubo(Q) + + save_test("setpacking_to_qubo", { + "problem": "SetPacking", + "source": { + "sets": sets, + "num_elements": n_elements, + "weights": weights, + "penalty": penalty, + }, + "qubo_matrix": Q.tolist(), + "qubo_num_vars": int(Q.shape[0]), + "qubo_optimal": qubo_result, + }, outdir) + + +def generate_max2sat(outdir: Path): + """Max 2-SAT: maximize satisfied clauses.""" + # 3 variables, 4 clauses: + # (x0 OR x1), (NOT x0 OR x2), (x1 OR NOT x2), (NOT x1 OR NOT x2) + literals = np.array([[0, 1], [0, 2], [1, 2], [1, 2]]) + signs = np.array( + [[True, True], [False, True], [True, False], [False, False]] + ) + + c = qubogen.Clauses(literals=literals, signs=signs) + Q = qubogen.qubo_max2sat(c) + + qubo_result = brute_force_qubo(Q) + + # Convert to list-of-clauses format matching our KSatisfiability model + clauses = [] + for i in range(len(literals)): + clause = [] + for j in range(2): + var = int(literals[i][j]) + negated = not bool(signs[i][j]) + clause.append({"variable": var, "negated": negated}) + clauses.append(clause) + + save_test("ksatisfiability_to_qubo", { + "problem": "KSatisfiability", + "source": {"num_variables": 3, "clauses": clauses}, + "qubo_matrix": Q.tolist(), + "qubo_num_vars": int(Q.shape[0]), + "qubo_optimal": qubo_result, + }, outdir) + + +def generate_ilp(outdir: Path): + """Binary ILP (General 0/1 Programming): min c^T x, s.t. Ax <= b.""" + # 3 variables + # minimize: x0 + 2*x1 + 3*x2 + # s.t.: x0 + x1 <= 1 + # x1 + x2 <= 1 + cost = np.array([1.0, 2.0, 3.0]) + A = np.array([[1.0, 1.0, 0.0], [0.0, 1.0, 1.0]]) + b = np.array([1.0, 1.0]) + sign = np.array([-1, -1]) # -1 means <= + penalty = 10.0 + + Q = qubogen.qubo_general01(cost, A, b, sign, penalty=penalty) + + qubo_result = brute_force_qubo(Q) + + save_test("ilp_to_qubo", { + "problem": "ILP", + "source": { + "num_variables": 3, + "objective": cost.tolist(), + "constraints_lhs": A.tolist(), + "constraints_rhs": b.tolist(), + "constraint_signs": sign.tolist(), + "penalty": penalty, + }, + "qubo_matrix": Q.tolist(), + "qubo_num_vars": int(Q.shape[0]), + "qubo_optimal": qubo_result, + }, outdir) + + +def main(): + outdir = Path(__file__).resolve().parent.parent / "tests" / "data" / "qubo" + outdir.mkdir(parents=True, exist_ok=True) + + print("Generating QUBO test datasets...") + generate_vertex_covering(outdir) + generate_independent_set(outdir) + generate_graph_coloring(outdir) + generate_set_packing(outdir) + generate_max2sat(outdir) + generate_ilp(outdir) + print("Done.") + + +if __name__ == "__main__": + main() diff --git a/scripts/pyproject.toml b/scripts/pyproject.toml new file mode 100644 index 00000000..a61d0e94 --- /dev/null +++ b/scripts/pyproject.toml @@ -0,0 +1,9 @@ +[project] +name = "scripts" +version = "0.1.0" +description = "Test data generation scripts for problem-reductions" +requires-python = ">=3.12" +dependencies = [ + "numpy>=1.26,<2", + "qubogen>=0.1.1", +] diff --git a/scripts/uv.lock b/scripts/uv.lock new file mode 100644 index 00000000..58b3004f --- /dev/null +++ b/scripts/uv.lock @@ -0,0 +1,55 @@ +version = 1 +revision = 2 +requires-python = ">=3.12" + +[[package]] +name = "networkx" +version = "3.6.1" +source = { registry = "https://pypi.org/simple" } +sdist = { url = "https://files.pythonhosted.org/packages/6a/51/63fe664f3908c97be9d2e4f1158eb633317598cfa6e1fc14af5383f17512/networkx-3.6.1.tar.gz", hash = "sha256:26b7c357accc0c8cde558ad486283728b65b6a95d85ee1cd66bafab4c8168509", size = 2517025, upload-time = "2025-12-08T17:02:39.908Z" } +wheels = [ + { url = "https://files.pythonhosted.org/packages/9e/c9/b2622292ea83fbb4ec318f5b9ab867d0a28ab43c5717bb85b0a5f6b3b0a4/networkx-3.6.1-py3-none-any.whl", hash = "sha256:d47fbf302e7d9cbbb9e2555a0d267983d2aa476bac30e90dfbe5669bd57f3762", size = 2068504, upload-time = "2025-12-08T17:02:38.159Z" }, +] + +[[package]] +name = "numpy" +version = "1.26.4" +source = { registry = "https://pypi.org/simple" } +sdist = { url = "https://files.pythonhosted.org/packages/65/6e/09db70a523a96d25e115e71cc56a6f9031e7b8cd166c1ac8438307c14058/numpy-1.26.4.tar.gz", hash = "sha256:2a02aba9ed12e4ac4eb3ea9421c420301a0c6460d9830d74a9df87efa4912010", size = 15786129, upload-time = "2024-02-06T00:26:44.495Z" } +wheels = [ + { url = "https://files.pythonhosted.org/packages/95/12/8f2020a8e8b8383ac0177dc9570aad031a3beb12e38847f7129bacd96228/numpy-1.26.4-cp312-cp312-macosx_10_9_x86_64.whl", hash = "sha256:b3ce300f3644fb06443ee2222c2201dd3a89ea6040541412b8fa189341847218", size = 20335901, upload-time = "2024-02-05T23:55:32.801Z" }, + { url = "https://files.pythonhosted.org/packages/75/5b/ca6c8bd14007e5ca171c7c03102d17b4f4e0ceb53957e8c44343a9546dcc/numpy-1.26.4-cp312-cp312-macosx_11_0_arm64.whl", hash = "sha256:03a8c78d01d9781b28a6989f6fa1bb2c4f2d51201cf99d3dd875df6fbd96b23b", size = 13685868, upload-time = "2024-02-05T23:55:56.28Z" }, + { url = "https://files.pythonhosted.org/packages/79/f8/97f10e6755e2a7d027ca783f63044d5b1bc1ae7acb12afe6a9b4286eac17/numpy-1.26.4-cp312-cp312-manylinux_2_17_aarch64.manylinux2014_aarch64.whl", hash = "sha256:9fad7dcb1aac3c7f0584a5a8133e3a43eeb2fe127f47e3632d43d677c66c102b", size = 13925109, upload-time = "2024-02-05T23:56:20.368Z" }, + { url = "https://files.pythonhosted.org/packages/0f/50/de23fde84e45f5c4fda2488c759b69990fd4512387a8632860f3ac9cd225/numpy-1.26.4-cp312-cp312-manylinux_2_17_x86_64.manylinux2014_x86_64.whl", hash = "sha256:675d61ffbfa78604709862923189bad94014bef562cc35cf61d3a07bba02a7ed", size = 17950613, upload-time = "2024-02-05T23:56:56.054Z" }, + { url = "https://files.pythonhosted.org/packages/4c/0c/9c603826b6465e82591e05ca230dfc13376da512b25ccd0894709b054ed0/numpy-1.26.4-cp312-cp312-musllinux_1_1_aarch64.whl", hash = "sha256:ab47dbe5cc8210f55aa58e4805fe224dac469cde56b9f731a4c098b91917159a", size = 13572172, upload-time = "2024-02-05T23:57:21.56Z" }, + { url = "https://files.pythonhosted.org/packages/76/8c/2ba3902e1a0fc1c74962ea9bb33a534bb05984ad7ff9515bf8d07527cadd/numpy-1.26.4-cp312-cp312-musllinux_1_1_x86_64.whl", hash = "sha256:1dda2e7b4ec9dd512f84935c5f126c8bd8b9f2fc001e9f54af255e8c5f16b0e0", size = 17786643, upload-time = "2024-02-05T23:57:56.585Z" }, + { url = "https://files.pythonhosted.org/packages/28/4a/46d9e65106879492374999e76eb85f87b15328e06bd1550668f79f7b18c6/numpy-1.26.4-cp312-cp312-win32.whl", hash = "sha256:50193e430acfc1346175fcbdaa28ffec49947a06918b7b92130744e81e640110", size = 5677803, upload-time = "2024-02-05T23:58:08.963Z" }, + { url = "https://files.pythonhosted.org/packages/16/2e/86f24451c2d530c88daf997cb8d6ac622c1d40d19f5a031ed68a4b73a374/numpy-1.26.4-cp312-cp312-win_amd64.whl", hash = "sha256:08beddf13648eb95f8d867350f6a018a4be2e5ad54c8d8caed89ebca558b2818", size = 15517754, upload-time = "2024-02-05T23:58:36.364Z" }, +] + +[[package]] +name = "qubogen" +version = "0.1.1" +source = { registry = "https://pypi.org/simple" } +dependencies = [ + { name = "networkx" }, +] +sdist = { url = "https://files.pythonhosted.org/packages/b1/c2/f8fad31da8e0029e867522ca94adb174ee4614e5c287a693b725961e4eff/qubogen-0.1.1.tar.gz", hash = "sha256:e54f4b454ded6deb292a3168ce0790bdabecc2fc464e7a225a585504fd87aa75", size = 4377, upload-time = "2019-02-27T14:14:47.914Z" } +wheels = [ + { url = "https://files.pythonhosted.org/packages/e1/d1/3b7177184c953a6a4a173d7838cc4291cd43a5b35dab27ca17277d8576f6/qubogen-0.1.1-py3-none-any.whl", hash = "sha256:f5a7ce03374c7e75977cf88c5e7bc6f0da4a4998ed7989893869728ec0e413fe", size = 5847, upload-time = "2019-02-27T14:14:46.446Z" }, +] + +[[package]] +name = "scripts" +version = "0.1.0" +source = { virtual = "." } +dependencies = [ + { name = "numpy" }, + { name = "qubogen" }, +] + +[package.metadata] +requires-dist = [ + { name = "numpy", specifier = ">=1.26,<2" }, + { name = "qubogen", specifier = ">=0.1.1" }, +] diff --git a/src/rules/coloring_qubo.rs b/src/rules/coloring_qubo.rs new file mode 100644 index 00000000..6e7f0952 --- /dev/null +++ b/src/rules/coloring_qubo.rs @@ -0,0 +1,120 @@ +//! Reduction from KColoring to QUBO. +//! +//! One-hot encoding: x_{v,c} = 1 iff vertex v gets color c. +//! QUBO variable index: v * K + c. +//! +//! One-hot penalty: P₁·Σ_v (1 - Σ_c x_{v,c})² +//! Edge penalty: P₂·Σ_{(u,v)∈E} Σ_c x_{u,c}·x_{v,c} +//! +//! QUBO has n·K variables. + +use crate::models::graph::KColoring; +use crate::models::optimization::QUBO; +use crate::poly; +use crate::reduction; +use crate::rules::registry::ReductionOverhead; +use crate::rules::traits::{ReduceTo, ReductionResult}; +use crate::topology::SimpleGraph; +use crate::traits::Problem; +use crate::types::ProblemSize; + +/// Result of reducing KColoring to QUBO. +#[derive(Debug, Clone)] +pub struct ReductionKColoringToQUBO<const K: usize> { + target: QUBO<f64>, + source_size: ProblemSize, + num_vertices: usize, +} + +impl<const K: usize> ReductionResult for ReductionKColoringToQUBO<K> { + type Source = KColoring<K, SimpleGraph, i32>; + type Target = QUBO<f64>; + + fn target_problem(&self) -> &Self::Target { + &self.target + } + + /// Decode one-hot: for each vertex, find which color bit is 1. + fn extract_solution(&self, target_solution: &[usize]) -> Vec<usize> { + (0..self.num_vertices) + .map(|v| { + (0..K) + .find(|&c| target_solution[v * K + c] == 1) + .unwrap_or(0) + }) + .collect() + } + + fn source_size(&self) -> ProblemSize { + self.source_size.clone() + } + + fn target_size(&self) -> ProblemSize { + self.target.problem_size() + } +} + +#[reduction( + source_graph = "SimpleGraph", + overhead = { ReductionOverhead::new(vec![("num_vars", poly!(num_vertices))]) } +)] +impl<const K: usize> ReduceTo<QUBO<f64>> for KColoring<K, SimpleGraph, i32> { + type Result = ReductionKColoringToQUBO<K>; + + fn reduce_to(&self) -> Self::Result { + let n = self.num_vertices(); + let edges = self.edges(); + let nq = n * K; + + // Penalty must be large enough to enforce one-hot constraints + // P1 for one-hot, P2 for edge conflicts; use same penalty + let penalty = 1.0 + n as f64; + + let mut matrix = vec![vec![0.0; nq]; nq]; + + // One-hot penalty: P₁·Σ_v (1 - Σ_c x_{v,c})² + // Expanding: (1 - Σ_c x_{v,c})² = 1 - 2·Σ_c x_{v,c} + (Σ_c x_{v,c})² + // = 1 - 2·Σ_c x_{v,c} + Σ_c x_{v,c}² + 2·Σ_{c<c'} x_{v,c}·x_{v,c'} + // Since x² = x for binary: = 1 - Σ_c x_{v,c} + 2·Σ_{c<c'} x_{v,c}·x_{v,c'} + for v in 0..n { + for c in 0..K { + let idx = v * K + c; + // Diagonal: -P₁ (from the linear term -Σ_c x_{v,c}) + matrix[idx][idx] -= penalty; + } + // Off-diagonal within same vertex: 2·P₁ for each pair of colors + for c1 in 0..K { + for c2 in (c1 + 1)..K { + let idx1 = v * K + c1; + let idx2 = v * K + c2; + matrix[idx1][idx2] += 2.0 * penalty; + } + } + } + + // Edge penalty: P₂·Σ_{(u,v)∈E} Σ_c x_{u,c}·x_{v,c} + let edge_penalty = penalty / 2.0; + for (u, v) in &edges { + for c in 0..K { + let idx_u = u * K + c; + let idx_v = v * K + c; + let (i, j) = if idx_u < idx_v { + (idx_u, idx_v) + } else { + (idx_v, idx_u) + }; + matrix[i][j] += edge_penalty; + } + } + + ReductionKColoringToQUBO { + target: QUBO::from_matrix(matrix), + source_size: self.problem_size(), + num_vertices: n, + } + } +} + +#[cfg(test)] +#[path = "../unit_tests/rules/coloring_qubo.rs"] +mod tests; diff --git a/src/rules/ilp_qubo.rs b/src/rules/ilp_qubo.rs new file mode 100644 index 00000000..bfb77a3e --- /dev/null +++ b/src/rules/ilp_qubo.rs @@ -0,0 +1,194 @@ +//! Reduction from binary ILP to QUBO. +//! +//! Binary ILP: optimize c^T x s.t. Ax {<=,>=,=} b, x ∈ {0,1}^n. +//! +//! Formulation (following qubogen): +//! 1. Normalize constraints to Ax = b by adding slack variables +//! 2. QUBO = -diag(c + 2·P·b·A) + P·A^T·A +//! +//! For Minimize sense, c is negated (convert to maximization). +//! Slack variables: ceil(log2(slack_range)) bits per inequality constraint. + +use crate::models::optimization::{Comparison, ObjectiveSense, ILP, QUBO}; +use crate::poly; +use crate::reduction; +use crate::rules::registry::ReductionOverhead; +use crate::rules::traits::{ReduceTo, ReductionResult}; +use crate::traits::Problem; +use crate::types::ProblemSize; + +/// Result of reducing binary ILP to QUBO. +#[derive(Debug, Clone)] +pub struct ReductionILPToQUBO { + target: QUBO<f64>, + source_size: ProblemSize, + num_original_vars: usize, +} + +impl ReductionResult for ReductionILPToQUBO { + type Source = ILP; + type Target = QUBO<f64>; + + fn target_problem(&self) -> &Self::Target { + &self.target + } + + /// Extract only the original variables (discard slack). + fn extract_solution(&self, target_solution: &[usize]) -> Vec<usize> { + target_solution[..self.num_original_vars].to_vec() + } + + fn source_size(&self) -> ProblemSize { + self.source_size.clone() + } + + fn target_size(&self) -> ProblemSize { + self.target.problem_size() + } +} + +#[reduction( + overhead = { ReductionOverhead::new(vec![("num_vars", poly!(num_vars))]) } +)] +impl ReduceTo<QUBO<f64>> for ILP { + type Result = ReductionILPToQUBO; + + fn reduce_to(&self) -> Self::Result { + let n = self.num_vars; + + // Verify all variables are binary + for (i, b) in self.bounds.iter().enumerate() { + assert!( + b.lower == Some(0) && b.upper == Some(1), + "ILP→QUBO requires binary variables (var {} has bounds {:?})", + i, + b + ); + } + + // Build dense constraint matrix A and rhs vector b + // Also compute slack sizes for inequality constraints + let num_constraints = self.constraints.len(); + let mut a_dense = vec![vec![0.0; n]; num_constraints]; + let mut b_vec = vec![0.0; num_constraints]; + let mut slack_sizes = vec![0usize; num_constraints]; + + for (k, constraint) in self.constraints.iter().enumerate() { + for &(var, coef) in &constraint.terms { + a_dense[k][var] += coef; + } + b_vec[k] = constraint.rhs; + + // Compute slack variable count: ceil(log2(slack_range + 1)) bits + // to represent integer values 0..slack_range with binary encoding. + // For binary variables, min_lhs = Σ min(0, a_i), max_lhs = Σ max(0, a_i). + match constraint.cmp { + Comparison::Eq => {} // no slack needed + Comparison::Le => { + // Ax <= b → Ax + s = b, s ∈ {0, ..., b - min_lhs} + let min_lhs: f64 = a_dense[k].iter().map(|&c| c.min(0.0)).sum(); + let slack_range = constraint.rhs - min_lhs; + if slack_range > 0.0 { + slack_sizes[k] = (slack_range + 1.0).log2().ceil() as usize; + } + } + Comparison::Ge => { + // Ax >= b → Ax - s = b, s ∈ {0, ..., max_lhs - b} + let max_lhs: f64 = a_dense[k].iter().map(|&c| c.max(0.0)).sum(); + let slack_range = max_lhs - constraint.rhs; + if slack_range > 0.0 { + slack_sizes[k] = (slack_range + 1.0).log2().ceil() as usize; + } + } + } + } + + let total_slack: usize = slack_sizes.iter().sum(); + let nq = n + total_slack; + + // Extend A with slack columns + let mut a_ext = vec![vec![0.0; nq]; num_constraints]; + for k in 0..num_constraints { + for j in 0..n { + a_ext[k][j] = a_dense[k][j]; + } + } + + // Add slack variable columns + let mut slack_col = n; + for (k, &ns) in slack_sizes.iter().enumerate() { + if ns > 0 { + let sign = match self.constraints[k].cmp { + Comparison::Le => 1.0, // Ax + s = b + Comparison::Ge => -1.0, // Ax - s = b + Comparison::Eq => 0.0, + }; + for s in 0..ns { + a_ext[k][slack_col + s] = sign * 2.0_f64.powi(s as i32); + } + slack_col += ns; + } + } + + // Build dense cost vector (nq elements) + let mut c_vec = vec![0.0; nq]; + for &(var, coef) in &self.objective { + c_vec[var] = coef; + } + + // For Minimize sense, negate the cost (formula assumes maximization) + if self.sense == ObjectiveSense::Minimize { + for c in c_vec.iter_mut() { + *c = -*c; + } + } + + // Penalty: must be large enough to enforce constraints + let penalty = 1.0 + + c_vec.iter().map(|c| c.abs()).sum::<f64>() + + b_vec.iter().map(|b| b.abs()).sum::<f64>(); + + // QUBO = -diag(c + 2·P·b·A) + P·A^T·A + let mut matrix = vec![vec![0.0; nq]; nq]; + + // Compute b·A (b_vec dot each column of a_ext) + let mut ba = vec![0.0; nq]; + for (j, ba_j) in ba.iter_mut().enumerate() { + for (k, &b_k) in b_vec.iter().enumerate() { + *ba_j += b_k * a_ext[k][j]; + } + } + + // Diagonal: -(c_j + 2·P·(b·A)_j) + for j in 0..nq { + matrix[j][j] = -(c_vec[j] + 2.0 * penalty * ba[j]); + } + + // A^T·A contribution (upper-triangular convention) + // Diagonal: P · Σ_k a_{ki}² + // Off-diagonal (i<j): 2·P · Σ_k a_{ki}·a_{kj} + for row in &a_ext { + for (i, row_i) in matrix.iter_mut().enumerate() { + if row[i].abs() < 1e-15 { + continue; + } + // Diagonal + row_i[i] += penalty * row[i] * row[i]; + // Off-diagonal + for j in (i + 1)..nq { + row_i[j] += 2.0 * penalty * row[i] * row[j]; + } + } + } + + ReductionILPToQUBO { + target: QUBO::from_matrix(matrix), + source_size: self.problem_size(), + num_original_vars: n, + } + } +} + +#[cfg(test)] +#[path = "../unit_tests/rules/ilp_qubo.rs"] +mod tests; diff --git a/src/rules/independentset_qubo.rs b/src/rules/independentset_qubo.rs new file mode 100644 index 00000000..ea6e3293 --- /dev/null +++ b/src/rules/independentset_qubo.rs @@ -0,0 +1,78 @@ +//! Reduction from IndependentSet to QUBO. +//! +//! Maximize Σ w_i·x_i s.t. x_i·x_j = 0 for (i,j) ∈ E +//! = Minimize -Σ w_i·x_i + P·Σ_{(i,j)∈E} x_i·x_j +//! +//! Q[i][i] = -w_i, Q[i][j] = P for edges. P = 1 + Σ w_i. + +use crate::models::graph::IndependentSet; +use crate::models::optimization::QUBO; +use crate::poly; +use crate::reduction; +use crate::rules::registry::ReductionOverhead; +use crate::rules::traits::{ReduceTo, ReductionResult}; +use crate::topology::SimpleGraph; +use crate::traits::Problem; +use crate::types::ProblemSize; + +/// Result of reducing IndependentSet to QUBO. +#[derive(Debug, Clone)] +pub struct ReductionISToQUBO { + target: QUBO<f64>, + source_size: ProblemSize, +} + +impl ReductionResult for ReductionISToQUBO { + type Source = IndependentSet<SimpleGraph, i32>; + type Target = QUBO<f64>; + + fn target_problem(&self) -> &Self::Target { + &self.target + } + + fn extract_solution(&self, target_solution: &[usize]) -> Vec<usize> { + target_solution.to_vec() + } + + fn source_size(&self) -> ProblemSize { + self.source_size.clone() + } + + fn target_size(&self) -> ProblemSize { + self.target.problem_size() + } +} + +#[reduction( + source_graph = "SimpleGraph", + overhead = { ReductionOverhead::new(vec![("num_vars", poly!(num_vertices))]) } +)] +impl ReduceTo<QUBO<f64>> for IndependentSet<SimpleGraph, i32> { + type Result = ReductionISToQUBO; + + fn reduce_to(&self) -> Self::Result { + let n = self.num_vertices(); + let edges = self.edges(); + let weights = self.weights_ref(); + let total_weight: f64 = weights.iter().map(|&w| w as f64).sum(); + let penalty = 1.0 + total_weight; + + let mut matrix = vec![vec![0.0; n]; n]; + for i in 0..n { + matrix[i][i] = -(weights[i] as f64); + } + for (u, v) in &edges { + let (i, j) = if u < v { (*u, *v) } else { (*v, *u) }; + matrix[i][j] += penalty; + } + + ReductionISToQUBO { + target: QUBO::from_matrix(matrix), + source_size: self.problem_size(), + } + } +} + +#[cfg(test)] +#[path = "../unit_tests/rules/independentset_qubo.rs"] +mod tests; diff --git a/src/rules/ksatisfiability_qubo.rs b/src/rules/ksatisfiability_qubo.rs new file mode 100644 index 00000000..5667b73b --- /dev/null +++ b/src/rules/ksatisfiability_qubo.rs @@ -0,0 +1,109 @@ +//! Reduction from KSatisfiability (K=2) to QUBO (Max-2-SAT). +//! +//! Each clause contributes to Q based on literal signs: +//! - (x_i ∨ x_j): penalty (1-x_i)(1-x_j) → Q[i][i]-=1, Q[j][j]-=1, Q[i][j]+=1, const+=1 +//! - (¬x_i ∨ x_j): penalty x_i(1-x_j) → Q[i][i]+=1, Q[i][j]-=1 +//! - (x_i ∨ ¬x_j): penalty (1-x_i)x_j → Q[j][j]+=1, Q[i][j]-=1 +//! - (¬x_i ∨ ¬x_j): penalty x_i·x_j → Q[i][j]+=1 +//! +//! CNFClause uses 1-indexed signed integers: positive = variable, negative = negated. + +use crate::models::optimization::QUBO; +use crate::models::satisfiability::KSatisfiability; +use crate::poly; +use crate::reduction; +use crate::rules::registry::ReductionOverhead; +use crate::rules::traits::{ReduceTo, ReductionResult}; +use crate::traits::Problem; +use crate::types::ProblemSize; + +/// Result of reducing KSatisfiability<2> to QUBO. +#[derive(Debug, Clone)] +pub struct ReductionKSatToQUBO { + target: QUBO<f64>, + source_size: ProblemSize, +} + +impl ReductionResult for ReductionKSatToQUBO { + type Source = KSatisfiability<2, i32>; + type Target = QUBO<f64>; + + fn target_problem(&self) -> &Self::Target { + &self.target + } + + fn extract_solution(&self, target_solution: &[usize]) -> Vec<usize> { + target_solution.to_vec() + } + + fn source_size(&self) -> ProblemSize { + self.source_size.clone() + } + + fn target_size(&self) -> ProblemSize { + self.target.problem_size() + } +} + +#[reduction( + overhead = { ReductionOverhead::new(vec![("num_vars", poly!(num_vars))]) } +)] +impl ReduceTo<QUBO<f64>> for KSatisfiability<2, i32> { + type Result = ReductionKSatToQUBO; + + fn reduce_to(&self) -> Self::Result { + let n = self.num_vars(); + let mut matrix = vec![vec![0.0; n]; n]; + + for clause in self.clauses() { + let lits = &clause.literals; + assert_eq!(lits.len(), 2, "Expected 2-SAT clause"); + + // Convert 1-indexed signed literals to 0-indexed (var, negated) pairs + let var_i = (lits[0].unsigned_abs() as usize) - 1; + let neg_i = lits[0] < 0; + let var_j = (lits[1].unsigned_abs() as usize) - 1; + let neg_j = lits[1] < 0; + + let (i, j, ni, nj) = if var_i <= var_j { + (var_i, var_j, neg_i, neg_j) + } else { + (var_j, var_i, neg_j, neg_i) + }; + + // Penalty for unsatisfied clause: minimize penalty + // (l_i ∨ l_j) unsatisfied when both literals false + match (ni, nj) { + (false, false) => { + // (x_i ∨ x_j): penalty = (1-x_i)(1-x_j) = 1 - x_i - x_j + x_i·x_j + matrix[i][i] -= 1.0; + matrix[j][j] -= 1.0; + matrix[i][j] += 1.0; + } + (true, false) => { + // (¬x_i ∨ x_j): penalty = x_i(1-x_j) = x_i - x_i·x_j + matrix[i][i] += 1.0; + matrix[i][j] -= 1.0; + } + (false, true) => { + // (x_i ∨ ¬x_j): penalty = (1-x_i)x_j = x_j - x_i·x_j + matrix[j][j] += 1.0; + matrix[i][j] -= 1.0; + } + (true, true) => { + // (¬x_i ∨ ¬x_j): penalty = x_i·x_j + matrix[i][j] += 1.0; + } + } + } + + ReductionKSatToQUBO { + target: QUBO::from_matrix(matrix), + source_size: self.problem_size(), + } + } +} + +#[cfg(test)] +#[path = "../unit_tests/rules/ksatisfiability_qubo.rs"] +mod tests; diff --git a/src/rules/mod.rs b/src/rules/mod.rs index f4518914..137a7c2e 100644 --- a/src/rules/mod.rs +++ b/src/rules/mod.rs @@ -9,18 +9,23 @@ pub use cost::{ pub use registry::{ReductionEntry, ReductionOverhead}; mod circuit_spinglass; +mod coloring_qubo; mod factoring_circuit; mod graph; +mod independentset_qubo; mod independentset_setpacking; +mod ksatisfiability_qubo; mod matching_setpacking; mod sat_coloring; mod sat_dominatingset; mod sat_independentset; mod sat_ksat; +mod setpacking_qubo; mod spinglass_maxcut; mod spinglass_qubo; mod traits; mod vertexcovering_independentset; +mod vertexcovering_qubo; mod vertexcovering_setcovering; pub mod unitdiskmapping; @@ -28,6 +33,8 @@ pub mod unitdiskmapping; #[cfg(feature = "ilp")] mod clique_ilp; #[cfg(feature = "ilp")] +mod ilp_qubo; +#[cfg(feature = "ilp")] mod coloring_ilp; #[cfg(feature = "ilp")] mod dominatingset_ilp; @@ -52,9 +59,13 @@ pub use factoring_circuit::ReductionFactoringToCircuit; pub use graph::{ EdgeJson, NodeJson, ReductionEdge, ReductionGraph, ReductionGraphJson, ReductionPath, }; +pub use coloring_qubo::ReductionKColoringToQUBO; +pub use independentset_qubo::ReductionISToQUBO; pub use independentset_setpacking::{ReductionISToSP, ReductionSPToIS}; +pub use ksatisfiability_qubo::ReductionKSatToQUBO; pub use matching_setpacking::ReductionMatchingToSP; pub use sat_coloring::ReductionSATToColoring; +pub use setpacking_qubo::ReductionSPToQUBO; pub use sat_dominatingset::ReductionSATToDS; pub use sat_independentset::{BoolVar, ReductionSATToIS}; pub use sat_ksat::{ReductionKSATToSAT, ReductionSATToKSAT}; @@ -62,6 +73,7 @@ pub use spinglass_maxcut::{ReductionMaxCutToSG, ReductionSGToMaxCut}; pub use spinglass_qubo::{ReductionQUBOToSG, ReductionSGToQUBO}; pub use traits::{ReduceTo, ReductionResult}; pub use vertexcovering_independentset::{ReductionISToVC, ReductionVCToIS}; +pub use vertexcovering_qubo::ReductionVCToQUBO; pub use vertexcovering_setcovering::ReductionVCToSC; #[cfg(feature = "ilp")] @@ -71,6 +83,8 @@ pub use coloring_ilp::{ReductionColoringToILP, ReductionKColoringToILP}; #[cfg(feature = "ilp")] pub use dominatingset_ilp::ReductionDSToILP; #[cfg(feature = "ilp")] +pub use ilp_qubo::ReductionILPToQUBO; +#[cfg(feature = "ilp")] pub use factoring_ilp::ReductionFactoringToILP; #[cfg(feature = "ilp")] pub use independentset_ilp::ReductionISToILP; diff --git a/src/rules/setpacking_qubo.rs b/src/rules/setpacking_qubo.rs new file mode 100644 index 00000000..e3924e5c --- /dev/null +++ b/src/rules/setpacking_qubo.rs @@ -0,0 +1,86 @@ +//! Reduction from SetPacking to QUBO. +//! +//! Same structure as IndependentSet on the intersection graph: +//! Maximize Σ w_i·x_i s.t. x_i·x_j = 0 for overlapping pairs (i,j). +//! = Minimize -Σ w_i·x_i + P·Σ_{overlapping (i,j)} x_i·x_j +//! +//! Q[i][i] = -w_i, Q[i][j] = P for overlapping pairs. P = 1 + Σ w_i. + +use crate::models::optimization::QUBO; +use crate::models::set::SetPacking; +use crate::poly; +use crate::reduction; +use crate::rules::registry::ReductionOverhead; +use crate::rules::traits::{ReduceTo, ReductionResult}; +use crate::traits::Problem; +use crate::types::{NumericWeight, ProblemSize}; + +use std::marker::PhantomData; + +/// Result of reducing SetPacking to QUBO. +#[derive(Debug, Clone)] +pub struct ReductionSPToQUBO<W> { + target: QUBO<f64>, + source_size: ProblemSize, + _phantom: PhantomData<W>, +} + +impl<W: NumericWeight + Into<f64>> ReductionResult for ReductionSPToQUBO<W> { + type Source = SetPacking<W>; + type Target = QUBO<f64>; + + fn target_problem(&self) -> &Self::Target { + &self.target + } + + fn extract_solution(&self, target_solution: &[usize]) -> Vec<usize> { + target_solution.to_vec() + } + + fn source_size(&self) -> ProblemSize { + self.source_size.clone() + } + + fn target_size(&self) -> ProblemSize { + self.target.problem_size() + } +} + +#[reduction( + source_weighted = true, + overhead = { ReductionOverhead::new(vec![("num_vars", poly!(num_sets))]) } +)] +impl<W: NumericWeight + Into<f64>> ReduceTo<QUBO<f64>> for SetPacking<W> { + type Result = ReductionSPToQUBO<W>; + + fn reduce_to(&self) -> Self::Result { + let n = self.num_sets(); + let weights = self.weights_ref(); + let total_weight: f64 = weights.iter().map(|w| w.clone().into()).sum(); + let penalty = 1.0 + total_weight; + + let mut matrix = vec![vec![0.0; n]; n]; + + // Diagonal: -w_i + for i in 0..n { + let w: f64 = weights[i].clone().into(); + matrix[i][i] = -w; + } + + // Off-diagonal: P for overlapping pairs + for (i, j) in self.overlapping_pairs() { + let (a, b) = if i < j { (i, j) } else { (j, i) }; + matrix[a][b] += penalty; + } + + ReductionSPToQUBO { + target: QUBO::from_matrix(matrix), + source_size: self.problem_size(), + _phantom: PhantomData, + } + } +} + +#[cfg(test)] +#[path = "../unit_tests/rules/setpacking_qubo.rs"] +mod tests; diff --git a/src/rules/vertexcovering_qubo.rs b/src/rules/vertexcovering_qubo.rs new file mode 100644 index 00000000..58824680 --- /dev/null +++ b/src/rules/vertexcovering_qubo.rs @@ -0,0 +1,90 @@ +//! Reduction from VertexCovering to QUBO. +//! +//! Minimize Σ w_i·x_i s.t. x_i + x_j ≥ 1 for (i,j) ∈ E +//! = Minimize Σ w_i·x_i + P·Σ_{(i,j)∈E} (1-x_i)(1-x_j) +//! +//! Expanding: Q[i][i] = w_i - P·deg(i), Q[i][j] = P for edges. +//! P = 1 + Σ w_i. + +use crate::models::graph::VertexCovering; +use crate::models::optimization::QUBO; +use crate::poly; +use crate::reduction; +use crate::rules::registry::ReductionOverhead; +use crate::rules::traits::{ReduceTo, ReductionResult}; +use crate::topology::SimpleGraph; +use crate::traits::Problem; +use crate::types::ProblemSize; + +/// Result of reducing VertexCovering to QUBO. +#[derive(Debug, Clone)] +pub struct ReductionVCToQUBO { + target: QUBO<f64>, + source_size: ProblemSize, +} + +impl ReductionResult for ReductionVCToQUBO { + type Source = VertexCovering<SimpleGraph, i32>; + type Target = QUBO<f64>; + + fn target_problem(&self) -> &Self::Target { + &self.target + } + + fn extract_solution(&self, target_solution: &[usize]) -> Vec<usize> { + target_solution.to_vec() + } + + fn source_size(&self) -> ProblemSize { + self.source_size.clone() + } + + fn target_size(&self) -> ProblemSize { + self.target.problem_size() + } +} + +#[reduction( + source_graph = "SimpleGraph", + overhead = { ReductionOverhead::new(vec![("num_vars", poly!(num_vertices))]) } +)] +impl ReduceTo<QUBO<f64>> for VertexCovering<SimpleGraph, i32> { + type Result = ReductionVCToQUBO; + + fn reduce_to(&self) -> Self::Result { + let n = self.num_vertices(); + let edges = self.edges(); + let weights = self.weights_ref(); + let total_weight: f64 = weights.iter().map(|&w| w as f64).sum(); + let penalty = 1.0 + total_weight; + + let mut matrix = vec![vec![0.0; n]; n]; + + // Compute degree of each vertex + let mut degree = vec![0usize; n]; + for (u, v) in &edges { + degree[*u] += 1; + degree[*v] += 1; + } + + // Diagonal: w_i - P * deg(i) + for i in 0..n { + matrix[i][i] = weights[i] as f64 - penalty * degree[i] as f64; + } + + // Off-diagonal: P for each edge + for (u, v) in &edges { + let (i, j) = if u < v { (*u, *v) } else { (*v, *u) }; + matrix[i][j] += penalty; + } + + ReductionVCToQUBO { + target: QUBO::from_matrix(matrix), + source_size: self.problem_size(), + } + } +} + +#[cfg(test)] +#[path = "../unit_tests/rules/vertexcovering_qubo.rs"] +mod tests; diff --git a/src/unit_tests/rules/coloring_qubo.rs b/src/unit_tests/rules/coloring_qubo.rs new file mode 100644 index 00000000..dfe2b089 --- /dev/null +++ b/src/unit_tests/rules/coloring_qubo.rs @@ -0,0 +1,75 @@ +use super::*; +use crate::solvers::{BruteForce, Solver}; + +#[test] +fn test_kcoloring_to_qubo_closed_loop() { + // Triangle K3, 3 colors → exactly 6 valid colorings (3! permutations) + let kc = KColoring::<3, SimpleGraph, i32>::new(3, vec![(0, 1), (1, 2), (0, 2)]); + let reduction = ReduceTo::<QUBO<f64>>::reduce_to(&kc); + let qubo = reduction.target_problem(); + + let solver = BruteForce::new(); + let qubo_solutions = solver.find_best(qubo); + + // All solutions should extract to valid colorings + for sol in &qubo_solutions { + let extracted = reduction.extract_solution(sol); + assert!(kc.solution_size(&extracted).is_valid); + } + + // Exactly 6 valid 3-colorings of K3 + assert_eq!(qubo_solutions.len(), 6); +} + +#[test] +fn test_kcoloring_to_qubo_path() { + // Path graph: 0-1-2, 2 colors + let kc = KColoring::<2, SimpleGraph, i32>::new(3, vec![(0, 1), (1, 2)]); + let reduction = ReduceTo::<QUBO<f64>>::reduce_to(&kc); + let qubo = reduction.target_problem(); + + let solver = BruteForce::new(); + let qubo_solutions = solver.find_best(qubo); + + for sol in &qubo_solutions { + let extracted = reduction.extract_solution(sol); + assert!(kc.solution_size(&extracted).is_valid); + } + + // 2-coloring of path: 0,1,0 or 1,0,1 → 2 solutions + assert_eq!(qubo_solutions.len(), 2); +} + +#[test] +fn test_kcoloring_to_qubo_reversed_edges() { + // Edge (2, 0) triggers the idx_v < idx_u swap branch (line 104). + // Path: 2-0-1 with reversed edge ordering + let kc = KColoring::<2, SimpleGraph, i32>::new(3, vec![(2, 0), (0, 1)]); + let reduction = ReduceTo::<QUBO<f64>>::reduce_to(&kc); + let qubo = reduction.target_problem(); + + let solver = BruteForce::new(); + let qubo_solutions = solver.find_best(qubo); + + for sol in &qubo_solutions { + let extracted = reduction.extract_solution(sol); + assert!(kc.solution_size(&extracted).is_valid); + } + + // Same as path graph: 2 valid 2-colorings + assert_eq!(qubo_solutions.len(), 2); +} + +#[test] +fn test_kcoloring_to_qubo_sizes() { + let kc = KColoring::<3, SimpleGraph, i32>::new(3, vec![(0, 1), (1, 2), (0, 2)]); + let reduction = ReduceTo::<QUBO<f64>>::reduce_to(&kc); + + let source_size = reduction.source_size(); + let target_size = reduction.target_size(); + assert!(!source_size.components.is_empty()); + assert!(!target_size.components.is_empty()); + + // QUBO should have n*K = 3*3 = 9 variables + assert_eq!(reduction.target_problem().num_variables(), 9); +} diff --git a/src/unit_tests/rules/graph.rs b/src/unit_tests/rules/graph.rs index 938851fa..fd3ebe34 100644 --- a/src/unit_tests/rules/graph.rs +++ b/src/unit_tests/rules/graph.rs @@ -38,12 +38,13 @@ fn test_has_direct_reduction() { } #[test] -fn test_no_path() { +fn test_is_to_qubo_path() { let graph = ReductionGraph::new(); - // No path between IndependentSet and QUBO (disconnected in graph topology) - let paths = - graph.find_paths::<IndependentSet<SimpleGraph, i32>, crate::models::optimization::QUBO<f64>>(); - assert!(paths.is_empty()); + // IS -> QUBO should now have a direct path + let path = + graph.find_shortest_path::<IndependentSet<SimpleGraph, i32>, crate::models::optimization::QUBO<f64>>(); + assert!(path.is_some()); + assert_eq!(path.unwrap().len(), 1); // Direct path } #[test] @@ -70,13 +71,13 @@ fn test_type_erased_paths() { fn test_find_paths_by_name() { let graph = ReductionGraph::new(); - let paths = graph.find_paths_by_name("MaxCut", "SpinGlass"); - assert!(!paths.is_empty()); - assert_eq!(paths[0].len(), 1); // Direct path + let shortest = graph.find_shortest_path_by_name("MaxCut", "SpinGlass"); + assert!(shortest.is_some()); + assert_eq!(shortest.unwrap().len(), 1); // Direct path - let paths = graph.find_paths_by_name("Factoring", "SpinGlass"); - assert!(!paths.is_empty()); - assert_eq!(paths[0].len(), 2); // Factoring -> CircuitSAT -> SpinGlass + let shortest = graph.find_shortest_path_by_name("Factoring", "SpinGlass"); + assert!(shortest.is_some()); + assert_eq!(shortest.unwrap().len(), 2); // Factoring -> CircuitSAT -> SpinGlass } #[test] @@ -582,12 +583,12 @@ fn test_find_cheapest_path_multi_step() { } #[test] -fn test_find_cheapest_path_no_path() { +fn test_find_cheapest_path_is_to_qubo() { let graph = ReductionGraph::new(); let cost_fn = MinimizeSteps; let input_size = ProblemSize::new(vec![("n", 10)]); - // No path from IndependentSet to QUBO + // Direct path from IndependentSet to QUBO let path = graph.find_cheapest_path( ("IndependentSet", "SimpleGraph"), ("QUBO", "SimpleGraph"), @@ -595,7 +596,8 @@ fn test_find_cheapest_path_no_path() { &cost_fn, ); - assert!(path.is_none()); + assert!(path.is_some()); + assert_eq!(path.unwrap().len(), 1); // Direct path } #[test] diff --git a/src/unit_tests/rules/ilp_qubo.rs b/src/unit_tests/rules/ilp_qubo.rs new file mode 100644 index 00000000..d53593d0 --- /dev/null +++ b/src/unit_tests/rules/ilp_qubo.rs @@ -0,0 +1,171 @@ +use super::*; +use crate::models::optimization::{LinearConstraint, ObjectiveSense}; +use crate::solvers::{BruteForce, Solver}; + +#[test] +fn test_ilp_to_qubo_closed_loop() { + // Binary ILP: maximize x0 + 2*x1 + 3*x2 + // s.t. x0 + x1 <= 1, x1 + x2 <= 1 + // Optimal: x = [1, 0, 1] with obj = 4 + let ilp = ILP::binary( + 3, + vec![ + LinearConstraint::le(vec![(0, 1.0), (1, 1.0)], 1.0), + LinearConstraint::le(vec![(1, 1.0), (2, 1.0)], 1.0), + ], + vec![(0, 1.0), (1, 2.0), (2, 3.0)], + ObjectiveSense::Maximize, + ); + let reduction = ReduceTo::<QUBO<f64>>::reduce_to(&ilp); + let qubo = reduction.target_problem(); + + let solver = BruteForce::new(); + let qubo_solutions = solver.find_best(qubo); + + for sol in &qubo_solutions { + let extracted = reduction.extract_solution(sol); + assert!(ilp.solution_size(&extracted).is_valid); + } + + // Optimal should be [1, 0, 1] + let best = reduction.extract_solution(&qubo_solutions[0]); + assert_eq!(best, vec![1, 0, 1]); +} + +#[test] +fn test_ilp_to_qubo_minimize() { + // Binary ILP: minimize x0 + 2*x1 + 3*x2 + // s.t. x0 + x1 >= 1 (at least one of x0, x1 selected) + // Optimal: x = [1, 0, 0] with obj = 1 + let ilp = ILP::binary( + 3, + vec![LinearConstraint::ge(vec![(0, 1.0), (1, 1.0)], 1.0)], + vec![(0, 1.0), (1, 2.0), (2, 3.0)], + ObjectiveSense::Minimize, + ); + let reduction = ReduceTo::<QUBO<f64>>::reduce_to(&ilp); + let qubo = reduction.target_problem(); + + let solver = BruteForce::new(); + let qubo_solutions = solver.find_best(qubo); + + for sol in &qubo_solutions { + let extracted = reduction.extract_solution(sol); + assert!(ilp.solution_size(&extracted).is_valid); + } + + let best = reduction.extract_solution(&qubo_solutions[0]); + assert_eq!(best, vec![1, 0, 0]); +} + +#[test] +fn test_ilp_to_qubo_equality() { + // Binary ILP: maximize x0 + x1 + x2 + // s.t. x0 + x1 + x2 = 2 + // Optimal: any 2 of 3 variables = 1 + let ilp = ILP::binary( + 3, + vec![LinearConstraint::eq( + vec![(0, 1.0), (1, 1.0), (2, 1.0)], + 2.0, + )], + vec![(0, 1.0), (1, 1.0), (2, 1.0)], + ObjectiveSense::Maximize, + ); + let reduction = ReduceTo::<QUBO<f64>>::reduce_to(&ilp); + let qubo = reduction.target_problem(); + + let solver = BruteForce::new(); + let qubo_solutions = solver.find_best(qubo); + + // Should have exactly 3 optimal solutions (C(3,2)) + assert_eq!(qubo_solutions.len(), 3); + + for sol in &qubo_solutions { + let extracted = reduction.extract_solution(sol); + assert!(ilp.solution_size(&extracted).is_valid); + assert_eq!(extracted.iter().filter(|&&x| x == 1).count(), 2); + } +} + +#[test] +fn test_ilp_to_qubo_ge_with_slack() { + // Ge constraint with slack_range > 1 to exercise slack variable code path. + // 3 vars: minimize x0 + x1 + x2 + // s.t. x0 + x1 + x2 >= 1 (max_lhs=3, b=1, slack_range=2, ns=ceil(log2(3))=2) + let ilp = ILP::binary( + 3, + vec![LinearConstraint::ge( + vec![(0, 1.0), (1, 1.0), (2, 1.0)], + 1.0, + )], + vec![(0, 1.0), (1, 1.0), (2, 1.0)], + ObjectiveSense::Minimize, + ); + let reduction = ReduceTo::<QUBO<f64>>::reduce_to(&ilp); + let qubo = reduction.target_problem(); + + // 3 original + ceil(log2(3))=2 slack = 5 QUBO variables + assert_eq!(qubo.num_variables(), 5); + + let solver = BruteForce::new(); + let qubo_solutions = solver.find_best(qubo); + + for sol in &qubo_solutions { + let extracted = reduction.extract_solution(sol); + assert!(ilp.solution_size(&extracted).is_valid); + } + + // Optimal: exactly one variable = 1 + let best = reduction.extract_solution(&qubo_solutions[0]); + assert_eq!(best.iter().sum::<usize>(), 1); +} + +#[test] +fn test_ilp_to_qubo_le_with_slack() { + // Le constraint with rhs > 1 to exercise Le slack variable code path. + // 3 vars: maximize x0 + x1 + x2 + // s.t. x0 + x1 + x2 <= 2 (min_lhs=0, b=2, slack_range=2, ns=ceil(log2(3))=2) + let ilp = ILP::binary( + 3, + vec![LinearConstraint::le( + vec![(0, 1.0), (1, 1.0), (2, 1.0)], + 2.0, + )], + vec![(0, 1.0), (1, 1.0), (2, 1.0)], + ObjectiveSense::Maximize, + ); + let reduction = ReduceTo::<QUBO<f64>>::reduce_to(&ilp); + let qubo = reduction.target_problem(); + + // 3 original + ceil(log2(3))=2 slack = 5 QUBO variables + assert_eq!(qubo.num_variables(), 5); + + let solver = BruteForce::new(); + let qubo_solutions = solver.find_best(qubo); + + for sol in &qubo_solutions { + let extracted = reduction.extract_solution(sol); + assert!(ilp.solution_size(&extracted).is_valid); + } + + // Optimal: exactly 2 of 3 variables = 1 (3 solutions) + let best = reduction.extract_solution(&qubo_solutions[0]); + assert_eq!(best.iter().sum::<usize>(), 2); +} + +#[test] +fn test_ilp_to_qubo_sizes() { + let ilp = ILP::binary( + 3, + vec![LinearConstraint::le(vec![(0, 1.0), (1, 1.0)], 1.0)], + vec![(0, 1.0), (1, 2.0), (2, 3.0)], + ObjectiveSense::Maximize, + ); + let reduction = ReduceTo::<QUBO<f64>>::reduce_to(&ilp); + + let source_size = reduction.source_size(); + let target_size = reduction.target_size(); + assert!(!source_size.components.is_empty()); + assert!(!target_size.components.is_empty()); +} diff --git a/src/unit_tests/rules/independentset_qubo.rs b/src/unit_tests/rules/independentset_qubo.rs new file mode 100644 index 00000000..f1ed1c11 --- /dev/null +++ b/src/unit_tests/rules/independentset_qubo.rs @@ -0,0 +1,66 @@ +use super::*; +use crate::solvers::{BruteForce, Solver}; + +#[test] +fn test_independentset_to_qubo_closed_loop() { + // Path graph: 0-1-2-3 (4 vertices, 3 edges) + // Maximum IS = {0, 2} or {1, 3} (size 2) + let is = IndependentSet::<SimpleGraph, i32>::new(4, vec![(0, 1), (1, 2), (2, 3)]); + let reduction = ReduceTo::<QUBO<f64>>::reduce_to(&is); + let qubo = reduction.target_problem(); + + let solver = BruteForce::new(); + let qubo_solutions = solver.find_best(qubo); + + for sol in &qubo_solutions { + let extracted = reduction.extract_solution(sol); + assert!(is.solution_size(&extracted).is_valid); + assert_eq!(extracted.iter().filter(|&&x| x == 1).count(), 2); + } +} + +#[test] +fn test_independentset_to_qubo_triangle() { + // Triangle: 0-1-2 (complete graph K3) + // Maximum IS = any single vertex (size 1) + let is = IndependentSet::<SimpleGraph, i32>::new(3, vec![(0, 1), (1, 2), (0, 2)]); + let reduction = ReduceTo::<QUBO<f64>>::reduce_to(&is); + let qubo = reduction.target_problem(); + + let solver = BruteForce::new(); + let qubo_solutions = solver.find_best(qubo); + + for sol in &qubo_solutions { + let extracted = reduction.extract_solution(sol); + assert!(is.solution_size(&extracted).is_valid); + assert_eq!(extracted.iter().filter(|&&x| x == 1).count(), 1); + } +} + +#[test] +fn test_independentset_to_qubo_empty_graph() { + // No edges: all vertices form the IS + let is = IndependentSet::<SimpleGraph, i32>::new(3, vec![]); + let reduction = ReduceTo::<QUBO<f64>>::reduce_to(&is); + let qubo = reduction.target_problem(); + + let solver = BruteForce::new(); + let qubo_solutions = solver.find_best(qubo); + + for sol in &qubo_solutions { + let extracted = reduction.extract_solution(sol); + assert!(is.solution_size(&extracted).is_valid); + assert_eq!(extracted.iter().filter(|&&x| x == 1).count(), 3); + } +} + +#[test] +fn test_independentset_to_qubo_sizes() { + let is = IndependentSet::<SimpleGraph, i32>::new(4, vec![(0, 1), (1, 2), (2, 3)]); + let reduction = ReduceTo::<QUBO<f64>>::reduce_to(&is); + + let source_size = reduction.source_size(); + let target_size = reduction.target_size(); + assert!(!source_size.components.is_empty()); + assert!(!target_size.components.is_empty()); +} diff --git a/src/unit_tests/rules/ksatisfiability_qubo.rs b/src/unit_tests/rules/ksatisfiability_qubo.rs new file mode 100644 index 00000000..641e9216 --- /dev/null +++ b/src/unit_tests/rules/ksatisfiability_qubo.rs @@ -0,0 +1,104 @@ +use super::*; +use crate::models::satisfiability::CNFClause; +use crate::solvers::{BruteForce, Solver}; + +#[test] +fn test_ksatisfiability_to_qubo_closed_loop() { + // 3 vars, 4 clauses (matches ground truth): + // (x1 ∨ x2), (¬x1 ∨ x3), (x2 ∨ ¬x3), (¬x2 ∨ ¬x3) + let ksat = KSatisfiability::<2, i32>::new( + 3, + vec![ + CNFClause::new(vec![1, 2]), // x1 ∨ x2 + CNFClause::new(vec![-1, 3]), // ¬x1 ∨ x3 + CNFClause::new(vec![2, -3]), // x2 ∨ ¬x3 + CNFClause::new(vec![-2, -3]), // ¬x2 ∨ ¬x3 + ], + ); + let reduction = ReduceTo::<QUBO<f64>>::reduce_to(&ksat); + let qubo = reduction.target_problem(); + + let solver = BruteForce::new(); + let qubo_solutions = solver.find_best(qubo); + + // Verify all solutions satisfy all clauses + for sol in &qubo_solutions { + let extracted = reduction.extract_solution(sol); + assert!(ksat.solution_size(&extracted).is_valid); + } +} + +#[test] +fn test_ksatisfiability_to_qubo_simple() { + // 2 vars, 1 clause: (x1 ∨ x2) → 3 satisfying assignments + let ksat = KSatisfiability::<2, i32>::new(2, vec![CNFClause::new(vec![1, 2])]); + let reduction = ReduceTo::<QUBO<f64>>::reduce_to(&ksat); + let qubo = reduction.target_problem(); + + let solver = BruteForce::new(); + let qubo_solutions = solver.find_best(qubo); + + for sol in &qubo_solutions { + let extracted = reduction.extract_solution(sol); + assert!(ksat.solution_size(&extracted).is_valid); + } +} + +#[test] +fn test_ksatisfiability_to_qubo_contradiction() { + // 1 var, 2 clauses: (x1 ∨ x1) and (¬x1 ∨ ¬x1) — can't satisfy both + // Actually, this is (x1) and (¬x1), which is a contradiction + // Max-2-SAT will satisfy 1 of 2 clauses + let ksat = KSatisfiability::<2, i32>::new( + 1, + vec![ + CNFClause::new(vec![1, 1]), // x1 ∨ x1 = x1 + CNFClause::new(vec![-1, -1]), // ¬x1 ∨ ¬x1 = ¬x1 + ], + ); + let reduction = ReduceTo::<QUBO<f64>>::reduce_to(&ksat); + let qubo = reduction.target_problem(); + + let solver = BruteForce::new(); + let qubo_solutions = solver.find_best(qubo); + + // Both x=0 and x=1 satisfy exactly 1 clause + assert_eq!(qubo_solutions.len(), 2); +} + +#[test] +fn test_ksatisfiability_to_qubo_reversed_vars() { + // Clause (3, -1) has var_i=2 > var_j=0, triggering the swap branch (line 71). + // 3 vars, clauses: (x3 ∨ ¬x1), (x1 ∨ x2) + let ksat = KSatisfiability::<2, i32>::new( + 3, + vec![ + CNFClause::new(vec![3, -1]), // var 2 > var 0 → swap + CNFClause::new(vec![1, 2]), + ], + ); + let reduction = ReduceTo::<QUBO<f64>>::reduce_to(&ksat); + let qubo = reduction.target_problem(); + + let solver = BruteForce::new(); + let qubo_solutions = solver.find_best(qubo); + + for sol in &qubo_solutions { + let extracted = reduction.extract_solution(sol); + assert!(ksat.solution_size(&extracted).is_valid); + } +} + +#[test] +fn test_ksatisfiability_to_qubo_sizes() { + let ksat = KSatisfiability::<2, i32>::new( + 3, + vec![CNFClause::new(vec![1, 2]), CNFClause::new(vec![-1, 3])], + ); + let reduction = ReduceTo::<QUBO<f64>>::reduce_to(&ksat); + + let source_size = reduction.source_size(); + let target_size = reduction.target_size(); + assert!(!source_size.components.is_empty()); + assert!(!target_size.components.is_empty()); +} diff --git a/src/unit_tests/rules/setpacking_qubo.rs b/src/unit_tests/rules/setpacking_qubo.rs new file mode 100644 index 00000000..e34d75bf --- /dev/null +++ b/src/unit_tests/rules/setpacking_qubo.rs @@ -0,0 +1,67 @@ +use super::*; +use crate::solvers::{BruteForce, Solver}; + +#[test] +fn test_setpacking_to_qubo_closed_loop() { + // 3 sets: {0,2}, {1,2}, {0,3} + // Overlaps: (0,1) share element 2, (0,2) share element 0 + // Max packing: sets 1 and 2 → {1,2} and {0,3} (no overlap) + let sp = SetPacking::<i32>::new(vec![vec![0, 2], vec![1, 2], vec![0, 3]]); + let reduction = ReduceTo::<QUBO<f64>>::reduce_to(&sp); + let qubo = reduction.target_problem(); + + let solver = BruteForce::new(); + let qubo_solutions = solver.find_best(qubo); + + for sol in &qubo_solutions { + let extracted = reduction.extract_solution(sol); + assert!(sp.solution_size(&extracted).is_valid); + assert_eq!(extracted.iter().filter(|&&x| x == 1).count(), 2); + } +} + +#[test] +fn test_setpacking_to_qubo_disjoint() { + // Disjoint sets: all can be packed + let sp = SetPacking::<i32>::new(vec![vec![0, 1], vec![2, 3], vec![4]]); + let reduction = ReduceTo::<QUBO<f64>>::reduce_to(&sp); + let qubo = reduction.target_problem(); + + let solver = BruteForce::new(); + let qubo_solutions = solver.find_best(qubo); + + for sol in &qubo_solutions { + let extracted = reduction.extract_solution(sol); + assert!(sp.solution_size(&extracted).is_valid); + // All 3 sets should be selected + assert_eq!(extracted.iter().filter(|&&x| x == 1).count(), 3); + } +} + +#[test] +fn test_setpacking_to_qubo_all_overlap() { + // All sets overlap: only 1 can be selected + let sp = SetPacking::<i32>::new(vec![vec![0, 1], vec![0, 2], vec![0, 3]]); + let reduction = ReduceTo::<QUBO<f64>>::reduce_to(&sp); + let qubo = reduction.target_problem(); + + let solver = BruteForce::new(); + let qubo_solutions = solver.find_best(qubo); + + for sol in &qubo_solutions { + let extracted = reduction.extract_solution(sol); + assert!(sp.solution_size(&extracted).is_valid); + assert_eq!(extracted.iter().filter(|&&x| x == 1).count(), 1); + } +} + +#[test] +fn test_setpacking_to_qubo_sizes() { + let sp = SetPacking::<i32>::new(vec![vec![0, 2], vec![1, 2], vec![0, 3]]); + let reduction = ReduceTo::<QUBO<f64>>::reduce_to(&sp); + + let source_size = reduction.source_size(); + let target_size = reduction.target_size(); + assert!(!source_size.components.is_empty()); + assert!(!target_size.components.is_empty()); +} diff --git a/src/unit_tests/rules/vertexcovering_qubo.rs b/src/unit_tests/rules/vertexcovering_qubo.rs new file mode 100644 index 00000000..83416b1c --- /dev/null +++ b/src/unit_tests/rules/vertexcovering_qubo.rs @@ -0,0 +1,66 @@ +use super::*; +use crate::solvers::{BruteForce, Solver}; + +#[test] +fn test_vertexcovering_to_qubo_closed_loop() { + // Cycle C4: 0-1-2-3-0 (4 vertices, 4 edges) + // Minimum VC = 2 vertices (e.g., {0, 2} or {1, 3}) + let vc = VertexCovering::<SimpleGraph, i32>::new(4, vec![(0, 1), (1, 2), (2, 3), (0, 3)]); + let reduction = ReduceTo::<QUBO<f64>>::reduce_to(&vc); + let qubo = reduction.target_problem(); + + let solver = BruteForce::new(); + let qubo_solutions = solver.find_best(qubo); + + for sol in &qubo_solutions { + let extracted = reduction.extract_solution(sol); + assert!(vc.solution_size(&extracted).is_valid); + assert_eq!(extracted.iter().filter(|&&x| x == 1).count(), 2); + } +} + +#[test] +fn test_vertexcovering_to_qubo_triangle() { + // Triangle K3: minimum VC = 2 (any two vertices) + let vc = VertexCovering::<SimpleGraph, i32>::new(3, vec![(0, 1), (1, 2), (0, 2)]); + let reduction = ReduceTo::<QUBO<f64>>::reduce_to(&vc); + let qubo = reduction.target_problem(); + + let solver = BruteForce::new(); + let qubo_solutions = solver.find_best(qubo); + + for sol in &qubo_solutions { + let extracted = reduction.extract_solution(sol); + assert!(vc.solution_size(&extracted).is_valid); + assert_eq!(extracted.iter().filter(|&&x| x == 1).count(), 2); + } +} + +#[test] +fn test_vertexcovering_to_qubo_star() { + // Star graph: center vertex 0 connected to 1, 2, 3 + // Minimum VC = {0} (just the center) + let vc = VertexCovering::<SimpleGraph, i32>::new(4, vec![(0, 1), (0, 2), (0, 3)]); + let reduction = ReduceTo::<QUBO<f64>>::reduce_to(&vc); + let qubo = reduction.target_problem(); + + let solver = BruteForce::new(); + let qubo_solutions = solver.find_best(qubo); + + for sol in &qubo_solutions { + let extracted = reduction.extract_solution(sol); + assert!(vc.solution_size(&extracted).is_valid); + assert_eq!(extracted.iter().filter(|&&x| x == 1).count(), 1); + } +} + +#[test] +fn test_vertexcovering_to_qubo_sizes() { + let vc = VertexCovering::<SimpleGraph, i32>::new(4, vec![(0, 1), (1, 2), (2, 3), (0, 3)]); + let reduction = ReduceTo::<QUBO<f64>>::reduce_to(&vc); + + let source_size = reduction.source_size(); + let target_size = reduction.target_size(); + assert!(!source_size.components.is_empty()); + assert!(!target_size.components.is_empty()); +} diff --git a/tests/data/qubo/coloring_to_qubo.json b/tests/data/qubo/coloring_to_qubo.json new file mode 100644 index 00000000..f10b02e8 --- /dev/null +++ b/tests/data/qubo/coloring_to_qubo.json @@ -0,0 +1 @@ +{"problem":"Coloring","source":{"num_vertices":3,"edges":[[0,1],[1,2],[0,2]],"num_colors":3,"penalty":10.0},"qubo_matrix":[[-10.0,10.0,10.0,5.0,0.0,0.0,5.0,0.0,0.0],[10.0,-10.0,10.0,0.0,5.0,0.0,0.0,5.0,0.0],[10.0,10.0,-10.0,0.0,0.0,5.0,0.0,0.0,5.0],[5.0,0.0,0.0,-10.0,10.0,10.0,5.0,0.0,0.0],[0.0,5.0,0.0,10.0,-10.0,10.0,0.0,5.0,0.0],[0.0,0.0,5.0,10.0,10.0,-10.0,0.0,0.0,5.0],[5.0,0.0,0.0,5.0,0.0,0.0,-10.0,10.0,10.0],[0.0,5.0,0.0,0.0,5.0,0.0,10.0,-10.0,10.0],[0.0,0.0,5.0,0.0,0.0,5.0,10.0,10.0,-10.0]],"qubo_num_vars":9,"qubo_optimal":{"value":-30.0,"configs":[[0,0,1,0,1,0,1,0,0],[0,0,1,1,0,0,0,1,0],[0,1,0,0,0,1,1,0,0],[0,1,0,1,0,0,0,0,1],[1,0,0,0,0,1,0,1,0],[1,0,0,0,1,0,0,0,1]]}} \ No newline at end of file diff --git a/tests/data/qubo/ilp_to_qubo.json b/tests/data/qubo/ilp_to_qubo.json new file mode 100644 index 00000000..a7408600 --- /dev/null +++ b/tests/data/qubo/ilp_to_qubo.json @@ -0,0 +1 @@ +{"problem":"ILP","source":{"num_variables":3,"objective":[1.0,2.0,3.0],"constraints_lhs":[[1.0,1.0,0.0],[0.0,1.0,1.0]],"constraints_rhs":[1.0,1.0],"constraint_signs":[-1,-1],"penalty":10.0},"qubo_matrix":[[-11.0,10.0,0.0],[10.0,-22.0,10.0],[0.0,10.0,-13.0]],"qubo_num_vars":3,"qubo_optimal":{"value":-24.0,"configs":[[1,0,1]]}} \ No newline at end of file diff --git a/tests/data/qubo/independentset_to_qubo.json b/tests/data/qubo/independentset_to_qubo.json new file mode 100644 index 00000000..571dba17 --- /dev/null +++ b/tests/data/qubo/independentset_to_qubo.json @@ -0,0 +1 @@ +{"problem":"IndependentSet","source":{"num_vertices":4,"edges":[[0,1],[1,2],[2,3],[0,3]],"penalty":8.0},"qubo_matrix":[[-1.0,8.0,0.0,8.0],[0.0,-1.0,8.0,0.0],[0.0,0.0,-1.0,8.0],[0.0,0.0,0.0,-1.0]],"qubo_num_vars":4,"qubo_optimal":{"value":-2.0,"configs":[[0,1,0,1],[1,0,1,0]]}} \ No newline at end of file diff --git a/tests/data/qubo/ksatisfiability_to_qubo.json b/tests/data/qubo/ksatisfiability_to_qubo.json new file mode 100644 index 00000000..585a784c --- /dev/null +++ b/tests/data/qubo/ksatisfiability_to_qubo.json @@ -0,0 +1 @@ +{"problem":"KSatisfiability","source":{"num_variables":3,"clauses":[[{"variable":0,"negated":false},{"variable":1,"negated":false}],[{"variable":0,"negated":true},{"variable":2,"negated":false}],[{"variable":1,"negated":false},{"variable":2,"negated":true}],[{"variable":1,"negated":true},{"variable":2,"negated":true}]]},"qubo_matrix":[[0.0,0.5,-0.5],[0.5,-1.0,0.0],[-0.5,0.0,1.0]],"qubo_num_vars":3,"qubo_optimal":{"value":-1.0,"configs":[[0,1,0]]}} \ No newline at end of file diff --git a/tests/data/qubo/setpacking_to_qubo.json b/tests/data/qubo/setpacking_to_qubo.json new file mode 100644 index 00000000..b9b86e2e --- /dev/null +++ b/tests/data/qubo/setpacking_to_qubo.json @@ -0,0 +1 @@ +{"problem":"SetPacking","source":{"sets":[[0,2],[1,2],[0,3]],"num_elements":4,"weights":[1.0,2.0,1.5],"penalty":8.0},"qubo_matrix":[[-1.0,4.0,4.0],[4.0,-2.0,0.0],[4.0,0.0,-1.5]],"qubo_num_vars":3,"qubo_optimal":{"value":-3.5,"configs":[[0,1,1]]}} \ No newline at end of file diff --git a/tests/data/qubo/vertexcovering_to_qubo.json b/tests/data/qubo/vertexcovering_to_qubo.json new file mode 100644 index 00000000..06479d1e --- /dev/null +++ b/tests/data/qubo/vertexcovering_to_qubo.json @@ -0,0 +1 @@ +{"problem":"VertexCovering","source":{"num_vertices":4,"edges":[[0,1],[1,2],[2,3],[0,3],[0,2]],"penalty":8.0},"qubo_matrix":[[-23.0,4.0,4.0,4.0],[4.0,-15.0,4.0,0.0],[4.0,4.0,-23.0,4.0],[4.0,0.0,4.0,-15.0]],"qubo_num_vars":4,"qubo_optimal":{"value":-38.0,"configs":[[1,0,1,0]]}} \ No newline at end of file diff --git a/tests/suites/reductions.rs b/tests/suites/reductions.rs index f02ebf2d..28fccf84 100644 --- a/tests/suites/reductions.rs +++ b/tests/suites/reductions.rs @@ -441,6 +441,329 @@ mod truth_table_tests { } } +/// Tests for QUBO reductions against ground truth JSON. +mod qubo_reductions { + use super::*; + use serde::Deserialize; + + #[derive(Deserialize)] + #[allow(dead_code)] + struct QuboOptimal { + value: f64, + configs: Vec<Vec<usize>>, + } + + #[derive(Deserialize)] + struct ISToQuboData { + source: ISSource, + qubo_num_vars: usize, + qubo_optimal: QuboOptimal, + } + + #[derive(Deserialize)] + struct ISSource { + num_vertices: usize, + edges: Vec<(usize, usize)>, + } + + #[test] + fn test_is_to_qubo_ground_truth() { + let json = std::fs::read_to_string("tests/data/qubo/independentset_to_qubo.json").unwrap(); + let data: ISToQuboData = serde_json::from_str(&json).unwrap(); + + let is = IndependentSet::<SimpleGraph, i32>::new( + data.source.num_vertices, + data.source.edges, + ); + let reduction = ReduceTo::<QUBO>::reduce_to(&is); + let qubo = reduction.target_problem(); + + assert_eq!(qubo.num_variables(), data.qubo_num_vars); + + let solver = BruteForce::new(); + let solutions = solver.find_best(qubo); + + // All QUBO optimal solutions should extract to valid IS solutions + for sol in &solutions { + let extracted = reduction.extract_solution(sol); + assert!(is.solution_size(&extracted).is_valid); + } + + // Optimal IS size should match ground truth + let gt_is_size: usize = data.qubo_optimal.configs[0].iter().sum(); + let our_is_size: usize = reduction.extract_solution(&solutions[0]).iter().sum(); + assert_eq!(our_is_size, gt_is_size); + } + + #[derive(Deserialize)] + struct VCToQuboData { + source: VCSource, + qubo_num_vars: usize, + qubo_optimal: QuboOptimal, + } + + #[derive(Deserialize)] + struct VCSource { + num_vertices: usize, + edges: Vec<(usize, usize)>, + } + + #[test] + fn test_vc_to_qubo_ground_truth() { + let json = + std::fs::read_to_string("tests/data/qubo/vertexcovering_to_qubo.json").unwrap(); + let data: VCToQuboData = serde_json::from_str(&json).unwrap(); + + let vc = VertexCovering::<SimpleGraph, i32>::new( + data.source.num_vertices, + data.source.edges, + ); + let reduction = ReduceTo::<QUBO>::reduce_to(&vc); + let qubo = reduction.target_problem(); + + assert_eq!(qubo.num_variables(), data.qubo_num_vars); + + let solver = BruteForce::new(); + let solutions = solver.find_best(qubo); + + for sol in &solutions { + let extracted = reduction.extract_solution(sol); + assert!(vc.solution_size(&extracted).is_valid); + } + + // Optimal VC size should match ground truth + let gt_vc_size: usize = data.qubo_optimal.configs[0].iter().sum(); + let our_vc_size: usize = reduction.extract_solution(&solutions[0]).iter().sum(); + assert_eq!(our_vc_size, gt_vc_size); + } + + #[derive(Deserialize)] + struct ColoringToQuboData { + source: ColoringSource, + qubo_num_vars: usize, + qubo_optimal: QuboOptimal, + } + + #[derive(Deserialize)] + struct ColoringSource { + num_vertices: usize, + edges: Vec<(usize, usize)>, + num_colors: usize, + } + + #[test] + fn test_coloring_to_qubo_ground_truth() { + let json = std::fs::read_to_string("tests/data/qubo/coloring_to_qubo.json").unwrap(); + let data: ColoringToQuboData = serde_json::from_str(&json).unwrap(); + + assert_eq!(data.source.num_colors, 3); + + let kc = KColoring::<3, SimpleGraph, i32>::new( + data.source.num_vertices, + data.source.edges, + ); + let reduction = ReduceTo::<QUBO>::reduce_to(&kc); + let qubo = reduction.target_problem(); + + assert_eq!(qubo.num_variables(), data.qubo_num_vars); + + let solver = BruteForce::new(); + let solutions = solver.find_best(qubo); + + for sol in &solutions { + let extracted = reduction.extract_solution(sol); + assert!(kc.solution_size(&extracted).is_valid); + } + + // Same number of optimal colorings as ground truth + assert_eq!(solutions.len(), data.qubo_optimal.configs.len()); + } + + #[derive(Deserialize)] + struct SPToQuboData { + source: SPSource, + qubo_num_vars: usize, + qubo_optimal: QuboOptimal, + } + + #[derive(Deserialize)] + struct SPSource { + sets: Vec<Vec<usize>>, + weights: Vec<f64>, + } + + #[test] + fn test_setpacking_to_qubo_ground_truth() { + let json = std::fs::read_to_string("tests/data/qubo/setpacking_to_qubo.json").unwrap(); + let data: SPToQuboData = serde_json::from_str(&json).unwrap(); + + let sp = SetPacking::with_weights(data.source.sets, data.source.weights); + let reduction = ReduceTo::<QUBO>::reduce_to(&sp); + let qubo = reduction.target_problem(); + + assert_eq!(qubo.num_variables(), data.qubo_num_vars); + + let solver = BruteForce::new(); + let solutions = solver.find_best(qubo); + + for sol in &solutions { + let extracted = reduction.extract_solution(sol); + assert!(sp.solution_size(&extracted).is_valid); + } + + // Optimal packing should match ground truth + let gt_selected: usize = data.qubo_optimal.configs[0].iter().sum(); + let our_selected: usize = reduction.extract_solution(&solutions[0]).iter().sum(); + assert_eq!(our_selected, gt_selected); + } + + #[derive(Deserialize)] + struct KSatToQuboData { + source: KSatSource, + qubo_num_vars: usize, + qubo_optimal: QuboOptimal, + } + + #[derive(Deserialize)] + struct KSatSource { + num_variables: usize, + clauses: Vec<Vec<KSatLiteral>>, + } + + #[derive(Deserialize)] + struct KSatLiteral { + variable: usize, + negated: bool, + } + + #[test] + fn test_ksat_to_qubo_ground_truth() { + let json = + std::fs::read_to_string("tests/data/qubo/ksatisfiability_to_qubo.json").unwrap(); + let data: KSatToQuboData = serde_json::from_str(&json).unwrap(); + + // Convert JSON clauses to CNFClause (1-indexed signed literals) + let clauses: Vec<CNFClause> = data + .source + .clauses + .iter() + .map(|lits| { + let signed: Vec<i32> = lits + .iter() + .map(|l| { + let var = (l.variable + 1) as i32; // 0-indexed to 1-indexed + if l.negated { -var } else { var } + }) + .collect(); + CNFClause::new(signed) + }) + .collect(); + + let ksat = KSatisfiability::<2, i32>::new(data.source.num_variables, clauses); + let reduction = ReduceTo::<QUBO>::reduce_to(&ksat); + let qubo = reduction.target_problem(); + + assert_eq!(qubo.num_variables(), data.qubo_num_vars); + + let solver = BruteForce::new(); + let solutions = solver.find_best(qubo); + + for sol in &solutions { + let extracted = reduction.extract_solution(sol); + assert!(ksat.solution_size(&extracted).is_valid); + } + + // Verify extracted solution matches ground truth assignment + let gt_config = &data.qubo_optimal.configs[0]; + let our_config = reduction.extract_solution(&solutions[0]); + assert_eq!(&our_config, gt_config); + } + + #[cfg(feature = "ilp")] + #[derive(Deserialize)] + struct ILPToQuboData { + source: ILPSource, + qubo_num_vars: usize, + qubo_optimal: QuboOptimal, + } + + #[cfg(feature = "ilp")] + #[derive(Deserialize)] + struct ILPSource { + num_variables: usize, + objective: Vec<f64>, + constraints_lhs: Vec<Vec<f64>>, + constraints_rhs: Vec<f64>, + constraint_signs: Vec<i32>, + } + + #[cfg(feature = "ilp")] + #[test] + fn test_ilp_to_qubo_ground_truth() { + let json = std::fs::read_to_string("tests/data/qubo/ilp_to_qubo.json").unwrap(); + let data: ILPToQuboData = serde_json::from_str(&json).unwrap(); + + // Build constraints from dense matrix + let constraints: Vec<LinearConstraint> = data + .source + .constraints_lhs + .iter() + .zip(data.source.constraints_rhs.iter()) + .zip(data.source.constraint_signs.iter()) + .map(|((row, &rhs), &sign)| { + let terms: Vec<(usize, f64)> = row + .iter() + .enumerate() + .filter(|(_, &c)| c.abs() > 1e-15) + .map(|(i, &c)| (i, c)) + .collect(); + match sign { + -1 => LinearConstraint::le(terms, rhs), + 0 => LinearConstraint::eq(terms, rhs), + 1 => LinearConstraint::ge(terms, rhs), + _ => panic!("Invalid constraint sign: {}", sign), + } + }) + .collect(); + + // Build objective (dense to sparse) + let objective: Vec<(usize, f64)> = data + .source + .objective + .iter() + .enumerate() + .filter(|(_, &c)| c.abs() > 1e-15) + .map(|(i, &c)| (i, c)) + .collect(); + + // The qubogen formula maximizes, so this is a Maximize ILP + let ilp = ILP::binary( + data.source.num_variables, + constraints, + objective, + ObjectiveSense::Maximize, + ); + let reduction = ReduceTo::<QUBO>::reduce_to(&ilp); + let qubo = reduction.target_problem(); + + // QUBO may have more variables (slack), but original count matches + assert!(qubo.num_variables() >= data.qubo_num_vars); + + let solver = BruteForce::new(); + let solutions = solver.find_best(qubo); + + for sol in &solutions { + let extracted = reduction.extract_solution(sol); + assert!(ilp.solution_size(&extracted).is_valid); + } + + // Optimal assignment should match ground truth + let gt_config = &data.qubo_optimal.configs[0]; + let our_config = reduction.extract_solution(&solutions[0]); + assert_eq!(&our_config, gt_config); + } +} + /// Tests for File I/O with reductions. mod io_tests { use super::*;