5528Count a b

题意：考虑一块$N*N$的板，$a[i][j]=i*j\%N$,记f(n)为这块板上非0的数目，求$g(n)=\sum_{i|n}f(i)$,$n<=10^9,T<=20000$

分析：考虑$i*j\%N==0的(i,j)$有什么性质；
不妨枚举$i$,那么合法的j应该满足$gcd(i,n)*j\%N==0$恰好有$gcd(i,n)个$，记$h(n)=\sum_{i=1}^{n}gcd(i,n)$那么$f(n)=n^2-h(n)$,因而$g(n)=(\sum_{i|n}i^2)-(\sum_{i|n}h(i))$
这个和式的两个部分都是积性函数，因此只需要考虑$h(p^k)$,这是一个比较简单的问题

#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;
}