1
Fork 0
mirror of https://github.com/RGBCube/uutils-coreutils synced 2025-07-31 04:57:45 +00:00

factor: Move the Miller-Rabin primality test to its own module.

This commit is contained in:
nicoo 2020-05-24 15:41:23 +02:00
parent 74054feb94
commit e3ecc81d97
5 changed files with 61 additions and 55 deletions

View file

@ -26,7 +26,7 @@ use std::path::Path;
use std::u64::MAX as MAX_U64;
#[cfg(test)]
use numeric::is_prime;
use miller_rabin::is_prime;
#[cfg(test)]
#[path = "src/numeric.rs"]

View file

@ -23,6 +23,7 @@ use std::fmt;
use std::io::{stdin, BufRead};
use std::ops;
mod miller_rabin;
mod numeric;
mod rho;
mod table;

View file

@ -0,0 +1,52 @@
use crate::numeric::*;
fn witness(mut a: u64, exponent: u64, m: u64) -> bool {
if a == 0 {
return false;
}
let mul = if m < 1 << 63 {
sm_mul as fn(u64, u64, u64) -> u64
} else {
big_mul as fn(u64, u64, u64) -> u64
};
if pow(a, m - 1, m, mul) != 1 {
return true;
}
a = pow(a, exponent, m, mul);
if a == 1 {
return false;
}
loop {
if a == 1 {
return true;
}
if a == m - 1 {
return false;
}
a = mul(a, a, m);
}
}
// Deterministic Miller-Rabin primality-checking algorithm, adapted to extract
// (some) dividers; it will fail to factor strong pseudoprimes.
pub(crate) fn is_prime(num: u64) -> bool {
if num < 2 {
return false;
}
if num % 2 == 0 {
return num == 2;
}
let mut exponent = num - 1;
while exponent & 1 == 0 {
exponent >>= 1;
}
// These witnesses detect all composites up to at least 2^64.
// Discovered by Jim Sinclair, according to http://miller-rabin.appspot.com
let witnesses = [2, 325, 9_375, 28_178, 450_775, 9_780_504, 1_795_265_022];
!witnesses
.iter()
.any(|&wit| witness(wit % num, exponent, num))
}

View file

@ -78,7 +78,7 @@ pub fn big_mul(mut a: u64, mut b: u64, m: u64) -> u64 {
}
// computes a.pow(b) % m
fn pow(mut a: u64, mut b: u64, m: u64, mul: fn(u64, u64, u64) -> u64) -> u64 {
pub(crate) fn pow(mut a: u64, mut b: u64, m: u64, mul: fn(u64, u64, u64) -> u64) -> u64 {
let mut result = 1;
while b > 0 {
if b & 1 != 0 {
@ -89,53 +89,3 @@ fn pow(mut a: u64, mut b: u64, m: u64, mul: fn(u64, u64, u64) -> u64) -> u64 {
}
result
}
fn witness(mut a: u64, exponent: u64, m: u64) -> bool {
if a == 0 {
return false;
}
let mul = if m < 1 << 63 {
sm_mul as fn(u64, u64, u64) -> u64
} else {
big_mul as fn(u64, u64, u64) -> u64
};
if pow(a, m - 1, m, mul) != 1 {
return true;
}
a = pow(a, exponent, m, mul);
if a == 1 {
return false;
}
loop {
if a == 1 {
return true;
}
if a == m - 1 {
return false;
}
a = mul(a, a, m);
}
}
// uses deterministic (i.e., fixed witness set) Miller-Rabin test
pub fn is_prime(num: u64) -> bool {
if num < 2 {
return false;
}
if num % 2 == 0 {
return num == 2;
}
let mut exponent = num - 1;
while exponent & 1 == 0 {
exponent >>= 1;
}
// These witnesses detect all composites up to at least 2^64.
// Discovered by Jim Sinclair, according to http://miller-rabin.appspot.com
let witnesses = [2, 325, 9_375, 28_178, 450_775, 9_780_504, 1_795_265_022];
!witnesses
.iter()
.any(|&wit| witness(wit % num, exponent, num))
}

View file

@ -1,4 +1,4 @@
use crate::Factors;
use crate::{miller_rabin, Factors};
use numeric::*;
use rand::distributions::{Distribution, Uniform};
use rand::rngs::SmallRng;
@ -41,8 +41,11 @@ fn find_divisor(num: u64) -> u64 {
pub(crate) fn factor(num: u64) -> Factors {
let mut factors = Factors::new();
if num == 1 { return factors; }
if is_prime(num) {
if num == 1 {
return factors;
}
if miller_rabin::is_prime(num) {
factors.push(num);
return factors;
}