解：
神奇的一批......参观yyb巨神的博客。
大致思路就是第一步枚举gcd，发现后面有个限制是gcd=1，用反演，得到的F(x)是两个等差数列求积。
然后发现有个地方我们除法的除数是乘积，于是换元枚举那个乘积。提到最前面。
稍微化一下，发现后面有个Id * miu，这个东西化成phi。
然后得到一个式子，前半部分是s2(n/i)这个整除分块，后面就要相应的求这个东西i2phi[i]的前缀和来迎合整除分块。
然后就是杜教筛，先设个g，把h(n)写出来发现要消掉一个d2，于是g(x) = x2。
没了。
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 ", ans); 88 return 0; 89 }