开篇，CSDN

以poj的大素数判断算法，作为我的第一篇。前几天就想弄的，可是写的东西都成一堆了，断不了行。不知道什么原因。今天写了下，又只可以 了，郁闷 。。。

题出处:http://acm.pku.edu.cn/JudgeOnline/problem?id=1811

 

 

 

资料：
http://book.51cto.com/art/200812/102576.htm
http://book.51cto.com/art/200812/102577.htm



问题：
        给定一个大整数(2<N<2^54)，判断它是不是素数，如果是，则输出"Prime"，否则输出它的最小素因子。
算法：
         Miller--Rabin素性测试法+ Pollard rho因子分解法。
算法介绍
Miller--Rabin：
         设n>2是一个奇数，n-2=2^(s)m,其中s是非负整数，m>0是奇数，则：
1.随机选取b=1,2,.....,n-1;
2.r=0；z=b^(m)mod n。如果z=1或者z=n-1,则n可能是素数，结束。
3.如果r=s-1,则n是合数，结束。
4.r=r+1;z=z^2 mod n;如果z=n-1，则n可能是素数。

注：每次测试得到是合数的概率<1/4，当测试次数很大时，其1/(4^k)是个很小的数。本题中
k=10;

Pollard rho:
原理：
1.假设有两个整数x1,x2，使得p可以整除x1-x2,但是n不能整除x1-x2;
2.可以证明gcd(x1-x2,n)=p。因为p可以整除x1-x2,可以写成x1-x2=q*p;由于n不能整除x1-x2，所以q不能整除n,这就表明gcd(x1-x2,n)既可以是1，也可以是n的一个因数。
算法：
1.选择x1,一个小的随机整数为种子。
2.运用函数算出x2,使得n不能整除x1-x2，所以函数x2=x1^2+a(在mod p)的条件下。
3.计算gcd(x1-x2,n)，如果不为1,则是n的因子，继续迭代。


注：x,a一开始随机获得。

注：为了减少反复的次数，算法做了一些改进。该算法用数对( , )开始，并且用 ，迭代计算 。在每一次迭代中，我们都应用上述函数式运算(从第二步)第一次计算数对中的第一个元素，第二次计算数对中的第二个元素。


Pollard p-1因子分解法：(一开始用了这个，死活不让过,到网上才查到还有一个Pollard rho,换了那种一下就过了。。。ft......)

算法：
1。输入奇整数n,输入整数B.
2。a=2;
3。对j=2到B,计算a=a^j mod n
4.计算d=gcd(a-1,n);
5,如果1<d<n,则d是n 的因子。



其他方法：
1。位运算：当n&1=1时，n为奇数。n=n>>1相当于n=n/2;
2。 a*b=c( mod m)时，而a*b超出所能表示的精度的处理。
3。a^b=c(mod m)；
4.最小公约数gcd(n,m);




Source Code
Problem: 1811  User: liuyuquan100
Memory: 308K  Time: 766MS
Language: G++  Result: Accepted
Source Code
#include<stdio.h>
#include<stdlib.h>
#include<math.h>
__int64 g(__int64 a,__int64 b,__int64 p)//a*b%p
{
__int64 R,D;
R=0;
D=a;
while(b>0)
{
  if(b&1) R=(R+D)%p;b为奇数b&1=1;偶数时为b&1=0;
  D=(D+D)%p;
  b>>=1;
}
return R;
}
__int64 f(__int64 b,__int64 m,__int64 n)//a^b % n
{
__int64 c;
c=1;
while(m)
{
  if(m&1) c=g(c,b,n);
  b=g(b,b,n);
  m>>=1;
}
return c;
}
//Miller-Rabin素数测试法。
bool check(__int64 b,__int64 n,__int64 m,__int64 s)//n-1=2^s * m
{
__int64 r,z;
z=f(b,m,n);
if(z==1||z==n-1) return true;
for(r=1;r<s;r++)
{
  z=g(z,z,n);
  if(z==n-1) return true;
}
return false;
}
//求最小公约数
__int64 gcd(__int64 n,__int64 m)
{
if(m==0)  return n;
return gcd(m,n%m);
}

//pollard的rho因子分解算法。
__int64 findminprime(__int64 n)
{
__int64 k,x,y,p,i,a;
x=rand()%n;
a=rand()%(n-1)+1;
y=x;
i=1;
k=2;
do
{
  i++;
  p=gcd(n+y-x,n);
  if(i==k) y=x,k=k*2;
  x=(g(x,x,n)+a)%n;
  if(p>1&&p<n) return p;
}while(x!=y);
return n;
}
//Miller-Rabin素性测试法。
bool checkprime(__int64 n)
{
__int64 k,s,i,b;
k=n-1;
s=0;
while(!(k&1))
{
  k=k>>1;
  s++;
}
for(i=0;i<10;i++)
{
  b=rand()%(n-1)+1;
  if(!check(b,n,k,s)) return false;
  
}
return true;
}
__int64 prim[]={2,3,5,7,11,13,17,19};
//找出最小的素因子。
__int64 findmin(__int64 n,__int64 t)
{
if(t==1&&checkprime(n)) return n;
else
{
  __int64 t,x,y,i;
  t=0;
  for(i=0;i<8;i++)
  {
   if(n%prim ==0) {t=1;x=prim;y=n/prim;break;}
  }
  if(t==0)
  {
   x=findminprime(n);
   y=n/x;
   x=findmin(x,1);
   y=findmin(y,1);
  }
  else
  {
   y=findmin(y,1);
  }
  if(x>y) x=y;
  return x;
}
}

int main()
{
__int64 m,n;
scanf("%I64d",&m);
while(m--)
{
  scanf("%I64d",&n);
  if(n==2) printf("Prime/n");
  else if(!(n&1)) printf("2/n");
  else
  {
   if(checkprime(n)) printf("Prime/n");
   else
   {
    printf("%I64d/n",findmin(n,0));
   }
  }
}
return 0;
}