From 824103769028240098c0f30b76fb3e067fedea3d Mon Sep 17 00:00:00 2001 From: nicoo Date: Sun, 24 May 2020 15:41:23 +0200 Subject: [PATCH] factor::miller_rabin: Extract dividers from the primality test Another 36% improvement. --- src/uu/factor/src/miller_rabin.rs | 33 +++++++++++++++++++++++++------ src/uu/factor/src/rho.rs | 20 ++++++++++++++----- 2 files changed, 42 insertions(+), 11 deletions(-) diff --git a/src/uu/factor/src/miller_rabin.rs b/src/uu/factor/src/miller_rabin.rs index 997cd1cd6..42c19b24f 100644 --- a/src/uu/factor/src/miller_rabin.rs +++ b/src/uu/factor/src/miller_rabin.rs @@ -4,14 +4,29 @@ use crate::numeric::*; // discovered by Jim Sinclair on 2011-04-20, see miller-rabin.appspot.com const BASIS: [u64; 7] = [2, 325, 9375, 28178, 450775, 9780504, 1795265022]; +#[derive(Eq, PartialEq)] +pub(crate) enum Result { + Prime, + Pseudoprime, + Composite(u64), +} + +impl Result { + pub(crate) fn is_prime(&self) -> bool { + return *self == Result::Prime; + } +} + // Deterministic Miller-Rabin primality-checking algorithm, adapted to extract // (some) dividers; it will fail to factor strong pseudoprimes. -pub(crate) fn is_prime(n: u64) -> bool { +pub(crate) fn test(n: u64) -> Result { + use self::Result::*; + if n < 2 { - return false; + return Pseudoprime; } if n % 2 == 0 { - return n == 2; + return if n == 2 { Prime } else { Composite(2) }; } let d = (n - 1).trailing_zeros(); @@ -30,7 +45,7 @@ pub(crate) fn is_prime(n: u64) -> bool { } if pow(x, n - 1, n, mul) != 1 { - return false; + return Pseudoprime; } x = pow(x, r, n, mul); if x == 1 || x == n - 1 { @@ -40,7 +55,7 @@ pub(crate) fn is_prime(n: u64) -> bool { loop { let y = mul(x, x, n); if y == 1 { - return false; + return Composite(gcd(x - 1, n)); } if y == n - 1 { // This basis element is not a witness of `n` being composite. @@ -51,5 +66,11 @@ pub(crate) fn is_prime(n: u64) -> bool { } } - true + Prime +} + +// Used by build.rs' tests +#[allow(dead_code)] +pub(crate) fn is_prime(n: u64) -> bool { + test(n).is_prime() } diff --git a/src/uu/factor/src/rho.rs b/src/uu/factor/src/rho.rs index 947eac478..e1a009810 100644 --- a/src/uu/factor/src/rho.rs +++ b/src/uu/factor/src/rho.rs @@ -1,3 +1,4 @@ +use crate::miller_rabin::Result::*; use crate::{miller_rabin, Factors}; use numeric::*; use rand::distributions::{Distribution, Uniform}; @@ -39,16 +40,25 @@ fn find_divisor(num: u64) -> u64 { } } -pub(crate) fn factor(num: u64) -> Factors { +pub(crate) fn factor(mut num: u64) -> Factors { let mut factors = Factors::new(); if num == 1 { return factors; } - if miller_rabin::is_prime(num) { - factors.push(num); - return factors; - } + match miller_rabin::test(num) { + Prime => { + factors.push(num); + return factors; + } + + Composite(d) => { + num = num / d; + factors *= factor(d); + } + + Pseudoprime => {} + }; let divisor = find_divisor(num); factors *= factor(divisor);