diff --git a/src/uu/factor/build.rs b/src/uu/factor/build.rs index 642b646c6..4c3125472 100644 --- a/src/uu/factor/build.rs +++ b/src/uu/factor/build.rs @@ -27,9 +27,11 @@ use self::sieve::Sieve; #[cfg(test)] use miller_rabin::is_prime; -#[path = "src/numeric.rs"] -mod numeric; -use numeric::modular_inverse; +#[path = "src/numeric/modular_inverse.rs"] +mod modular_inverse; +#[path = "src/numeric/traits.rs"] +mod traits; +use modular_inverse::modular_inverse; mod sieve; diff --git a/src/uu/factor/src/numeric/gcd.rs b/src/uu/factor/src/numeric/gcd.rs new file mode 100644 index 000000000..54067c56e --- /dev/null +++ b/src/uu/factor/src/numeric/gcd.rs @@ -0,0 +1,67 @@ +// * This file is part of the uutils coreutils package. +// * +// * (c) 2015 Wiktor Kuropatwa +// * (c) 2020 nicoo +// * +// * For the full copyright and license information, please view the LICENSE file +// * that was distributed with this source code. + +use std::cmp::min; +use std::mem::swap; + +pub fn gcd(mut n: u64, mut m: u64) -> u64 { + // Stein's binary GCD algorithm + // Base cases: gcd(n, 0) = gcd(0, n) = n + if n == 0 { + return m; + } else if m == 0 { + return n; + } + + // Extract common factor-2: gcd(2ⁱ n, 2ⁱ m) = 2ⁱ gcd(n, m) + // and reducing until odd gcd(2ⁱ n, m) = gcd(n, m) if m is odd + let k = { + let k_n = n.trailing_zeros(); + let k_m = m.trailing_zeros(); + n >>= k_n; + m >>= k_m; + min(k_n, k_m) + }; + + loop { + // Invariant: n odd + debug_assert!(n % 2 == 1, "n = {} is even", n); + + if n > m { + swap(&mut n, &mut m); + } + m -= n; + + if m == 0 { + return n << k; + } + + m >>= m.trailing_zeros(); + } +} + +#[cfg(test)] +mod tests { + use super::*; + use quickcheck::quickcheck; + + quickcheck! { + fn gcd(a: u64, b: u64) -> bool { + // Test against the Euclidean algorithm + let g = { + let (mut a, mut b) = (a, b); + while b > 0 { + a %= b; + swap(&mut a, &mut b); + } + a + }; + super::gcd(a, b) == g + } + } +} diff --git a/src/uu/factor/src/numeric/mod.rs b/src/uu/factor/src/numeric/mod.rs new file mode 100644 index 000000000..d73195b81 --- /dev/null +++ b/src/uu/factor/src/numeric/mod.rs @@ -0,0 +1,17 @@ +// * This file is part of the uutils coreutils package. +// * +// * (c) 2020 nicoo +// * +// * For the full copyright and license information, please view the LICENSE file +// * that was distributed with this source code. + +mod gcd; +pub use gcd::gcd; + +mod traits; + +mod modular_inverse; +pub(crate) use modular_inverse::modular_inverse; + +mod montgomery; +pub(crate) use montgomery::{Arithmetic, Montgomery}; diff --git a/src/uu/factor/src/numeric/modular_inverse.rs b/src/uu/factor/src/numeric/modular_inverse.rs new file mode 100644 index 000000000..5310a1886 --- /dev/null +++ b/src/uu/factor/src/numeric/modular_inverse.rs @@ -0,0 +1,62 @@ +// * This file is part of the uutils coreutils package. +// * +// * (c) 2015 Wiktor Kuropatwa +// * (c) 2020 nicoo +// * +// * For the full copyright and license information, please view the LICENSE file +// * that was distributed with this source code. + +use super::traits::Int; + +// extended Euclid algorithm +// precondition: a is odd +pub(crate) fn modular_inverse(a: T) -> T { + let zero = T::zero(); + let one = T::one(); + debug_assert!(a % (one + one) == one, "{:?} is not odd", a); + + let mut t = zero; + let mut newt = one; + let mut r = zero; + let mut newr = a; + + while newr != zero { + let quot = if r == zero { + // special case when we're just starting out + // This works because we know that + // a does not divide 2^64, so floor(2^64 / a) == floor((2^64-1) / a); + T::max_value() + } else { + r + } / newr; + + let newtp = t.wrapping_sub(".wrapping_mul(&newt)); + t = newt; + newt = newtp; + + let newrp = r.wrapping_sub(".wrapping_mul(&newr)); + r = newr; + newr = newrp; + } + + debug_assert_eq!(r, one); + t +} + +#[cfg(test)] +mod tests { + use super::{super::traits::Int, *}; + use crate::parametrized_check; + + fn test_inverter() { + // All odd integers from 1 to 20 000 + let one = T::from(1).unwrap(); + let two = T::from(2).unwrap(); + let mut test_values = (0..10_000) + .map(|i| T::from(i).unwrap()) + .map(|i| two * i + one); + + assert!(test_values.all(|x| x.wrapping_mul(&modular_inverse(x)) == one)); + } + parametrized_check!(test_inverter); +} diff --git a/src/uu/factor/src/numeric.rs b/src/uu/factor/src/numeric/montgomery.rs similarity index 57% rename from src/uu/factor/src/numeric.rs rename to src/uu/factor/src/numeric/montgomery.rs index 814d72612..730477bb5 100644 --- a/src/uu/factor/src/numeric.rs +++ b/src/uu/factor/src/numeric/montgomery.rs @@ -1,58 +1,14 @@ // * This file is part of the uutils coreutils package. // * -// * (c) 2015 Wiktor Kuropatwa -// * (c) 2020 nicoo +// * (c) 2020 Alex Lyon +// * (c) 2020 nicoo // * // * For the full copyright and license information, please view the LICENSE file // * that was distributed with this source code. -use num_traits::{ - identities::{One, Zero}, - int::PrimInt, - ops::wrapping::{WrappingMul, WrappingNeg, WrappingSub}, -}; -use std::cmp::min; -use std::fmt::{Debug, Display}; -use std::mem::swap; - -// This is incorrectly reported as dead code, -// presumably when included in build.rs. -#[allow(dead_code)] -pub fn gcd(mut n: u64, mut m: u64) -> u64 { - // Stein's binary GCD algorithm - // Base cases: gcd(n, 0) = gcd(0, n) = n - if n == 0 { - return m; - } else if m == 0 { - return n; - } - - // Extract common factor-2: gcd(2ⁱ n, 2ⁱ m) = 2ⁱ gcd(n, m) - // and reducing until odd gcd(2ⁱ n, m) = gcd(n, m) if m is odd - let k = { - let k_n = n.trailing_zeros(); - let k_m = m.trailing_zeros(); - n >>= k_n; - m >>= k_m; - min(k_n, k_m) - }; - - loop { - // Invariant: n odd - debug_assert!(n % 2 == 1, "n = {} is even", n); - - if n > m { - swap(&mut n, &mut m); - } - m -= n; - - if m == 0 { - return n << k; - } - - m >>= m.trailing_zeros(); - } -} +use super::traits::{DoubleInt, Int, OverflowingAdd}; +use super::*; +use num_traits::identities::{One, Zero}; pub(crate) trait Arithmetic: Copy + Sized { // The type of integers mod m, in some opaque representation @@ -223,161 +179,10 @@ impl Arithmetic for Montgomery { } } -// NOTE: Trait can be removed once num-traits adds a similar one; -// see https://github.com/rust-num/num-traits/issues/168 -pub(crate) trait OverflowingAdd: Sized { - fn overflowing_add_(self, n: Self) -> (Self, bool); -} -macro_rules! overflowing { - ($x:ty) => { - impl OverflowingAdd for $x { - fn overflowing_add_(self, n: Self) -> (Self, bool) { - self.overflowing_add(n) - } - } - }; -} -overflowing!(u32); -overflowing!(u64); -overflowing!(u128); - -pub(crate) trait Int: - Display + Debug + PrimInt + OverflowingAdd + WrappingNeg + WrappingSub + WrappingMul -{ - fn as_u64(&self) -> u64; - fn from_u64(n: u64) -> Self; - - #[cfg(debug_assertions)] - fn as_u128(&self) -> u128; -} - -pub(crate) trait DoubleInt: Int { - /// An integer type with twice the width of `Self`. - /// In particular, multiplications (of `Int` values) can be performed in - /// `Self::DoubleWidth` without possibility of overflow. - type DoubleWidth: Int; - - fn as_double_width(self) -> Self::DoubleWidth; - fn from_double_width(n: Self::DoubleWidth) -> Self; -} - -macro_rules! int { - ( $x:ty ) => { - impl Int for $x { - fn as_u64(&self) -> u64 { - *self as u64 - } - fn from_u64(n: u64) -> Self { - n as _ - } - #[cfg(debug_assertions)] - fn as_u128(&self) -> u128 { - *self as u128 - } - } - }; -} -macro_rules! double_int { - ( $x:ty, $y:ty ) => { - int!($x); - impl DoubleInt for $x { - type DoubleWidth = $y; - - fn as_double_width(self) -> $y { - self as _ - } - fn from_double_width(n: $y) -> $x { - n as _ - } - } - }; -} - -double_int!(u32, u64); -double_int!(u64, u128); -int!(u128); - -// extended Euclid algorithm -// precondition: a is odd -pub(crate) fn modular_inverse(a: T) -> T { - let zero = T::zero(); - let one = T::one(); - debug_assert!(a % (one + one) == one, "{:?} is not odd", a); - - let mut t = zero; - let mut newt = one; - let mut r = zero; - let mut newr = a; - - while newr != zero { - let quot = if r == zero { - // special case when we're just starting out - // This works because we know that - // a does not divide 2^64, so floor(2^64 / a) == floor((2^64-1) / a); - T::max_value() - } else { - r - } / newr; - - let newtp = t.wrapping_sub(".wrapping_mul(&newt)); - t = newt; - newt = newtp; - - let newrp = r.wrapping_sub(".wrapping_mul(&newr)); - r = newr; - newr = newrp; - } - - debug_assert_eq!(r, one); - t -} - #[cfg(test)] mod tests { use super::*; - use quickcheck::quickcheck; - - quickcheck! { - fn gcd(a: u64, b: u64) -> bool { - // Test against the Euclidean algorithm - let g = { - let (mut a, mut b) = (a, b); - while b > 0 { - a %= b; - swap(&mut a, &mut b); - } - a - }; - super::gcd(a, b) == g - } - } - - macro_rules! parametrized_check { - ( $f:ident ) => { - paste::item! { - #[test] - fn [< $f _ u32 >]() { - $f::() - } - #[test] - fn [< $f _ u64 >]() { - $f::() - } - } - }; - } - - fn test_inverter() { - // All odd integers from 1 to 20 000 - let one = T::from(1).unwrap(); - let two = T::from(2).unwrap(); - let mut test_values = (0..10_000) - .map(|i| T::from(i).unwrap()) - .map(|i| two * i + one); - - assert!(test_values.all(|x| x.wrapping_mul(&modular_inverse(x)) == one)); - } - parametrized_check!(test_inverter); + use crate::parametrized_check; fn test_add() { for n in 0..100 { diff --git a/src/uu/factor/src/numeric/traits.rs b/src/uu/factor/src/numeric/traits.rs new file mode 100644 index 000000000..50f5ab5c1 --- /dev/null +++ b/src/uu/factor/src/numeric/traits.rs @@ -0,0 +1,105 @@ +// * This file is part of the uutils coreutils package. +// * +// * (c) 2020 Alex Lyon +// * (c) 2020 nicoo +// * +// * For the full copyright and license information, please view the LICENSE file +// * that was distributed with this source code. + +use num_traits::{ + int::PrimInt, + ops::wrapping::{WrappingMul, WrappingNeg, WrappingSub}, +}; +use std::fmt::{Debug, Display}; + +// NOTE: Trait can be removed once num-traits adds a similar one; +// see https://github.com/rust-num/num-traits/issues/168 +pub(crate) trait OverflowingAdd: Sized { + fn overflowing_add_(self, n: Self) -> (Self, bool); +} + +macro_rules! overflowing { + ($x:ty) => { + impl OverflowingAdd for $x { + fn overflowing_add_(self, n: Self) -> (Self, bool) { + self.overflowing_add(n) + } + } + }; +} +overflowing!(u32); +overflowing!(u64); +overflowing!(u128); + +pub(crate) trait Int: + Display + Debug + PrimInt + OverflowingAdd + WrappingNeg + WrappingSub + WrappingMul +{ + fn as_u64(&self) -> u64; + fn from_u64(n: u64) -> Self; + + #[cfg(debug_assertions)] + fn as_u128(&self) -> u128; +} + +pub(crate) trait DoubleInt: Int { + /// An integer type with twice the width of `Self`. + /// In particular, multiplications (of `Int` values) can be performed in + /// `Self::DoubleWidth` without possibility of overflow. + type DoubleWidth: Int; + + fn as_double_width(self) -> Self::DoubleWidth; + fn from_double_width(n: Self::DoubleWidth) -> Self; +} + +macro_rules! int { + ( $x:ty ) => { + impl Int for $x { + fn as_u64(&self) -> u64 { + *self as u64 + } + fn from_u64(n: u64) -> Self { + n as _ + } + #[cfg(debug_assertions)] + fn as_u128(&self) -> u128 { + *self as u128 + } + } + }; +} +macro_rules! double_int { + ( $x:ty, $y:ty ) => { + int!($x); + impl DoubleInt for $x { + type DoubleWidth = $y; + + fn as_double_width(self) -> $y { + self as _ + } + fn from_double_width(n: $y) -> $x { + n as _ + } + } + }; +} +double_int!(u32, u64); +double_int!(u64, u128); +int!(u128); + +/// Helper macro for instantiating tests over u32 and u64 +#[cfg(test)] +#[macro_export] +macro_rules! parametrized_check { + ( $f:ident ) => { + paste::item! { + #[test] + fn [< $f _ u32 >]() { + $f::() + } + #[test] + fn [< $f _ u64 >]() { + $f::() + } + } + }; +}