Problem
Description
设 (d(x)) 为 (x) 的约数个数,给定 (N)、(M),求
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_{pin prime}p^{x_p}) 次方这样的形式,其中 (x) 数列为 (i) 的每个质因子的次数序列。
所以,我们先考虑怎么求 (d(i)) ,我们可以枚举 (i) 的每个质因子选多少个乘起来,由于每一个不同的数分解成的质因数乘法必定是唯一的,所以我们只需要枚举有多少个本质不同的质因数乘法序列即可。由于还要再加上不选这个因数的选择,所以 (d(i)) 的计算方法即为:
从这个角度入手,我们就可以考虑如何计算 (d(ij)) 了,设 (y) 数列为 (j) 的每个质因子次数序列,则也有:
那么如何计算 (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)) 的表达式了:
接下来就是很简单的套路了!
先设 (nle m)
接着,根据常见套路,我们改为枚举 (x,y) :
再把莫比乌斯函数套进去:
再次根据常见套路,我们改为枚举 (d) :
那么又一次根据常见套路,继续减少枚举次数:
为了更直观地体现,我们需要根据乘法分配率将式子分成两个部分分别计算,然后乘在一起:
那么我们就这样得到了最终的式子。
接下来,我们只需要正常整除分块就好了!因为对于式子中的 (lfloorfrac{n}{xd} floor) 来讲呢,我们可以拆成 ((n/d)/x) ,所以就可以预处理一个函数 (f(x)=sum_{i=1}^x lfloorfrac{x}{i} floor) ,则:
再预处理一下 (mu) 函数的前缀和,那么一切都只是整除分块而已,分块后面两个函数的值即可,吸溜。
Code
#include <cstdio>
#include <iostream>
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() {
freopen("3327.in", "r", stdin);
freopen("3327.out", "w", stdout);
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
", ans);
}
}