Fix #116: Add Knapsack to QUBO reduction (#172)

* Add plan for #116: Knapsack to QUBO reduction

Co-Authored-By: Claude Opus 4.6 <noreply@anthropic.com>

* Add plan for #114: Knapsack

* feat: add Knapsack to QUBO reduction rule

Co-Authored-By: Claude Opus 4.6 <noreply@anthropic.com>

* feat: add Knapsack to QUBO example program

Co-Authored-By: Claude Opus 4.6 <noreply@anthropic.com>

* chore: regenerate reduction graph after Knapsack->QUBO rule

Co-Authored-By: Claude Opus 4.6 <noreply@anthropic.com>

* docs: add Knapsack to QUBO reduction in paper

Co-Authored-By: Claude Opus 4.6 <noreply@anthropic.com>

* style: fix formatting (cargo fmt)

Co-Authored-By: Claude Opus 4.6 <noreply@anthropic.com>

* fix: address PR #172 review comments

- Delete plan files per user request
- Remove KS alias per user request (keep only lowercase "knapsack")
- Fix num_slack_bits() to use integer bit-length instead of float log2
  (avoids precision loss for large i64 capacities above 2^53)
- Fix 1i64 << j to 1u64 << j in QUBO slack coefficients
  (avoids negative coefficient from signed overflow)

Co-Authored-By: Claude Opus 4.6 <noreply@anthropic.com>

* Fix Knapsack to QUBO review findings

* Fix final review findings: remove README example, use nonempty check, fmt

- Remove specific Knapsack→QUBO example from README (keep README high-level)
- Replace hardcoded example count (43) with nonempty check in rule_builders test
- Fix rustfmt formatting in example and test files

Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>

* Remove legacy example binary and its integration test

The per-reduction example file examples/reduction_knapsack_to_qubo.rs
follows a deprecated pattern. Canonical examples belong in
src/example_db/rule_builders.rs (already registered there).

Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>

* Derive paper example values from loaded data instead of hardcoding

Replace hardcoded formulas and item selections in the Knapsack→QUBO
paper extra section with data-driven computations from load-example.
Constraint equation, penalty P, objective H, selected items, and
weight/value sums are now all derived from the canonical example data.

Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>

---------

Co-authored-by: Claude Opus 4.6 <noreply@anthropic.com>
Co-authored-by: GiggleLiu <cacate0129@gmail.com>

Author: 3 people

docs/paper/reductions.typ

@@ -1626,6 +1626,47 @@ where $P$ is a penalty weight large enough that any constraint violation costs m
   _Solution extraction._ Discard slack variables: return $bold(x)' [0..n]$.
 ]
 
