【2017多校訓練08 1002】【HDOJ 6134】Battlestation Operational

典型的數列反演題。

運用莫比烏斯反演的一個結論 $[n = 1] = \sum_{d | n} \mu(d)$，將表達式作以下轉化：

$$ ans = \sum_{i=1}^n \sum_{j=1}^i (\lfloor \frac{i-1}{j} \rfloor + 1) \sum_{d | i \land d | j} \mu(d) \\ = \sum_{d=1}^n \mu(d) \sum_{i=1}^{\lfloor \frac{n}{d} \rfloor} \sum_{j=1}^i (\lfloor \frac{i-1}{j} \rfloor + 1) $$

令$$F_n = \sum_{i=1}^n \sum_{j=1}^i (\lfloor \frac{i-1}{j} \rfloor + 1)$$

則有$$ans = \sum_{d=1}^n \mu(d) F(\lfloor \frac{n}{d} \rfloor)$$

先考慮如何計算Fn.

觀察得知，內層求和與n無直接關聯，不妨直接對F相鄰兩項作差：

$$dF_n = F_n - F_{n-1} \\= \sum_{j=1}^n (\lfloor \frac{n-1}{j} \rfloor + 1) $$

考慮每一個j對每一個n的貢獻。

對於一個給定的j，咱們能夠枚舉$\lfloor \frac{n-1}{j} \rfloor$的取值t，此時j和t將對$[j*t+1, j*(t+1)+1)$這一範圍內全部的n對應的$F_n$產生t+1的貢獻。

由調和級數可知，對全部j枚舉它們在N之內的倍數只須要O(Nlog(N))級別的時間複雜度。咱們只要在枚舉j和t的同時維護一下$dF_n$的相鄰兩項差，最後作兩次前綴和就能夠獲得$F_n$數列了。

再來考慮如何由$F_n$計算$ans_n$。

由$ans_n = \sum_{d=1}^n \mu(d) F(\lfloor \frac{n}{d} \rfloor)$，與上一步相似，一樣能夠考慮每一個d對每一個n的答案的貢獻。對於每一個d，枚舉$\lfloor \frac{n}{d} \rfloor$的取值t，此時d和t對$[t*d, (t+1)*d)$範圍內全部的n對應的$ans_n$產生$\mu(d) * F_t$ 的貢獻。枚舉結束後再作一遍前綴和便可。

 1 #include <iostream>
 2 #include <cmath>
 3 #include <cstdio>
 4 #include <algorithm>
 5 using namespace std;
 6 const int maxn = 1000005, mod = 1000000007;
 7 typedef long long LL;
 8 int mu[maxn], muS[maxn], F[maxn], ans[maxn], P[maxn], pcnt, N;
 9 bool not_p[maxn];
10 
11 void sieve()
12 {
13     mu[1] = 1;
14     for(int i = 2;i < maxn;++i)
15     {
16         if(!not_p[i]) P[pcnt++] = i, mu[i] = -1;
17         for(int j = 0;j < pcnt;++j)
18         {
19             if(i * P[j] >= maxn) break;
20             not_p[i * P[j]] = true;
21             if(i % P[j] == 0)
22             {
23                 mu[i * P[j]] = 0;
24                 break;
25             }
26             else mu[i * P[j]] = -mu[i];
27         }
28     }
29 }
30 
31 void init()
32 {
33     sieve();
34     for(int i = 2;i <= N;++i) muS[i] = muS[i-1] + mu[i];
35     N = 1000000;
36     int L, R;
37     for(int k = 1;k < N;++k)
38     {
39         L = k, R = k+k;
40         for(int t = 1;L < N;++t, L += k, R += k)
41         {
42             F[L+1] = (F[L+1] + t) % mod;
43             if(R < N) F[R+1] = (F[R+1] - t) % mod;
44             //if(L < 3) printf("k = %d, (%d, %d] = %d\n", k, L, N, t);
45         }
46     }
47     for(int i = 2;i <= N;++i) F[i] = (F[i] + F[i-1]) % mod;
48     for(int i = 2;i <= N;++i) F[i] = (F[i] + F[i-1]) % mod;
49     for(int i = 1;i <= N;++i)
50         F[i] = (F[i] + (LL) i * (i+1) / 2) % mod;
51     for(int i = 1;i <= N;++i) //F[d]
52     {
53         for(int j = 1, k = i;k <= N;++j, k += i)
54         {
55             int tmp = (mod + F[j] * mu[i]) % mod;
56             ans[k] = (ans[k] + tmp) % mod;
57             if(k+i <= N) ans[k+i] = (ans[k+i] + mod - tmp) % mod;
58         }
59     }
60     for(int i = 2;i <= N;++i) ans[i] = (ans[i] + ans[i-1]) % mod;
61 }
62 
63 void work()
64 {
65     while(~scanf("%d", &N))
66     {
67         printf("%d\n", ans[N]);
68     }
69 }
70 int main() {
71     // your code goes here
72     int T;
73     init();
74     work();
75     return 0;
76 }

莫比烏斯函數，差分數列