Add BoxedUint::inv_mod

Cargo.toml

@@ -31,7 +31,7 @@ zeroize = { version = "1", optional = true,  default-features = false }
 bincode = "1"
 criterion = { version = "0.5", features = ["html_reports"] }
 hex-literal = "0.4"
-num-bigint = "0.4"
+num-bigint = { package = "num-bigint-dig", version = "0.8" }
 num-integer = "0.1"
 num-traits = "0.2"
 proptest = "1"

src/uint/boxed.rs

@@ -10,6 +10,7 @@ mod div;
 pub(crate) mod encoding;
 mod inv_mod;
 mod mul;
+mod neg;
 mod shl;
 mod shr;
 mod sub;
@@ -18,7 +19,7 @@ mod sub_mod;
 use crate::{Limb, NonZero, Uint, Word, Zero, U128, U64};
 use alloc::{boxed::Box, vec, vec::Vec};
 use core::fmt;
-use subtle::{Choice, ConditionallySelectable};
+use subtle::{Choice, ConditionallySelectable, ConstantTimeEq};
 
 #[cfg(feature = "zeroize")]
 use zeroize::Zeroize;
@@ -85,6 +86,13 @@ impl BoxedUint {
             .fold(Choice::from(1), |acc, limb| acc & limb.is_zero())
     }
 
+    /// Is this [`BoxedUint`] equal to one?
+    pub fn is_one(&self) -> Choice {
+        let mut iter = self.limbs.iter();
+        let choice = iter.next().copied().unwrap_or(Limb::ZERO).ct_eq(&Limb::ONE);
+        iter.fold(choice, |acc, limb| acc & limb.is_zero())
+    }
+
     /// Is this integer value an odd number?
     ///
     /// # Returns
@@ -197,6 +205,25 @@ impl BoxedUint {
         Self { limbs }
     }
 
+    /// Conditionally assign `other` to `self`, according to `choice`.
+    ///
+    /// This function should execute in constant time.
+    #[inline]
+    pub fn conditional_assign(&mut self, other: &Self, choice: Choice) {
+        *self = Self::conditional_select(self, other, choice);
+    }
+
+    /// Conditionally swap `self` and `other` if `choice == 1`; otherwise,
+    /// reassign both unto themselves.
+    ///
+    /// This function should execute in constant time.
+    #[inline]
+    fn conditional_swap(a: &mut Self, b: &mut Self, choice: Choice) {
+        let t = a.clone();
+        a.conditional_assign(b, choice);
+        b.conditional_assign(&t, choice);
+    }
+
     /// Widen this type's precision to the given number of bits.
     ///
     /// Panics if `bits_precision` is not a multiple of `Limb::BITS` or smaller than the current

src/uint/boxed/add.rs

@@ -1,7 +1,7 @@
 //! [`BoxedUint`] addition operations.
 
 use crate::{BoxedUint, CheckedAdd, Limb, Zero};