+#let ks_qubo = load-example("Knapsack", "QUBO")
+#let ks_qubo_sol = ks_qubo.solutions.at(0)
+#let ks_qubo_num_items = ks_qubo.source.instance.weights.len()
+#let ks_qubo_num_slack = ks_qubo.target.instance.num_vars - ks_qubo_num_items
+#let ks_qubo_penalty = 1 + ks_qubo.source.instance.values.fold(0, (a, b) => a + b)
+#let ks_qubo_selected = ks_qubo_sol.source_config.enumerate().filter(((i, x)) => x == 1).map(((i, x)) => i)
+#let ks_qubo_sel_weight = ks_qubo_selected.fold(0, (a, i) => a + ks_qubo.source.instance.weights.at(i))
+#let ks_qubo_sel_value = ks_qubo_selected.fold(0, (a, i) => a + ks_qubo.source.instance.values.at(i))
+#reduction-rule("Knapsack", "QUBO",
+  example: true,
+  example-caption: [$n = #ks_qubo_num_items$ items, capacity $C = #ks_qubo.source.instance.capacity$],
+  extra: [
+    *Step 1 -- Source instance.* The canonical knapsack instance has weights $(#ks_qubo.source.instance.weights.map(str).join(", "))$, values $(#ks_qubo.source.instance.values.map(str).join(", "))$, and capacity $C = #ks_qubo.source.instance.capacity$.
+
+    *Step 2 -- Introduce slack variables.* The inequality $sum_i w_i x_i lt.eq C$ becomes an equality by adding $B = #ks_qubo_num_slack$ binary slack bits that encode unused capacity:
+    $ #ks_qubo.source.instance.weights.enumerate().map(((i, w)) => $#w x_#i$).join($+$) + #range(ks_qubo_num_slack).map(j => $#calc.pow(2, j) s_#j$).join($+$) = #ks_qubo.source.instance.capacity $
+    This gives $n + B = #ks_qubo_num_items + #ks_qubo_num_slack = #ks_qubo.target.instance.num_vars$ QUBO variables.
+
+    *Step 3 -- Add the penalty objective.* With penalty $P = 1 + sum_i v_i = #ks_qubo_penalty$, the QUBO minimizes
+    $ H = -(#ks_qubo.source.instance.values.enumerate().map(((i, v)) => $#v x_#i$).join($+$)) + #ks_qubo_penalty (#ks_qubo.source.instance.weights.enumerate().map(((i, w)) => $#w x_#i$).join($+$) + #range(ks_qubo_num_slack).map(j => $#calc.pow(2, j) s_#j$).join($+$) - #ks_qubo.source.instance.capacity)^2 $
+    so any violation of the equality is more expensive than the entire knapsack value range.
+
+    *Step 4 -- Verify a solution.* The QUBO ground state $bold(z) = (#ks_qubo_sol.target_config.map(str).join(", "))$ extracts to the knapsack choice $bold(x) = (#ks_qubo_sol.source_config.map(str).join(", "))$. This selects items $\{#ks_qubo_selected.map(str).join(", ")\}$ with total weight $#ks_qubo_selected.map(i => str(ks_qubo.source.instance.weights.at(i))).join(" + ") = #ks_qubo_sel_weight$ and total value $#ks_qubo_selected.map(i => str(ks_qubo.source.instance.values.at(i))).join(" + ") = #ks_qubo_sel_value$, so the slack bits are all zero and the penalty term vanishes #sym.checkmark.
+
+    *Count:* #ks_qubo.solutions.len() optimal QUBO solution. The source optimum is unique because items $\{#ks_qubo_selected.map(str).join(", ")\}$ are the only feasible selection achieving value #ks_qubo_sel_value.
+  ],
+)[
+  For a standard 0-1 Knapsack instance with nonnegative weights, nonnegative values, and nonnegative capacity, the inequality $sum_i w_i x_i lt.eq C$ is converted to equality using binary slack variables that encode the unused capacity. When $C > 0$, one can take $B = floor(log_2 C) + 1$ slack bits; when $C = 0$, a single slack bit also suffices. The penalty method (@sec:penalty-method) combines the negated value objective with a quadratic constraint penalty, producing a QUBO with $n + B$ binary variables.
+][
+  _Construction._ Given $n$ items with nonnegative weights $w_0, dots, w_(n-1)$, nonnegative values $v_0, dots, v_(n-1)$, and nonnegative capacity $C$, introduce $B = floor(log_2 C) + 1$ binary slack variables $s_0, dots, s_(B-1)$ when $C > 0$ (or one slack bit when $C = 0$) to convert the capacity inequality to equality:
+  $ sum_(i=0)^(n-1) w_i x_i + sum_(j=0)^(B-1) 2^j s_j = C $
+  Let $a_k$ denote the constraint coefficient of the $k$-th binary variable ($a_k = w_k$ for $k < n$, $a_(n+j) = 2^j$ for $j < B$). The QUBO objective is:
+  $ f(bold(z)) = -sum_(i=0)^(n-1) v_i x_i + P (sum_k a_k z_k - C)^2 $
+  where $bold(z) = (x_0, dots, x_(n-1), s_0, dots, s_(B-1))$ and $P = 1 + sum_i v_i$. Expanding the quadratic penalty using $z_k^2 = z_k$ (binary):
+  $ Q_(k k) = P a_k^2 - 2 P C a_k - [k < n] v_k, quad Q_(i j) = 2 P a_i a_j quad (i < j) $
+
+  _Correctness._ ($arrow.r.double$) If $bold(x)^*$ is a feasible knapsack solution with value $V^*$, then there exist slack values $bold(s)^*$ satisfying the equality constraint (encoding $C - sum w_i x_i^*$ in binary), so $f(bold(z)^*) = -V^*$. ($arrow.l.double$) If the equality constraint is violated, the penalty $(sum a_k z_k - C)^2 gt.eq 1$ contributes at least $P > sum_i v_i$ to the objective. Since all values are nonnegative, every feasible assignment has objective in the range $[-sum_i v_i, 0]$, so that penalty exceeds the entire feasible value range. Among feasible assignments (penalty zero), $f$ reduces to $-sum v_i x_i$, minimized at the knapsack optimum.
+
+  _Solution extraction._ Discard slack variables: return $bold(z)[0..n]$.
+]
+
 #let qubo_ilp = load-example("QUBO", "ILP")
 #let qubo_ilp_sol = qubo_ilp.solutions.at(0)
 #reduction-rule("QUBO", "ILP",

docs/src/reductions/problem_schemas.json

@@ -229,17 +229,17 @@
       {
         "name": "weights",
         "type_name": "Vec<i64>",
-        "description": "Item weights w_i"
+        "description": "Nonnegative item weights w_i"
       },
       {
         "name": "values",
         "type_name": "Vec<i64>",
-        "description": "Item values v_i"
+        "description": "Nonnegative item values v_i"
       },
       {
         "name": "capacity",
         "type_name": "i64",
-        "description": "Knapsack capacity C"
+        "description": "Nonnegative knapsack capacity C"
       }
     ]
   },

docs/src/reductions/reduction_graph.json

@@ -727,6 +727,17 @@
       ],
       "doc_path": "rules/sat_ksat/index.html"
     },
