洛谷P3768 简单的数学题

解：

神奇的一批......参观yyb巨神的博客。

大致思路就是第一步枚举gcd，发现后面有个限制是gcd=1，用反演，得到的F(x)是两个等差数列求积。

然后发现有个地方我们除法的除数是乘积，于是换元枚举那个乘积。提到最前面。

稍微化一下，发现后面有个Id * miu，这个东西化成phi。

然后得到一个式子，前半部分是s²(n/i)这个整除分块，后面就要相应的求这个东西i²phi[i]的前缀和来迎合整除分块。

然后就是杜教筛，先设个g，把h(n)写出来发现要消掉一个d²，于是g(x) = x²。

没了。

 1 #include <cstdio>
 2 #include <map>
 3 
 4 typedef long long LL;
 5 const int N = 5000010, T = 5000008;
 6 
 7 LL MO;
 8 
 9 inline LL qpow(LL a, LL b) {
10     LL ans = 1;
11     while(b) {
12         if(b & 1) ans = ans * a % MO;
13         a = a * a % MO;
14         b = b >> 1;
15     }
16     return ans;
17 }
18 
19 std::map<LL, LL> mp;
20 int p[N], top, phi[N];
21 LL F[N], inv2, inv6;
22 bool vis[N];
23 
24 inline void getp(int n) {
25     phi[1] = 1;
26     for(int i = 2; i <= n; i++) {
27         if(!vis[i]) {
28             p[++top] = i;
29             phi[i] = i - 1;
30         }
31         for(int j = 1; j <= top && i * p[j] <= n; j++) {
32             vis[i * p[j]] = 1;
33             if(i % p[j] == 0) {
34                 phi[i * p[j]] = phi[i] * p[j];
35                 break;
36             }
37             phi[i * p[j]] = phi[i] * (p[j] - 1);
38         }
39     }
40     for(int i = 1; i <= n; i++) {
41        F[i] = (F[i - 1] + (1ll * i * i % MO * phi[i] % MO)) % MO;
42     }
43     return;
44 }
45 
46 inline LL s2(LL x) { /// sum[1~n] ^ 2
47     x %= MO;
48     LL temp = (x + 1) * x / 2 % MO;
49     return temp * temp % MO;
50 }
51 
52 inline LL H(LL x) { /// sum of n^3
53     return s2(x);
54 }
55 
56 inline LL G(LL x) { /// sum of n^2
57     x %= MO;
58     return (x << 1 | 1) % MO * (x + 1) % MO * x % MO * inv6 % MO;
59 }
60 
61 LL getF(LL x) {
62     if(x <= 0) return 0;
63     if(x <= T) return F[x];
64     if(mp.count(x)) return mp[x];
65     LL ans = H(x);
66     for(LL i = 2, j; i <= x; i = j + 1) {
67         j = x / (x / i);
68         ans -= (G(j) - G(i - 1) + MO) % MO * getF(x / i) % MO;
69         ans = (ans % MO + MO) % MO;
70     }
71     return mp[x] = ans;
72 }
73 
74 int main() {
75     LL n;
76     scanf("%lld%lld", &MO, &n);
77     getp(T);
78     inv6 = qpow(6, MO - 2);
79     inv2 = (MO + 1) / 2;
80 
81     LL ans = 0;
82     for(LL i = 1, j; i <= n; i = j + 1) {
83         j = n / (n / i);
84         ans += s2(n / i) * (getF(j) - getF(i - 1) + MO) % MO;
85         ans = (ans % MO + MO) % MO;
86     }
87     printf("%lld\n", ans);
88     return 0;
89 }