題目連結:http://poj.org/problem?id=1811
題目解析:2<=n<2^54,如果n是素數直接輸出,否則求N的最小質因數。
求大整數最小質因數的算法沒看懂,不打算看了,直接貼代碼,以後當模版用。
資料比較大,隻能先用Miller_Rabin算法進行素數判斷。
在用Pollard_rho分解因子。
#include <iostream>
#include <stdio.h>
#include <string.h>
#include <stdlib.h>
#include <math.h>
#include <time.h>
#include <algorithm>
typedef long long ll;
#define Time 15 //随機算法判定次數,Time越大,判錯機率越小
using namespace std;
ll n,ans,factor[10001];//質因數分解結果(剛傳回時是無序的)
ll tol;//質因數的個數,數組下标從0開始
//****************************************************************
// Miller_Rabin 算法進行素數測試
//速度快,而且可以判斷 <2^63的數
//****************************************************************
long long mult_mod(ll a,ll b,ll c)//計算 (a*b)%c. a,b都是ll的數,直接相乘可能溢出的
{
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;
}
ll pow_mod(ll x,ll n,ll mod)//x^n%n
{
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||n==3||n==5||n==7)return true;
if(n==1||(n%2==0)||(n%3==0)||(n%5==0)||(n%7==0)) return false;//偶數
ll x=n-1;
ll t=0;
while((x&1)==0)
{
x>>=1;
t++;
}
for(int i=0; i<Time; i++)
{
ll a=rand()%(n-1)+1;//rand()需要stdlib.h頭檔案
if(check(a,n,x,t))
return false;//合數
}
return true;
}
//************************************************
//pollard_rho 算法進行質因數分解
//************************************************
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;
}
ll Pollard_rho(ll x,ll c)
{
ll i=1,k=2;
ll x0=rand()%x;
ll y=x0;
while(1)
{
i++;
x0=(mult_mod(x0,x0,x)+c)%x;
long long 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()
{
int T;
//srand(time(NULL));加上RE不懂
scanf("%d",&T);
while(T--)
{
scanf("%lld",&n);//(n>=2)
/*if(n==1)
{
printf("1\n");
continue;
}*/
if(Miller_Rabin(n))
{
printf("Prime\n");
continue;
}
tol=0;
findfac(n);//對n分解質因子
ll ans=factor[0];
for(int i=1; i<tol; i++)
if(factor[i]<ans)
ans=factor[i];
/*for(int i=0;i<tol;i++)
{
printf("%lld\n",factor[i]);
}*/
printf("%lld\n",ans);
}
return 0;
}
算法解析:
由費馬小定理可以知道,若p是素數且a是整數,則滿足a^p==a(mod p)。若存在正整數a不滿足a^p==a(mod p),那麼p是合數。
定義:令a是一個正整數,若p是合數且滿足a^p==a(mod p),則p稱為以a為基的僞素數。
Miller-Rabin素數測試算法原理: 假如p是素數,且(a,p)==1,(a為任意小于p的正整數),那麼a^p-1==1(mod p)。如果a^p-1==1(mod p),
則可認為n是素數,取多個底進行試驗,次數越多,n為素數機率越大。(我的個人了解多次試驗為p換基,使之成為僞素數的可能性大大減小)。
(Miller-Rabin測試:不斷選取不超過n-1的基b(s次),計算是否每次都有bn-1 ≡ 1(mod n),若每次都成立則n是素數,否則為合數。)
轉載:說Miller-Rabin測試以前先說兩個比較高效的求a*b% n 和 ab %n 的函數,這裡都是用到二進制思想,将b拆分成二進制,然後與a相加(相乘)
// a * b % n
//例如: b = 1011101那麼a * b mod n = (a * 1000000 mod n + a * 10000 mod n + a * 1000 mod n + a * 100 mod n + a * 1 mod n) mod n
ll mod_mul(ll a, ll b, ll n) {
ll res = 0;
while(b) {
if(b&1) res = (res + a) % n;
a = (a + a) % n;//a=(a<<1)%n
b >>= 1;
}
return res;
}
這代碼很棒,以後計算a*b時,如果裡面有一個數很大,則可以選擇上面的算法,(nlogn)的時間複雜度。
//a^b % n
//同理
ll mod_exp(ll a, ll b, ll n) {
ll res = 1;
while(b) {
if(b&1) res = mod_mul(res, a, n);
a = mod_mul(a, a, n);
b >>= 1;
}
return res;
}
快速幂,沒什麼好說的。
核心代碼:
開始程式時需加srand(time(NULL));
bool miller_rabin(ll n)
{
for(int i=1; i<=N; i++) //N為你打算測試的次數,N(10~20)
{
ll a=random(n-2)+1;//需頭檔案stdlib.h,random(X)産生0~X的随機數,+1産生1~n-1
if(mod_exp(a,n-1,mod)!=1)
{
"合數";
}
}
}
注意,MIller-Rabin測試是機率型的,不是确定型的,不過由于多次運作後出錯的機率非常小,是以實際應用還是可行的。(一次Miller-Rabin測試其成功的機率為3/4)
二次探測定理:(改進)
一個合數n,若對所有滿足(b,n)=1的正整數b都有b^n-1==1(mod n)成立,(上面的反例,但出現這種數的幾率不大),則稱之為卡邁克爾數。
二次探測 如果p是奇素數,則 x2 ≡ 1(mod p)的解為 x = 1 || x = p - 1(mod p);
可以利用二次探測定理在實作Miller-Rabin上添加一些細節,具體實作如下:
bool miller_rabin(ll n) {
if(n == 2 || n == 3 || n == 5 || n == 7 || n == 11) return true;
if(n == 1 || !(n%2) || !(n%3) || !(n%5) || !(n%7) || !(n%11)) return false;
ll x, pre, u;
int i, j, k = 0;
u = n - 1; //要求x^u % n
while(!(u&1)) { //如果u為偶數則u右移,用k記錄移位數
k++; u >>= 1;
}
srand((ll)time(0));
for(i = 0; i < S; ++i) { //進行S次測試
x = rand()%(n-2) + 2; //在[2, n)中取随機數
if((x%n) == 0) continue;
x = mod_exp(x, u, n); //先計算(x^u) % n,
pre = x;
for(j = 0; j < k; ++j) { //把移位減掉的量補上,并在這地方加上二次探測
x = mod_mul(x, x, n);
if(x == 1 && pre != 1 && pre != n-1) return false; //二次探測定理,這裡如果x = 1則pre 必須等于 1,或則 n-1否則可以判斷不是素數
pre = x;
}
if(x != 1) return false; //費馬小定理
}
return true;
}
前邊這個算法經過測試還是比較靠譜的,可以用作模闆。