lib-modulo 0.6.0

Fast modular arithmetic
Documentation
use lib_modulo::Modulus64;

fn main() {
    assert!(primality_test(2));
    assert!(!primality_test(2 * 3));

    // 9-th Mersenne number
    assert!(primality_test((1 << 61) - 1));
    // 2^62 - 1 = (2^32 + 1)(2^32 - 1)
    assert!(!primality_test(u64::MAX));
}

/// Miller-Rabin primality test
pub fn primality_test(x: u64) -> bool {
    let (s, d) = {
        let x = x - 1;
        let s = x.trailing_zeros();
        (s - 1, x >> s)
    };

    let modulus = Modulus64::new(x);
    let one = modulus.residue(1);

    'test: for witness in [2, 325, 9375, 28178, 450775, 9780504, 1795265022] {
        let mut witness = modulus.residue(witness);

        if witness.is_zero() {
            continue;
        }

        witness = witness.pow(d);
        if witness == one || witness == -one {
            continue;
        }

        let mut s = s;
        while s > 0 {
            s -= 1;

            witness = witness * witness;
            if witness == -one {
                continue 'test;
            }
        }

        return false;
    }

    true
}