cp_library_rs/number_theory/
miller_rabin_test.rs

1//! ミラー・ラビン素数判定法
2
3/// 余りをとる累乗
4fn powmod(a: usize, b: usize, m: usize) -> usize {
5    let (mut a, mut b, m) = (a as u128, b as u128, m as u128);
6    let mut res = 1;
7    while b > 0 {
8        if b & 1 == 1 {
9            res = (res * a) % m;
10        }
11        a = (a * a) % m;
12        b >>= 1;
13    }
14    res as usize
15}
16
17/// ## ミラーラビン素数判定法
18/// 参考: <https://zenn.dev/kaki_xxx/articles/40a92b43200215>
19pub fn is_prime_MR(N: usize) -> bool {
20    if N <= 2 {
21        return N == 2;
22    }
23    if N % 2 == 0 {
24        return false;
25    }
26
27    let (mut s, mut d) = (0, N - 1);
28    while d % 2 == 0 {
29        s += 1;
30        d >>= 1;
31    }
32
33    // n < 2^64 の場合、以下を調べれば十分
34    let A = [2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37];
35    for &a in &A {
36        if a % N == 0 {
37            break;
38        }
39        let mut t = 0;
40        let mut x = powmod(a, d, N);
41        if x != 1 {
42            while t < s {
43                if x == N - 1 {
44                    break;
45                }
46                x = ((x as u128).pow(2) % (N as u128)) as usize;
47                t += 1;
48            }
49            if t == s {
50                return false;
51            }
52        }
53    }
54    true
55}