From 032886ac829f4b086bc571b6634cabf2d34860bc Mon Sep 17 00:00:00 2001 From: Evan Hines Date: Mon, 6 Dec 2021 20:32:35 -0500 Subject: [PATCH 1/4] Bump dependencies and move allocations outside of the hotpath Clean up some old code comments that were distracting Use a new factorization algorithm so that the tests actually finish. This will make working much less tedious. --- Cargo.toml | 6 +-- README.md | 12 ------ src/lib.rs | 115 +++++++++++++++++++++++------------------------------ 3 files changed, 53 insertions(+), 80 deletions(-) diff --git a/Cargo.toml b/Cargo.toml index ddbdc50..ed35a3e 100644 --- a/Cargo.toml +++ b/Cargo.toml @@ -20,10 +20,10 @@ travis-ci = { repository = "AtropineTears/num-primes", branch = "master" } [dependencies] # num num = { version = "0.4.0", default-features = false } -num-traits = "0.2.11" -num-bigint = { version = "0.2.6", features = ["rand"] } +num-traits = "0.2.14" +num-bigint = { version = "0.4.3", features = ["rand"] } # random -rand = "0.5.6" +rand = "0.8.4" # Logging log = "0.4.14" diff --git a/README.md b/README.md index 0519e12..99a66eb 100644 --- a/README.md +++ b/README.md @@ -11,10 +11,6 @@ It takes full advantage of the [num](https://crates.io/crates/num) crate on **st * Read the [License](#license) * Read the [Contribution](#contribution) -## Notice - -Please note there is a critical bug in this program that I cannot seem to fix where it marks some prime numbers as not prime. It is in the miller-rabin implementation and I cannot seem to fix it. If anyone is up to it, feel free to look through the issues tab for information about the bug and submit a PR if you find a fix. - ## Usage Add this to your `Cargo.toml`: @@ -24,10 +20,6 @@ Add this to your `Cargo.toml`: num-primes = "0.2.0" ``` -## Warning - -There is currently a major bug in `is_prime()` and `is_composite()` that makes some values return wrong. For example, a prime can sometimes be marked as composite unless it was generated as they use the same tests to test for primality. - ## How To Use There are three main structs that are included in this library @@ -98,10 +90,6 @@ fn main(){ } ``` -## Verification - -WARNING: There is currently a bug that makes verification of certain prime numbers fail. Be careful when using this feature. - ### Verify Composite Number This function will verify whether a `BigUint` type is a **composite** by returning a boolean value. diff --git a/src/lib.rs b/src/lib.rs index 7de7569..028ebf8 100644 --- a/src/lib.rs +++ b/src/lib.rs @@ -8,8 +8,10 @@ extern crate num_bigint as bigint; use core::ops::Sub; +use core::convert::TryInto; +use num::integer::gcd; use num::Integer; -pub use bigint::{BigUint,RandBigInt}; +pub use bigint::{BigUint, BigInt, RandBigInt}; use num_traits::{Zero, One}; use num_traits::*; @@ -19,7 +21,7 @@ use log::info; // Settings // NIST recomends 5 rounds for miller rabin. This implementation does 8. Apple uses 16. Three iterations has a probability of 2^80 of failing -const MILLER_RABIN_ROUNDS: usize = 8usize; +const MILLER_RABIN_ROUNDS: usize = 64usize; /// # Generator @@ -97,11 +99,13 @@ impl Generator { /// ``` pub fn new_composite(n: usize) -> BigUint { let mut rng = rand::thread_rng(); + let n = n.try_into().unwrap(); + loop { // Make mutable and set LSB and MSB - let candidate: BigUint = rng.gen_biguint(n); - //candidate.set_bit(0, true); - //candidate.set_bit((n-1) as u32, true); + let mut candidate: BigUint = rng.gen_biguint(n); + candidate.set_bit(0, true); + candidate.set_bit((n-1) as u64, true); if is_prime(&candidate) == false { return candidate; } @@ -123,7 +127,7 @@ impl Generator { /// ``` pub fn new_uint(n: usize) -> BigUint { let mut rng = rand::thread_rng(); - return rng.gen_biguint(n); + return rng.gen_biguint(n.try_into().unwrap()); } /// # Generate Prime Number @@ -146,13 +150,14 @@ impl Generator { /// ``` pub fn new_prime(n: usize) -> BigUint { let mut rng = rand::thread_rng(); - + let n = n.try_into().unwrap(); + loop { // Make mutable and set LSB and MSB - let candidate: BigUint = rng.gen_biguint(n); - - //candidate.set_bit(0, true); - //candidate.set_bit((n-1) as u32, true); + let mut candidate: BigUint = rng.gen_biguint(n); + + candidate.set_bit(0, true); + candidate.set_bit((n-1) as u64, true); if is_prime(&candidate) == true { return candidate; @@ -172,11 +177,15 @@ impl Generator { /// ``` pub fn safe_prime(n: usize) -> BigUint { let mut rng = rand::thread_rng(); + let n = n.try_into().unwrap(); + loop { // Make mutable and set LSB and MSB - let candidate: BigUint = rng.gen_biguint(n); - //candidate.set_bit(0, true); - //candidate.set_bit((n-1) as u32, true); + let mut candidate: BigUint = rng.gen_biguint(n); + + candidate.set_bit(0, true); + candidate.set_bit((n-1) as u64, true); + if is_prime(&candidate) == true { if is_safe_prime(&candidate) == true { // checks with (p-1/n) @@ -219,7 +228,7 @@ impl Verification { /// let x: BigUint = BigUint::from_u64(7u64).unwrap(); /// /// // Verify Its A Smooth Number - /// let result: bool = Verification::is_smooth_number(&x,31.0,5); + /// let result: bool = Verification::is_very_smooth_number(&x,31.0,5); /// /// println!("Is A {} Smooth Number: {}",x,result); /// } @@ -255,34 +264,32 @@ impl Factorization { if is_prime(&n) { return Some(n) } + + let one = BigInt::one(); + let two: BigInt = &one + &one; - let one = BigUint::one(); - let two = &one + &one; - - // STEP 1 | n divided by 2 - while n.is_even() { - n = n / &two; - } - - // STEP 2 | 3..sqrt(n) | Divide i by n. On failure, add 2 to i - let n_sqrt = n.sqrt().to_usize().unwrap(); - - for mut i in 3..n_sqrt { - while n.divides(&BigUint::from(i)) { - n = n / BigUint::from(i); - } - i = i + 2usize; + let mut x = two.clone(); + let mut y = two.clone(); + let mut d = one.clone(); + let n_int = n.clone().try_into().unwrap(); + + let g = |x: &BigInt| x.modpow(&two, &n_int); + + while d == one { + x = g(&x); + y = g(&g(&y)); + d = gcd((&x - &y).abs(), &n_int - 0); } - // Step 3 - if n > two { - return Some(n) + let d = d.try_into().unwrap(); + + if d == n { + None } else { - return None + Some(d) } - } } @@ -318,7 +325,6 @@ fn fermat(candidate: &BigUint) -> bool { let result = random.modpow(&(exponent), candidate); - //let result = random.pow_mod(&(candidate - One::one()), candidate); result == One::one() } @@ -327,61 +333,41 @@ fn miller_rabin(candidate: &BigUint, limit: usize) -> bool { let zero: BigUint = Zero::zero(); let one = BigUint::one(); let two: BigUint = &one + &one; - let two_2 = num_bigint::ToBigUint::to_biguint(&2).unwrap(); // Check Whether Candidate Is 2 (which is prime) - if candidate == &two_2 { + if candidate == &two { return true } - - // Check Whether Candidate Is Even - /* - if candidate.mod(two) { - return false - } - */ let (d,s) = rewrite(&candidate); let step = s.sub(&one).to_usize().unwrap(); + let one_usize = one.to_usize().unwrap(); + let zero_usize = zero.to_usize().unwrap(); + let c_minus_1 = candidate-&one; let mut rng = rand::thread_rng(); for _i in 0..limit { // Generate Random Number between [2,n-1) | Exclusive End Range; Uses (n-1), not (n-2) - let a = rng.gen_biguint_range(&two, &(candidate-&one)); - - // Reference Implementation - // Pretty sure `sample_range()` has an inclusive end - //let basis = Int::sample_range(&two, &(candidate-&two)); + let a = rng.gen_biguint_range(&two, &c_minus_1); // (a^d mod n) let mut x = a.modpow(&d, &candidate); - // Reference Implementation - //let mut y = Int::modpow(&basis, &d, candidate); - - if x == one || x == (candidate - &one) { + if x == one || x == c_minus_1 { continue // return true } else { - // Convert To Usizes For Loop - // step = (s - 1) - let one_usize = one.to_usize().unwrap(); - let zero_usize = zero.to_usize().unwrap(); - let mut break_early = false; // Issue #1 | Changed one_usize to zero_usize; step (s-1) was equal to iterations-1 and therefore needed an extra iteration for _ in zero_usize..step { x = x.modpow(&two,candidate); - // Reference Implementation - //y = Int::modpow(&y, &two, candidate); - if x == one { return false } - else if x == (candidate - BigUint::one()) { + else if x == c_minus_1 { break_early = true; break; } @@ -450,7 +436,6 @@ fn is_prime(candidate: &BigUint) -> bool { else { return true } - return true } // p = 2q + 1 From 92ae30f0858ec0c73188ebd64fa3ab573fd9bb0e Mon Sep 17 00:00:00 2001 From: Evan Hines Date: Mon, 6 Dec 2021 23:27:00 -0500 Subject: [PATCH 2/4] Fix bug in small prime factor code. For some reason it would check modulus the prime and return maybe prime although we can just as easily check if it actually that prime right here. This was causing small primes to go through the more expensive primality tests although we didn't need to. Add a test for factoring two prime numbers --- src/lib.rs | 86 +++++++++++++++++++++++++----------------- tests/factorization.rs | 19 +++++----- 2 files changed, 61 insertions(+), 44 deletions(-) diff --git a/src/lib.rs b/src/lib.rs index 028ebf8..f614b1a 100644 --- a/src/lib.rs +++ b/src/lib.rs @@ -1,4 +1,4 @@ -#![no_std] +//#![no_std] #![forbid(unsafe_code)] extern crate num; @@ -102,9 +102,7 @@ impl Generator { let n = n.try_into().unwrap(); loop { - // Make mutable and set LSB and MSB let mut candidate: BigUint = rng.gen_biguint(n); - candidate.set_bit(0, true); candidate.set_bit((n-1) as u64, true); if is_prime(&candidate) == false { return candidate; @@ -153,7 +151,6 @@ impl Generator { let n = n.try_into().unwrap(); loop { - // Make mutable and set LSB and MSB let mut candidate: BigUint = rng.gen_biguint(n); candidate.set_bit(0, true); @@ -180,7 +177,6 @@ impl Generator { let n = n.try_into().unwrap(); loop { - // Make mutable and set LSB and MSB let mut candidate: BigUint = rng.gen_biguint(n); candidate.set_bit(0, true); @@ -264,38 +260,55 @@ impl Factorization { if is_prime(&n) { return Some(n) } - + + let mut n: BigInt = n.try_into().unwrap(); + + let mut f = BigInt::zero(); + let one = BigInt::one(); let two: BigInt = &one + &one; - let mut x = two.clone(); - let mut y = two.clone(); - let mut d = one.clone(); - let n_int = n.clone().try_into().unwrap(); + // This while loop will run until n is completely factorized, + // thus giving the largest factor. If you only care about A factor, + // this can be rewritten as below for much better performance. + while n > one { + // while f == BigInt::zero() { - let g = |x: &BigInt| x.modpow(&two, &n_int); - while d == one { - x = g(&x); - y = g(&g(&y)); - d = gcd((&x - &y).abs(), &n_int - 0); - } + let mut x = two.clone(); + let mut y = two.clone(); + let mut d = one.clone(); + let n_int = n.clone().try_into().unwrap(); - let d = d.try_into().unwrap(); + let g = |x: &BigInt| (x.modpow(&two, &n_int) + 1) % &n_int; - if d == n { - None - } - else { - Some(d) + while d == one { + x = g(&x); + y = g(&g(&y)); + d = gcd((&x - &y).abs(), &n_int - 0); + } + + n /= &d; + + if &d > &f { + f = d; + } } + return Some(f.try_into().unwrap()); } } +#[derive(PartialEq)] +enum SmallPrimeResult { + Prime, + NotPrime, + MaybePrime +} + // if true, then is not prime // if false, then maybe prime -fn div_small_primes(numb: &BigUint) -> bool { +fn div_small_primes(numb: &BigUint) -> SmallPrimeResult { let zero: BigUint = Zero::zero(); let one: BigUint = One::one(); @@ -303,15 +316,19 @@ fn div_small_primes(numb: &BigUint) -> bool { static SMALL_PRIMES: [u32; 2048] = [2,3,5,7,11,13,17,19,23,29,31,37,41,43,47,53,59,61,67,71,73,79,83,89,97,101,103,107,109,113,127,131,137,139,149,151,157,163,167,173,179,181,191,193,197,199,211,223,227,229,233,239,241,251,257,263,269,271,277,281,283,293,307,311,313,317,331,337,347,349,353,359,367,373,379,383,389,397,401,409,419,421,431,433,439,443,449,457,461,463,467,479,487,491,499,503,509,521,523,541,547,557,563,569,571,577,587,593,599,601,607,613,617,619,631,641,643,647,653,659,661,673,677,683,691,701,709,719,727,733,739,743,751,757,761,769,773,787,797,809,811,821,823,827,829,839,853,857,859,863,877,881,883,887,907,911,919,929,937,941,947,953,967,971,977,983,991,997,1009,1013,1019,1021,1031,1033,1039,1049,1051,1061,1063,1069,1087,1091,1093,1097,1103,1109,1117,1123,1129,1151,1153,1163,1171,1181,1187,1193,1201,1213,1217,1223,1229,1231,1237,1249,1259,1277,1279,1283,1289,1291,1297,1301,1303,1307,1319,1321,1327,1361,1367,1373,1381,1399,1409,1423,1427,1429,1433,1439,1447,1451,1453,1459,1471,1481,1483,1487,1489,1493,1499,1511,1523,1531,1543,1549,1553,1559,1567,1571,1579,1583,1597,1601,1607,1609,1613,1619,1621,1627,1637,1657,1663,1667,1669,1693,1697,1699,1709,1721,1723,1733,1741,1747,1753,1759,1777,1783,1787,1789,1801,1811,1823,1831,1847,1861,1867,1871,1873,1877,1879,1889,1901,1907,1913,1931,1933,1949,1951,1973,1979,1987,1993,1997,1999,2003,2011,2017,2027,2029,2039,2053,2063,2069,2081,2083,2087,2089,2099,2111,2113,2129,2131,2137,2141,2143,2153,2161,2179,2203,2207,2213,2221,2237,2239,2243,2251,2267,2269,2273,2281,2287,2293,2297,2309,2311,2333,2339,2341,2347,2351,2357,2371,2377,2381,2383,2389,2393,2399,2411,2417,2423,2437,2441,2447,2459,2467,2473,2477,2503,2521,2531,2539,2543,2549,2551,2557,2579,2591,2593,2609,2617,2621,2633,2647,2657,2659,2663,2671,2677,2683,2687,2689,2693,2699,2707,2711,2713,2719,2729,2731,2741,2749,2753,2767,2777,2789,2791,2797,2801,2803,2819,2833,2837,2843,2851,2857,2861,2879,2887,2897,2903,2909,2917,2927,2939,2953,2957,2963,2969,2971,2999,3001,3011,3019,3023,3037,3041,3049,3061,3067,3079,3083,3089,3109,3119,3121,3137,3163,3167,3169,3181,3187,3191,3203,3209,3217,3221,3229,3251,3253,3257,3259,3271,3299,3301,3307,3313,3319,3323,3329,3331,3343,3347,3359,3361,3371,3373,3389,3391,3407,3413,3433,3449,3457,3461,3463,3467,3469,3491,3499,3511,3517,3527,3529,3533,3539,3541,3547,3557,3559,3571,3581,3583,3593,3607,3613,3617,3623,3631,3637,3643,3659,3671,3673,3677,3691,3697,3701,3709,3719,3727,3733,3739,3761,3767,3769,3779,3793,3797,3803,3821,3823,3833,3847,3851,3853,3863,3877,3881,3889,3907,3911,3917,3919,3923,3929,3931,3943,3947,3967,3989,4001,4003,4007,4013,4019,4021,4027,4049,4051,4057,4073,4079,4091,4093,4099,4111,4127,4129,4133,4139,4153,4157,4159,4177,4201,4211,4217,4219,4229,4231,4241,4243,4253,4259,4261,4271,4273,4283,4289,4297,4327,4337,4339,4349,4357,4363,4373,4391,4397,4409,4421,4423,4441,4447,4451,4457,4463,4481,4483,4493,4507,4513,4517,4519,4523,4547,4549,4561,4567,4583,4591,4597,4603,4621,4637,4639,4643,4649,4651,4657,4663,4673,4679,4691,4703,4721,4723,4729,4733,4751,4759,4783,4787,4789,4793,4799,4801,4813,4817,4831,4861,4871,4877,4889,4903,4909,4919,4931,4933,4937,4943,4951,4957,4967,4969,4973,4987,4993,4999,5003,5009,5011,5021,5023,5039,5051,5059,5077,5081,5087,5099,5101,5107,5113,5119,5147,5153,5167,5171,5179,5189,5197,5209,5227,5231,5233,5237,5261,5273,5279,5281,5297,5303,5309,5323,5333,5347,5351,5381,5387,5393,5399,5407,5413,5417,5419,5431,5437,5441,5443,5449,5471,5477,5479,5483,5501,5503,5507,5519,5521,5527,5531,5557,5563,5569,5573,5581,5591,5623,5639,5641,5647,5651,5653,5657,5659,5669,5683,5689,5693,5701,5711,5717,5737,5741,5743,5749,5779,5783,5791,5801,5807,5813,5821,5827,5839,5843,5849,5851,5857,5861,5867,5869,5879,5881,5897,5903,5923,5927,5939,5953,5981,5987,6007,6011,6029,6037,6043,6047,6053,6067,6073,6079,6089,6091,6101,6113,6121,6131,6133,6143,6151,6163,6173,6197,6199,6203,6211,6217,6221,6229,6247,6257,6263,6269,6271,6277,6287,6299,6301,6311,6317,6323,6329,6337,6343,6353,6359,6361,6367,6373,6379,6389,6397,6421,6427,6449,6451,6469,6473,6481,6491,6521,6529,6547,6551,6553,6563,6569,6571,6577,6581,6599,6607,6619,6637,6653,6659,6661,6673,6679,6689,6691,6701,6703,6709,6719,6733,6737,6761,6763,6779,6781,6791,6793,6803,6823,6827,6829,6833,6841,6857,6863,6869,6871,6883,6899,6907,6911,6917,6947,6949,6959,6961,6967,6971,6977,6983,6991,6997,7001,7013,7019,7027,7039,7043,7057,7069,7079,7103,7109,7121,7127,7129,7151,7159,7177,7187,7193,7207,7211,7213,7219,7229,7237,7243,7247,7253,7283,7297,7307,7309,7321,7331,7333,7349,7351,7369,7393,7411,7417,7433,7451,7457,7459,7477,7481,7487,7489,7499,7507,7517,7523,7529,7537,7541,7547,7549,7559,7561,7573,7577,7583,7589,7591,7603,7607,7621,7639,7643,7649,7669,7673,7681,7687,7691,7699,7703,7717,7723,7727,7741,7753,7757,7759,7789,7793,7817,7823,7829,7841,7853,7867,7873,7877,7879,7883,7901,7907,7919,7927,7933,7937,7949,7951,7963,7993,8009,8011,8017,8039,8053,8059,8069,8081,8087,8089,8093,8101,8111,8117,8123,8147,8161,8167,8171,8179,8191,8209,8219,8221,8231,8233,8237,8243,8263,8269,8273,8287,8291,8293,8297,8311,8317,8329,8353,8363,8369,8377,8387,8389,8419,8423,8429,8431,8443,8447,8461,8467,8501,8513,8521,8527,8537,8539,8543,8563,8573,8581,8597,8599,8609,8623,8627,8629,8641,8647,8663,8669,8677,8681,8689,8693,8699,8707,8713,8719,8731,8737,8741,8747,8753,8761,8779,8783,8803,8807,8819,8821,8831,8837,8839,8849,8861,8863,8867,8887,8893,8923,8929,8933,8941,8951,8963,8969,8971,8999,9001,9007,9011,9013,9029,9041,9043,9049,9059,9067,9091,9103,9109,9127,9133,9137,9151,9157,9161,9173,9181,9187,9199,9203,9209,9221,9227,9239,9241,9257,9277,9281,9283,9293,9311,9319,9323,9337,9341,9343,9349,9371,9377,9391,9397,9403,9413,9419,9421,9431,9433,9437,9439,9461,9463,9467,9473,9479,9491,9497,9511,9521,9533,9539,9547,9551,9587,9601,9613,9619,9623,9629,9631,9643,9649,9661,9677,9679,9689,9697,9719,9721,9733,9739,9743,9749,9767,9769,9781,9787,9791,9803,9811,9817,9829,9833,9839,9851,9857,9859,9871,9883,9887,9901,9907,9923,9929,9931,9941,9949,9967,9973,10007,10009,10037,10039,10061,10067,10069,10079,10091,10093,10099,10103,10111,10133,10139,10141,10151,10159,10163,10169,10177,10181,10193,10211,10223,10243,10247,10253,10259,10267,10271,10273,10289,10301,10303,10313,10321,10331,10333,10337,10343,10357,10369,10391,10399,10427,10429,10433,10453,10457,10459,10463,10477,10487,10499,10501,10513,10529,10531,10559,10567,10589,10597,10601,10607,10613,10627,10631,10639,10651,10657,10663,10667,10687,10691,10709,10711,10723,10729,10733,10739,10753,10771,10781,10789,10799,10831,10837,10847,10853,10859,10861,10867,10883,10889,10891,10903,10909,10937,10939,10949,10957,10973,10979,10987,10993,11003,11027,11047,11057,11059,11069,11071,11083,11087,11093,11113,11117,11119,11131,11149,11159,11161,11171,11173,11177,11197,11213,11239,11243,11251,11257,11261,11273,11279,11287,11299,11311,11317,11321,11329,11351,11353,11369,11383,11393,11399,11411,11423,11437,11443,11447,11467,11471,11483,11489,11491,11497,11503,11519,11527,11549,11551,11579,11587,11593,11597,11617,11621,11633,11657,11677,11681,11689,11699,11701,11717,11719,11731,11743,11777,11779,11783,11789,11801,11807,11813,11821,11827,11831,11833,11839,11863,11867,11887,11897,11903,11909,11923,11927,11933,11939,11941,11953,11959,11969,11971,11981,11987,12007,12011,12037,12041,12043,12049,12071,12073,12097,12101,12107,12109,12113,12119,12143,12149,12157,12161,12163,12197,12203,12211,12227,12239,12241,12251,12253,12263,12269,12277,12281,12289,12301,12323,12329,12343,12347,12373,12377,12379,12391,12401,12409,12413,12421,12433,12437,12451,12457,12473,12479,12487,12491,12497,12503,12511,12517,12527,12539,12541,12547,12553,12569,12577,12583,12589,12601,12611,12613,12619,12637,12641,12647,12653,12659,12671,12689,12697,12703,12713,12721,12739,12743,12757,12763,12781,12791,12799,12809,12821,12823,12829,12841,12853,12889,12893,12899,12907,12911,12917,12919,12923,12941,12953,12959,12967,12973,12979,12983,13001,13003,13007,13009,13033,13037,13043,13049,13063,13093,13099,13103,13109,13121,13127,13147,13151,13159,13163,13171,13177,13183,13187,13217,13219,13229,13241,13249,13259,13267,13291,13297,13309,13313,13327,13331,13337,13339,13367,13381,13397,13399,13411,13417,13421,13441,13451,13457,13463,13469,13477,13487,13499,13513,13523,13537,13553,13567,13577,13591,13597,13613,13619,13627,13633,13649,13669,13679,13681,13687,13691,13693,13697,13709,13711,13721,13723,13729,13751,13757,13759,13763,13781,13789,13799,13807,13829,13831,13841,13859,13873,13877,13879,13883,13901,13903,13907,13913,13921,13931,13933,13963,13967,13997,13999,14009,14011,14029,14033,14051,14057,14071,14081,14083,14087,14107,14143,14149,14153,14159,14173,14177,14197,14207,14221,14243,14249,14251,14281,14293,14303,14321,14323,14327,14341,14347,14369,14387,14389,14401,14407,14411,14419,14423,14431,14437,14447,14449,14461,14479,14489,14503,14519,14533,14537,14543,14549,14551,14557,14561,14563,14591,14593,14621,14627,14629,14633,14639,14653,14657,14669,14683,14699,14713,14717,14723,14731,14737,14741,14747,14753,14759,14767,14771,14779,14783,14797,14813,14821,14827,14831,14843,14851,14867,14869,14879,14887,14891,14897,14923,14929,14939,14947,14951,14957,14969,14983,15013,15017,15031,15053,15061,15073,15077,15083,15091,15101,15107,15121,15131,15137,15139,15149,15161,15173,15187,15193,15199,15217,15227,15233,15241,15259,15263,15269,15271,15277,15287,15289,15299,15307,15313,15319,15329,15331,15349,15359,15361,15373,15377,15383,15391,15401,15413,15427,15439,15443,15451,15461,15467,15473,15493,15497,15511,15527,15541,15551,15559,15569,15581,15583,15601,15607,15619,15629,15641,15643,15647,15649,15661,15667,15671,15679,15683,15727,15731,15733,15737,15739,15749,15761,15767,15773,15787,15791,15797,15803,15809,15817,15823,15859,15877,15881,15887,15889,15901,15907,15913,15919,15923,15937,15959,15971,15973,15991,16001,16007,16033,16057,16061,16063,16067,16069,16073,16087,16091,16097,16103,16111,16127,16139,16141,16183,16187,16189,16193,16217,16223,16229,16231,16249,16253,16267,16273,16301,16319,16333,16339,16349,16361,16363,16369,16381,16411,16417,16421,16427,16433,16447,16451,16453,16477,16481,16487,16493,16519,16529,16547,16553,16561,16567,16573,16603,16607,16619,16631,16633,16649,16651,16657,16661,16673,16691,16693,16699,16703,16729,16741,16747,16759,16763,16787,16811,16823,16829,16831,16843,16871,16879,16883,16889,16901,16903,16921,16927,16931,16937,16943,16963,16979,16981,16987,16993,17011,17021,17027,17029,17033,17041,17047,17053,17077,17093,17099,17107,17117,17123,17137,17159,17167,17183,17189,17191,17203,17207,17209,17231,17239,17257,17291,17293,17299,17317,17321,17327,17333,17341,17351,17359,17377,17383,17387,17389,17393,17401,17417,17419,17431,17443,17449,17467,17471,17477,17483,17489,17491,17497,17509,17519,17539,17551,17569,17573,17579,17581,17597,17599,17609,17623,17627,17657,17659,17669,17681,17683,17707,17713,17729,17737,17747,17749,17761,17783,17789,17791,17807,17827,17837,17839,17851,17863]; for p in SMALL_PRIMES.iter() { - if numb % &BigUint::from(*p) == zero { - return false + let prime = &BigUint::from(*p); + + if numb == prime { + return SmallPrimeResult::Prime; } - // Fixes part of Issue 1 but may slow down generation | https://github.com/AtropineTearz/num-primes/issues/1 - if numb / &BigUint::from(*p) == one { - return true + + // numb % prime == 0 and numb != prime => numb is not prime + if numb % prime == zero { + return SmallPrimeResult::NotPrime; } } - return true + + return SmallPrimeResult::MaybePrime; } @@ -323,7 +340,7 @@ fn fermat(candidate: &BigUint) -> bool { // p - 1 let exponent: BigUint = candidate.sub(BigUint::one()); - let result = random.modpow(&(exponent), candidate); + let result = random.modpow(&exponent, candidate); result == One::one() } @@ -356,7 +373,6 @@ fn miller_rabin(candidate: &BigUint, limit: usize) -> bool { if x == one || x == c_minus_1 { continue - // return true } else { let mut break_early = false; @@ -419,9 +435,11 @@ fn is_prime(candidate: &BigUint) -> bool { return false } - // First, simple trial divide - if div_small_primes(candidate) == false { - return false + + let p = div_small_primes(candidate); + + if p != SmallPrimeResult::MaybePrime { + return p == SmallPrimeResult::Prime; } // Second, Fermat's little theo test on the candidate diff --git a/tests/factorization.rs b/tests/factorization.rs index 59dcd30..ca2eba3 100644 --- a/tests/factorization.rs +++ b/tests/factorization.rs @@ -2,18 +2,17 @@ use num_primes::{Generator,Factorization}; #[test] fn factor_uint(){ - // Generate a Large Unsigned Integer of 32 bits - let x = Generator::new_uint(32); - - println!("Number: {}",x); + // Generate a Unsigned Integer of 16 bits + let p = Generator::new_prime(16); + let q = Generator::new_prime(16); + let f = &p * &q; // Factor The Largest Prime of x - let prime_factor = Factorization::prime_factor(x); + let prime_factor = Factorization::prime_factor(f).unwrap(); - // Print Out The Statements - match prime_factor { - Some(prime_factor) => println!("Prime Factor: {}",prime_factor), - None => println!("There are no prime factors") + if &p > &q { + assert_eq!(prime_factor, p); + } else { + assert_eq!(prime_factor, q); } - } \ No newline at end of file From 1bd4baa2eb2e625e718c3399b508f75aab3315f9 Mon Sep 17 00:00:00 2001 From: Evan Hines Date: Mon, 6 Dec 2021 23:30:01 -0500 Subject: [PATCH 3/4] Remove commented out no_std --- src/lib.rs | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/lib.rs b/src/lib.rs index f614b1a..8a105e0 100644 --- a/src/lib.rs +++ b/src/lib.rs @@ -1,4 +1,4 @@ -//#![no_std] +#![no_std] #![forbid(unsafe_code)] extern crate num; From a83a76cb584db400c9616bbdaa0a4959a914f66e Mon Sep 17 00:00:00 2001 From: Evan Hines Date: Tue, 7 Dec 2021 00:06:52 -0500 Subject: [PATCH 4/4] Fall back to brute force factorization in cases where rho fails to find a solution. --- src/lib.rs | 33 ++++++++++++++++++++++++++++++++- 1 file changed, 32 insertions(+), 1 deletion(-) diff --git a/src/lib.rs b/src/lib.rs index 8a105e0..dbca7f1 100644 --- a/src/lib.rs +++ b/src/lib.rs @@ -260,7 +260,8 @@ impl Factorization { if is_prime(&n) { return Some(n) } - + + let org = n.clone().try_into().unwrap(); let mut n: BigInt = n.try_into().unwrap(); let mut f = BigInt::zero(); @@ -271,6 +272,7 @@ impl Factorization { // This while loop will run until n is completely factorized, // thus giving the largest factor. If you only care about A factor, // this can be rewritten as below for much better performance. + while n > one { // while f == BigInt::zero() { @@ -288,6 +290,10 @@ impl Factorization { d = gcd((&x - &y).abs(), &n_int - 0); } + if &d == &org { + return Some(brute_force_prime_factor(d.try_into().unwrap())); + } + n /= &d; if &d > &f { @@ -299,6 +305,31 @@ impl Factorization { } } +fn brute_force_prime_factor(mut n: BigUint) -> BigUint { + + let one = BigUint::one(); + let two = &one + &one; + + // STEP 1 | n divided by 2 + while n.is_even() { + n = n / &two; + } + + // STEP 2 | 3..sqrt(n) | Divide i by n. On failure, add 2 to i + let n_sqrt = n.sqrt().to_usize().unwrap(); + + for mut i in 3..n_sqrt + 1 { + let divisor = BigUint::from(i); + while n.divides(&divisor) { + n /= &divisor; + } + i = i + 1usize; + } + + return if n == one { BigUint::from(n_sqrt) } else { n }; +} + + #[derive(PartialEq)] enum SmallPrimeResult { Prime,