真 TM 毒瘤……
典型的式子好推不可实现题……
一个显然的想法是先枚举 \(r = \lfloor\sqrt[3]i\rfloor\)，然后再考虑 \(i\) 的贡献。
即 \(\sum\limits_{i=1}^{\lfloor\sqrt[3]i\rfloor}\gcd(\lfloor\sqrt[3]i\rfloor,i) = \sum\limits_{r = 1}^{\lfloor\sqrt[3]i\rfloor}\sum\limits_{i=1}^n[\lfloor\sqrt[3]i\rfloor=r]\gcd(r,i)\)。
又因为 \(\lfloor\sqrt[3]i\rfloor = r\) 等价于 \(r^3 \le i < (r+1)^3\)，相当于把这些数分成了 \(\lfloor\sqrt[3]i\rfloor\) 块。
得 \(\sum\limits_{i=1}^{\lfloor\sqrt[3]i\rfloor}\gcd(\lfloor\sqrt[3]i\rfloor,i) = \sum\limits_{i={\lfloor\sqrt[3]n\rfloor}^3}^n\gcd(\lfloor\sqrt[3]i\rfloor,i) + \sum\limits_{r=1}^{\lfloor\sqrt[3]i\rfloor-1}\sum\limits_{i=r^3}^{(r+1)^3-1}\gcd(r,i)\)。
根据基本反演套路易得 \(\sum\limits_{i=1}^n\gcd(r,i) = \sum\limits_{d|r}\lfloor\frac nd\rfloor\varphi(d)\)。
于是两个和式分别算即可。
第一个和式需要 \(\lfloor\sqrt[3]n\rfloor\) 的约数的欧拉函数值，可以考虑类似杜教筛的思想，预处理 \(n^{2/9}\) 之内的值，然后别的值分解质因数硬算（只需要 \(n^{1/6}\) 以下的质数）。
复杂度比较难算，姑且认为是 \(O(n^{2/9})\)（实际上跑不满）。
（实际上可以通过 DFS 做到 \(O(n^{1/6})\)，因为都是 \(\lfloor\sqrt[3]n\rfloor\) 的约数，但是我懒得写）
第二个和式，根据上式得 \[ \begin{align*} &\sum\limits_{r=1}^{\lfloor\sqrt[3]i\rfloor-1}\sum\limits_{i=r^3}^{(r+1)^3-1}\gcd(r,i)\\ =&\sum\limits_{r=1}^{\lfloor\sqrt[3]i\rfloor-1}\sum\limits_{d|r}\left(\left\lfloor\frac{(r+1)^3-1}{d}\right\rfloor - \left\lfloor\frac{r^3-1}{d}\right\rfloor\right)\varphi(d) \\ =&\sum\limits_{d=1}^{\lfloor\sqrt[3]i\rfloor-1}\varphi(d)\sum\limits_{r=1}^{\left\lfloor\frac{\lfloor\sqrt[3]i\rfloor-1}d\right\rfloor}\left(\left\lfloor\frac{(dr+1)^3-1}{d}\right\rfloor - \left\lfloor\frac{(dr)^3-1}{d}\right\rfloor\right) \end{align*} \]
注意到 \(\left\lfloor\frac{(dr+1)^3-1}{d}\right\rfloor = \left\lfloor\frac{d^3r^3 + 3d^2r^2 + 3dr}{d}\right\rfloor = d^2r^3 + 3dr^2 + 3r\)，且 \(\left\lfloor\frac{(dr)^3-1}{d}\right\rfloor = \left\lfloor\frac{d^3r^3-1}{d}\right\rfloor = d^2r^3 - 1\)。
故 \[
\begin{align*}
&\sum\limits_{d=1}^{\lfloor\sqrt[3]i\rfloor-1}\varphi(d)\sum\limits_{r=1}^{\left\lfloor\frac{\lfloor\sqrt[3]i\rfloor-1}d\right\rfloor}\left(\left\lfloor\frac{(dr+1)^3-1}{d}\right\rfloor - \left\lfloor\frac{(dr)^3-1}{d}\right\rfloor\right)\\
=&\sum\limits_{d=1}^{\lfloor\sqrt[3]i\rfloor-1}\varphi(d)\sum\limits_{r=1}^{\left\lfloor\frac{\lfloor\sqrt[3]i\rfloor-1}d\right\rfloor}(3dr^2 + 3r + 1)
\end{align*}
\]
然后稍微杜教筛一下就行了。
代码： 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
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
using namespace std;
template<class T>
inline void read(T &x)
{
char ch = 0,w = 0;
for(;ch < '0' || ch > '9';w |= ch == '-',ch = getchar());
for(x = 0;ch >= '0' && ch <= '9';x = (x << 3) + (x << 1) + (ch ^ '0'),ch = getchar());
w && (x = -x);
}
const int MX = 5e6;
const int mod = 998244353;
const int inv2 = 499122177;
const int inv6 = 166374059;
__int128 n;
long long rt;
int vis[MX + 5],cnt,prime[MX + 5],phi[MX + 5],sum[MX + 5];
inline long long cubert(__int128 n)
{
__int128 l = 1,r = n,mid,ret;
while(l <= r)
mid = l + r >> 1,mid <= n / mid / mid ? (l = mid + 1,ret = mid) : (r = mid - 1);
return ret;
}
int w[3][MX + 5];
inline int &mem_calc(long long x)
{
return w[0][rt / x];
}
inline int &mem_sphi(long long x)
{
return w[1][rt / x];
}
inline int &mem_sidphi(long long x)
{
return w[2][rt / x];
}
int calc(long long n)
{
if(n <= MX)
return phi[n] - phi[n - 1];
if(~mem_calc(n))
return mem_calc(n);
long long ret = n,t = n;
for(register int i = 1;(long long)prime[i] * prime[i] <= n;++i)
{
!(t % prime[i]) && (ret -= ret / prime[i]);
for(;!(t % prime[i]);t /= prime[i]);
}
t > 1 && (ret -= ret / t);
return mem_calc(n) = ret % mod;
}
int sphi(long long n)
{
if(n <= MX)
return phi[n];
if(~mem_sphi(n))
return mem_sphi(n);
long long ret = (n % mod) * ((n + 1) % mod) % mod * inv2 % mod;
for(register long long l = 2,r;l <= n;l = r + 1)
{
r = n / (n / l);
ret = (ret - (r - l + 1) % mod * sphi(n / l) % mod + mod) % mod;
}
return mem_sphi(n) = ret;
}
int sidphi(long long n)
{
if(n <= MX)
return sum[n];
if(~mem_sidphi(n))
return mem_sidphi(n);
long long ret = (n % mod) * ((n + 1) % mod) % mod * ((2 * n + 1) % mod) % mod * inv6 % mod;
for(register long long l = 2,r;l <= n;l = r + 1)
{
r = n / (n / l);
ret = (ret - ((l + r) % mod) * ((r - l + 1) % mod) % mod * inv2 % mod * sidphi(n / l) % mod + mod) % mod;
}
return mem_sidphi(n) = ret;
}
int ans;
int main()
{
memset(w,-1,sizeof w),sum[1] = phi[1] = 1;
for(register int i = 2;i <= MX;++i)
{
if(!vis[i])
phi[prime[++cnt] = i] = i - 1;
for(register int j = 1;j <= cnt && i * prime[j] <= MX;++j)
{
vis[i * prime[j]] = 1;
if(!(i % prime[j]))
{
phi[i * prime[j]] = phi[i] * prime[j];
break;
}
phi[i * prime[j]] = phi[i] * (prime[j] - 1);
}
sum[i] = (sum[i - 1] + (long long)i * phi[i] % mod) % mod,phi[i] = (phi[i] + phi[i - 1]) % mod;
}
read(n),rt = cubert(n);
for(register int i = 1;(long long)i * i <= rt;++i)
if(!(rt % i))
{
ans = (ans + (long long)((int)(n / i % mod) - (int)(((__int128)rt * rt * rt - 1) / i % mod) + mod) * calc(i) % mod) % mod;
if(i ^ rt / i)
ans = (ans + (long long)((int)(n / (rt / i) % mod) - (int)(((__int128)rt * rt * rt - 1) / (rt / i) % mod) + mod) * calc(rt / i) % mod) % mod;
}
--rt;
for(register long long l = 1,r;l <= rt;l = r + 1)
{
r = rt / (rt / l);
ans = (ans + 3LL * (sidphi(r) - sidphi(l - 1) + mod) % mod * (rt / l % mod) % mod * ((rt / l + 1) % mod) % mod * ((rt / l * 2 + 1) % mod) % mod * inv6 % mod) % mod;
ans = (ans + (long long)(sphi(r) - sphi(l - 1) + mod) * (3LL * (rt / l % mod) % mod * ((rt / l + 1) % mod) % mod * inv2 % mod + (rt / l) % mod) % mod) % mod;
}
printf("%d\n",ans);
}