zoukankan      html  css  js  c++  java
  • Miller-Rabin素数测试

    费马小定理:对于素数p和任意整数a,有ap ≡ a(mod p)(同余)。反过来,满足ap ≡ a(mod p),p也几乎一定是素数。

    伪素数:如果n是一个正整数,如果存在和n互素的正整数a满足 an-1 ≡ 1(mod n),我们说n是基于a的伪素数。如果一个数是伪素数,那么它几乎肯定是素数。

    Miller-Rabin测试:不断选取不超过n-1的基b(s次),计算是否每次都有bn-1 ≡ 1(mod n),若每次都成立则n是素数,否则为合数。 

    二次探测定理:如果n是一个素数,且0<x<p,则方程x^2%p=1的解为:x=1或    x=p-1.

      1 //Achen
      2 #include<algorithm>
      3 #include<iostream>
      4 #include<cstring>
      5 #include<cstdlib>
      6 #include<cstdio>
      7 #include<vector>
      8 #include<queue>
      9 #include<cmath>
     10 #include<ctime>
     11 const int N=1e5+7;
     12 typedef long long LL;
     13 using namespace std;
     14 int T;
     15 LL x,p[N]; 
     16 
     17 template<typename T> void read(T &x) {
     18     T f=1; x=0; char ch=getchar();
     19     while(ch!='-'&&(ch<'0'||ch>'9')) ch=getchar();
     20     if(ch=='-') f=-1,ch=getchar();
     21     for(;ch>='0'&&ch<='9';ch=getchar()) x=x*10+ch-'0'; x*=f;
     22 }
     23 
     24 LL ksc(LL a,LL b,LL mod) {
     25     LL base=a%mod,res=0;
     26     while(b) {
     27         if(b&1) res=(res+base)%mod;
     28         base=(base+base)%mod;
     29         b>>=1;
     30     } 
     31     return res;
     32 }
     33 
     34 LL ksm(LL a,LL b,LL mod) {
     35     LL base=a,res=1;
     36     while(b) {
     37         if(b&1) res=ksc(res,base,mod);
     38         base=ksc(base,base,mod);
     39         b>>=1;
     40     }
     41     return res;
     42 }
     43 
     44 int miller_rabin(LL n) {
     45     LL u=n-1;
     46     int k=0;
     47     if(n==2||n==3||n==5||n==7||n==11) return 1; 
     48     if(!(n%2)||!(n%3)||!(n%5)||!(n%7)||!(n%11)) return 0;
     49     while(!(u&1)) {
     50         u>>=1; k++;
     51     }
     52     for(int i=1;i<=20;i++) {
     53         LL x=rand()%(n-2)+2;
     54         LL tp=ksm(x,u,n);
     55         for(int j=1;j<=k;j++) {
     56             LL tpp=ksc(tp,tp,n);
     57             if(tpp==1&&tp!=1&&tp!=n-1) return 0;
     58             tp=tpp;
     59         }
     60         if(tp!=1) return 0;
     61     }
     62     return 1;
     63 } 
     64 
     65 LL gcd(LL a,LL b) {return (!b)?a:gcd(b,a%b);}
     66 
     67 LL pallord_rho(LL n,int c) {
     68     LL x=rand()%n,y=x;
     69     int k=2,i=1;
     70     for(;;) {
     71         i++;
     72         x=(ksc(x,x,n)+c)%n;
     73         LL tp=gcd((x-y+n)%n,n);
     74         if(tp>1&&tp<n) return tp;
     75         if(x==y) return n; 
     76         if(i==k) y=x,k+=k;
     77     }
     78 }
     79 
     80 void find(LL n) {
     81     if(miller_rabin(n)) {
     82         p[++p[0]]=n;
     83         return ;
     84     } 
     85     LL p=n;
     86     for(int c=13;;c++) {
     87         p=pallord_rho(n,c);
     88         if(p>1&&p<n) break;
     89     }
     90     find(p); find(n/p);
     91 }
     92 
     93 int main() {
     94     read(T);
     95     while(T--) {
     96         p[0]=0;
     97         read(x);
     98         if(x==1||miller_rabin(x)) puts("Prime");
     99         else {    
    100             find(x);
    101             LL ans=p[1];
    102             for(int i=2;i<=p[0];i++) ans=min(ans,p[i]);  
    103             printf("%lld
    ",ans);
    104         }
    105     }
    106     return 0;
    107 }
    View Code

     

  • 相关阅读:
    数据结构前言
    Linux---远程连接、命令行基础、文件及目录管理
    HTTP协议
    Docker---dockerfile
    Docker---指令
    Docker---介绍
    进程模块的使用
    numpy---(精简)
    OpenJudge/Poj 2105 IP Address
    OpenJudge 2786 Pell数列
  • 原文地址:https://www.cnblogs.com/Achenchen/p/7501304.html
Copyright © 2011-2022 走看看