【题意】
求 (displaystyle sum_{i=1}^n gcd(lfloorsqrt[3] i floor, i)mod 998244353, nleq 10^{21})
【分析】
正面求解很难,考虑枚举 (lfloorsqrt[3] i floor)
当 (lfloorsqrt[3] i floor)=t$ 时 (iin [t^3, (t+1)^3-1])
故考虑枚举开三方的数字:
记 (m=lfloorsqrt[3] n floor)
(displaystyle sum_{i=1}^n gcd(lfloorsqrt[3] i floor, i)=sum_{t=1}^{m-1} sum_{i=t^3}^{(t+1)^3-1} gcd(t, i)+sum_{i=m^3}^n gcd(m, i))
考虑到 (gcd(t, i+t)=gcd(t, i+t-t)=gcd(i, t))
因此每 (t) 个连续的 (gcd(t, i)) 的求和都等于 (displaystyle sum_{i=1}^t gcd(i, t))
记 (displaystyle G_t(n)=sum_{i=1}^n gcd(i, t)) 则 (G_t(n)=(n/t)cdot G_t(t)+G_t(nmod t))
因此 (displaystyle sum_{i=t^3}^{(t+1)^3-1} gcd(t, i)=G_t(t^3+3t^2+3t)-G_t(t^3)+gcd(t^3, t)=(t^2+3t+3-t^2)G_t(t)+t=(3t+3)G_t(t)+t)
而同理, (displaystyle sum_{i=m^3}^n gcd(m, i)=(n/m-m^2)G_m(m)+G_m(nmod m)+m)
(egin{aligned} herefore& sum_{i=1}^n gcd(lfloorsqrt[3] i floor, i) \\=&sum_{t=1}^{m-1} sum_{i=t^3}^{(t+1)^3-1} gcd(t, i)+sum_{i=m^3}^n gcd(m, i) \\=&sum_{t=1}^{m-1} [(3t+3)G_t(t)+t]+(n/m-m^2)G_m(m)+m-G_m(nmod m) end{aligned})
预处理 (G_t(t)) 极其前缀和,每次查询的时候就获取 (m=lfloorsqrt[3] n floor) ,然后 (O(1)) 算出除 (G_m(nmod m)) 的部分
当然,由于 (nleq 10^{21}) 可能涉及 __int128
,可能无法使用 pow()
函数。可以考虑二分 (m)
考虑剩下的部分 (displaystyle G_m(nmod m)=sum_{i=1}^{nmod m}gcd(i, m))
方便期间,记 (displaystyle g(n, m)=sum_{i=1}^{nmod m}gcd(i, m))
由莫比乌斯反演得:
(egin{aligned} g(n, m)&=sum_{i=1}^ngcd(i, m) \\&=sum_{i=1}^n sum_{dmid iwedge dmid m}oldsymbol varphi(d)& ext{(欧拉反演)} \\&=sum_{dmid m}oldsymbol varphi(d)sum_{i=1}^n [dmid i] \\&=sum_{dmid m}oldsymbol varphi(d)(n/d) end{aligned})
答案为 (G_m(nmod m)=g(nmod m, m)) ,直接预处理欧拉函数,然后 (O(sqrt m)) 枚举因数即可
最后剩下 (G_t(t)) 没解决,等价于 (displaystyle g(t, t)=sum_{dmid t} oldsymbol varphi(d) ({tover d})=(oldsymbol varphi*oldsymbol {id})(t))
这个为两个积性函数的狄利克雷卷积,也为积性函数,为方便,记 (oldsymbol f=oldsymbol varphi*oldsymbol {id})
考虑如何线筛得到:
(egin{aligned} oldsymbol f(p^k)&=sum_{i=1}^{k-1} oldsymbol varphi(p^i)oldsymbol {id}(p^{k-i})+oldsymbol varphi(p^k)+oldsymbol {id}(p^k)&(k>1) \\&=sum_{i=1}^{k-1} (p^{i-1}(p-1)cdot p^{k-i})+p^{k-1}(p-1)+p^k \\&=p^{k-1}(p-1)cdot k+p^k end{aligned})
则 (oldsymbol f(p^{k-1})=p^{k-2}(p-1)cdot (k-1)+p^{k-1})
因而 (pcdot oldsymbol f(p^{k-1})=p^{k-1}(p-1)cdot (k-1)+p^k=oldsymbol f(p^k)-p^{k-1}(p-1)=oldsymbol f(p^k)-oldsymbol varphi(p^k))
因此 (oldsymbol f(p^k)=pcdot oldsymbol f(p^{k-1})+oldsymbol varphi(p^k)=pcdot left[oldsymbol f(p^{k-1})+oldsymbol varphi(p^{k-1}) ight])
记 (fc_n) 表示 (n) 的最小质因数, (fck_n) 表示 (n) 的最小质因子的指数,例如 (fck_{12}=2^2)
故当 (p eq fc_n) 时 (oldsymbol f(ncdot p)=oldsymbol f(n)cdot oldsymbol f(p)=oldsymbol f(n)cdot (2p-1))
当 (p=fc_n) 时
(egin{aligned} &oldsymbol f(ncdot p) \\=&oldsymbol f({nover fck_n})cdot oldsymbol f(fck_ncdot p) \\=&oldsymbol f({nover fck_n}) cdot pcdot left[oldsymbol f(fck_n)+oldsymbol varphi(fck_n) ight] \\=&pcdot [oldsymbol f(n)+oldsymbol f({nover fck_n})cdot oldsymbol varphi(fck_n)] end{aligned})
考虑线筛的时候,同时维护最小质因数 fc[i]
,最小质因子的指数 fck[i]
,欧拉函数 phi[i]
,以及该函数值 f[i]
即可直接转移
总复杂度 (O(n+T(log n+sqrt n) )=O(n+Tsqrt n))
【代码】
只跑了 1.8s
#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
const int Lim=1e7, MAXN=Lim+10, P=998244353;
#define de(x) cout << #x << " = " << x << endl
#define dd(x) cout << #x << " = " << x << " "
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 fc[MAXN], prime[MAXN/10], cntprime, fck[MAXN], phi[MAXN], f[MAXN], sumit[MAXN];
inline void init() {
phi[1]=f[1]=1;
for(int i=2; i<=Lim; ++i) {
if(!fc[i]) {
fc[i]=prime[++cntprime]=fck[i]=i;
phi[i]=i-1;
f[i]=i+i-1;
}
for(int j=1; j<=cntprime; ++j)
if(prime[j]*i>Lim) break;
else if(prime[j]<fc[i]) {
fc[prime[j]*i]=prime[j];
fck[prime[j]*i]=prime[j];
phi[prime[j]*i]=phi[i]*(prime[j]-1);
f[prime[j]*i]=f[i]*(2ll*prime[j]-1)%P;
}
else{
fc[prime[j]*i]=prime[j];
fck[prime[j]*i]=prime[j]*fck[i];
phi[prime[j]*i]=phi[i]*prime[j];
f[prime[j]*i]=(f[i]+(ll)f[i/fck[i]]*phi[fck[i]])%P*prime[j]%P;
}
}
for(int i=1; i<=Lim; ++i)
sumit[i]=((3ll*i+3)*f[i]+i+sumit[i-1])%P;
}
inline int g(int n, int m) {
int res=0;
for(int i=1, j; i*i<=m; ++i) if(m%i==0) {
j=m/i;
res=(res+(ll)phi[i]*(n/i))%P;
if(i!=j) res=(res+(ll)phi[j]*(n/j))%P;
}
return res;
}
inline int findit(__int128 m) {
__int128 l=1, r=1e7, mid, ans=0;
while(l<=r) {
mid=l+r>>1;
if(mid*mid*mid<=m) {
ans=mid;
l=mid+1;
}
else r=mid-1;
}
return ans;
}
inline int ans(__int128 n) {
int m=findit(n);
int res=((n/m-(ll)m*m)*f[m]+m)%P+g(n%m, m);
res=(res%P+sumit[m-1])%P;
return res<0?res+P:res;
}
int main() {
init();
int t; read(t);
__int128 n;
while(t--) read(n), printf("%d
", ans(n));
return 0;
}