【BZOJ3529】数表（SDOI2014）-莫比乌斯反演+树状数组

测试地址：数表
做法：本题需要用到莫比乌斯反演+树状数组。
首先忽略a$a$的限制，我们要求的是：
ans=∑ni=1∑mj=1D(gcd(i,j))$ans=\sum _{i=1}^n\sum _{j=1}^{m}D(gcd(i,j))$
其中D(i)$D(i)$指i$i$的约数和。
我们把上式改成枚举公因数d$d$，不妨设n≤m$n\le m$，则：
ans=∑nd=1D(d)∑ni=1∑mj=1[gcd(i,j)=d]$ans=\sum _{d=1}^{n}D(d)\sum _{i=1}^n\sum _{j=1}^m[gcd(i,j)=d]$
令f(d)=∑ni=1∑mj=1[gcd(i,j)=d]$f(d)=\sum _{i=1}^n\sum _{j=1}^m[gcd(i,j)=d]$，我们发现这就是求gcd$gcd$等于d$d$的数对个数（这简直是句废话）。
又令F(d)=∑d|if(i)$F(d)=\sum _{d|i}f(i)$，显然F(d)$F(d)$为gcd$gcd$等于d$d$倍数的数对个数，那么也有F(d)=⌊nd⌋⌊md⌋$F(d)=\lfloor \frac{n}{d}\rfloor \lfloor \frac{m}{d}\rfloor$。根据莫比乌斯反演定理的第二种形式，有：
f(d)=∑d|iμ(id)F(i)=∑d|iμ(id)⌊ni⌋⌊mi⌋$f(d)=\sum _{d|i}\mu (\frac{i}{d})F(i)=\sum _{d|i}\mu (\frac{i}{d})\lfloor \frac{n}{i}\rfloor \lfloor \frac{m}{i}\rfloor$
把这个式子带进ans$ans$那个式子，有：
ans=∑nd=1D(d)∑d|iμ(id)⌊ni⌋⌊mi⌋$ans=\sum _{d=1}^{n}D(d)\sum _{d|i}\mu (\frac{i}{d})\lfloor \frac{n}{i}\rfloor \lfloor \frac{m}{i}\rfloor$
交换d,i$d,i$的位置，有：
ans=∑ni=1⌊ni⌋⌊mi⌋∑d|iD(d)μ(id)$ans=\sum _{i=1}^n\lfloor \frac{n}{i}\rfloor \lfloor \frac{m}{i}\rfloor \sum _{d|i}D(d)\mu (\frac{i}{d})$
显然我们如果预处理出g(i)=∑d|iD(d)μ(id)$g(i)=\sum _{d|i}D(d)\mu (\frac{i}{d})$的前缀和，就可以用数论分块处理每个询问了。
那么现在我们考虑a$a$的限制，实际上a$a$是在限制只有D(i)≤a$D(i)\le a$的D(i)$D(i)$对g(i)$g(i)$有贡献，因此我们把所有询问按a$a$从小到大排序，对新产生的贡献暴力在g$g$中单点修改，因为D(i)$D(i)$对g(j)$g(j)$产生贡献当且仅当i|j$i|j$，那么修改次数显然是nlogn$n\log n$级别的，然后要在数论分块中维护g$g$的前缀和查询，这个显然可以用常数小又好写的树状数组解决。
那么我们就解决了此题，时间复杂度为O(Qn−−√logn+nlog2n)$O(Q\sqrt{n}\log n+n{\log }^{2}n)$。注意到这题模数很特殊，直接用int自然溢出，最后输出时对231−1$2^{31}-1$取个按位与就行了。
以下是本人代码：

#include <bits/stdc++.h>
using namespace std;
int T,p[100010];
int maxn,D[100010],mu[100010],s[100010],ans[20010];
int prime[100010];
bool vis[100010]={0};
struct query
{
    int id,n,m,a;
}q[20010];

bool cmp(query a,query b)
{
    return a.a<b.a;
}

bool cmpp(int a,int b)
{
    return D[a]<D[b];
}

int lowbit(int x)
{
    return x&(-x);
}

void add(int x,int c)
{
    for(int i=x;i<=maxn;i+=lowbit(i))
        s[i]+=c;
}

int sum(int x)
{
    int ans=0;
    for(int i=x;i;i-=lowbit(i))
        ans+=s[i];
    return ans;
}

void init()
{
    scanf("%d",&T);
    maxn=0;
    for(int i=1;i<=T;i++)
    {
        scanf("%lld%lld%lld",&q[i].n,&q[i].m,&q[i].a);
        if (q[i].n>q[i].m) swap(q[i].n,q[i].m);
        maxn=max(maxn,q[i].n);
        q[i].id=i;
    }
    sort(q+1,q+T+1,cmp);

    for(int i=1;i<=maxn;i++)
        for(int j=1;i*j<=maxn;j++)
            D[i*j]=D[i*j]+i;
    for(int i=1;i<=maxn;i++)
        p[i]=i;
    sort(p+1,p+maxn+1,cmpp);
}

void calc_mu()
{
    mu[1]=1;
    prime[0]=0;
    for(int i=2;i<=maxn;i++)
    {
        if (!vis[i])
        {
            prime[++prime[0]]=i;
            mu[i]=-1;
        }
        for(int j=1;j<=prime[0]&&i*prime[j]<=maxn;j++)
        {
            vis[i*prime[j]]=1;
            if (i%prime[j]==0)
            {
                mu[i*prime[j]]=0;
                break;
            }
            mu[i*prime[j]]=-mu[i];
        }
    }
}

int main()
{
    init();
    calc_mu();

    int nowp=1;
    for(int i=1;i<=T;i++)
    {
        while(D[p[nowp]]<=q[i].a)
        {
            for(int j=1;j*p[nowp]<=maxn;j++)
                add(j*p[nowp],D[p[nowp]]*mu[j]);
            nowp++;
        }
        int n=q[i].n,m=q[i].m;
        ans[q[i].id]=0;
        for(int j=n;j>=1;j=max(n/(n/j+1),m/(m/j+1)))
        {
            int l=max(n/(n/j+1)+1,m/(m/j+1)+1),r=j;
            ans[q[i].id]=ans[q[i].id]+(sum(r)-sum(l-1))*(n/j)*(m/j);
        }
    }

    for(int i=1;i<=T;i++)
        printf("%d
",ans[i]&2147483647);

    return 0;
}