Skip to main content

problemreductions/models/algebraic/
closest_vector_problem.rs

1//! Closest Vector Problem (CVP).
2//!
3//! Given an integer lattice basis `B` and a target vector `t`, find integer
4//! coefficients `x` minimizing the squared distance `||Bx - t||_2^2`.
5
6use crate::registry::{ConstructionError, CreateSpec, ProblemSchemaEntry, VariantDimension};
7use crate::traits::{EvaluationError, Problem};
8use crate::types::Min;
9use num_bigint::BigInt;
10use num_traits::Zero;
11use serde::{Deserialize, Serialize};
12
13#[derive(Debug, Deserialize, crate::CreateSpec)]
14struct ClosestVectorProblemCreateSpec {
15    /// Integer basis matrix as semicolon-separated column vectors.
16    #[create(codec = "semicolon-separated")]
17    basis: Vec<Vec<i64>>,
18    /// Integer target vector.
19    #[create(name = "target_vec", codec = "comma-separated")]
20    target: Vec<i64>,
21}
22
23impl TryFrom<ClosestVectorProblemCreateSpec> for ClosestVectorProblem {
24    type Error = ConstructionError;
25
26    fn try_from(spec: ClosestVectorProblemCreateSpec) -> Result<Self, Self::Error> {
27        Self::new(spec.basis, spec.target)
28    }
29}
30
31inventory::submit! {
32    ProblemSchemaEntry {
33        name: "ClosestVectorProblem",
34        display_name: "Closest Vector Problem",
35        aliases: &["CVP"],
36        dimensions: &[VariantDimension::new("coefficient", "i64", &["i64"])],
37        category: crate::registry::ProblemCategory::Algebraic,
38        module_path: module_path!(),
39        description: "Find the closest point in an integer lattice to a target vector",
40        fields: ClosestVectorProblemCreateSpec::FIELDS,
41    }
42}
43
44/// Euclidean Closest Vector Problem over an integer lattice basis.
45#[derive(Debug, Clone, Serialize)]
46pub struct ClosestVectorProblem {
47    /// Basis matrix stored as column vectors.
48    basis: Vec<Vec<i64>>,
49    /// Target vector in the ambient space.
50    target: Vec<i64>,
51    #[serde(skip)]
52    coefficient_box_bits: u64,
53}
54
55impl ClosestVectorProblem {
56    /// Construct a CVP instance with a full-column-rank integer basis.
57    pub fn new(basis: Vec<Vec<i64>>, target: Vec<i64>) -> Result<Self, ConstructionError> {
58        let ambient_dimension = target.len();
59        for (index, column) in basis.iter().enumerate() {
60            if column.len() != ambient_dimension {
61                return Err(ConstructionError::Conversion(format!(
62                    "basis vector {index} has length {}, expected {ambient_dimension}",
63                    column.len()
64                )));
65            }
66        }
67        if basis.len() > ambient_dimension {
68            return Err(ConstructionError::Conversion(format!(
69                "{} basis vectors cannot be independent in ambient dimension {ambient_dimension}",
70                basis.len()
71            )));
72        }
73        let (_, determinant) = independent_rows(&basis, ambient_dimension).ok_or_else(|| {
74            ConstructionError::Conversion(
75                "closest-vector basis columns must be linearly independent".into(),
76            )
77        })?;
78        let coefficient_box_bits = if basis.is_empty() {
79            0
80        } else {
81            let norms: Vec<BigInt> = basis
82                .iter()
83                .map(|column| column.iter().map(|&entry| BigInt::from(entry).pow(2)).sum())
84                .collect();
85            // Cauchy-Binet and Hadamard bound every adjugate row's squared
86            // norm by the product of the other column norms. The determinant
87            // is already available from the rank check; no source solve is needed.
88            let norm_product = norms.iter().product::<BigInt>();
89            let mut squared_radius: BigInt =
90                target.iter().map(|&entry| BigInt::from(entry).pow(2)).sum();
91            if basis.len() == ambient_dimension {
92                // With full rank, rounding A^-1 t leaves coefficient errors
93                // <=1/2. Cauchy-Schwarz gives ||B error||^2 <=r sum(norms)/4.
94                squared_radius = squared_radius
95                    .min(BigInt::from(basis.len()) * norms.iter().sum::<BigInt>() / 4);
96            }
97            // Every coefficient interval has integer width at most
98            // floor(sqrt(4 P R^2 / det(A)^2)); width's bit length is the
99            // exact-range binary bit count. BigInt keeps large cancellations exact.
100            let scaled_radius = norm_product * squared_radius * 4u32;
101            let squared_determinant = determinant.pow(2);
102            norms
103                .iter()
104                .map(|norm| {
105                    (&scaled_radius / norm / &squared_determinant)
106                        .to_biguint()
107                        .expect("squared width is nonnegative")
108                        .sqrt()
109                        .bits()
110                })
111                .fold(0u64, u64::saturating_add)
112            // Constructible target counts fit u64, so capping an upper bound
113            // at u64::MAX remains sound even if the geometric relaxation exceeds it.
114        };
115        Ok(Self {
116            basis,
117            target,
118            coefficient_box_bits,
119        })
120    }
121
122    /// Number of basis vectors.
123    pub fn num_basis_vectors(&self) -> usize {
124        self.basis.len()
125    }
126
127    /// Dimension of the ambient space.
128    pub fn ambient_dimension(&self) -> usize {
129        self.target.len()
130    }
131
132    /// Maximum bit length of an absolute basis or target entry, at least one.
133    pub fn max_numeric_magnitude_bits(&self) -> u64 {
134        crate::types::max_numeric_magnitude_bits(
135            self.basis.iter().flatten().chain(&self.target).copied(),
136        )
137    }
138
139    /// Cached upper bound on total bits in a constructible closest-vector box.
140    ///
141    /// Uses the selected coordinate determinant, Hadamard bounds on inverse
142    /// row norms, and the zero or full-rank rounded candidate's distance bound.
143    /// Computing it adds norm sums and integer arithmetic to the existing rank
144    /// check, without computing coefficient centers or solving CVP.
145    pub fn coefficient_box_bits(&self) -> u64 {
146        self.coefficient_box_bits
147    }
148
149    /// Integer basis columns.
150    pub fn basis(&self) -> &[Vec<i64>] {
151        &self.basis
152    }
153
154    /// Target coordinates.
155    pub fn target(&self) -> &[i64] {
156        &self.target
157    }
158
159    pub(crate) fn independent_rows(&self) -> Result<Vec<usize>, ConstructionError> {
160        independent_rows(&self.basis, self.ambient_dimension())
161            .map(|(rows, _)| rows)
162            .ok_or_else(|| {
163                ConstructionError::Conversion(
164                    "closest-vector basis columns must be linearly independent".into(),
165                )
166            })
167    }
168}
169
170fn independent_rows(basis: &[Vec<i64>], ambient_dimension: usize) -> Option<(Vec<usize>, BigInt)> {
171    let num_columns = basis.len();
172    if num_columns == 0 {
173        return Some((Vec::new(), BigInt::from(1)));
174    }
175
176    let mut matrix = (0..ambient_dimension)
177        .map(|row| {
178            basis
179                .iter()
180                .map(|column| BigInt::from(column[row]))
181                .collect::<Vec<_>>()
182        })
183        .collect::<Vec<_>>();
184    // Rank is an exact predicate; elimination intermediates are not model fields.
185    let mut previous_pivot = BigInt::from(1);
186    let mut row_indices = (0..ambient_dimension).collect::<Vec<_>>();
187
188    for column in 0..num_columns {
189        let pivot_row = (column..ambient_dimension).find(|&row| !matrix[row][column].is_zero())?;
190        matrix.swap(column, pivot_row);
191        row_indices.swap(column, pivot_row);
192        let pivot = matrix[column][column].clone();
193
194        for row in (column + 1)..ambient_dimension {
195            for next_column in (column + 1)..num_columns {
196                matrix[row][next_column] = (&matrix[row][next_column] * &pivot
197                    - &matrix[row][column] * &matrix[column][next_column])
198                    / &previous_pivot;
199            }
200            matrix[row][column] = BigInt::zero();
201        }
202        previous_pivot = pivot;
203    }
204    row_indices.truncate(num_columns);
205    Some((row_indices, previous_pivot))
206}
207
208impl<'de> Deserialize<'de> for ClosestVectorProblem {
209    fn deserialize<D>(deserializer: D) -> Result<Self, D::Error>
210    where
211        D: serde::Deserializer<'de>,
212    {
213        #[derive(Deserialize)]
214        struct Raw {
215            basis: Vec<Vec<i64>>,
216            target: Vec<i64>,
217        }
218
219        let raw = Raw::deserialize(deserializer)?;
220        Self::new(raw.basis, raw.target).map_err(serde::de::Error::custom)
221    }
222}
223
224impl Problem for ClosestVectorProblem {
225    const NAME: &'static str = "ClosestVectorProblem";
226    type Solution = Vec<i64>;
227    type Value = Min<i64>;
228
229    crate::problem_parameters![
230        ("ambient_dimension", ambient_dimension),
231        ("num_basis_vectors", num_basis_vectors),
232        ("max_numeric_magnitude_bits", max_numeric_magnitude_bits),
233        ("coefficient_box_bits", coefficient_box_bits),
234    ];
235
236    fn evaluate(&self, solution: &Self::Solution) -> Result<Min<i64>, EvaluationError> {
237        if solution.len() != self.num_basis_vectors() {
238            return Err(EvaluationError::InvalidConfiguration(format!(
239                "expected {} closest-vector coefficients, got {}",
240                self.num_basis_vectors(),
241                solution.len()
242            )));
243        }
244
245        let overflow = || EvaluationError::IntegerOverflow("computing CVP squared distance".into());
246        let mut squared = 0_i64;
247        for (row, &target) in self.target.iter().enumerate() {
248            let mut coordinate = 0_i64;
249            for (&coefficient, column) in solution.iter().zip(&self.basis) {
250                coordinate = coordinate
251                    .checked_add(coefficient.checked_mul(column[row]).ok_or_else(overflow)?)
252                    .ok_or_else(overflow)?;
253            }
254            let difference = coordinate.checked_sub(target).ok_or_else(overflow)?;
255            squared = squared
256                .checked_add(difference.checked_mul(difference).ok_or_else(overflow)?)
257                .ok_or_else(overflow)?;
258        }
259        Ok(Min(Some(squared)))
260    }
261
262    fn variant() -> Vec<(&'static str, &'static str)> {
263        vec![("coefficient", "i64")]
264    }
265}
266
267crate::declare_variants! {
268    default ClosestVectorProblem => "2^(num_basis_vectors * log(num_basis_vectors))" create ClosestVectorProblemCreateSpec,
269}
270
271#[cfg(feature = "example-db")]
272pub(crate) fn canonical_model_example_specs() -> Vec<crate::example_db::specs::ModelExampleSpec> {
273    vec![crate::example_db::specs::ModelExampleSpec {
274        id: "closest_vector_problem",
275        instance: Box::new(
276            ClosestVectorProblem::new(vec![vec![2, 0], vec![1, 2]], vec![3_i64, 2])
277                .expect("canonical closest-vector instance must be valid"),
278        ),
279        optimal_config: serde_json::json!(vec![1, 1]),
280        optimal_value: serde_json::json!(0),
281    }]
282}
283
284#[cfg(test)]
285#[path = "../../unit_tests/models/algebraic/closest_vector_problem.rs"]
286mod tests;
287
288crate::decision_problem_meta!(ClosestVectorProblem, "DecisionClosestVectorProblem");
289
290inventory::submit! {
291    crate::registry::ProblemSchemaEntry {
292        name: "DecisionClosestVectorProblem", display_name: "Decision ClosestVectorProblem", aliases: &[],
293        dimensions: &[VariantDimension::new("coefficient", "i64", &["i64"])], category: crate::registry::ProblemCategory::Algebraic, module_path: module_path!(),
294        description: "Does a feasible solution have objective value <= the bound?",
295        fields: &[
296        crate::registry::FieldInfo { name: "basis", type_name: "Vec<Vec<i64>>", description: "Basis matrix as semicolon-separated column vectors." },
297        crate::registry::FieldInfo { name: "target_vec", type_name: "Vec<i64>", description: "Target vector." },
298        crate::registry::FieldInfo { name: "bound", type_name: "i64", description: "Accept objective values <= this bound" },
299    ],
300    }
301}
302crate::declare_variants! {
303    default crate::models::decision::Decision<ClosestVectorProblem> => "2^(num_basis_vectors * log(num_basis_vectors))" create crate::models::decision::DecisionCreateSpec<ClosestVectorProblem>,
304}
305crate::register_decision_variant!(@edges ClosestVectorProblem, "DecisionClosestVectorProblem");
306
307#[cfg(feature = "example-db")]
308pub(crate) fn decision_canonical_rule_example_specs(
309) -> Vec<crate::example_db::specs::RuleExampleSpec> {
310    vec![crate::example_db::specs::RuleExampleSpec {
311        id: "decision_closest_vector_problem_to_closest_vector_problem",
312        build: || {
313            let source = crate::models::decision::Decision::new(
314                ClosestVectorProblem::new(vec![vec![2, 0], vec![1, 2]], vec![3_i64, 2])
315                    .expect("canonical closest-vector instance must be valid"),
316                0,
317            );
318            let witness = serde_json::json!(vec![1, 1]);
319            crate::example_db::specs::rule_example_with_witness::<_, ClosestVectorProblem>(
320                source,
321                crate::export::SolutionPair {
322                    source_config: witness.clone(),
323                    target_config: witness,
324                },
325            )
326        },
327    }]
328}