浅谈Miller Robin

· · 个人记录

1. 二次探测定理

可解得 $P|(x-1)$ 或 $P|(x+1)

于是 x = 1 或 P-1

2. 费马小定理

又可表示为 $a^{2^s\cdot t}\equiv1\pmod{P}

3. 算法思路

随机数字做 a，从 a^t 开始自乘，根据二次探测定理若 a^{2^i*t}\equiv1\pmod{P} 则a^{2^{i-1}*t}=1或(P-1)，由此检验 P；注意最后需要验证 a^{2^s*t}\equiv1\pmod{P}

建议 ll 20，int 10。（来源）

CODE

#include <cstdio>
#include <cstdlib>
#define ll long long

ll Pow(ll x, ll k, const ll M){
    ll ret =1;
    for(; k; x =x*x%M, k >>=1) if(k&1) ret =ret*x%M;
    return ret;
}

bool mr(ll p){
    if(p < 2) return 0;
    if(p == 2) return 1;
    if(p == 3) return 1;
    ll d =p-1, r =0;
    while(!(d&1)) ++r, d >>=1;
    for(int k =0; k < 20; ++k){
        ll a =1ll*rand()*rand()*rand()*rand()%(p-2-1)+2;/*[2, p-2]*/
        ll x =Pow(a, d, p);
        if(x == 1 || x == p-1) continue;
        for(int i =0; i < r-1; ++i){
            x =x*x%p;
            if(x == p-1) break;
        }
        if(x != p-1) return 0;
    }
    return 1;
}

以下是之前的版本，仅供参考，不保证正确性

参考：AAA, BBB, CCC, DDD

1. 二次探测定理

可解得$P|(x-1)$或$P|(x+1)

于是(x\pm1) = 0或P

因此只要枚举比P小的数即可确定P是素数

2. 费马小定理

又可表示为 $a^{2^s\cdot t}\equiv1\pmod{P}

可以通过证明得到这种方法(?且选取质数做a?)得到的x的优越性

(关于此条没有找到完整的证明，但通过检验(10^7，时间10s左右；可以自行验证更大数据；但据经验分析取100个质数(建议筛法)非高精度情况下不可能出错，int范围20个即可)接下来的算法思路是可以应付OI竞赛 )

3. 算法思路

用质数做a，从a^t开始自乘，根据二次探测定理若a^{2^i*t}\equiv1\pmod{P}则a^{2^{i-1}*t}=1或(P-1)，由此检验P；注意最后a^{2^s*t}\equiv1\pmod{P}

寻找500个质数a，检验的错误率可忽略不计

4. 一些性质

CODE

#include <cstdio>
#define ll long long

inline ll read(){
    ll x = 0; bool f =0; char c =getchar();
    while(c < '0' || c > '9') (c == '-') ? f =1, c =getchar() : c =getchar();
    while(c >= '0' && c <= '9') x = (x<<3) + (x<<1) + (48^c), c =getchar();
    return (f) ? -x : x;
}

ll qPow(ll x, ll k, const ll &M){
    ll ans =1;
    do if(k&1) ans =ans*x%M; while((x =x*x%M) && (k >>=1));
    return ans;
}

const int N =20;
ll TRY[] ={2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47, 53, 59, 61, 67, 71};

bool Miller_Rabin(ll P){
    for(int i =0; i < N; ++i)
        if(P == TRY[i]) return 1;
    for(int i =0; i < N; ++i){
        /*事实上对于一般取500样本，这条可以检测的不会出错，质数太大的又没用*/
        //if(TRY[i]*TRY[i]%P == 1) return 0;//TRY一定奇数，不可能P-1
        ll t =P-1, k =0; while(!t&1) ++k, t >>=1;
        ll res =P-1/*!避免cur第一次就为1*/, cur =qPow(TRY[i], t, P)%P;
        while(k--){
            if(cur == 1){
                if(res == P-1) break;
                else return 0;/*by 20.7 : 这个 else 基本就不可能用到*/
            }
            res =cur, cur =cur*cur%P;
        }
        if(cur != 1) return 0;
    }
    return 1;
}

int main(){
    while(1) printf("%s", (Miller_Rabin(read())) ? "Yes\n" : "No\n");
}
同步于剪贴板