• 微信公众号:美女很有趣。 工作之余,放松一下,关注即送10G+美女照片!

Miiler Rabin 算法

开发技术 开发技术 6小时前 3次浏览

Miller Rabin 算法

参考:【朝夕的ACM笔记】数论-Miller Rabin素数判定 – ~(o°ω°o)

作用在于判断素数,它有前置两个定理

费马小定理

Theorem(p) 是一个素数,(ain Z^+) 且不是 (p) 的倍数,那么有 (a^{p-1}equiv1(bmod p))

Proof 我们知道,若有 (i)(j) 在模 (p) 意义下不同余,那么 (itimes a)(jtimes a) 在模 (p) 意义下也不同余

那么有: (1times 2times 3timesdotstimes(p-1)equiv(1times a)times(2times a)timesdotstimes(p-1)*a(bmod p))

(W=Pi_{i=1}^{p-1}) ,那么 (Wequiv Wtimes a^{p-1}(bmod p))

又知:(gcd(W,p)=1), 所以,(a^{p-1}equiv 1pmod p)

证毕。

我们可以猜想,费马小定理的逆定理是否成立

也就是说,如果 (ain Z^+) 且不是 (p) 的倍数,而且 (a^{p-1}equiv1pmod p) ,是否能说 (p) 是一个素数?

答案是不能,反例:对 (a = 2),满足 (2^{n−1} bmod n = 1) 的非素数 (n) 是存在的,比如 (n = 341 = 11 × 31)

伪素数:对于整数 (a),称满足 (a^{n−1} bmod n = 1) 的合数为以 (a) 为底的伪素数。

经测试,前 10 亿的自然数中,同时以 2 和 3 为底的伪素数有 1272 个。

我们用费马小定理验证素数的话,出错的概率大概只有 0.000025。

于是我们初步得到了一个判断素数的算法

Solution 随机选取若干个小于待测整数的正整数作为底 a,然后用费马小定理来测试

但是,有一些坑爹的合数能通过所有的测试,也就是 (text{Carmichael})

二次探测定理

Theorem(p) 为素数, $a^2≡1pmod p $,那么 (a≡±1pmod p)

Proof 其实感觉是显然的。由题,(p|(a^2-1)) ,即 (p|(a-1)(a+1))

那么,因为 (p) 是素数,根据唯一分解定理,(p|(a-1)) 或者 (p|(a+1)) 。得证。

所以,如果 $a^2≡1pmod p $ ,但是不满足 (a≡±1pmod p) ,那么 (p) 不是素数。

MR

二次探测定理和素数判定有什么用呢?

首先,根据 (Miller Rabin) 算法的过程

假设需要判断的数是 (p) (假设 (p) 是奇数,是偶数时显然不是素数)

我们把 (p−1) 分解为 (2^k+t) 的形式(当 (p) 是素数,(a^{2^k+t}equiv 1pmod p)

我们随机一个 (a) ,计算出 (a^t) ,让它不断的自乘,

如果在某时,(a^{2k}equiv1mod p) ,但是不满足 (a^k≡±1pmod p) ,那么 (p) 不是素数

如果最后算到 (a^{2^k+t}bmod pne1)(p) 不是素数。

  • 正确性:首先是底数的选取:Miller Rabin检测依旧有错误的概率,虽然微乎其微,但我们显然不想接受。但通过选取适合的底数,可以避免这一情况的发生。

    (int) 范围内,选取 (2,7,61) 三个数作为底数可以保证 (100%) 正确。

    (long long) 范围内,选取 (2,325,9375,28178,450775,9780504,1795265022) 七个数作为底数可保证 (100%) 正确。

    故正常情况下时间复杂度最多为 (O(n log^3 n)) ,常数为 (7)

代码实现:

#include<cstdio>
#define LL long long 
inline int read() {
    char c = getchar(); int x = 0, f = 1;
    while(c < '0' || c > '9') {if(c == '-') f = -1; c = getchar();}
    while(c >= '0' && c <= '9') x = x * 10 + c - '0', c = getchar();
    return x * f;
}
int N, M, Test[10] = {2, 3, 5, 7, 11, 13, 17};
int pow(int a, int p, int mod) {
    int base = 1;
    for(; p; p >>= 1, a = (1ll * a * a) % mod) 
        if(p & 1) base = (1ll * base * a) % mod;
    return base % mod;
}
bool Query(int P) {
    if(P == 1) return 0;
    int t = P - 1, k = 0;
    while(!(t & 1)) k++, t >>= 1;
    for(int i = 0; i < 4; i++) {
        if(P == Test[i]) return 1;
        LL a = pow(Test[i], t, P), nxt = a;
        for(int j = 1; j <= k; j++) {
            nxt = (a * a) % P;
            if(nxt == 1 && a != 1 && a != P - 1) return 0;
            a = nxt;
        }
        if(a != 1) return 0;
    }
    return 1;
}
main() { 
    N = read(); M = read();    
    while(M--) puts(Query(read()) ? "Yes" : "No");
}

程序员灯塔
转载请注明原文链接:Miiler Rabin 算法
喜欢 (0)