BZOJ3994: [SDOI2015]约数个数和（洛谷P3327）

莫比乌斯反演

BZOJ题目传送门
洛谷题目传送门

这道题有一个我不知道的结论：

$$d(nm)=\sum_{i|n}\sum_{j|m}[gcd(i,j)=1]$$
证明过程：

设 $nm = p_1^{x_1}\cdot p_2^{x_2}\cdot p_3^{x_3}\cdots p_k^{x_k}$且 $n = p_1^{y_1}\cdot p_2^{y_2}\cdot p_3^{y_3}\cdots p_k^{y_k}$，则 $m = p_1^{x_1 - y_1}\cdot p_2^{x_2 - y_2}\cdot p_3^{x_3 - y_3}\cdots p_k^{x_k - y_k}$。
设 $i = p_1^{a_1}\cdot p_2^{a_2}\cdot p_3^{a_3}\cdots p_k^{a_k}$，$j = p_1^{b_1}\cdot p_2^{b_2}\cdot p_3^{b_3}\cdots p_k^{b_k}$，要使 $\gcd(i, j) = 1$，必然有 $a_1$ 或 $b_1$ 为 0，如果$a_1$ 为 0，那么 $b_1$有 $x_i - y_1 + 1$种取值；如果 $b_1$为 0，那么$a_1$有$y_1 + 1$种取值，一共有 $x_1 + 1$种取值，对于 $x_i$也如此。
所以满足条件的$i, j$对数为 $\prod(x_i + 1)$，与约数定理给出的形式一样，得证。

——By sengxian

现在我们来推式子：

$$\sum_{i=1}^n\sum_{j=1}^md(ij)$$

$$=\sum_{i=1}^n\sum_{j=1}^m\sum_{a|i}\sum_{b|j}[gcd(a,b)=1]$$

$$=\sum_{i=1}^n\sum_{j=1}^m\sum_{a|i}\sum_{b|j}\sum_{x|gcd(a,b)}\mu(x)$$

$$=\sum_{i=1}^n\sum_{j=1}^m\sum_{x|i,x|j}\mu(x)d(\frac ix)d(\frac jx)$$

$$=\sum_{x=1}^n\mu(x)\sum_{x|i}d(\frac ix)\sum_{x|j}d(\frac jx)$$

$$=\sum_{x=1}^n\mu(x)\sum_{i=1}^{\lfloor\frac nx\rfloor}d(i)\sum_{j=1}^{\lfloor\frac mx\rfloor}d(j)$$

因为$d(x)$是积性函数，所以可以线性筛预处理并算出前缀和，然后分块统计。

代码：

#include<cctype>
#include<cstdio>
#include<cstring>
#include<algorithm>
#define N 50005
#define F inline
using namespace std;
typedef long long LL;
int t,n,m,c,mu[N],p[N],a[N],b[N];
bool f[N];
F char readc(){
    static char buf[100000],*l=buf,*r=buf;
    if (l==r) r=(l=buf)+fread(buf,1,100000,stdin);
    if (l==r) return EOF; return *l++;
}
F int _read(){
    int x=0; char ch=readc();
    while (!isdigit(ch)) ch=readc();
    while (isdigit(ch)) x=(x<<3)+(x<<1)+(ch^48),ch=readc();
    return x;
}
F void mkp(){
    mu[1]=a[1]=1,f[1]=true;
    for (int i=2;i<N;i++){
        if (!f[i]) p[++c]=i,mu[i]=-1,a[i]=2,b[i]=1;
        for (int j=1,k;j<=c&&(k=p[j]*i)<N;j++){
            f[k]=true;
            if (i%p[j]==0){
                mu[k]=0,a[k]=a[i]/(b[i]+1)*(b[i]+2);
                b[k]=b[i]+1; break;
            }
            else mu[k]=-mu[i],a[k]=a[i]*2,b[k]=1;
        }
    }
    for (int i=2;i<N;i++) mu[i]+=mu[i-1],a[i]+=a[i-1];
}
int main(){
    for (mkp(),t=_read();t;t--){
        n=_read(),m=_read(); if (n>m) swap(n,m);
        LL ans=0; int p;
        for (int i=1;i<=n;i=p+1){
            p=min(n/(n/i),m/(m/i));
            ans+=1ll*(mu[p]-mu[i-1])*a[n/i]*a[m/i];
        }
        printf("%lld\n",ans);
    }
    return 0;
}