[BZOJ2693]jzptab
试题描述
给出 (n, m)，求 (sum_{i=1}^n sum_{j=1}^m mathrm{lcm}(i, j)) 对 (100000009) 取模后的结果。
多组询问。
输入
一个正整数 (T) 表示数据组数
接下来 (T) 行，每行两个正整数，表示 (n)、(m)
输出
(T) 行，每行一个整数，表示第 (i) 组数据的结果
输入示例
1
4 5
输出示例
122
数据规模及约定
(T le 10000)
(n, m le 10000000)
题解
这道题是那道题的加强版，强化的地方是多组数据而不是 (n, m) 的范围，所以那篇文章后面的优化做法是没有用的。
那篇博文的公式有点丑，所以在这里重推一下……
[sum_{i=1}^n sum_{j=1}^m mathrm{lcm}(i, j) \
= sum_{i=1}^n sum_{j=1}^m frac{ij}{gcd(i, j)}
]
现在考虑枚举最大公约数
[= sum_{g=1}^{min {n, m} } frac{1}{g} sum_{i=1}^n i sum_{j=1}^m j cdot [gcd(i, j) = g] \
= sum_{g=1}^{min {n, m} } frac{1}{g} sum_{i=1}^{lfloor frac{n}{g}
floor} ig sum_{j=1}^{lfloor frac{m}{g}
floor} jg sum_{d|i, d|j} mu(d) \
= sum_{g=1}^{min {n, m} } g sum_{d=1}^{lfloor frac{min {n, m}}{g}
floor} mu(d) sum_{i=1}^{lfloor frac{n}{gd}
floor} id sum_{j=1}^{lfloor frac{m}{gd}
floor} jd \
= sum_{g=1}^{min {n, m} } g sum_{d=1}^{lfloor frac{min {n, m}}{g}
floor} mu(d) cdot d^2 frac{lfloor frac{n}{gd}
floor left( lfloor frac{n}{gd}
floor + 1
ight)}{2} cdot frac{lfloor frac{m}{gd}
floor left( lfloor frac{m}{gd}
floor + 1
ight)}{2} \
= frac{1}{4} sum_{g=1}^{min {n, m} } g sum_{d=1}^{lfloor frac{min {n, m}}{g}
floor} mu(d) cdot d^2 lfloor frac{n}{gd}
floor left( lfloor frac{n}{gd}
floor + 1
ight) lfloor frac{m}{gd}
floor left( lfloor frac{m}{gd}
floor + 1
ight)
]
从这里开始将和那道题不一样了，现在用 (T) 替换 (gd)
[= frac{1}{4} sum_{T=1}^{min {n, m} } lfloor frac{n}{T}
floor left( lfloor frac{n}{T}
floor + 1
ight) lfloor frac{m}{T}
floor left( lfloor frac{m}{T}
floor + 1
ight) sum_{d|T} mu(d) cdot d^2 cdot frac{T}{d} \
= frac{1}{4} sum_{T=1}^{min {n, m} } lfloor frac{n}{T}
floor left( lfloor frac{n}{T}
floor + 1
ight) lfloor frac{m}{T}
floor left( lfloor frac{m}{T}
floor + 1
ight) T sum_{d|T} mu(d) cdot d
]
令 (f(T) = T sum_{d|T} mu(d) cdot d)，则 (f(T)) 是积性函数（两个积性函数的乘积是积性函数，两个积性函数的狄利克雷卷积也是积性函数），我们可以线性筛出来。
如何线性筛？假设 (p) 是质数，则 (f(p) = p(1 - p))，我们要解决的就是已知 (f(T)) 求 (f(Tp))，分两种情况：
- (p|T)，令 (T = T' p^k)，且 (T') 中不含因子 (p)，那么 (f(Tp) = f(T') cdot f(p^{k+1}))，在线性筛的过程 (f(T')) 已知，(f(p^k)) 已知，那么 (f(p^{k+1}) = f(p^k) cdot p)，这是因为 (f(p^k) = p^k sum_{i=0}^k mu(p^i) p^i)，后面 (mu(p^i)(i > 1)) 都是 (0)，所以只有前面那个 (p^k) 会乘上一个 (p)，所以最后得到 (f(Tp) = f(T) cdot p)；
- (p mid T)，显然 (f(Tp) = f(T) cdot f(p))。
这样我们可以 (O(n)) 的时间求出 (f(T)) 的前缀和，后面的询问就变成 (O(sqrt n)) 一次啦。
注意 (100000009) 不是素数，但这并不影响我们求 (4) 的逆元。
#include <iostream>
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <cctype>
#include <algorithm>
using namespace std;
#define rep(i, s, t) for(int i = (s), mi = (t); i <= mi; i++)
#define dwn(i, s, t) for(int i = (s), mi = (t); i >= mi; i--)
int read() {
int x = 0, f = 1; char c = getchar();
while(!isdigit(c)){ if(c == '-') f = -1; c = getchar(); }
while(isdigit(c)){ x = x * 10 + c - '0'; c = getchar(); }
return x * f;
}
#define maxn 10000010
#define MOD 100000009
#define inv4 75000007
#define LL long long
bool vis[maxn];
int prime[maxn], cp, f[maxn], sumf[maxn];
void init() {
int n = maxn - 1;
f[1] = sumf[1] = 1;
rep(i, 2, n) {
if(!vis[i]) prime[++cp] = i, f[i] = (LL)i * (1 - i + MOD) % MOD;
for(int j = 1; j <= cp && i * prime[j] <= n; j++) {
vis[i*prime[j]] = 1;
if(i % prime[j] == 0) {
f[i*prime[j]] = (LL)f[i] * prime[j] % MOD;
break;
}
f[i*prime[j]] = (LL)f[i] * f[prime[j]] % MOD;
}
sumf[i] = sumf[i-1] + f[i];
if(sumf[i] >= MOD) sumf[i] -= MOD;
}
return ;
}
void work() {
int n = read(), m = read(), ans = 0;
for(int T = 1; T <= min(n, m); ) {
int r = min(min(n / (n / T), m / (m / T)), min(n, m));
ans += (LL)(n / T) * (n / T + 1) % MOD * (m / T) % MOD * (m / T + 1) % MOD * (sumf[r] - sumf[T-1] + MOD) % MOD;
if(ans >= MOD) ans -= MOD;
T = r + 1;
}
ans = (LL)ans * inv4 % MOD;
printf("%d
", ans);
return ;
}
int main() {
init();
int T = read();
while(T--) work();
return 0;
}