qfall_math/rational/q/
rounding.rs

1// Copyright © 2023 Marvin Beckmann
2//
3// This file is part of qFALL-math.
4//
5// qFALL-math is free software: you can redistribute it and/or modify it under
6// the terms of the Mozilla Public License Version 2.0 as published by the
7// Mozilla Foundation. See <https://mozilla.org/en-US/MPL/2.0/>.
8
9//! This module includes functionality for rounding instances of [`Q`].
10
11use super::Q;
12use crate::{error::MathError, integer::Z, traits::Distance};
13use flint_sys::{
14    fmpq::{fmpq_sgn, fmpq_simplest_between},
15    fmpz::{fmpz_cdiv_q, fmpz_fdiv_q},
16};
17
18impl Q {
19    /// Rounds the given rational [`Q`] down to the next integer [`Z`].
20    ///
21    /// # Examples
22    /// ```
23    /// use qfall_math::rational::Q;
24    /// use qfall_math::integer::Z;
25    ///
26    /// let value = Q::from((5, 2));
27    /// assert_eq!(Z::from(2), value.floor());
28    ///
29    /// let value = Q::from((-5, 2));
30    /// assert_eq!(Z::from(-3), value.floor());
31    ///
32    /// let value = Q::from(2);
33    /// assert_eq!(Z::from(2), value.floor());
34    /// ```
35    pub fn floor(&self) -> Z {
36        let mut out = Z::default();
37        unsafe { fmpz_fdiv_q(&mut out.value, &self.value.num, &self.value.den) };
38        out
39    }
40
41    /// Rounds the given rational [`Q`] up to the next integer [`Z`].
42    ///
43    /// # Examples
44    /// ```
45    /// use qfall_math::rational::Q;
46    /// use qfall_math::integer::Z;
47    ///
48    /// let value = Q::from((5, 2));
49    /// assert_eq!(Z::from(3), value.ceil());
50    ///
51    /// let value = Q::from((-5, 2));
52    /// assert_eq!(Z::from(-2), value.ceil());
53    ///
54    /// let value = Q::from(2);
55    /// assert_eq!(Z::from(2), value.ceil());
56    /// ```
57    pub fn ceil(&self) -> Z {
58        let mut out = Z::default();
59        unsafe { fmpz_cdiv_q(&mut out.value, &self.value.num, &self.value.den) };
60        out
61    }
62
63    /// Rounds the given rational [`Q`] to the closest integer [`Z`].
64    /// If the distance is equal, it rounds up.
65    ///
66    /// # Examples
67    /// ```
68    /// use qfall_math::rational::Q;
69    /// use qfall_math::integer::Z;
70    ///
71    /// let value = Q::from((5, 2));
72    /// assert_eq!(Z::from(3), value.round());
73    ///
74    /// let value = Q::from((-5, 2));
75    /// assert_eq!(Z::from(-2), value.round());
76    ///
77    /// let value = Q::from(2);
78    /// assert_eq!(Z::from(2), value.round());
79    /// ```
80    pub fn round(&self) -> Z {
81        if Q::from(self.floor()).distance(self) < Q::from(0.5) {
82            self.floor()
83        } else {
84            self.ceil()
85        }
86    }
87
88    /// Returns the smallest rational with the smallest denominator in the range
89    /// `\[self - |precision|, self + |precision|\]`.
90    ///
91    /// This function allows to free memory in exchange for the specified loss of
92    /// precision (see Example 3). Be aware that this loss of precision is propagated by
93    /// arithmetic operations and can be significantly increased depending on the
94    /// performed operations.
95    ///
96    /// This function ensures that simplifying does not change the sign of `self`.
97    ///
98    /// Parameters:
99    /// - `precision`: the precision the new value can differ from `self`.
100    ///   Note that the absolute value is relevant, not the sign.
101    ///
102    /// Returns the simplest [`Q`] within the defined range.
103    ///
104    /// # Examples
105    /// ```
106    /// use qfall_math::rational::Q;
107    ///
108    /// let value = Q::from((17, 20));
109    /// let precision = Q::from((1, 20));
110    ///
111    /// let simplified = Q::from((4, 5));
112    /// assert_eq!(simplified, value.simplify(&precision));
113    /// ```
114    ///
115    /// ```
116    /// use qfall_math::rational::Q;
117    ///
118    /// let value = Q::from((3, 2));
119    ///
120    /// assert_eq!(Q::ONE, value.simplify(0.5));
121    /// ```
122    ///
123    /// ## Simplify with reasonable precision loss
124    /// This example uses [`Q::INV_MAX32`], i.e. a loss of precision of at most `1 / 2^31 - 2` behind the decimal point.
125    /// If you require higher precision, [`Q::INV_MAX62`] is available.
126    /// ```
127    /// use qfall_math::rational::Q;
128    /// let value = Q::PI;
129    ///
130    /// let simplified = value.simplify(Q::INV_MAX32);
131    ///
132    /// assert_ne!(&Q::PI, &simplified);
133    /// assert!(&simplified >= &(Q::PI - Q::INV_MAX32));
134    /// assert!(&simplified <= &(Q::PI + Q::INV_MAX32));
135    /// ```
136    pub fn simplify(&self, precision: impl Into<Q>) -> Self {
137        let precision = precision.into();
138
139        let lower = self - &precision;
140        let upper = self + &precision;
141        let mut out = Q::default();
142        unsafe { fmpq_simplest_between(&mut out.value, &lower.value, &upper.value) };
143
144        if unsafe { fmpq_sgn(&self.value) != fmpq_sgn(&out.value) } {
145            return Q::MINUS_ONE * out;
146        }
147
148        out
149    }
150
151    /// Performs the randomized rounding algorithm
152    /// by sampling from a discrete Gaussian over the integers shifted
153    /// by `self` with gaussian parameter `r`.
154    ///
155    /// Parameters:
156    /// - `r`: specifies the Gaussian parameter, which is proportional
157    ///   to the standard deviation `sigma * sqrt(2 * pi) = r`
158    ///
159    /// Returns the rounded value as an [`Z`] or an error if `r < 0`.
160    ///
161    /// # Examples
162    /// ```
163    /// use qfall_math::rational::Q;
164    ///
165    /// let value = Q::from((5, 2));
166    /// let rounded = value.randomized_rounding(3).unwrap();
167    /// ```
168    ///
169    /// # Errors and Failures
170    /// - Returns a [`MathError`] of type [`InvalidIntegerInput`](MathError::InvalidIntegerInput)
171    ///   if `r < 0`.
172    ///
173    /// This function implements randomized rounding according to:
174    /// - Peikert, C. (2010, August).
175    ///   An efficient and parallel Gaussian sampler for lattices.
176    ///   In: Annual Cryptology Conference (pp. 80-97).
177    ///   <https://link.springer.com/chapter/10.1007/978-3-642-14623-7_5>
178    pub fn randomized_rounding(&self, r: impl Into<Q>) -> Result<Z, MathError> {
179        Z::sample_discrete_gauss(self, r)
180    }
181}
182
183#[cfg(test)]
184mod test_floor {
185    use crate::{integer::Z, rational::Q};
186
187    // Ensure that positive rationals are rounded correctly
188    #[test]
189    fn positive() {
190        let val_1 = Q::from((i64::MAX, 2));
191        let val_2 = Q::from((1, u64::MAX));
192
193        assert_eq!(Z::from((i64::MAX - 1) / 2), val_1.floor());
194        assert_eq!(Z::ZERO, val_2.floor());
195    }
196
197    // Ensure that negative rationals are rounded correctly
198    #[test]
199    fn negative() {
200        let val_1 = Q::from((-i64::MAX, 2));
201        let val_2 = Q::from((-1, u64::MAX));
202
203        assert_eq!(Z::from((-i64::MAX - 1) / 2), val_1.floor());
204        assert_eq!(Z::MINUS_ONE, val_2.floor());
205    }
206}
207
208#[cfg(test)]
209mod test_ceil {
210    use crate::{integer::Z, rational::Q};
211
212    // Ensure that positive rationals are rounded correctly
213    #[test]
214    fn positive() {
215        let val_1 = Q::from((i64::MAX, 2));
216        let val_2 = Q::from((1, u64::MAX));
217
218        assert_eq!(Z::from((i64::MAX - 1) / 2 + 1), val_1.ceil());
219        assert_eq!(Z::ONE, val_2.ceil());
220    }
221
222    // Ensure that negative rationals are rounded correctly
223    #[test]
224    fn negative() {
225        let val_1 = Q::from((-i64::MAX, 2));
226        let val_2 = Q::from((-1, u64::MAX));
227
228        assert_eq!(Z::from((-i64::MAX - 1) / 2 + 1), val_1.ceil());
229        assert_eq!(Z::ZERO, val_2.ceil());
230    }
231}
232
233#[cfg(test)]
234mod test_round {
235    use crate::{integer::Z, rational::Q};
236
237    // Ensure that positive rationals are rounded correctly
238    #[test]
239    fn positive() {
240        let val_1 = Q::from((i64::MAX, 2));
241        let val_2 = Q::from((1, u64::MAX));
242
243        assert_eq!(Z::from((i64::MAX - 1) / 2 + 1), val_1.round());
244        assert_eq!(Z::ZERO, val_2.round());
245    }
246
247    // Ensure that negative rationals are rounded correctly
248    #[test]
249    fn negative() {
250        let val_1 = Q::from((-i64::MAX, 2));
251        let val_2 = Q::from((-1, u64::MAX));
252
253        assert_eq!(Z::from((-i64::MAX - 1) / 2 + 1), val_1.round());
254        assert_eq!(Z::ZERO, val_2.round());
255    }
256}
257
258#[cfg(test)]
259mod test_simplify {
260    use crate::{integer::Z, rational::Q, traits::Distance};
261
262    /// Ensure that negative precision works as expected
263    #[test]
264    fn precision_absolute_value() {
265        let value_1 = Q::from((17, 20));
266        let value_2 = Q::from((-17, 20));
267        let precision = Q::from((-1, 20));
268
269        let simplified_1 = Q::from((4, 5));
270        let simplified_2 = Q::from((-4, 5));
271        assert_eq!(simplified_1, value_1.simplify(&precision));
272        assert_eq!(simplified_2, value_2.simplify(precision));
273    }
274
275    /// Ensure that large values with pointer representations are reduced
276    #[test]
277    fn large_pointer_representation() {
278        let value = Q::from((u64::MAX - 1, u64::MAX));
279        let precision = Q::from((1, u64::MAX));
280
281        assert_eq!(Q::ONE, value.simplify(&precision));
282    }
283
284    /// Ensure that the simplified value stays in range
285    #[test]
286    fn stay_in_precision() {
287        let value = Q::from((i64::MAX - 1, i64::MAX));
288        let precision = Q::from((1, u64::MAX - 1));
289
290        let simplified = value.simplify(&precision);
291        assert!(&value - &precision <= simplified && simplified <= &value + &precision);
292        assert!(Q::from((i64::MAX - 2, i64::MAX)) <= simplified && simplified <= 1);
293    }
294
295    /// Ensure max_bits of denominator are not bigger than 1/2 * precision
296    #[test]
297    fn max_bits_denominator() {
298        let value = Q::PI;
299        let precisions = [Q::INV_MAX8, Q::INV_MAX16, Q::INV_MAX32];
300
301        for precision in precisions {
302            let inv_precision = precision.get_denominator();
303            let inv_precision = inv_precision.div_ceil(2);
304
305            let simplified = value.simplify(&precision);
306            let denominator = simplified.get_denominator();
307
308            assert!(denominator.distance(Z::ZERO) < inv_precision);
309        }
310    }
311
312    /// Ensure that a value which can not be simplified is not changed
313    #[test]
314    fn no_change() {
315        let precision = Q::from((1, u64::MAX - 1));
316
317        assert_eq!(Q::ONE, Q::ONE.simplify(&precision));
318        assert_eq!(Q::MINUS_ONE, Q::MINUS_ONE.simplify(&precision));
319        assert_eq!(Q::ZERO, Q::ZERO.simplify(precision));
320    }
321
322    /// Ensure that no sign change can occurr.
323    #[test]
324    fn no_change_of_sign() {
325        assert!(Q::ZERO < Q::ONE.simplify(2));
326        assert!(Q::ZERO > Q::MINUS_ONE.simplify(2));
327    }
328}
329
330#[cfg(test)]
331mod test_randomized_rounding {
332    use crate::rational::Q;
333
334    /// Ensure that a `r < 0` throws an error
335    #[test]
336    fn negative_r() {
337        let value = Q::from((2, 3));
338        assert!(value.randomized_rounding(-1).is_err());
339    }
340}
qfall_math/rational/q/rounding.rs

qfall_math/rational/q/
rounding.rs
