[BZOJ2820]YY的GCD（莫比乌斯反演）

题目描述

传送门

题解

这题和上一题差不多的…
令$p(i)$表示第i个质数，假设$n<m$

$$\sum\limits_{i=1}^n\sum\limits_{j=1}^m\sum\limits_{d=1}^k[(i,j)=p(d)]$$
$$=\sum\limits_{d=1}^k\sum\limits_{i=1}^{n\over p(d)}\sum\limits_{j=1}^{m\over p(d)}[(i,j)=1]$$
利用反演公式 $[n=1]=\sum\limits_{d|n}\mu(d)$
$$=\sum\limits_{d=1}^k\sum\limits_{i=1}^{n\over p(d)}\sum\limits_{j=1}^{m\over p(d)}\sum\limits_{t|(i,j)}\mu(t)$$
$$=\sum\limits_{d=1}^k\sum\limits_{t=1}^{n\over p(d)}\sum\limits_{i=1}^{n\over p(d)}[t|i]\sum\limits_{j=1}^{m\over p(d)}[t|j]\mu(t)$$
$$=\sum\limits_{d=1}^k\sum\limits_{t=1}^{n\over p(d)}\lfloor{n\over p(d)t}\rfloor\lfloor{m\over p(d)t}\rfloor\mu(t)$$
令 $T=p(d)t$，那么现在 $p(d)$变成了 $T$的约数
$$=\sum\limits_{T=1}^n\lfloor{n\over T}\rfloor\lfloor{m\over T}\rfloor\sum\limits_{p(d)|T}\mu({T\over p(d)})$$
令 $f(n)=\sum\limits_{p(d)|n}\mu({n\over p(d)})$，我们如果想要 $O(\sqrt n)$求解的话就必须做到 $O(1)$求 $f$
直接求$f$只能是 $O(\sqrt n)$的
如何预处理 $f$？
考虑线性筛，但是求$f(n)$的过程中枚举的是n的每一个质因子，并不是n的约数，所以 $f$并不一定是一个积性函数
但是我们仍然可以用类似线性筛的方法求出$f$，关键是要利用莫比乌斯函数的性质
观察 $f$的表达式，由$\mu$的性质可以知道， $f(n)$不为0，当且仅当将n至多含有一个指数大于1的质因子
咦这和 $\mu$的定义很像啊…
显然当n为质数时 $f(n)=1$
考虑用 $f(a)$和质数p推出 $f(a*p)$
1° 当 $\mu(a)=0$时， $f(a*p)=0$
2° 当 $(a,p)=1$时，相当于是在 $f(a)$的每一项上都乘上了一个 $\mu(p)$，然后再加入了一个 $\mu(a)$，即 $f(a*p)=-f(a)+\mu(a)$
3° 当 $(a,p)\ne 1$时， $f(a)$中所有含有 $\mu(...p)$的项都变为0，只有不含p的一项不变，也就是说 $f(a*p)=\mu(a)$

这样求出$f$之后原式变成

$$=\sum\limits_{T=1}^n\lfloor{n\over T}\rfloor\lfloor{m\over T}\rfloor f(T)$$
分块之后$O(n+T(\sqrt n+\sqrt m))$求解

代码

#include<algorithm>
#include<iostream>
#include<cstring>
#include<cstdio>
#include<cmath>
using namespace std;
#define N 10000005
#define LL long long
#define Mod 1000000007

int T,n,m;
int p[N],prime[N],mu[N],f[N];
LL ans;

void get(int n)
{
    mu[1]=1;f[1]=0;
    for (int i=2;i<=n;++i)
    {
        if (!p[i])
        {
            prime[++prime[0]]=i;
            f[i]=1;
            mu[i]=-1;
        }
        for (int j=1;j<=prime[0]&&i*prime[j]<=n;++j)
        {
            p[i*prime[j]]=1;
            if (i%prime[j]==0)
            {
                mu[i*prime[j]]=0;
                if (mu[i]) f[i*prime[j]]=mu[i];
                else f[i*prime[j]]=0;
                break;
            }
            else
            {
                mu[i*prime[j]]=-mu[i];
                f[i*prime[j]]=-f[i]+mu[i];
            }
        }
    }
    for (int i=1;i<=n;++i) f[i]+=f[i-1];
}
int main()
{
    get(10000000);
    scanf("%d",&T);
    while (T--)
    {
        scanf("%d%d",&n,&m);
        if (n>m) swap(n,m);
        ans=0;
        for (int i=1,j=0;i<=n;i=j+1)
        {
            j=min(n/(n/i),m/(m/i));
            ans+=(LL)(n/i)*(m/i)*(f[j]-f[i-1]);
        }
        printf("%lld\n",ans);
    }
}