-use subtle::CtOption;
+use subtle::{Choice, CtOption};
 
 impl BoxedUint {
     /// Computes `a + b + carry`, returning the result along with the new carry.
@@ -14,6 +14,18 @@ impl BoxedUint {
     pub fn wrapping_add(&self, rhs: &Self) -> Self {
         self.adc(rhs, Limb::ZERO).0
     }
+
+    /// Perform wrapping addition, returning the truthy value as the second element of the tuple
+    /// if an overflow has occurred.
+    pub(crate) fn conditional_wrapping_add(&self, rhs: &Self, choice: Choice) -> (Self, Choice) {
+        let actual_rhs = Self::conditional_select(
+            &Self::zero_with_precision(rhs.bits_precision()),
+            rhs,
+            choice,
+        );
+        let (sum, carry) = self.adc(&actual_rhs, Limb::ZERO);
+        (sum, Choice::from((carry.0 & 1) as u8))
+    }
 }
 
 impl CheckedAdd<&BoxedUint> for BoxedUint {

src/uint/boxed/bits.rs

@@ -29,6 +29,20 @@ impl BoxedUint {
         self.limbs.len() * Limb::BITS
     }
 
+    /// Calculate the number of trailing zeros in the binary representation of this number.
+    pub fn trailing_zeros(&self) -> usize {
+        let mut count: Word = 0;
+        let mut nonzero_limb_not_encountered = Choice::from(1u8);
+
+        for l in &*self.limbs {
+            let z = l.trailing_zeros() as Word;
+            count += Word::conditional_select(&0, &z, nonzero_limb_not_encountered);
+            nonzero_limb_not_encountered &= l.is_zero();
+        }
+
+        count as usize
+    }
+
     /// Sets the bit at `index` to 0 or 1 depending on the value of `bit_value`.
     pub(crate) fn set_bit(&mut self, index: usize, bit_value: Choice) {
         let limb_num = index / Limb::BITS;

src/uint/boxed/encoding.rs

@@ -24,6 +24,10 @@ impl BoxedUint {
     /// If the length of `bytes` (when interpreted as bits) is larger than `bits_precision`, this
     /// function will return [`DecodeError::InputSize`].
     pub fn from_be_slice(bytes: &[u8], bits_precision: usize) -> Result<Self, DecodeError> {
+        if bytes.is_empty() && bits_precision == 0 {
+            return Ok(Self::zero());
+        }
+
         if bits_precision % Limb::BITS != 0 {
             return Err(DecodeError::Precision);
         }
@@ -50,6 +54,10 @@ impl BoxedUint {
     /// If the length of `bytes` (when interpreted as bits) is larger than `bits_precision`, this
     /// function will return [`DecodeError::InputSize`].
     pub fn from_le_slice(bytes: &[u8], bits_precision: usize) -> Result<Self, DecodeError> {
+        if bytes.is_empty() && bits_precision == 0 {
+            return Ok(Self::zero());
+        }
+
         if bits_precision % Limb::BITS != 0 {
             return Err(DecodeError::Precision);
         }

src/uint/boxed/inv_mod.rs

@@ -1,13 +1,46 @@
 //! [`BoxedUint`] modular inverse (i.e. reciprocal) operations.
 
 use crate::{BoxedUint, Word};
-use subtle::{Choice, ConstantTimeLess};
+use subtle::{Choice, ConstantTimeEq, ConstantTimeLess, CtOption};
 
 impl BoxedUint {
+    /// Computes the multiplicative inverse of `self` mod `modulus`.
+    /// Returns `None` if an inverse does not exist.
+    pub fn inv_mod(&self, modulus: &Self) -> CtOption<Self> {
+        debug_assert_eq!(self.bits_precision(), modulus.bits_precision());
+
+        // Decompose `modulus = s * 2^k` where `s` is odd
+        let k = modulus.trailing_zeros();
+        let s = modulus.shr(k);
+
+        // Decompose `self` into RNS with moduli `2^k` and `s` and calculate the inverses.
+        // Using the fact that `(z^{-1} mod (m1 * m2)) mod m1 == z^{-1} mod m1`
+        let (a, a_is_some) = self.inv_odd_mod(&s);
+        let b = self.inv_mod2k(k);
+        // inverse modulo 2^k exists either if `k` is 0 or if `self` is odd.
+        let b_is_some = (k as Word).ct_eq(&0) | self.is_odd();
+
+        // Restore from RNS:
+        // self^{-1} = a mod s = b mod 2^k
+        // => self^{-1} = a + s * ((b - a) * s^(-1) mod 2^k)
+        // (essentially one step of the Garner's algorithm for recovery from RNS).
+
+        let m_odd_inv = s.inv_mod2k(k); // `s` is odd, so this always exists
+
+        // This part is mod 2^k
+        let mask = Self::one().shl(k).wrapping_sub(&Self::one());
+        let t = (b.wrapping_sub(&a).wrapping_mul(&m_odd_inv)).bitand(&mask);
+
+        // Will not overflow since `a <= s - 1`, `t <= 2^k - 1`,
+        // so `a + s * t <= s * 2^k - 1 == modulus - 1`.
+        let result = a.wrapping_add(&s.wrapping_mul(&t));
+        CtOption::new(result, a_is_some & b_is_some)
+    }
+
     /// Computes 1/`self` mod `2^k`.
     ///
     /// Conditions: `self` < 2^k and `self` must be odd
-    pub fn inv_mod2k(&self, k: usize) -> Self {
+    pub(crate) fn inv_mod2k(&self, k: usize) -> Self {
         // This is the same algorithm as in `inv_mod2k_vartime()`,
         // but made constant-time w.r.t `k` as well.
 
@@ -33,6 +66,83 @@ impl BoxedUint {
 
         x
     }
+
+    /// Computes the multiplicative inverse of `self` mod `modulus`, where `modulus` is odd.
+    /// Returns `None` if an inverse does not exist.
+    fn inv_odd_mod(&self, modulus: &Self) -> (Self, Choice) {
+        self.inv_odd_mod_bounded(modulus, self.bits_precision(), modulus.bits_precision())
+    }
+
+    /// Computes the multiplicative inverse of `self` mod `modulus`, where `modulus` is odd.
+    /// In other words `self^-1 mod modulus`.
+    ///
+    /// `bits` and `modulus_bits` are the bounds on the bit size
+    /// of `self` and `modulus`, respectively.
+    ///
+    /// (the inversion speed will be proportional to `bits + modulus_bits`).
+    /// The second element of the tuple is the truthy value if an inverse exists,
+    /// otherwise it is a falsy value.
+    ///
+    /// **Note:** variable time in `bits` and `modulus_bits`.
+    ///
+    /// The algorithm is the same as in GMP 6.2.1's `mpn_sec_invert`.
+    fn inv_odd_mod_bounded(
+        &self,
+        modulus: &Self,
+        bits: usize,
+        modulus_bits: usize,
+    ) -> (Self, Choice) {
+        debug_assert_eq!(self.bits_precision(), modulus.bits_precision());
+
+        let bits_precision = self.bits_precision();
+        debug_assert!(bool::from(modulus.is_odd()));
+
+        let mut a = self.clone();
+        let mut u = Self::one_with_precision(bits_precision);
+        let mut v = Self::zero_with_precision(bits_precision);
+        let mut b = modulus.clone();
+
+        // `bit_size` can be anything >= `self.bits()` + `modulus.bits()`, setting to the minimum.
+        let bit_size = bits + modulus_bits;
+
+        let mut m1hp = modulus.clone();
+        let (m1hp_new, carry) = m1hp.shr_1();
+        debug_assert!(bool::from(carry));
+        m1hp = m1hp_new.wrapping_add(&Self::one_with_precision(bits_precision));
+
+        for _ in 0..bit_size {
+            debug_assert!(bool::from(b.is_odd()));
+
+            let self_odd = a.is_odd();
+
+            // Set `self -= b` if `self` is odd.
+            let (new_a, swap) = a.conditional_wrapping_sub(&b, self_odd);
+            // Set `b += self` if `swap` is true.
+            b = Self::conditional_select(&b, &b.wrapping_add(&new_a), swap);
+            // Negate `self` if `swap` is true.
+            a = new_a.conditional_wrapping_neg(swap);
+
+            let mut new_u = u.clone();
+            let mut new_v = v.clone();
+            Self::conditional_swap(&mut new_u, &mut new_v, swap);
+            let (new_u, cy) = new_u.conditional_wrapping_sub(&new_v, self_odd);
+            let (new_u, cyy) = new_u.conditional_wrapping_add(modulus, cy);
+            debug_assert!(bool::from(cy.ct_eq(&cyy)));
+
+            let (new_a, overflow) = a.shr_1();
+            debug_assert!(!bool::from(overflow));
+            let (new_u, cy) = new_u.shr_1();
+            let (new_u, cy) = new_u.conditional_wrapping_add(&m1hp, cy);
+            debug_assert!(!bool::from(cy));
+
+            a = new_a;
+            u = new_u;
+            v = new_v;
+        }
+
+        debug_assert!(bool::from(a.is_zero()));
+        (v, b.is_one())
+    }
 }
 
 #[cfg(test)]

src/uint/boxed/mul.rs

@@ -27,6 +27,11 @@ impl BoxedUint {
         ret
     }
 
+    /// Perform wrapping multiplication, wrapping to the width of `self`.
+    pub fn wrapping_mul(&self, rhs: &Self) -> Self {
+        self.mul_wide(rhs).shorten(self.bits_precision())
+    }
+
     /// Multiply `self` by itself.
     pub fn square(&self) -> Self {
         // TODO(tarcieri): more optimized implementation