hdu 3864 D_num Pollard_rho算法和Miller_Rabin算法

D_num

Time Limit: 2000/1000 MS (Java/Others)    Memory Limit: 65536/32768 K (Java/Others)


Problem Description
Oregon Maple was waiting for Bob When Bob go back home. Oregon Maple asks Bob a problem that as a Positive number N, if there are only four Positive number M makes Gcd(N, M) == M then we called N is a D_num. now, Oregon Maple has some Positive numbers, and if a Positive number N is a D_num , he want to know the four numbers M. But Bob have something to do, so can you help Oregon Maple? 
Gcd is Greatest common divisor.
 
Input
Some cases (case < 100);
Each line have a numeral N(1<=N<10^18)
 
Output
For each N, if N is a D_NUM, then output the four M (if M > 1) which makes Gcd(N, M) = M. output must be Small to large, else output “is not a D_num”.
 
Sample Input
6 10 9
 
Sample Output
2 3 6 2 5 10 is not a D_num
 
Source

题意:一个数的因数是否为4个;是,输出>1的因数;

思路:Pollard_rho算法和Miller_Rabin算法,前一个大质数分解,后面判断大数是否为质数;

   acdream的模板;

#include<iostream>
#include<cstdio>
#pragma comment(linker, "/STACK:1024000000,1024000000")
#include<iostream>
#include<cstdio>
#include<cmath>
#include<string>
#include<queue>
#include<algorithm>
#include<stack>
#include<cstring>
#include<vector>
#include<list>
#include<set>
#include<map>
#include<bitset>
using namespace std;
#define LL unsigned long long
#define pi (4*atan(1.0))
#define eps 1e-4
#define bug(x)  cout<<"bug"<<x<<endl;
const int N=5500,M=1e5+10,inf=1e9+10;
const LL INF=1e18+10,mod=2147493647;

const int Times = 10;
LL ct, cnt;
LL fac[N], num[N];

LL gcd(LL a, LL b)
{
    return b? gcd(b, a % b) : a;
}

LL multi(LL a, LL b, LL m)
{
    LL ans = 0;
    a %= m;
    while(b)
    {
        if(b & 1)
        {
            ans = (ans + a) % m;
            b--;
        }
        b >>= 1;
        a = (a + a) % m;
    }
    return ans;
}

LL quick_mod(LL a, LL b, LL m)
{
    LL ans = 1;
    a %= m;
    while(b)
    {
        if(b & 1)
        {
            ans = multi(ans, a, m);
            b--;
        }
        b >>= 1;
        a = multi(a, a, m);
    }
    return ans;
}

bool Miller_Rabin(LL n)
{
    if(n == 2) return true;
    if(n < 2 || !(n & 1)) return false;
    LL m = n - 1;
    int k = 0;
    while((m & 1) == 0)
    {
        k++;
        m >>= 1;
    }
    for(int i=0; i<Times; i++)
    {
        LL a = rand() % (n - 1) + 1;
        LL x = quick_mod(a, m, n);
        LL y = 0;
        for(int j=0; j<k; j++)
        {
            y = multi(x, x, n);
            if(y == 1 && x != 1 && x != n - 1) return false;
            x = y;
        }
        if(y != 1) return false;
    }
    return true;
}

LL pollard_rho(LL n, LL c)
{
    LL i = 1, k = 2;
    LL x = rand() % (n - 1) + 1;
    LL y = x;
    while(true)
    {
        i++;
        x = (multi(x, x, n) + c) % n;
        LL d = gcd((y - x + n) % n, n);
        if(1 < d && d < n) return d;
        if(y == x) return n;
        if(i == k)
        {
            y = x;
            k <<= 1;
        }
    }
}

void Find(LL n, int c)
{
    if(n == 1) return;
    if(Miller_Rabin(n))
    {
        fac[ct++] = n;
        return ;
    }
    LL p = n;
    LL k = c;
    while(p >= n) p = pollard_rho(p, c--);
    Find(p, k);
    Find(n / p, k);
}
int main()
{
    LL n;
    while(~scanf("%I64d",&n))
    {
        ct = 0;
        Find(n, 120);
        sort(fac, fac + ct);
        num[0] = 1;
        int k = 1;
        for(int i=1; i<ct; i++)
        {
            if(fac[i] == fac[i-1])
                ++num[k-1];
            else
            {
                num[k] = 1;
                fac[k++] = fac[i];
            }
        }
        cnt = k;
        LL ans = 0;
        if(cnt==1)
        {
            if(num[0]==3)printf("%lld %lld %lld
",fac[0],fac[0]*fac[0],n);
            else printf("is not a D_num
");
        }
        else if(cnt==2)
        {
            if(num[0]==1&&num[1]==1)printf("%lld %lld %lld
",fac[0],fac[1],n);
            else printf("is not a D_num
");
        }
        else printf("is not a D_num
");
        //for(int i=0;i<cnt;i++)
        //cout<<num[i]<<" "<<fac[i]<<endl;
    }
    return 0;
}
原文地址:https://www.cnblogs.com/jhz033/p/7107047.html