目录
- @[email protected]
- @[email protected]
- @accepted [email protected]
- @[email protected]
@[email protected]
给定 n,求:
\[\sum_{i=1}^{n}gcd(\lfloor^3\sqrt{i}\rfloor, i)\mod 998244353\]
Input
第一行包含一个整数 T(1≤T≤11) 描述数据组数。
接下来 T 行,每行一个整数 n (1≤n≤10^21) 描述询问。
Output
对于每组询问,输出答案 mod 998244353。
Sample Input
10
64
180
526
267
775
649
749
908
300
255
Sample Output
103
327
1069
522
1693
1379
1631
1998
606
492
@[email protected]
考虑 \(\lfloor^3\sqrt{n}\rfloor\) 其实只有 10^7 种取值,令 \(x =\lfloor^3\sqrt{n}\rfloor\),我们的问题其实是在求:
\[ans=\sum_{i=1}^{x-1}\sum_{j=i^3}^{(i+1)^3-1}gcd(i, j) + \sum_{j=x^3}^{n}gcd(x, j)\]
我们现在考虑怎么求 \(f(n, i)=\sum_{j=1}^{n}gcd(i, j)\),然后用前缀和相减的方法搞一搞。
考虑一波反演变成 \(f(n, i)=\sum_{p|i}\sum_{d|p}\lfloor\frac{n}{p}\rfloor*\mu(\frac{p}{d})*d\)。
注意到 \(i|i^3, i|(i+1)^3-1\),放在上式相当于 \(n=i^3 或 (i+1)^3-1\),于是可以得到 \(p|n\)。
由此,考虑对式子进行进一步地变形,得到 \(f(n, i)=n*\sum_{p|i}\sum_{d|p}\frac{d}{p}*\mu(\frac{p}{d})\)。
令整数 \(a = \frac{p}{d}\)。考虑对于某个 i,一共有 \(\sigma(\frac{i}{a})\) 对 (p, d) 满足 \(a = \frac{p}{d}\)。
所以得到 \(f(n, i)=n*\sum_{a|i}\frac{\mu(a)}{a}*\sigma(\frac{i}{a})\)。
然后后面那个 \(\sum_{a|i}\frac{\mu(a)}{a}*\sigma(\frac{i}{a})\) 显然是个积性函数,筛一筛就好了,不妨令为 g(i)。
现在答案的表达式变为:
\[ans=\sum_{i=1}^{x-1}(((i+1)^3-1)*g(i) - i^3*g(i) + gcd(i, i^3)) + \sum_{j=1}^{n}gcd(x, j) - g(x)*x^3 + gcd(x, x^3)\]
因为 gcd(i, i^3) = i,所以只需要一个等差数列求和公式套上去就可以解决这些项了。那些与 n 无关的项可以直接预处理。
现在只剩下 \(\sum_{j=1}^{n}gcd(x, j)\) 这一项需要求解。因为 \(x | n\) 不一定成立,所以上面的推导不一定成立。
考虑回到最初的反演式:\(\sum_{j=1}^{n}gcd(x, j) = \sum_{p|x}\sum_{d|p}\lfloor\frac{n}{p}\rfloor*\mu(\frac{p}{d})*d\)
假如枚举出 p,那么剩下的 \(\sum_{d|p}\mu(\frac{p}{d})*d\) 就又是个积性函数了,筛就完事儿了。
最后只需要枚举因数,因此复杂度是 \(O(10^7 + \sqrt{10^7}*T)\)。好像比正解要优秀些?
@accepted [email protected]
#include<cstdio>
const int MAXN = int(1E7);
const int MOD = 998244353;
template <class T>
void read(T &x) {
static char ch;static bool neg;
for(ch=neg=0;ch<'0' || '9'<ch;neg|=ch=='-',ch=getchar());
for(x=0;'0'<=ch && ch<='9';(x*=10)+=ch-'0',ch=getchar());
x=neg?-x:x;
}
int pow_mod(int b, int p) {
int ret = 1;
while( p ) {
if( p & 1 ) ret = 1LL*ret*b%MOD;
b = 1LL*b*b%MOD;
p >>= 1;
}
return ret;
}
int club(__int128 x) {
int l = 1, r = int(1E7);
while( l < r ) {
int mid = (l + r + 1) >>1;
__int128 y = mid;
if( y*y*y <= x ) l = mid;
else r = mid - 1;
}
return l;
}
int inv[MAXN + 5], f[MAXN + 5], g[MAXN + 5], h[MAXN + 5];
int prm[MAXN + 5], low[MAXN + 5], pcnt = 0;
int gcd(int a, int b) {
return (b == 0) ? a : gcd(b, a%b);
}
void init() {
inv[1] = 1;
for(int i=2;i<=MAXN;i++)
inv[i] = MOD - 1LL*inv[MOD%i]*(MOD/i)%MOD;
low[1] = f[1] = g[1] = 1;
for(int i=2;i<=MAXN;i++) {
if( !low[i] ) {
low[i] = i, prm[++pcnt] = i;
f[i] = (2 + MOD - inv[i])%MOD;
g[i] = (i - 1)%MOD;
long long p = 1LL*i*i;
for(int j=2;p<=MAXN;p*=i,j++)
low[p] = p, f[p] = (j + 1 + (MOD - 1LL*inv[i]*j%MOD))%MOD, g[p] = p/i*(i-1);
}
for(int j=1;1LL*i*prm[j]<=MAXN;j++) {
if( i % prm[j] == 0 ) {
if( i != low[i] ) {
low[i*prm[j]] = low[i]*prm[j];
f[i*prm[j]] = 1LL*f[i/low[i]]*f[prm[j]*low[i]]%MOD;
g[i*prm[j]] = 1LL*g[i/low[i]]*g[prm[j]*low[i]]%MOD;
}
break;
}
else {
low[i*prm[j]] = prm[j];
f[i*prm[j]] = 1LL*f[i]*f[prm[j]]%MOD;
g[i*prm[j]] = 1LL*g[i]*g[prm[j]]%MOD;
}
}
}
for(int i=1;i<=MAXN;i++) {
int tmp = (1LL*(i+1)*(i+1)%MOD*(i+1)%MOD + MOD - 1)%MOD;
h[i] = 1LL*f[i]*tmp%MOD;
h[i] = (h[i-1] + h[i])%MOD;
tmp = 1LL*i*i%MOD*i%MOD;
f[i] = 1LL*f[i]*tmp%MOD;
f[i] = (f[i-1] + f[i])%MOD;
}
}
int main() {
init();
int T; scanf("%d", &T);
while( T-- ) {
__int128 n; read(n);
int p = club(n), ans = 1LL*p*(p+1)%MOD*inv[2]%MOD;
for(int i=1;i*i<=p;i++) {
if( p % i == 0 ) {
ans = (ans + 1LL*g[i]*((n/i)%MOD)%MOD)%MOD;
if( i*i != p )
ans = (ans + 1LL*g[p/i]*((n/(p/i))%MOD)%MOD)%MOD;
}
}
ans = (ans + MOD - f[p])%MOD;
printf("%d\n", (ans + h[p-1])%MOD);
}
}
@[email protected]
因为怕 __int128 没办法用 cmath 里的 pow 函数,所以自己手写了个二分(居然没有问题)。
考场上因为某个地方的 n/i 写成了 n/(p/i),没有调出来然后 GG 了。考后 debug 一会儿就出来了 QAQ。。。
原文地址:https://www.cnblogs.com/Tiw-Air-OAO/p/11304899.html