prime-formula 0.3.1

//! Library for generating likely primes beyond 10^38 upwards
//!
//! The algorithm uses The Periodic Table of Primes to select candidates
//! for prime numbers after the number given, and then in parallel runs
//! a BigInt optimized version of the Miller-Rabin test to eliminate
//! composite numbers. Anything left has a good chance of being prime.
//!
//! For primes below 10^38 then recommended to use the other prime
//! formula as they will run faster in general due to the use of fixed
//! precision integers.
//
// Copyright: Adam Cottrell (cottrela@gmail.com)

use num_bigint::BigUint;

use crate::common::PRIME_ROOTS;
pub mod big_primality;

/// Generate huge prime numbers using a primality test
///
/// Using primality tests is quicker for very large prime
/// numbers, but considerably longer than a sieve for smaller
/// prime numbers.
///
/// Use the formula of primes to find large primes
pub struct BigPrimeFormula {
    k: u32, // Number of tests to use (limited to 48 in deterministic mode)
}

impl BigPrimeFormula {
    pub fn new(k: u32) -> Self {
        BigPrimeFormula { k }
    }

    /// Get next primes after start, for primes beyond 10^38
    ///
    /// Note that primes below 11 will not be found using this
    /// formula, and so must be manually added by the caller if
    /// needed.
    ///
    /// This function returns the primes in order of their roots,
    /// therefore may need sorting before use. Set get_primes.
    pub fn get_next_prime(&self, start: BigUint) -> BigUint {
        let mut i = &start / 210u8; // Align to floor of prime cycle

        // Iterate over `i` (cycles of 210) and roots in order
        loop {
            for root in PRIME_ROOTS.iter().copied() {
                let candidate = &i * 210u8 + BigUint::from(root);
                if candidate > start && big_primality::is_likely_big_prime(&candidate, self.k) {
                    return candidate;
                }
            }
            i += 1u8; // Increment the cycle counter
        }
    }

    /// Get next primes after start, for primes beyond 10^38
    ///
    /// Note that primes below 11 will not be found using this
    /// formula, and so must be manually added by the caller if
    /// needed.
    ///
    /// This function returns the primes in order of their roots,
    /// therefore may need sorting before use. Set get_primes.
    pub fn get_prev_prime(&self, start: BigUint) -> BigUint {
        let mut i = (&start + 209u8) / 210u8 + 1u8; // Align to ceiling of prime cycle
        // Iterate over `i` (cycles of 210) and roots in order
        loop {
            for root in PRIME_ROOTS.iter().rev().copied() {
                let candidate = &i * 210u8 + BigUint::from(root);
                if candidate < start && big_primality::is_likely_big_prime(&candidate, self.k) {
                    return candidate;
                }
            }
            i -= 1u8; // Increment the cycle counter
        }
    }
}

#[cfg(test)]
mod tests {
    use super::*;

    #[test]
    fn test_large_prime() {
        let mut next_prime;
        let inst = BigPrimeFormula::new(2);

        // Get the next prime from a sparse 128-bit number
        next_prime = inst.get_next_prime(BigUint::from(0x8000_8080_8888_8008_u128));
        assert_eq!(next_prime, 0x8000808088888027u128.into());

        // Get the next prime from a sparse 128-bit number
        next_prime = inst.get_next_prime(BigUint::from(0x8000_0000_0000_0000_u128 - 1));
        assert_eq!(next_prime, 0x800000000000001Du128.into());

        // Get the previous prime from 128-bit MAX
        next_prime = inst.get_prev_prime(BigUint::from(u128::MAX));
        assert_eq!(
            next_prime,
            340282366920938463463374607431768211297u128.into()
        );

        // Get the next prime after 128-bit MAX
        next_prime = inst.get_next_prime(340282366920938463463374607431768211297u128.into());
        assert_eq!(next_prime.clone(), BigUint::from(u128::MAX) + 52u8);

        // Get the previous-previous prime from 128-bit MAX
        next_prime =
            inst.get_prev_prime(BigUint::from(340282366920938463463374607431768211297u128));
        assert_eq!(
            next_prime,
            340282366920938463463374607431768211283u128.into()
        );
    }

    #[test]
    fn test_huge_prime() {
        // Very low chance of failure for low order primes
        let next_prime;
        let inst = BigPrimeFormula::new(4);

        // Check for a huge prime (bigger than a googol)
        next_prime = inst.get_next_prime(BigUint::from(2u8).pow(300));
        let huge_prime_str = "2037035976334486086268445688409378161051468393665936250636140449354381299763336706183397533";
        assert_eq!(format!("{next_prime}"), huge_prime_str);
    }

    #[test]
    fn test_rsa_prime() {
        // Run 1000 iterations of a very large prime
        let next_prime;
        let inst = BigPrimeFormula::new(1000); // Increase confidence

        // Find the next RSA prime (U1024)
        let rsa_init = BigUint::from(1u8) << 1023u32
            | BigUint::from(1u8) << 250u32
            | BigUint::from(1u8) << 64u32;
        next_prime = inst.get_next_prime(rsa_init);
        let huge_prime_str = "89884656743115795386465259539451236680898848947115328636715040578866337902750481566354238661203768010560056939935696678829394884407208311246423715319737062188883946712432742638151109800623047059726541476042502884419075341171231440738765806664746684135168551983053897680180966479640491965552649659009464271097";
        assert_eq!(format!("{next_prime}"), huge_prime_str);
    }
}