zoukankan      html  css  js  c++  java
  • Miller Rabin素数检测与Pollard Rho算法

    一些前置知识可以看一下我的联赛前数学知识

    如何判断一个数是否为质数

    方法一:试除法

    扫描(2sim sqrt{n})之间的所有整数,依次检查它们能否整除(n),若都不能整除,则(n)是质数,否则(n)是合数。

    代码

    bool is_prime(int n){
       if(n<2) return 0;
       int m=sqrt(n);
       for(int i=2;i<=m;i++){
       	if(n%i==0) return 0;
       }
       return 1;
    }
    
    

    方法二、线性筛

    (O(n)) 的复杂度筛出所有的质数,然后 (O(1)) 判断

    代码

    const int maxn=1e6+5;
    int pri[maxn]; 
    bool not_prime[maxn];
    void xxs(int n){
    	not_prime[0]=not_prime[1]=1;
    	for(int i=2;i<=n;++i){
    		if(!not_prime[i]){
    			pri[++pri[0]]=i;
    		}
    		for (int j=1;j<=pri[0] && i*pri[j]<=n;j++){
    			not_prime[i*pri[j]]=1;
    			if(i%pri[j]== 0) break;
    		}
    	}
    }
    

    但是,我们来看一下这道题 LOJ#143. 质数判定

    (5s) 的时间内判断 (10^5)(10^{18}) 级别的数是否是质数

    以上两种方法显然都不可行

    所以我们要用到一种更高效的算法 (Miller Rabin)

    Miller Rabin素数检测

    根据费马小定理

    (p) 为素数,对于任意整数 (a),都有

    (a^{p-1} equiv 1(mod p))

    它的逆定理在大多数情况下是成立的

    所以,可以每一次随机挑选一些数 (a)

    判断是否存在 (a^{p-1} equiv 1(mod p))

    只要有一个 (a) 不满足条件,就将其标记为合数

    否则标记为质数

    代码

    #include<cstdio>
    #include<ctime>
    #include<cstdlib>
    #define rg register
    typedef long long ll;
    ll gsc(rg ll ds,rg ll zs,rg ll mod){
    	return ((unsigned long long)(ds*zs)-(unsigned long long)((long double)ds/mod*zs)*mod+mod)%mod;
    }
    ll ksm(rg ll ds,rg ll zs,rg ll mod){
    	rg ll nans=1;
    	while(zs){
    		if(zs&1LL) nans=gsc(nans,ds,mod);
    		ds=gsc(ds,ds,mod);
    		zs>>=1LL;
    	}
    	return nans;
    }
    int main(){
    	rg ll aa,bb;
    	rg bool jud=1;
    	while(scanf("%lld",&aa)!=EOF){
    		jud=1;
    		if(aa==1){
    			printf("N
    ");
    			continue;
    		}
    		for (rg int i=1;i<=100;i++) {
    			bb=1LL*rand()*rand()%(aa-1)+1;
    			if(ksm(bb,aa-1,aa)!=1){
    				jud=0;
    				break;
    			}
    		}
    		if(jud) printf("Y
    ");
    		else printf("N
    ");
    	}
    	return 0;
    }
    

    但是也有例外,即存在一种极端反例卡迈克尔数(一种合数)

    对于任何卡迈克尔数,费马定理都成立

    虽然这种数很少,但是还是有可能会被卡(只有(33)分)

    所以我们需要借助另一个定理来优化它

    二次探测定理:若 (p) 为质数且 (x∈(0,p)),则方程 (x^2≡1(mod p)) 的解为 (x_1=1,x_2=p-1)

    对于任意一个奇素数 (p)(p-1) 一定可以写成 (r2^t) 的形式

    因此我们可以把费马小定理写成 (a^{r2^t}equiv 1(mod p)) 的形式

    一开始,我们先求出 (a^r mod p) 的值

    然后每一次给这个值平方,一共平方 (t)

    算一下每次得出来的结果是否满足二次探测定理

    如果不满足,说明这个数不是质数

    最后再看一下最终的值是否满足费马小定理即可

    时间复杂度 (klogn)

    其中 (k) 为每次检测的次数

    事实证明,这个算法出错的概率非常小

    (frac{1}{4^k})

    代码

    #include<cstdio>
    #include<algorithm>
    #include<cstring>
    #include<cstdlib>
    #include<iostream>
    #define rg register
    inline int read(){
    	rg int x=0,fh=1;
    	rg char ch=getchar();
    	while(ch<'0' || ch>'9'){
    		if(ch=='-') fh=-1;
    		ch=getchar();
    	}
    	while(ch>='0' && ch<='9'){
    		x=(x<<1)+(x<<3)+(ch^48);
    		ch=getchar();
    	}
    	return x*fh;
    }
    typedef long long ll;
    const int times=7;
    ll gsc(rg ll ds,rg ll zs,rg ll mod){
        return ((unsigned long long)(ds*zs)-(unsigned long long)((long double)ds/mod*zs)*mod+mod)%mod;
    }
    ll ksm(rg ll ds,rg ll zs,rg ll mod){
    	rg ll nans=1;
    	while(zs){
    		if(zs&1LL) nans=gsc(nans,ds,mod);
    		ds=gsc(ds,ds,mod);
    		zs>>=1LL;
    	}
    	return nans;
    }
    bool check(rg ll a,rg ll r,rg ll t,rg ll mod){
    	ll nans=ksm(a,r,mod),tmp=nans;
    	for(rg ll i=1;i<=t;i++){
    		nans=gsc(tmp,tmp,mod);
    		if(nans==1 && tmp!=1 && tmp!=mod-1) return 0;
    		tmp=nans;
    	}
    	if(tmp==1) return 1;
    	else return 0;
    }
    bool Miller_Robin(rg ll n){
    	if(n==2) return 1;
    	if(n<2 || (n&1LL)==0) return 0;
    	rg ll t=0,r=n-1;
    	while((r&1LL)==0){
    		r>>=1;
    		t++;
    	}
    	for(rg int i=1;i<=times;i++){
    		rg ll a=rand()%(n-1)+1;
    		if(!check(a,r,t,n)) return 0;
    	}
    	return 1;
    }
    ll x;
    int main(){
    	srand(time(0));
    	while(scanf("%lld",&x)!=EOF){
    		if(Miller_Robin(x)) printf("Y
    ");
    		else printf("N
    ");
    	}
    	return 0;
    }
    

    Pollard Rho算法

    模板题

    这道题不仅要判断是否是质数,还要求输出最大质因子

    判质数的话用 (Miller Rabin) 判一下就好了

    关键是怎么找出最大质因子

    其实还是随机的思想

    我们每一次随机一个数 (n),判断 (n) 是不是 (p) 的因数

    如果是,就分成两个子问题 (n)(frac{p}{n}) 递归求解

    如果当前的数为质数就更新答案

    但是这样随机概率非常小

    利用生日悖论,采用组合随机采样的方法,满足答案的组合比单个个体要多一些.这样可以提高正确率

    具体方法是随机两个整数 (n)(m)

    判断 (|n-m|) 是不是 (p) 的因数

    看起来没有什么不同,实际上效率提高了不少

    剩下的就在于怎么构造随机的 (n,m)

    构造的方法会影响到我们程序的效率

    实践证明,直接 (rand) 效率极低

    而构造一种形如 (f(x)=x^2+c) 的数列效率很高

    我们首先 (rand) 一个 (x_1)

    然后令 (f(x_2)=x_1^2+c)

    (f(x_3)=x_2^2+c)

    依次类推,每次把序列相邻两项作差判断是否是 (p) 的因数

    但是因为我们是在模意义下进行运算,所以这个序列一定有循环节

    所以当遇到循环节的时候我们就退出循环,重新构造一个数列

    因为每次求 (gcd) 的开销太大

    所以我们可以先把相邻两项的差连乘起来,这是不影响结果的

    当累乘的次数等于我们设定的一个值时再进行求 (gcd) 的运算

    一般把这个值设为 (2^k)

    代码

    #include<cstdio>
    #include<algorithm>
    #include<cstring>
    #include<cstdlib>
    #include<iostream>
    #include<cmath>
    #define rg register
    inline int read(){
    	rg int x=0,fh=1;
    	rg char ch=getchar();
    	while(ch<'0' || ch>'9'){
    		if(ch=='-') fh=-1;
    		ch=getchar();
    	}
    	while(ch>='0' && ch<='9'){
    		x=(x<<1)+(x<<3)+(ch^48);
    		ch=getchar();
    	}
    	return x*fh;
    }
    typedef long long ll;
    const int times=10;
    int nt;
    ll x,ans;
    ll gsc(rg ll ds,rg ll zs,rg ll mod){
        return ((unsigned long long)(ds*zs)-(unsigned long long)((long double)ds/mod*zs)*mod+mod)%mod;
    }
    ll ksm(rg ll ds,rg ll zs,rg ll mod){
    	rg ll nans=1;
    	while(zs){
    		if(zs&1LL) nans=gsc(nans,ds,mod);
    		ds=gsc(ds,ds,mod);
    		zs>>=1LL;
    	}
    	return nans;
    }
    bool check(rg ll a,rg ll r,rg ll t,rg ll mod){
    	ll nans=ksm(a,r,mod),tmp=nans;
    	for(rg ll i=1;i<=t;i++){
    		nans=gsc(tmp,tmp,mod);
    		if(nans==1 && tmp!=1 && tmp!=mod-1) return 0;
    		tmp=nans;
    	}
    	if(tmp==1) return 1;
    	else return 0;
    }
    bool Miller_Robin(rg ll n){
    	if(n==2) return 1;
    	if(n<2 || (n&1LL)==0) return 0;
    	rg ll t=0,r=n-1;
    	while((r&1LL)==0){
    		r>>=1;
    		t++;
    	}
    	for(rg int i=1;i<=times;i++){
    		rg ll a=rand()%(n-1)+1;
    		if(!check(a,r,t,n)) return 0;
    	}
    	return 1;
    }
    ll gcd(rg ll aa,rg ll bb){
    	if(bb==0) return aa;
    	return gcd(bb,aa%bb);
    }
    ll rp(rg ll x,rg ll y){
    	rg int t=0,k=1;
    	rg ll v0=rand()%(x-1)+1,v=v0,d,s=1;
    	while(1){
    		v=(gsc(v,v,x)+y)%x;
    		s=gsc(s,std::abs(v-v0),x);
    		if(v==v0 || !s) return x;
    		if(++t==k){
    			d=gcd(s,x);
    			if(d!=1) return d;
    			s=1,v0=v,k<<=1;
    		}
    	}
    
    }
    void dfs(rg ll n){
    	if(n<=ans || n==1) return;
    	if(Miller_Robin(n)){
    		ans=n;
    		return;
    	}
    	rg ll now=n;
    	while(now==n) now=rp(n,rand()%n+1);
    	while(n%now==0) n/=now;
    	dfs(n),dfs(now);
    }
    ll solve(rg ll n){
    	ans=1;
    	dfs(n);
    	return ans;
    }
    int main(){
    	srand(time(0));
    	scanf("%d",&nt);
    	while(nt--){
    		scanf("%lld",&x);
    		if(Miller_Robin(x)) printf("Prime
    ");
    		else printf("%lld
    ",solve(x));
    	}
    	return 0;
    }
    
  • 相关阅读:
    nginx-1.8.1的安装
    ElasticSearch 在3节点集群的启动
    The type java.lang.CharSequence cannot be resolved. It is indirectly referenced from required .class files
    sqoop导入导出对mysql再带数据库test能跑通用户自己建立的数据库则不行
    LeetCode 501. Find Mode in Binary Search Tree (找到二叉搜索树的众数)
    LeetCode 437. Path Sum III (路径之和之三)
    LeetCode 404. Sum of Left Leaves (左子叶之和)
    LeetCode 257. Binary Tree Paths (二叉树路径)
    LeetCode Questions List (LeetCode 问题列表)- Java Solutions
    LeetCode 561. Array Partition I (数组分隔之一)
  • 原文地址:https://www.cnblogs.com/liuchanglc/p/14245195.html
Copyright © 2011-2022 走看看