+    {
+      "source": 22,
+      "target": 47,
+      "overhead": [
+        {
+          "field": "num_vars",
+          "formula": "num_items + num_slack_bits"
+        }
+      ],
+      "doc_path": "rules/knapsack_qubo/index.html"
+    },
     {
       "source": 23,
       "target": 11,

src/example_db/rule_builders.rs

@@ -12,10 +12,10 @@ mod tests {
     use super::*;
 
     #[test]
-    fn builds_all_42_canonical_rule_examples() {
+    fn builds_all_canonical_rule_examples() {
         let examples = build_rule_examples();
 
-        assert_eq!(examples.len(), 42);
+        assert!(!examples.is_empty());
         assert!(examples
             .iter()
             .all(|example| !example.source.problem.is_empty()));

src/models/misc/knapsack.rs

@@ -17,16 +17,17 @@ inventory::submit! {
         module_path: module_path!(),
         description: "Select items to maximize total value subject to weight capacity constraint",
         fields: &[
-            FieldInfo { name: "weights", type_name: "Vec<i64>", description: "Item weights w_i" },
-            FieldInfo { name: "values", type_name: "Vec<i64>", description: "Item values v_i" },
-            FieldInfo { name: "capacity", type_name: "i64", description: "Knapsack capacity C" },
+            FieldInfo { name: "weights", type_name: "Vec<i64>", description: "Nonnegative item weights w_i" },
+            FieldInfo { name: "values", type_name: "Vec<i64>", description: "Nonnegative item values v_i" },
+            FieldInfo { name: "capacity", type_name: "i64", description: "Nonnegative knapsack capacity C" },
         ],
     }
 }
 
 /// The 0-1 Knapsack problem.
 ///
-/// Given `n` items, each with weight `w_i` and value `v_i`, and a capacity `C`,
+/// Given `n` items, each with nonnegative weight `w_i` and nonnegative value `v_i`,
+/// and a nonnegative capacity `C`,
 /// find a subset `S ⊆ {0, ..., n-1}` such that `∑_{i∈S} w_i ≤ C`,
 /// maximizing `∑_{i∈S} v_i`.
 ///
@@ -47,22 +48,35 @@ inventory::submit! {
 /// ```
 #[derive(Debug, Clone, Serialize, Deserialize)]
 pub struct Knapsack {
+    #[serde(deserialize_with = "nonnegative_i64_vec::deserialize")]
     weights: Vec<i64>,
+    #[serde(deserialize_with = "nonnegative_i64_vec::deserialize")]
     values: Vec<i64>,
+    #[serde(deserialize_with = "nonnegative_i64::deserialize")]
     capacity: i64,
 }
 
 impl Knapsack {
     /// Create a new Knapsack instance.
     ///
     /// # Panics
-    /// Panics if `weights` and `values` have different lengths.
+    /// Panics if `weights` and `values` have different lengths, or if any
+    /// weight, value, or the capacity is negative.
     pub fn new(weights: Vec<i64>, values: Vec<i64>, capacity: i64) -> Self {
         assert_eq!(
             weights.len(),
             values.len(),
             "weights and values must have the same length"
         );
+        assert!(
+            weights.iter().all(|&weight| weight >= 0),
+            "Knapsack weights must be nonnegative"
+        );
+        assert!(
+            values.iter().all(|&value| value >= 0),
+            "Knapsack values must be nonnegative"
+        );
+        assert!(capacity >= 0, "Knapsack capacity must be nonnegative");
         Self {
             weights,
             values,
@@ -89,6 +103,18 @@ impl Knapsack {
     pub fn num_items(&self) -> usize {
         self.weights.len()
     }
+
+    /// Returns the number of binary slack bits used by the QUBO encoding.
+    ///
+    /// For positive capacity this is `floor(log2(C)) + 1`; for zero capacity we
+    /// keep one slack bit so the encoding shape remains uniform.
+    pub fn num_slack_bits(&self) -> usize {
+        if self.capacity == 0 {
+            1
+        } else {
+            (u64::BITS - (self.capacity as u64).leading_zeros()) as usize
+        }
+    }
 }
 
 impl Problem for Knapsack {
@@ -141,6 +167,42 @@ crate::declare_variants! {
     default opt Knapsack => "2^(num_items / 2)",
 }
 
+mod nonnegative_i64 {
+    use serde::de::Error;
+    use serde::{Deserialize, Deserializer};
+
+    pub fn deserialize<'de, D>(deserializer: D) -> Result<i64, D::Error>
+    where
+        D: Deserializer<'de>,
+    {
+        let value = i64::deserialize(deserializer)?;
+        if value < 0 {
+            return Err(D::Error::custom(format!(
+                "expected nonnegative integer, got {value}"
+            )));
+        }
+        Ok(value)
+    }
+}
+
+mod nonnegative_i64_vec {
+    use serde::de::Error;
+    use serde::{Deserialize, Deserializer};
+
+    pub fn deserialize<'de, D>(deserializer: D) -> Result<Vec<i64>, D::Error>
+    where
+        D: Deserializer<'de>,
+    {
+        let values = Vec::<i64>::deserialize(deserializer)?;
+        if let Some(value) = values.iter().copied().find(|value| *value < 0) {
+            return Err(D::Error::custom(format!(
+                "expected nonnegative integers, got {value}"
+            )));
+        }
+        Ok(values)
+    }
+}
+
 #[cfg(test)]
 #[path = "../../unit_tests/models/misc/knapsack.rs"]
 mod tests;

src/rules/knapsack_qubo.rs

@@ -0,0 +1,116 @@
+//! Reduction from Knapsack to QUBO.
+//!
+//! Converts a nonnegative 0-1 Knapsack instance into QUBO by turning the
+//! capacity inequality sum(w_i * x_i) <= C into equality using binary slack
+//! variables, then constructing a QUBO that combines the objective
+//! -sum(v_i * x_i) with a quadratic penalty
+//! P * (sum(w_i * x_i) + sum(2^j * s_j) - C)^2.
+//! For nonnegative values, penalty P > sum(v_i) ensures any infeasible solution
+//! costs more than any feasible one.
+//!
+//! Reference: Lucas, 2014, "Ising formulations of many NP problems".
+
+use crate::models::algebraic::QUBO;
+use crate::models::misc::Knapsack;
+use crate::reduction;
+use crate::rules::traits::{ReduceTo, ReductionResult};
+
+/// Result of reducing Knapsack to QUBO.
+#[derive(Debug, Clone)]
+pub struct ReductionKnapsackToQUBO {
+    target: QUBO<f64>,
+    num_items: usize,
+}
+
+impl ReductionResult for ReductionKnapsackToQUBO {
+    type Source = Knapsack;
+    type Target = QUBO<f64>;
+
+    fn target_problem(&self) -> &Self::Target {
+        &self.target
+    }
+
+    fn extract_solution(&self, target_solution: &[usize]) -> Vec<usize> {
+        target_solution[..self.num_items].to_vec()
+    }
+}
+
+#[reduction(overhead = { num_vars = "num_items + num_slack_bits" })]
+impl ReduceTo<QUBO<f64>> for Knapsack {
+    type Result = ReductionKnapsackToQUBO;
+
+    fn reduce_to(&self) -> Self::Result {
+        let n = self.num_items();
+        let c = self.capacity();
+        let b = self.num_slack_bits();
+        let total = n + b;
+
+        // Penalty must exceed sum of all values
+        let sum_values: i64 = self.values().iter().sum();
+        let penalty = (sum_values + 1) as f64;
+
+        // Build QUBO matrix
+        // H = -sum(v_i * x_i) + P * (sum(w_i * x_i) + sum(2^j * s_j) - C)^2
+        //
+        // Let a_k be the coefficient of variable k in the constraint:
+        //   a_k = w_k for k < n (item variables)
+        //   a_{n+j} = 2^j for j < B (slack variables)
+        //
+        // Expanding the penalty:
+        //   P * (sum(a_k * z_k) - C)^2 = P * sum_i sum_j a_i * a_j * z_i * z_j
+        //                                 - 2P * C * sum(a_k * z_k) + P * C^2
+        // Since z_k is binary, z_k^2 = z_k, so diagonal terms become:
+        //   Q[k][k] = P * a_k^2 - 2P * C * a_k  (from penalty)
+        //   Q[k][k] -= v_k                       (from objective, item vars only)
+        // Off-diagonal terms (i < j):
+        //   Q[i][j] = 2P * a_i * a_j
+
+        let mut coeffs = vec![0.0f64; total];
+        for (i, coeff) in coeffs.iter_mut().enumerate().take(n) {
+            *coeff = self.weights()[i] as f64;
+        }
+        for j in 0..b {
+            coeffs[n + j] = (1u64 << j) as f64;
+        }
+
+        let c_f = c as f64;
+        let mut matrix = vec![vec![0.0f64; total]; total];
+
+        // Diagonal: P * a_k^2 - 2P * C * a_k - v_k (for items)
+        for k in 0..total {
+            matrix[k][k] = penalty * coeffs[k] * coeffs[k] - 2.0 * penalty * c_f * coeffs[k];
+            if k < n {
+                matrix[k][k] -= self.values()[k] as f64;
+            }
+        }
+
+        // Off-diagonal (upper triangular): 2P * a_i * a_j
+        for i in 0..total {
+            for j in (i + 1)..total {
+                matrix[i][j] = 2.0 * penalty * coeffs[i] * coeffs[j];
+            }
+        }
+
+        ReductionKnapsackToQUBO {
+            target: QUBO::from_matrix(matrix),
+            num_items: n,
+        }
+    }
+}
+
+#[cfg(feature = "example-db")]
+pub(crate) fn canonical_rule_example_specs() -> Vec<crate::example_db::specs::RuleExampleSpec> {
+    vec![crate::example_db::specs::RuleExampleSpec {
+        id: "knapsack_to_qubo",
+        build: || {
+            crate::example_db::specs::direct_best_example::<_, QUBO<f64>, _>(
+                Knapsack::new(vec![2, 3, 4, 5], vec![3, 4, 5, 7], 7),
+                |_, _| true,
+            )
+        },
+    }]
+}
+
+#[cfg(test)]
+#[path = "../unit_tests/rules/knapsack_qubo.rs"]
+mod tests;

src/rules/mod.rs

@@ -11,6 +11,7 @@ pub(crate) mod coloring_qubo;
 pub(crate) mod factoring_circuit;
 mod graph;
 mod kcoloring_casts;
+mod knapsack_qubo;
 mod ksatisfiability_casts;
 pub(crate) mod ksatisfiability_qubo;
 pub(crate) mod ksatisfiability_subsetsum;
@@ -79,6 +80,7 @@ pub(crate) fn canonical_rule_example_specs() -> Vec<crate::example_db::specs::Ru
     specs.extend(circuit_spinglass::canonical_rule_example_specs());
     specs.extend(coloring_qubo::canonical_rule_example_specs());
     specs.extend(factoring_circuit::canonical_rule_example_specs());
+    specs.extend(knapsack_qubo::canonical_rule_example_specs());
     specs.extend(ksatisfiability_qubo::canonical_rule_example_specs());
     specs.extend(ksatisfiability_subsetsum::canonical_rule_example_specs());
     specs.extend(maximumclique_maximumindependentset::canonical_rule_example_specs());