【BZOJ3529】数表，莫比乌斯反演+BIT

传送门
思路：
复习一下莫比乌斯反演
好像是数论入门题，然后我做了近一天……
式子是这么化没错
$n≤m$

$$\sum^n_{i=1}\sum^m_{j=1}[\sigma(gcd(i,j))≤a]\sigma(gcd(i,j))$$
$$=\sum^n_{i=1}\sum^m_{j=1}\sum^n_{d=1}[gcd(i,j)=d][\sigma(d)≤a]\sigma(d)$$
$i=di,j=dj$
$$=\sum^n_{d=1}[\sigma(d)≤a]\sigma(d)\sum^{\lfloor \frac n d\rfloor}_{i=1}\sum^{\lfloor\frac m d\rfloor}_{j=1}[gcd(i,j)=1]$$
$$=\sum^{n}_{d=1}[\sigma(d)≤a]\sigma(d)\sum^{\lfloor \frac n d\rfloor}_{k=1}\mu(k)\lfloor \frac n {kd}\rfloor\lfloor \frac m {kd}\rfloor$$
一开始我化到这里就停下开始求解了
写了一发离线+BIT交上去发现比暴力多跑了一个点(?)
分析了一下复杂度好像是 $O(Tn)$
被reflash婊了一通
悲伤过后重新来搞，在reflash的提(jiao)示(yu)下化成了下面这样
$P=kd$
$$=\sum^n_{d=1}[\sigma(d)≤a]\sigma(d)\sum^n_{P=1}[d|P]\mu(\frac P d)\lfloor \frac n P\rfloor\lfloor \frac m P\rfloor$$
$$=\sum^n_{P=1}\lfloor \frac n P\rfloor\lfloor \frac m P\rfloor\sum^P_{d=1}[d|P][\sigma(d)≤a]\sigma(d)\mu(\frac P d)$$
设 $F(n)=\sum^n_{i=1}[i|n][\sigma(i)≤a]\sigma(i)\mu(\frac n i)$
$$=\sum^n_{P=1}\lfloor \frac n P\rfloor\lfloor \frac m P\rfloor F(P)$$
用BIT维护 $F$的前缀和，然后处理每个询问的时候$O(\sqrt n)$前缀和搞就可以了
一开始把所有的 $\sigma(i)\mu(\frac n i)$保存下来，大约有 $10^6$个，按照 $\sigma(i)$排序
在BZ上跑了10s……rank倒数第5
又被reflash婊了一通
悲伤过后重新来搞，在reflash的提(jiao)示(yu)下开始现求 $\mu$里的东西（就是枚举d的倍数……）
然后跑了7s……
其实还可以线筛 $\sigma$，因为这是个积性函数（要记录最小质因子的若干次方）
改完以后终于跑到了4s左右……
不过还是很慢
复杂度 $O(n\log^2n+T\sqrt n\log n)$
写一下这篇博客，记住自己菜到家的数论知识，顺便%%%reflash

#include<cstdio>
#include<iostream>
#include<algorithm>
#define LL long long
#define mo 2147483647
#define N 100000
using namespace std;
int cnt,prime[N+5],mu[N+5],t[N+5];
LL d[N+5],ans[N+5];
typedef pair<int,int>xy;
xy sigma[N+5];
bool vis[N+5];
struct node{
    int id,n,m,A;
    bool operator <(const node other)const
    {
        return A<other.A;
    }
}q[20005];
void add(int x,int val)
{
    for (;x<=N;x+=x&-x) d[x]+=val;
}
LL get(int x)
{
    LL t=0;
    for (;x;x-=x&-x) t+=d[x];
    return t;
}
main()
{
    mu[1]=1;
    sigma[1]=make_pair(1,1);
    t[1]=1;
    for (int i=2;i<=N;++i)
    {
        if (!vis[i])
            prime[++prime[0]]=i,
            t[i]=i,
            sigma[i]=make_pair(i+1,i),
            mu[i]=-1;
        for (int j=1;j<=prime[0];++j)
        {
            if (prime[j]*i>N) break;
            vis[prime[j]*i]=1;
            if (i%prime[j])
                t[i*prime[j]]=prime[j],
                sigma[i*prime[j]]=make_pair(sigma[prime[j]].first*sigma[i].first,i*prime[j]),
                mu[i*prime[j]]=-mu[i];
            else
            {
                t[i*prime[j]]=t[i]*prime[j];
                sigma[i*prime[j]]=make_pair(sigma[i].first+sigma[i/t[i]].first*t[i]*prime[j],i*prime[j]);
                mu[i*prime[j]]=0;
                break;
            }
        }
    }
    sort(sigma+1,sigma+N+1);
    int T;
    scanf("%d",&T);
    for (int i=1;i<=T;++i)
    {
        scanf("%d%d%d",&q[i].n,&q[i].m,&q[i].A);
        if (q[i].n>q[i].m) swap(q[i].n,q[i].m);
        q[i].id=i;
    }
    sort(q+1,q+T+1);
    int n,m,last=1;
    for (int i=1;i<=T;++i)
    {
        n=q[i].n;m=q[i].m;
        for (;last<=N;++last)
            if (sigma[last].first<=q[i].A)
            {
                int fi=sigma[last].first,se=sigma[last].second;
                for (int j=1;se*j<=N;++j)
                    add(j*se,fi*mu[j]);
            }
            else
                break;
        for (int last,j=1;j<=n;j=last+1)
            last=min(n/(n/j),m/(m/j)),
            ans[q[i].id]+=1LL*(n/j)*(m/j)*(get(last)-get(j-1));
    }
    for (int i=1;i<=T;++i) printf("%d\n",ans[i]&mo); 
}