hdu 5297 Y sequence 容斥+迭代

hdu 5297 Y sequence 容斥+迭代

题目链接：hdu 5297
题意：定义数列A为不能被表示成pow(a,b)(a>0,2<=b<=r)的正整数的升序排列。求A(n)。(0<n<2*1e18，2<=r<=62。)
思路：

思路主要分为两个部分：
首先，先求A(n)的反函数，最基本的思想是把n开次方直接取整去掉，但是由于会产生重排，所以要用容斥去搞，我的做法是先利用莫比乌斯反演预处理出容斥的系数。
之后，外层的做法：由于不合法的数非常少，所以可以直接进行迭代求解。

上代码（虽然比较丑：

#include <cmath>
#include <iostream>
#include <cstring>
#include <cstdio>
#include <cstdlib>
using namespace std;
long long n;
int r;
const long long INF=0x7fffffffffffffff;
const int M=100000;
const int N=63;
const double eps=1e-7;
long long biao[64][M];


int prime[N];
int tt=0;
void get_prime(){
    bool is_prime[N];
    for(int i=1;i<N;i++)is_prime[i]=true;
    for(int i=2;i<N;i++)
        for(int j=i*i;j<N;j+=i)
            is_prime[j]=false;
    for(int i=1;i<N;i++)if(is_prime[i])
        prime[tt++]=i;
}
bool haha[63][63];
void get(int r){
    for(int i=0;prime[i]<=r;i++){
        haha[r][1]=1;
        for(int j=0;prime[j]<=r;j++){
            for(int k=1;k<N;k++){
                if(haha[r][k]&&k*prime[j]<N){
                    haha[r][k*prime[j]]=1;
                }
            }
        }
    }
}

void prepare(){
    get_prime();
    memset(haha,0,sizeof(haha));
    for(int i=0;i<N;i++){get(i);}
    for(int i=0;i<M;i++)biao[0][i]=1;
    for(int i=1;i<N;i++){
        for(int j=1;j<M;j++){
            if(biao[i-1][j]<INF/j)
            biao[i][j]=biao[i-1][j]*j;
            else biao[i][j]=INF;
          //  cout<<j<<" "<<i<<" "<<biao[i][j]<<endl;
        }
    }
}


int mu[N+10];


int hehe[100],tot=0;
void get_mu(){
    mu[1]=1;

    for(int i=1;i<=N;i++)
    for(int j=i+i;j<=N;j+=i){
        mu[j]-=mu[i];
    }
    for(int i=2;i<N;i++)
    {
        if(mu[i]!=0){
            hehe[tot++]=i;
        }
    }
}
int kgh(long long a,int b){
    if(b==2)return (int)(sqrt(a)+eps);
    else if(b==3) return (int)(pow(a,1.0/3)+eps);
    int l=1,h=14,mid;
    if(b==5) h=10000;
    if(b==6) h=5000;
    if(b==7) h=500;
    if(b==10) h=73;
    if(b==11) h=50;
    if(b==12) h=40;
    if(b==13) h=28;
    if(b==14) h=21;
    if(b==15) h=17;
    while(l<h){
        mid=(l+h)>>1;
        if(biao[b][mid]>a) h=mid;
        else l=mid+1;
    }
    return l-1;
}
int data[N][N];
long long ok(long long tp){
    int ans=0;
   for(int t=1;t<data[r][0];t++){
        int i=data[r][t];
        if(mu[i]==1) ans+=(kgh(tp,i)-1);
        else if(mu[i]==-1)ans-=(kgh(tp,i)-1);
    }
    return n-(tp-1+ans);
}
void pp(){
    for(int i=2;i<N;i++){
        data[i][0]=1;
        for(int j=2;j<N;j++){
            if(haha[i][j]&&mu[j]!=0){
                data[i][data[i][0]++]=j;
            }
        }
    }
}
int main(){
//    freopen("1010.in","r",stdin);
//    freopen("data1.out","w",stdout);
    int T;get_mu();
    prepare();
    pp();

   // cout<<kgh(INF/2,11)<<endl;
    scanf("%d",&T);
    long long l,h,mid;
    while(T--){
        scanf("%I64d%d",&n,&r);
        long long tp=ok(n),tn=n;
        while(tp){
            tn+=tp;
            tp=ok(tn);
        }
        printf("%I64d\n",tn);
    }
    return 0;
}