【BZOJ3529】【莫比乌斯反演 + 树状数组】[Sdoi2014]数表
Description
有一张N×m的数表，其第i行第j列（1 < =i < =礼，1 < =j < =m）的数值为
能同时整除i和j的所有自然数之和。给定a，计算数表中不大于a的数之和。
Input
输入包含多组数据。
输入的第一行一个整数Q表示测试点内的数据组数，接下来Q行，每行三个整数n，m，a(|a| < =10^9)描述一组数据。
Output
对每组数据，输出一行一个整数，表示答案模2^31的值。
Sample Input
2
4 4 3
10 10 5
4 4 3
10 10 5
Sample Output
20
148
148
HINT
1 < =N．m < =10^5 ， 1 < =Q < =2×10^4
Source
【分析】
其实是个比较水的题目。
$\sum\limits_{i < = n\ {\rm{ j < = m}}}^{} { F(gcd(i,j))} $其中gcd(i,j) <= a，F(i)，为i的因子和。
我们可以将公式变形,先忽略gcd(i, j) <= a这个条件。
$\sum\limits_{i = 1}^{\min (n,m)} {F(i) \times num(i)} $ num(i)为i在 gcd(i,j) (i <=n, j <= m)中出现的次数。
易得：
$num(i) = \sum\limits_{i|d} {\mu (\frac{d}{i}) \times \left\lfloor {\frac{n}{d}} \right\rfloor } \times \left\lfloor {\frac{m}{d}} \right\rfloor $
然后原式可以化为:
$\sum\limits_{i = 1}^{\min (n,m)} {F(i)} \times \sum\limits_{i|d} {\mu (\frac{d}{i}) \times \left\lfloor {\frac{n}{d}} \right\rfloor } \times \left\lfloor {\frac{m}{d}} \right\rfloor$
换元:
$\sum\limits_{d = 1}^{\min (n,m)} {\left\lfloor {\frac{m}{d}} \right\rfloor \left\lfloor {\frac{m}{d}} \right\rfloor } \sum\limits_{i|d} {F(i)\mu (\frac{d}{i})} $
然后再把$\sum\limits_{i|d} {F(i)\mu (\frac{d}{i})} $预处理出来,对询问以a为关键字排序,用树状数组记录一段的和，前面用分块，然后就可以做了。
1 /* 2 唐代杜牧的《遣怀》 3 落魄江南载酒行，楚腰纤细掌中轻。 4 十年一觉扬州梦，赢得青楼薄幸名。 5 */ 6 #include <cstdio> 7 #include <cstring> 8 #include <algorithm> 9 #include <cmath> 10 #include <queue> 11 #include <vector> 12 #include <iostream> 13 #include <string> 14 #include <ctime> 15 #include <map> 16 #include <set> 17 #define LOCAL 18 long long MOD = 1000000000 + 7; 19 const int MAXM = 1000 * 1000 + 10; 20 const int MAXN = 100000 + 10; 21 using namespace std; 22 //输入输出优化 23 int read(){ 24 int x = 0, flag = 1; 25 char ch = getchar(); 26 while (ch < '1' || ch >'9') {if (ch == '-') flag = -1; ch = getchar();} 27 while (ch >= '0' && ch <= '9') {x = x * 10 + (ch - '0'); ch = getchar();} 28 return x * flag; 29 } 30 struct FF{ 31 int order;//order表示该num对应的gcd值 32 int num; 33 bool operator < (const FF &b)const{ 34 return num < b.num; 35 } 36 }F[MAXN]; 37 struct QUERY{ 38 int l, r, a, order; 39 bool operator < (const QUERY &b)const{ 40 return a < b.a; 41 } 42 }q[MAXN]; 43 int mu[MAXN], prime[MAXN]; 44 int g[MAXN], C[MAXN], Q, Ans[MAXN]; 45 46 int lowbit(int x){return x & -x;} 47 void add(int x, int val){ 48 while (x <= 100000){ 49 C[x] += val; 50 x += lowbit(x); 51 } 52 return; 53 } 54 int sum(int x){ 55 int cnt = 0; 56 while (x > 0){ 57 cnt += C[x]; 58 x -= lowbit(x); 59 } 60 return cnt; 61 } 62 void prepare(){ 63 memset(prime, 0, sizeof(prime)); 64 memset(C, 0, sizeof(C)); 65 66 mu[1] = 1; 67 for (int i = 2;i <= 100000; i++){ 68 if (!prime[i]){ 69 prime[++prime[0]] = i; 70 mu[i] = -1; 71 } 72 for (int j = 1; j <= prime[0]; j++){ 73 if (i * prime[j] > 100000) break; 74 prime[i * prime[j]] = 1; 75 if (i % prime[j] == 0){ 76 mu[i * prime[j]] = 0; 77 break; 78 }else mu[i * prime[j]] = -mu[i]; 79 } 80 } 81 //F[i]代表i的因数和 82 F[1].num = F[1].order = 1; 83 for (int i = 2; i <= 100000; i++){ 84 int cnt = 0; 85 for (long long j = 1; j * (long long)j <= (long long)i; j++){ 86 if (j * j == i){cnt += j; break;} 87 if (i % j != 0) continue; 88 cnt += j + (i / j); 89 } 90 F[i].num = cnt; 91 F[i].order = i; 92 } 93 sort(F + 1, F + 1 + 100000); 94 //for (int i = 1; i <= 100; i++) printf("%d\n", mu[i]); 95 } 96 void init(){ 97 scanf("%d", &Q); 98 for (int i = 1; i <= Q; i++){ 99 int l = read(), r = read(), a = read(); 100 q[i].l = l; q[i].r = r; 101 q[i].a = a; q[i].order = i; 102 } 103 sort(q + 1, q + 1 + Q); 104 } 105 //直接回答第x个询问 106 int query(int x){ 107 int cnt = 0; 108 for (int i = 1; i <= min(q[x].l, q[x].r); i++){ 109 int t = min(q[x].l / (q[x].l / i), q[x].r / (q[x].r / i)); 110 cnt += (q[x].l / i) * (q[x].r / i) * (sum(t) - sum(i - 1)); 111 i = t; 112 } 113 return cnt; 114 } 115 void work(){ 116 int pos = 1;//表示a现在的大小 117 for (int i = 1; i <= Q; i++){ 118 while (F[pos].num <= q[i].a && pos <= 100000){ 119 for (int j = 1; j * F[pos].order <= 100000; j++) 120 add(j * F[pos].order, F[pos].num * mu[j]); 121 pos++; 122 } 123 Ans[q[i].order] = query(i); 124 } 125 for (int i = 1; i <= Q; i++) printf("%d\n", Ans[i] & 0x7fffffff); 126 } 127 128 int main(){ 129 int T; 130 131 prepare(); 132 init(); 133 work(); 134 return 0; 135 } 136
(∑k=1nakbk)2≤(∑k=1na2k)(∑k=1nb2k)
