BZOJ2082: [Poi2010]Divine divisor 数论

Description
Tz耍畸形，在寂寞的时候玩一个游戏，他随便找出一个n，然后算出n的所有因子，最后找出一个最大的k，即有一个因子d的k次方为n的因子，那个因子d就是非凡因子啦 。比如48他的非凡因子d就是2，k最大为4，因为16也是48的因子。一个整数的非凡因子可能不止一个，比如6就有3个：2，3，6（k最大是1）。
Input
有两行，第一行给出一个整数m(1<=m<=600)，第二行是由空格分开的m个数ai(2<=ai<=1018(10的18次方))，n就是ai的乘积（n=a1*a2*a3*….*am）。
Output
也有两行（Tz真畸形），第一行为最大的k，第二行为非凡因子d的个数。

分解质因子，取出现次数最大值即可。
然而pollard_pho被卡死了。。。观察最耗时的情况，实际上是两个非常大的质数相乘，这时pollard_pho想试探出其中一个就非常慢。不过这题我们只关心个数，不关心具体是什么，因此可以特殊处理这种情况。
先预处理出10^6以内所有质数，这样一个ai被这些数除完后最多就只剩两个质因子了，可以分情况讨论：
ai为1，可以不再理会。
ai还剩一个质因子，可以直接用miller_rabin检查，然后扔到哈希表中。
ai还剩两个质因子，若相同，可以开方判断，注意此时不要用math库中的取整函数否则会炸精度。若开方成功，将其扔到哈希表中。
否则ai只可能是两个不同质数的乘积了。对所有数去重后两两取gcd，可以得到所有可能在多个数中出现的质因子，全部加入哈希表中，注意这时不要统计出现次数。
再对每个数枚举哈希表中的质数来分解，有不能分解的因子直接当做无名质数，不关心它是多少，只知道它与其它数都无关了，它的出现次数也是可以知道的，就是重复次数。
就可以统计了，时刻注意精度问题，注意快速乘的写法先除后乘。
代码：

#include<cstdio>
#include<cmath>
#include<ext/pb_ds/assoc_container.hpp>
#include<tr1/random>
typedef long long ll;
struct Istream
{
    static const size_t str=1<<16;
    char buf[str],*s,*t;
    Istream():buf(),s(),t(){}
    inline char get()
    {
        return (s==t)?(t=buf+fread(s=buf,1,str,stdin),*s++):(*s++);
    }
    inline Istream& operator>>(int &x)
    {
        register char c;
        do c=get(); while(c<'0'||c>'9');
        x=0;
        while(c>='0'&&c<='9') x=(x<<3)+(x<<1)+c-'0',c=get();
        return *this;
    }
    inline Istream& operator>>(ll &x)
    {
        register char c;
        do c=get(); while(c<'0'||c>'9');
        x=0;
        while(c>='0'&&c<='9') x=(x<<3)+(x<<1)+c-'0',c=get();
        return *this;
    }
}cin;
typedef __gnu_pbds::gp_hash_table<ll,int> map;
inline ll mul(const ll &a,const ll &b,const ll &m)
{
    return (a*b-(ll)((long double)a/m*b)*m+m)%m;
}
inline ll qpow(ll a,ll x,const ll &m)
{
    ll ans=1;
    while(x)
    {
        if(x&1) ans=mul(ans,a,m);
        a=mul(a,a,m);
        x>>=1;
    }
    return ans;
}
inline bool miller_rabin(ll a,const ll &u,int c,const ll &p)
{
    a=qpow(a,u,p);
    while(c--)
    {
        ll b=mul(a,a,p);
        if(b==1&&(a!=1&&a!=p-1)) return 0;
        a=b;
    }
    return a==1;
}
std::tr1::mt19937 rnd(114514);
inline bool is_prime(const ll &p)
{
    if(p==2||p==3||p==5) return 1;
    if(!(p&1)||p<=1) return 0;
    int c=0;
    ll u=p-1;
    while(!(u&1)) ++c,u>>=1;
    if(!miller_rabin(2,u,c,p)||!miller_rabin(3,u,c,p)||!miller_rabin(5,u,c,p)) return 0;
    return 1;
}
inline ll gcd(ll a,ll b)
{
    ll c;
    while(b)
    c=a,a=b,b=c%b;
    return a;
}
struct Int
{
    int a[1000];
    int tot;
    static const int mob=1000000000;
    Int():a(),tot(){a[0]=1;}
    Int& operator <<= (int x)
    {
        while(x--)
        {
            for(int i=0;i<=tot;++i)
            a[i]<<=1;
            for(int i=0;i<=tot;++i)
            if(a[i]>=mob)
            {
                a[i+1]+=a[i]/mob;
                a[i]%=mob;
            }
            if(a[tot+1]) ++tot;
        }
        return *this;
    }
    void print()
    {
        --a[0];
        printf("%d",a[tot]);
        for(int i=tot-1;~i;--i)
        printf("%09d",a[i]);
    }
}fin;
#define gm 1000001
int stk[gm],top=0;
int cnt[gm];
bool np[gm];
void liner()
{
    for(int i=2;i<gm;++i)
    {
        if(!np[i])
        {
            stk[++top]=i;
        }
        for(int j=1;j<=top;++j)
        {
            int x=stk[j];
            if(ll(i)*x>=gm) break;
            np[i*x]=1;
            if(i%x==0) break;
        }
    }
}
map ta,cc,oy;
ll b[301],c[301];
int nb,nc;
int main()
{
    liner();
    int m;ll a;
    cin>>m;
    while(m--)
    {
        cin>>a;
        for(int i=1;i<=top&&a!=1;++i)
        {
            int x=stk[i];
            while(a%x==0)
            {
                ++cnt[i];
                a/=x;
            }
        }
        if(a!=1)
        {
            if(is_prime(a)) ta[a]+=1;
            else
            {
                ll w=(int)sqrt(a);
                if(w*w==a) ta[w]+=2;
                else ++cc[b[++nb]=a];
            }
        }
    }
    if(nb)
    {
        std::sort(b+1,b+nb+1);
        nc=std::unique_copy(b+1,b+nb+1,c+1)-c-1;
        for(int i=1;i<nc;++i)
        {
            for(int j=i+1;j<=nc;++j)
            {
                ll d=gcd(c[i],c[j]);
                if(d!=1) ta[d]+=0;
            }
        }
        for(int i=1;i<=nb;++i)
        {
            a=b[i];
            bool found=0;
            for(map::iterator i=ta.begin(),j=ta.end();i!=j;++i)
            {
                if(a%(i->first)==0)
                {
                    found=1;
                    ++i->second;
                    ++ta[a/(i->first)];
                    break;
                }
            }
            if(!found&&a!=b[i-1]) oy[cc[a]]+=2;
        }
    }
    int ans=1;
    for(int i=1;i<=top;++i)
    if(cnt[i]>ans) ans=cnt[i];
    for(map::iterator i=ta.begin(),j=ta.end();i!=j;++i)
    if(i->second>ans) ans=i->second;
    int sum=oy[ans];
    for(int i=1;i<=top;++i)
    if(cnt[i]==ans) ++sum;
    for(map::iterator i=ta.begin(),j=ta.end();i!=j;++i)
    if(i->second==ans) ++sum;
    fin<<=sum;
    printf("%d\n",ans);
    fin.print();putchar('\n');
    return 0;
}