【SDOI 2015】约数个数和
Problem
Description
设 \(d(x)\) 为 \(x\) 的约数个数,给定 \(N\)、\(M\),求
\[ \sum_{i=1}^N \sum_{j=1}^M d(ij) \]
Input Format
输入文件包含多组测试数据。
第一行,一个整数 \(T\),表示测试数据的组数。
接下来的 \(T\) 行,每行两个整数 \(N\)、\(M\)。
Output Format
\(T\) 行,每行一个整数,表示你所求的答案。
Sample
Input
2
7 4
5 6
Output
110
121
Range
对于所有的数据,\(1 \leq N, M \leq 50000,\ 1 \leq T \leq 50000\)。
Algorithm
莫比乌斯反演
Mentality
显然是莫比乌斯反演啊!
不过做这题的话,我们不难发现:\(d(ij)\) 非常的不好求!如果枚举因数的话,那肯定会出锅,因为 \(ij\) 最高有 \(10^{10}\) ,所以必定超时。那怎么办呢?
可以考虑从 \(i\) 与 \(j\) 的因数入手。
那么我们先考虑 \(i\) ,对于任意一个数 \(i\) ,它必定都能分解成 \(i=\prod_{p\in prime}p^{x_p}\) 次方这样的形式,其中 \(x\) 数列为 \(i\) 的每个质因子的次数序列。
所以,我们先考虑怎么求 \(d(i)\) ,我们可以枚举 \(i\) 的每个质因子选多少个乘起来,由于每一个不同的数分解成的质因数乘法必定是唯一的,所以我们只需要枚举有多少个本质不同的质因数乘法序列即可。由于还要再加上不选这个因数的选择,所以 \(d(i)\) 的计算方法即为:
\[ d(i)=\prod_{p\in prime}(x_p+1) \]
从这个角度入手,我们就可以考虑如何计算 \(d(ij)\) 了,设 \(y\) 数列为 \(j\) 的每个质因子次数序列,则也有:
\[ d(j)=\prod_{p\in prime}(y_p+1) \]
那么如何计算 \(d(ij)\) 呢?其实从枚举质因子的角度来看,就很简单了。我们可以枚举 \(i\) 和 \(j\) 的因数 \(a|i\) 与 \(b|j\) ,则 \(ab\) 亦为 \(ij\) 的因数之一。不难发现,我们的重点是在枚举两个因数的时候,避免两个因数的质因数分解序列的乘积与之前出现过的乘积有重合,那怎么办呢?
回忆一下:当我们思考怎么计算 \(d(i)\) 的时候,我们的第一想法是可以枚举每个质因子选多少个,这种时候,对于一个质因子 \(p\) 而言,它在枚举序列中会被枚举出 \(x_p\) 种不同的大小,即:\(p^1,p^2,\dots ,p^{x_p}\) 这几个不同的数。也即 \(p\) 这一位总共会被枚举 \(x_p\) 次。
那么我们发现,当在计算 \(d(ij)\) 时,只需要保证每个质因子被枚举相应次数即可。
换句话说,对于 \(p|ij\) ,它理应被枚举 \(x_p+y_p\) 次!
那么当我们枚举 \(a|i\) 和 \(b|j\) 时,应该怎么处理呢?
其实这时,你会发现一个很简单的事情,那就是我们只需要确保 \(gcd(a,b)=1\) ,那么 \(a,b\) 中就不会含有相同的质因子。那么当我们枚举某个质因子 \(p|a\) 的时候,它会被枚举 \(x_p\) 次,而当我们枚举 \(p|b\) 时,它会被枚举 \(y_p\) 次!
换句话说,只需要确保 \(gcd(a,b)=1\) ,就能使枚举过程中,不会出现重复的质因子序列了!
所以,我们可以得到 \(d(ij)\) 的表达式了:
\[ d(ij)=\sum_{x|i}\sum_{y|j}[gcd(x,y)=1] \]
接下来就是很简单的套路了!
先设 \(n\le m\)
\[ \sum_{i=1}^n \sum_{j=1}^m d(ij) \]
\[ \sum_{i=1}^n \sum_{j=1}^m \sum_{x|i}\sum_{y|j}[gcd(x,y)=1] \]
接着,根据常见套路,我们改为枚举 \(x,y\) :
\[ \sum_{x=1}^n \sum_{y=1}^m [gcd(x,y)=1]\sum_{i=1}^{n}\sum_{j=1}^m [x|n\&\&y|j] \]
\[ \sum_{x=1}^n \sum_{y=1}^m [gcd(x,y)=1]\lfloor\frac{n}{x}\rfloor\lfloor\frac{m}{y}\rfloor \]
再把莫比乌斯函数套进去:
\[ \sum_{x=1}^n \sum_{y=1}^m\sum_{d|gcd(x,y)}\mu(d)\lfloor\frac{n}{x}\rfloor\lfloor\frac{m}{y}\rfloor \]
再次根据常见套路,我们改为枚举 \(d\) :
\[ \sum_{d=1}^n\mu(d)\sum_{x=1}^n\sum_{y=1}^m[d|x\&\&d|y]\lfloor\frac{n}{x}\rfloor\lfloor\frac{m}{y}\rfloor \]
那么又一次根据常见套路,继续减少枚举次数:
\[ \sum_{d=1}^n\mu(d)\sum_{x=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{y=1}^{\lfloor\frac{m}{d}\rfloor}\lfloor\frac{n}{xd}\rfloor\lfloor\frac{m}{yd}\rfloor \]
为了更直观地体现,我们需要根据乘法分配率将式子分成两个部分分别计算,然后乘在一起:
\[ \sum_{d=1}^n\mu(d)(\sum_{x=1}^{\lfloor\frac{n}{d}\rfloor}\lfloor\frac{n}{xd}\rfloor)(\sum_{y=1}^{\lfloor\frac{m}{d}\rfloor}\lfloor\frac{m}{yd}\rfloor) \]
那么我们就这样得到了最终的式子。
接下来,我们只需要正常整除分块就好了!因为对于式子中的 \(\lfloor\frac{n}{xd}\rfloor\) 来讲呢,我们可以拆成 \((n/d)/x\) ,所以就可以预处理一个函数 \(f(x)=\sum_{i=1}^x \lfloor\frac{x}{i}\rfloor\) ,则:
\[ ans=\sum_{d=1}^n\mu(d)*f(\lfloor\frac{n}{d}\rfloor)*f(\lfloor\frac{m}{d}\rfloor) \]
再预处理一下 \(\mu\) 函数的前缀和,那么一切都只是整除分块而已,分块后面两个函数的值即可,吸溜。
Code
#include<iostream>
#include<cstdio>
using namespace std;
int T,n,m;
int cntp,pri[50001];
long long ans,mu[50001],f[50001];
bool vis[50001];
void init()
{
mu[1]=1;
for(int i=2;i<=50000;i++)
{
if(!vis[i])pri[++cntp]=i,mu[i]=-1;
for(int j=1;j<=cntp&&pri[j]*i<=50000;j++)
{
vis[pri[j]*i]=true;
if(!(i%pri[j]))break;
mu[pri[j]*i]=-mu[i];
}
}//求莫比乌斯函数
for(int i=1;i<=50000;i++)mu[i]+=mu[i-1];//前缀和
for(int i=1;i<=50000;i++)
for(int l=1,r;l<=i;l=r+1)
{
r=i/(i/l);
f[i]+=1ll*(r-l+1)*(i/l);
}//求 f 函数
}
int main()
{
cin>>T;
init();
while(T--)
{
scanf("%d%d",&n,&m);
if(n>m)swap(n,m);
ans=0;
for(int l=1,r;l<=n;l=r+1)
{
r=min(n/(n/l),m/(m/l));
ans+=1ll*(mu[r]-mu[l-1])*f[n/l]*f[m/l];//整除分块计算答案
}
printf("%lld\n",ans);
}
}