HDU 4344 随机法判素数(费马小定理

#include <cstdio>
#include <ctime>
#include <cmath>
#include <algorithm>
using namespace std;
typedef long long ll;
const int N = 108;
const int S = 10;

ll mult_mod(ll a, ll b, ll c) {
    a %= c;
    b %= c;
    ll ret = 0;
    while(b) {
        if(b&1) ret = (ret + a) % c;
        a = (a + a) % c;
        b >>= 1;
    }
    return ret;
}
ll pow_mod(ll x, ll n, ll mod) {
    if(n == 1) return x % mod;
    x %= mod;
    ll tmp = x, ret = 1;
    while(n > 0){
        if(n&1) ret = mult_mod(ret, tmp, mod);
        tmp = mult_mod(tmp, tmp, mod);
        n >>= 1;
    }
    return ret;
}

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;
    
}
bool Miller_Rabin(ll n) {
    if(n < 2) return false;
    if(n==2||n==3||n==5||n==7) return true;
    if(n%2==0||n%3==0||n%5==0||n%7==0) return false;
    
    ll x = n - 1, t = 0;
    while((x&1)==0) {
        x >>= 1;
        t ++;
    }
    for(int i = 0; i < S; i ++) {
        ll a = rand()%(n-1) +1;
        if(check(a, n, x, t)) return false;
    }
    return true;
}
ll gcd(ll a, ll b) {
    if(a < 0) return gcd(-a, b);
    if(b < 0) return gcd(a, -b);
    while(a > 0 && b > 0) {
        if(a > b) a %= b;
        else b %= a;
    }
    return a+b;
}
ll Pollard_rho(ll x, ll c) {
    ll i = 1, k = 2;
    ll x0 = ((rand() % x) + x) % 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;
        }
    }
    
}
ll P[N], tot;

void findfac(ll n) {
    if(Miller_Rabin(n)) {
        P[tot++] = 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;scanf("%d", &T);
    while(T-- > 0) {
        ll n;scanf("%I64d", &n);
        
        tot = 0;
        findfac(n);
        sort(P, P + tot);
        ll t = 0, ans = 0;
        for(int i = 0; i < tot; i ++) {
            if(!i || P[i] != P[i-1]) {
                ans += pow_mod(P[i], count(P, P+tot, P[i]), n);
                t ++;
            }
        }
        printf("%I64d %I64d\n", t, t==1?n/P[0]:ans);
    }
    return 0;
}

posted @ 2015-01-03 10:33  zfyouxi  阅读(161)  评论(0编辑  收藏  举报