代码如下:
#include<stdio.h>
#include<string.h>
#include<stdlib.h>
#include<time.h>
#include<iostream>
#include<algorithm>
#define LL long long
using namespace std;
//****************************************************************
// Miller_Rabin 算法进行素数测试
//速度快,而且可以判断 <2^63的数
//****************************************************************
const int S=20;//随机算法判定次数,S越大,判错概率越小
//计算 (a*b)%c. a,b都是long long的数,直接相乘可能溢出的
// a,b,c <2^63
LL mult_mod(LL a, LL b, LL c) {
a %= c;
b %= c;
LL ret=0;
while(b) {
if(b&1) {
ret += a;
ret %= c;
}
a <<= 1;
if(a >= c)
a %= c;
b >>= 1;
}
return ret;
}
//计算 x^n %c
LL pow_mod(LL x, LL n, LL mod) {//x^n%c
if(n==1)return x % mod;
x %= mod;
LL tmp=x;
LL ret=1;
while(n){
if(n&1)
ret=mult_mod(ret, tmp, mod);
tmp = mult_mod(tmp, tmp, mod);
n >>= 1;
}
return ret;
}
//以a为基,n-1=x*2^t a^(n-1)=1(mod n) 验证n是不是合数
//一定是合数返回true,不一定返回false
bool check(LL a, LL n, LL x, LL t){
LL ret = pow_mod(a,x,n);
LL last = ret;
for(int i=1; i<=t; i++) {
ret = mult_mod(ret, ret, n);
if(ret==1 && last!=1 && last!=n-1)
return true;//合数
last = ret;
}
if(ret != 1) return true;
return false;
}
// Miller_Rabin()算法素数判定
//是素数返回true.(可能是伪素数,但概率极小)
//合数返回false;
bool Miller_Rabin(LL n){
if(n < 2)return false;
if(n == 2)return true;
if((n&1) == 0) return false;//偶数
LL x = n-1;
LL t = 0;
while((x&1) == 0) {
x >>= 1;
t++;
}
for(int i=0; i<S; i++){
LL a = rand()%(n-1)+1;//rand()需要stdlib.h头文件
if(check(a, n, x, t))
return false;//合数
}
return true;
}
//************************************************
//pollard_rho 算法进行质因数分解
//************************************************
LL factor[100];//质因数分解结果(刚返回时是无序的)
int tol;//质因数的个数。数组小标从0开始
LL gcd(LL a,LL b) {
if(a==0)return 1;//???????
if(a<0) return gcd(-a,b);
while(b)
{
long long t=a%b;
a=b;
b=t;
}
return a;
}
long long Pollard_rho(LL x, LL c) {
LL i = 1,k = 2;
LL x0 = rand() % x;
LL y = x0;
while(true){
i++;
x0 = (mult_mod(x0, x0, x) + c) % x;
LL d=gcd(y-x0, x);
if(d!=1 && d!=x) return d;
if(y == x0) return x;
if(i == k) {
y = x0;
k += k;
}
}
}
//对n进行素因子分解
void findfac(LL n)
{
if(Miller_Rabin(n)) {//素数
factor[tol++] = n;
return;
}
LL p = n;
while(p >= n)
p=Pollard_rho(p, rand()%(n-1)+1);
findfac(p);
findfac(n/p);
}
int main(){
//srand(time(NULL));//需要time.h头文件//POJ上G++不能加这句话
LL n;
while(scanf("%I64d",&n) != EOF){
tol = 0;
findfac(n);
for(int i=0; i<tol; i++)
printf("%lld ", factor[i]);
if(Miller_Rabin(n))
printf("Yes\n");
else printf("No\n");
}
return 0;
}
再附上一份简短的,如果对于比表小的数,用这个模板 比较好,上面的会超时不知道什么原因
bool witness(__int64 a,__int64 n)
{
__int64 t,d,x;
d=1;
int i=ceil(log(n-1.0)/log(2.0)) - 1;
for(;i>=0;i--)
{
x=d; d=(d*d)%n;
if(d==1 && x!=1 && x!=n-1) return true;
if( ((n-1) & (1<<i)) > 0)
d=(d*a)%n;
}
return d==1? false : true;
}
bool miller_rabin(__int64 n)
{
if(n==2) return true;
if(n==1 || ((n&1)==0)) return false;
for(int i=0;i<50;i++){
__int64 a=rand()*(n-2)/RAND_MAX +1;
if(witness(a, n)) return false;
}
return true;
}