prime_formula/common/
mod.rs

1//! Common library for generating prime numbers
2//
3// Copyright: Adam Cottrell (cottrela@gmail.com)
4
5/// Base primes excluded from root calculations (2, 3, 5, 7)
6///
7/// These form the foundation of the prime period wheel factorization
8pub const NON_ROOT_PRIMES: [u8; 4] = [2, 3, 5, 7];
9
10/// The prime period (210) - product of base primes
11///
12/// All prime roots >7 will be coprime with this value
13pub const PRIME_PERIOD: u8 = {
14    let mut product = 1;
15    let mut i = 0;
16    while i < NON_ROOT_PRIMES.len() {
17        product *= NON_ROOT_PRIMES[i];
18        i += 1;
19    }
20    product
21};
22
23/// Prime roots coprime with the base primes (2,3,5,7)
24///
25/// These 48 roots form the basis for all prime number generation in this sieve.
26/// Any prime number >7 can be expressed as: root + k×210
27pub const PRIME_ROOTS: [u8; 48] = {
28    let mut primes = [0u8; 48];
29    let mut num = 2;
30    let mut idx = 0;
31
32    // Iterate through all candidates in range
33    while num < PRIME_PERIOD + 2 {
34        // Check if number is coprime with base primes
35        if num % 2 != 0 && num % 3 != 0 && num % 5 != 0 && num % 7 != 0 {
36            primes[idx] = num;
37            idx += 1;
38        }
39        num += 1;
40    }
41
42    // Compile-time validation of array size
43    if idx != 48 {
44        panic!("Prime roots count mismatch");
45    }
46
47    primes
48};
49
50/// The number of prime roots (48)
51pub const NUM_ROOTS: usize = PRIME_ROOTS.len();
52
53/// Public function to align start and end on 210 boundaries
54pub fn align_to_cycle(start: u128, end: u128) -> (u128, u128) {
55    // Find w/c aligned range for a given start/end window
56    //
57    // For any prime number (P):
58    //    P = r + b * 210
59    //    b = (P - r) / 210
60    //
61    // Where:
62    //   r is the prime root (one of 48)
63    //   b is the multiplier (1..N)
64    let r_min = PRIME_ROOTS[0] as u128; // 13
65    let r_max = PRIME_ROOTS[NUM_ROOTS - 1] as u128; // 211
66
67    // Get ceil(b_start) and floor(b_end) to find the tightest range
68    let b_start = (start.saturating_sub(r_max)).div_ceil(210); // Use integer ceil
69    let b_end = 1 + (end - r_min) / 210; // Use integer floor
70
71    // If b_start == b_end, then there is nothing to do
72    (b_start, b_end)
73}
74
75/// Helper function to find the inverse prime root values.
76/// The inverse of a prime root (q) is the congruent
77/// of q * q_hat - r = 0, solved for q_hat.
78/// Where r is a different prime root.
79///
80/// n.b. there is one solution per prime root that satisfies
81/// this equation.
82pub fn find_inv_prime_root_vec(root_idx: usize) -> Vec<u8> {
83    let r = PRIME_ROOTS[root_idx] as u32;
84    PRIME_ROOTS
85        .iter()
86        .copied()
87        .map(|q| {
88            let period = PRIME_PERIOD as u32;
89            let q = q as u32;
90            PRIME_ROOTS
91                .iter()
92                .copied()
93                .find(|&h_candidate| (q * h_candidate as u32) % period == r % period)
94                .expect("No inverse prime root found")
95        })
96        .collect()
97}
98
99/// Function to generate cyclic table of composites (CTC)
100/// These offsets are used as look-ups for each
101/// prime root.
102///
103/// Constructing the CTC is the first step in creating
104/// the Periodic Table of Primes (PTP).
105pub fn get_cyclic_composite_vec(root_index: usize, inv_prime_root_table_row: &[u8]) -> Vec<u8> {
106    let r = PRIME_ROOTS[root_index] as u32;
107    let period = PRIME_PERIOD as u32;
108    PRIME_ROOTS
109        .iter()
110        .enumerate()
111        .map(|(j, &q)| {
112            let q_hat = inv_prime_root_table_row[j] as u32;
113            let q = q as u32;
114            let val = if r > (q * q_hat) {
115                1 + (period + q * q_hat - r) / period
116            } else {
117                1 + (q * q_hat - r) / period
118            };
119            val as u8
120        })
121        .collect()
122}
123
124/// Get primes 2,3,5,7 within range [start, end]
125pub fn get_small_primes(start: u128, end: u128) -> Vec<u128> {
126    [2, 3, 5, 7]
127        .iter()
128        .filter(|&&p| p >= start && p <= end)
129        .cloned()
130        .collect()
131}
132
133/// High-level mode for the selecting algorithm to use
134#[derive(Debug, Copy, Clone, PartialEq)]
135pub enum ConfigMode {
136    /// Novel sieve method - best for medium ranges
137    Sieve,
138    /// Windowed search - optimized for small number ranges
139    Window,
140    /// Use Miller-Rabin test with deterministic witnesses
141    MillerRabinDeterministic,
142    /// Use Miller-Rabin test with random witnesses
143    MillerRabinProbabilistic,
144}
145
146/// Helper function to select ConfigMode based on the
147/// start/end range from big-O approximations. The goal
148/// being to select the fastest algorithm for the range.
149pub fn get_auto_mode(_start: u128, end: u128) -> ConfigMode {
150    let big_o_sieve = end.isqrt(); // Dominated by sqrt(n)
151    let big_o_miller = end.ilog10().pow(6) as u128;
152
153    if big_o_sieve > big_o_miller {
154        // For very large numbers, use best-effort
155        // to reduce the number of options - this
156        // may return pseudoprimes.
157        if end > u64::MAX as u128 {
158            ConfigMode::MillerRabinProbabilistic
159        } else {
160            ConfigMode::MillerRabinDeterministic
161        }
162    } else {
163        // For everything else, sieve mode should give
164        // the best overall performance. It is exhaustive
165        // and will quickly eliminate a lot of composites.
166        ConfigMode::Sieve
167    }
168}
169
170#[cfg(test)]
171mod tests {
172    use super::*;
173
174    #[test]
175    fn test_constants() {
176        // Verify base primes
177        assert_eq!(NON_ROOT_PRIMES, [2, 3, 5, 7]);
178
179        // Verify prime period calculation
180        assert_eq!(PRIME_PERIOD, 2 * 3 * 5 * 7);
181
182        // Validate prime roots count and properties
183        assert_eq!(PRIME_ROOTS.len(), 48);
184        for &root in &PRIME_ROOTS {
185            assert!(root >= 11);
186            assert!(root <= 211);
187            assert_ne!(root % 2, 0);
188            assert_ne!(root % 3, 0);
189            assert_ne!(root % 5, 0);
190            assert_ne!(root % 7, 0);
191        }
192    }
193
194    #[test]
195    fn test_align_to_cycle() {
196        // Test exact boundaries
197        assert_eq!(align_to_cycle(0, 210), (0, 1));
198        assert_eq!(align_to_cycle(210, 420), (0, 2)); // 210 < 211 so should include first cycle as well as second
199
200        // Test values within cycle
201        assert_eq!(align_to_cycle(100, 200), (0, 1));
202        assert_eq!(align_to_cycle(300, 500), (1, 3));
203
204        // Test edge cases
205        assert_eq!(align_to_cycle(211, 211), (0, 1)); // Upper boundary
206        assert_eq!(align_to_cycle(13, 13), (0, 1)); // Lower boundary
207    }
208
209    #[test]
210    fn test_inverse_prime_roots() {
211        // Test first and last roots
212        let first_inverses = find_inv_prime_root_vec(0);
213        assert_eq!(first_inverses.len(), 48);
214
215        let last_inverses = find_inv_prime_root_vec(47);
216        assert_eq!(last_inverses.len(), 48);
217
218        // Verify inverse property for sample root
219        let test_root = PRIME_ROOTS[10]; // Arbitrary index
220        let inverses = find_inv_prime_root_vec(10);
221        for (i, &h) in inverses.iter().enumerate() {
222            let q = PRIME_ROOTS[i];
223            assert_eq!(
224                (q as u32 * h as u32) % PRIME_PERIOD as u32,
225                test_root as u32 % PRIME_PERIOD as u32
226            );
227        }
228    }
229
230    #[test]
231    fn test_cyclic_composites() {
232        // Test known values
233        let test_index = PRIME_ROOTS.iter().position(|&r| r == 11).unwrap();
234        let inverses = find_inv_prime_root_vec(test_index);
235        let composites = get_cyclic_composite_vec(test_index, &inverses);
236
237        // Should generate predictable offsets
238        assert_eq!(composites.len(), 48);
239        assert!(composites.iter().all(|&v| v > 0));
240    }
241
242    #[test]
243    fn test_small_primes() {
244        // Full range
245        assert_eq!(get_small_primes(2, 10), vec![2, 3, 5, 7]);
246
247        // Partial range
248        assert_eq!(get_small_primes(5, 7), vec![5, 7]);
249
250        // Empty range
251        assert!(get_small_primes(8, 10).is_empty());
252    }
253
254    #[test]
255    fn test_auto_mode_selection() {
256        // Small numbers should use sieve
257        assert_eq!(get_auto_mode(0, 1_000), ConfigMode::Sieve);
258
259        // Very large numbers should use probabilistic MR
260        assert_eq!(
261            get_auto_mode(0, u128::MAX),
262            ConfigMode::MillerRabinProbabilistic
263        );
264
265        // Borderline case (sqrt(n) vs log^6(n))
266        let borderline = 10u128.pow(15);
267        let mode = get_auto_mode(0, borderline);
268        assert!(matches!(
269            mode,
270            ConfigMode::Sieve | ConfigMode::MillerRabinDeterministic
271        ));
272    }
273
274    #[test]
275    fn test_prime_root_validity() {
276        // Verify all roots are in correct range
277        assert!(PRIME_ROOTS.iter().all(|&r| (11..=211).contains(&r)));
278
279        // Verify no duplicates
280        let mut sorted = PRIME_ROOTS.to_vec();
281        sorted.sort_unstable();
282        sorted.dedup();
283        assert_eq!(sorted.len(), 48);
284    }
285
286    #[test]
287    fn test_boundary_alignment() {
288        // Test alignment near cycle boundaries
289        assert_eq!(align_to_cycle(209, 211), (0, 1));
290        assert_eq!(align_to_cycle(211, 421), (0, 2));
291
292        // Test with extremely large numbers
293        let big_num = 210u128.pow(10);
294        assert_eq!(
295            align_to_cycle(big_num, big_num + 100),
296            (big_num / 210 - 1, big_num / 210 + 1)
297        );
298    }
299}
prime_formula/common/mod.rs

prime_formula/common/
mod.rs
