【jzoj3623】【SDOI2014】【数表】【莫比乌斯反演】

题目大意

解题思路

$ans=\sum_{i=1}^n\sum_{j=1}^m[f(gcd(i,j))<=a]$其中$f[x]=\sum_{d|x}d$

经过一波莫比乌斯反演，转换主题之后，$ans=\sum_{t=1}^n \lfloor{n\over t}\rfloor \lfloor{m\over t}\rfloor g(t)$其中$g(t)=\sum_{d|t}f(d)*\mu( \lfloor{t\over d}\rfloor )$

可以用埃氏筛法求出f，用欧拉筛法求出mu，对f从小到大排序，对询问按a从小到大离线，每次询问之前把符合条件的f枚举另一个数组合成t，在树状数组中添加。询问时分快处理，在树状数组中查询即可。

code

#include<cstdio>
#include<cmath>
#include<cstring>
#include<algorithm>
#define LF double
#define LL long long
#define ULL unsigned int
#define fo(i,j,k) for(int i=j;i<=k;i++)
#define fd(i,j,k) for(int i=j;i>=k;i--)
#define fr(i,j) for(int i=begin[j];i;i=next[i])
using namespace std;
int const mn=1e5+9,mq=2*1e4+9,inf=1e9+7;
int T,lim,b[mn],ss[mn],mu[mn],tag[mn];
ULL f[mn],t[mn],anss[mn];
struct rec{
    int n,m,a,p;
};
rec a[mq];
bool cmp(int i,int j){
    return f[i]<f[j];
}
bool cm2(rec i,rec j){
    return i.a<j.a;
}
void add(int x,int y){
    while(x<=lim){
        t[x]+=y;
        x+=x&(-x);
    }
}
ULL qury(int x){
    LL ans=0;
    while(x>0){
        ans+=t[x];
        x-=x&(-x);
    }
    return ans;
}
int main(){
    //freopen("table.in","r",stdin);
    //freopen("table.out","w",stdout);
    freopen("d.in","r",stdin);
    freopen("d.out","w",stdout);
    scanf("%d",&T);
    lim=1e5;
    fo(i,1,lim)fo(j,1,lim/i)f[i*j]+=i;
    mu[1]=1;
    fo(i,2,lim){
        if(!tag[i])ss[++ss[0]]=i,mu[i]=-1;
        fo(j,1,ss[0]){
            if(i*ss[j]>lim)break;
            tag[i*ss[j]]=1;
            mu[i*ss[j]]=-mu[i];
            if(i%ss[j]==0){mu[i*ss[j]]=0;break;}
        }
    }
    fo(i,1,lim)b[i]=i;
    sort(b+1,b+lim,cmp);
    fo(i,1,T)scanf("%d%d%d",&a[i].n,&a[i].m,&a[i].a),a[i].p=i;
    sort(a+1,a+T+1,cm2);
    int tmp=1;
    fo(cas,1,T){
        int tt=a[cas].a;
        for(;(tmp<=lim)&&(f[b[tmp]]<=tt);tmp++){
            int t1=lim/b[tmp];
            fo(j,1,t1)add(b[tmp]*j,f[b[tmp]]*mu[j]);
        }
        int n=a[cas].n,m=a[cas].m,i=1;
        if(n>m)swap(n,m);
        ULL ans=0;
        while(i<=n){
            int t1=(n/i),t2(m/i),t3=n/t1,t4=m/t2,j=min(t3,t4);
            ans+=t1*t2*(qury(j)-qury(i-1));
            i=j+1;
        }
        anss[a[cas].p]=ans;
    }
    fo(i,1,T)printf("%llu\n",(unsigned long long)((anss[i]>2147483647)?anss[i]-2147483648:anss[i]));
    return 0;
}