5528Count a b

来源:互联网 发布:me域名多少钱 编辑:程序博客网 时间:2024/05/17 22:16

题意:考虑一块NN的板,a[i][j]=ij%N,记f(n)为这块板上非0的数目,求g(n)=i|nf(i),n<=109,T<=20000

分析:考虑ij%N==0(i,j)有什么性质;
不妨枚举i,那么合法的j应该满足gcd(i,n)j%N==0恰好有gcd(i,n),记h(n)=ni=1gcd(i,n)那么f(n)=n2h(n),因而g(n)=(i|ni2)(i|nh(i))
这个和式的两个部分都是积性函数,因此只需要考虑h(pk),这是一个比较简单的问题

#include<stdio.h>#include<algorithm>#include<math.h>#include<string.h>#include<string>#include<vector>#include<set>#include<queue>#include<time.h>#include<assert.h>#include<iostream>using namespace std;typedef long long LL;typedef unsigned long long ULL;const int Maxn=1000020;LL f[Maxn],g[Maxn];bool isp[Maxn];int least[Maxn];vector<int>pri;void getp(){    for(int i=2;i<Maxn;i++){        if(!isp[i])pri.push_back(i),least[i]=i;        for(int j=0;j<pri.size()&&pri[j]*i<Maxn;j++){            isp[pri[j]*i]=1;            least[pri[j]*i]=pri[j];            if(i%pri[j]==0){break;}        }    }}vector<int>yinzi;int mul(int a,int b,int n){return 1LL*a*b%n;}int powmod(int x,int y,int mod){    int ret=1;    x%=mod;    while(y){        if(y&1)ret=mul(ret,x,mod);        y>>=1;        x=mul(x,x,mod);    }    return ret;}bool check(int a,int n){    int m=n-1;    int s=0;    while(~m&1)m>>=1,s++;    int x=powmod(a,m,n);    for(int i=1;i<=s;i++){        int y=mul(x,x,n);        if(y==1&&x!=1&&x!=n-1)return 0;        x=y;    }    return x==1;}int Miller_Rabin(int x){    if(x==1)return 0;    if(x==2)return 1;    if(~x&1)return 0;    int times=3;//Pro=1/4^times    while(times--)if(!check(rand()%(x-2)+2,x))return 0;    return 1;}int Pollard_rho(int n,int c){    int factor=1;    int cirsize=1;    int x=rand()%n,xfixed=x;    while(factor==1){        for(int i=0;i<cirsize;i++){            x=(mul(x,x,n)+c)%n;            if(x==xfixed)return 1;            int d=__gcd(n,(x-xfixed+n)%n);//            if(d>1&&d<n)return d;        }        cirsize<<=1;        xfixed=x;    }    assert(0);}void findfac(int x,int c){    if(x==1)return;    if(x<Maxn){        yinzi.push_back(least[x]);        findfac(x/least[x],c);        return;    }    if(Miller_Rabin(x)){//ispri        yinzi.push_back(x);        return;    }    int tmp;    while((tmp=Pollard_rho(x,c--))==1);    findfac(x/tmp,c);findfac(tmp,c);}int main(){    /*       for(int i=1;i<=100;i++){       LL tot=0;       for(int j=1;j<=i;j++)tot+=__gcd(i,j);       f[i]=tot;       }    //for(int i=1;i<=100;i++)printf("%lld ",f[i]);puts("");    for(int i=1;i<=100;i++){    for(int j=1;j<=i;j++){    if(i%j==0)g[i]+=f[j];    }    }    for(int i=1;i<=100;i++)printf("%lld ",f[i]);puts("");    */    getp();    int _;scanf("%d",&_);    while(_--){        yinzi.clear();        int x;scanf("%d",&x);        findfac(x,2178);        sort(yinzi.begin(),yinzi.end());        ULL ans1=1,ans2=1;        for(int i=0,j;i<yinzi.size();i=j){            for(j=i;j<yinzi.size()&&yinzi[j]==yinzi[i];j++);            //printf("pri=%d k=%d\n",yinzi[i],j-i);            ULL tmp=1,tmp2=1;            ULL tmpans=1,tmpans2=1;            ULL mul=(ULL)yinzi[i]*(ULL)yinzi[i];            for(int it=1;it<=j-i;it++){                tmpans+=(it+1)*tmp*(ULL)yinzi[i]-it*tmp;                tmp*=(ULL)yinzi[i];                tmp2*=mul;                tmpans2+=tmp2;            }            ans2*=tmpans;            ans1*=tmpans2;        }        printf("%llu\n",ans1-ans2);    }    return 0;}
0 0
原创粉丝点击