Problem
DZY Loves Math
TimeLimit:20Sec
MemoryLimit:512MB
Description
对于正整数n,定义f(n)为n所含质因子的最大幂指数。例如f(1960)=f(23×51×7×2)=3, f(10007)=1, f(1)=0。
给定正整数a,b,求∑i=1a∑j=1bf(gcd(i,j))。
第一行一个数T,表示询问数。
接下来T行,每行两个数a,b,表示一个询问。
Output
对于每一个询问,输出一行一个非负整数作为回答。
1 2 3 4 5
| 4 7558588 9653114 6514903 4451211 7425644 1189442 6335198 4957
|
Sample Output
1 2 3 4
| 35793453939901 14225956593420 4332838845846 15400094813
|
HINT
T≤104
1≤a,b≤107
标签:莫比乌斯反演
Solution
好题,get线筛新姿势。
首先套路转化出莫比乌斯函数:
Ans=i=1∑aj=1∑bf(gcd(i,j))=d=1∑min(a,b)i=1∑aj=1∑bf(d)[gcd(i,j)=d]=d=1∑min(a,b)f(d)i=1∑⌊da⌋j=1∑⌊db⌋k∣gcd(i,j)∑μ(k)=d=1∑min(a,b)f(d)k=1∑⌊dmin(a,b)⌋μ(k)×⌊d×ka⌋×⌊d×kb⌋=t=1∑min(a,b)⌊ta⌋⌊tb⌋d∣t∑μ(dt)×f(d)
这时会发现前面用根号分块很好处理,而后面的部分需要O(1)计算,所以需要线性筛预处理。
令g(t)=∑d∣tμ(dt)×f(d),考虑通过μ与f的性质找到其积性关系。
设t=p1a1×p2a2×⋯×pkak,dt=p1a1′×p2a2′×⋯×pkak′,d=p1a1−a1′×p2a2−a2′×⋯×pkak−ak′
那么一定有0≤a1′,a2′,⋯,ak′≤1,否则μ(⌊dt⌋)=0,不计入总贡献。
- 若a1=a2=⋯=ak=max{a}
- 对于f(d)=max{a}−1的情况,只有一种,即a1′=a2′=⋯=ak′=0。而μ(dt)=(−1)k。故贡献为(a−1)×(−1)k=a×(−1)k−(−1)k;
- 对于f(d)=max{a}的情况,根据组合原理,有∑i=0k−1(−1)i(ki),而又由二项式基本定理知∑i=0k(−1)k(ki)=0,因而贡献为(−1)k(kk)=(−1)k。
- 故此情况g(d)=a×(−1)k−(−1)k+(−1)k=a×(−1)k。
- 若∃i,j使得i=j,ai=aj
- 不论f(d)=max{a}还是f(d)=max{a}−1,都存在至少一个质因数pr使得ar′不论取0还是1对f(d)的取值都没有影响。然而ar′取0或1会使得μ(dt)取到−1或1,此处的贡献为f(d)+(−f(d))=0,一定全部被抵消。
- 故此情况g(d)=0。
综上,线性筛预处理g(t)需要知道每个数最小的质因数的次数num和最小质因数的幂指数次幂sp,这样看是否有num[i×pri[j]]=num[i/sp[i]]即可知i×pri[j]的最小与次小质因数的次数是否相等,由此可判断是情况1还是情况2。这样先线性筛预处理后,对每次询问O(min(a,b))进行根号分块,即可达到O(Tmin(a,b))的复杂度。
Code
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36
| #include <bits/stdc++.h> #define MAX_N 10000000 using namespace std; typedef long long lnt; template <class T> inline void read(T &x) { x = 0; int c = getchar(), f = 1; for (; !isdigit(c); c = getchar()) if (c == 45) f = -1; for (; isdigit(c); c = getchar()) (x *= 10) += f*(c-'0'); } int g[MAX_N+5], sp[MAX_N+5], num[MAX_N+5]; int pri[MAX_N+5], cnt; lnt ans; bool NotPri[MAX_N+5]; void init(int n) { NotPri[0] = NotPri[1] = true; for (int i = 2; i <= n; i++) { if (!NotPri[i]) pri[cnt++] = sp[i] = i, g[i] = num[i] = 1; for (int j = 0; j < cnt; j++) { if (i*pri[j] > n) break; NotPri[i*pri[j]] = true; if (i%pri[j]) sp[i*pri[j]] = pri[j], num[i*pri[j]] = 1, g[i*pri[j]] = num[i] == 1 ? -g[i] : 0; else sp[i*pri[j]] = sp[i]*pri[j], num[i*pri[j]] = num[i]+1, g[i*pri[j]] = sp[i] == i ? 1 : (num[i/sp[i]] == num[i*pri[j]] ? -g[i/sp[i]] : 0); if (i%pri[j] == 0) break; } } for (int i = 1; i <= n; i++) g[i] += g[i-1]; } int main() { int T; read(T), init(MAX_N); while (T--) { int a, b; read(a), read(b), ans = 0; for (int l = 1, r; l <= min(a,b); l = r+1) r = min(a/(a/l), b/(b/l)), ans += 1LL*(a/l)*(b/l)*(g[r]-g[l-1]); printf("%lld\n", ans); } return 0; }
|