## 3667: Rabin-Miller算法

Time Limit: 60 Sec  Memory Limit: 512 MB
Submit: 1214  Solved: 368
[Submit][Status][Discuss]

6
2
13
134
8897
1234567654321
1000000000000

Prime
Prime
67
41
4649
5

## HINT

裸的rho。跑得好慢，30多秒，是不是我乘法那一块比别人慢啊？

#include<cstdio>
#include<iostream>
#include<cstdlib>
#include<algorithm>
#define rep(i,j,k) for(register int i = j; i <= k; i++)
#define dow(i,j,k) for(register int i = j; i >= k; i--)
#define ll long long
using namespace std;

ll s = 0, t = 1; char c = getchar();
while( !isdigit(c) ) { if( c == '-' ) t = -1; c = getchar(); }
while( isdigit(c) ) s = s * 10 + c - 48, c = getchar();
return s * t;
}

const int times = 15;
inline void M(ll &x,ll p) { if( x >= p ) x -= p; }
inline void Max(ll &x,ll v) { if( v > x ) x = v; }
inline ll mul(ll x,ll t,ll p) {
ll ret = 0;
while( t ) {
if( t & 1 ) M(ret += x,p);
M(x += x,p), t >>= 1;
} return ret;
}
inline ll pow(ll x,ll t,ll p) {
ll ret = 1;
while( t ) {
if( t & 1 ) ret = mul(ret,x,p);
x = mul(x,x,p), t >>= 1;
} return ret;
}
inline bool miller_rabin(ll n) {
if( n == 2 || n == 3 || n == 5 ) return 1;
if( n < 2 || !(n&1) || n % 3 == 0 || n % 5 == 0 ) return 0;
ll m = n - 1; int k = 0;
while( !(m & 1) ) m >>= 1, k++;
rep(i,1,times) {
ll x = pow(rand() % (n-1) + 1,m,n), y;
rep(j,1,k) {
y = mul(x,x,n);
if( y == 1 && x != 1 && x != n - 1 ) return 0;
x = y;
}
if( y != 1 ) return 0;
}
return 1;
}

inline ll gcd(ll x,ll y) {
ll t;
while( y ) t = x, x = y, y = t % y;
return x;
}
inline ll rho(ll n,int c) {
ll i = 1, k = 2, x = rand() % (n-1) + 1, y = x;
while( 1 ) {
i++;
x = (mul(x,x,n) + c) % n;
ll d = gcd((y-x+n)%n,n);
if( 1 < d && d < n ) return d;
if( x == y ) return n;
if( i == k ) y = x, k <<= 1;
}
}

ll ans = 0;
inline void find(ll x,int k) {
if( x == 1 || x <= ans ) return;
if( miller_rabin(x) ) { Max(ans,x); return; }
ll p = x; int c = k;
while( p >= x ) p = rho(x,k--);
find(p,c), find(x/p,c);
}

int main() {
while( T-- ) {
ll n = read(); ans = 0;
find(n,120);
if( ans == n ) puts("Prime");
else printf("%lld\n",ans);
}
return 0;
}


posted on 2016-12-15 13:04  83131  阅读(195)  评论(0编辑  收藏  举报