欧拉函数+反演——2019hdu多校6588

[求sum_{i=1}^{n}(sqrt[3]i,i)\ 首先转化一下这个式子，考虑对于iin[j^3,(j+1)^3-1]，sqrt[3]i=j\ 所以可以枚举所有j，然后对iin[j^3,(j+1)^3-1]区间的(i,j)求和即可 那么我们把n分成两部分，分别求和：\ sum_{i=1}^{n}(sqrt[3]i,i)=sum_{i={lfloor sqrt[3]n floor}^3}^{n}(sqrt[3]n,n) +sum_{j=1}^{sqrt[3]n-1}sum_{i=j^3}^{(j+1)^3-1}(i,j)\ ]

[先来看前面一个合式sum_{i={lfloor sqrt[3]n floor}^3}^{n}(sqrt[3]n,i)，直接用欧拉替换式\ 原式=sum_{i={lfloor sqrt[3]n floor}^3}^{n} sum_{d|i,d|{lfloor sqrt[3]n floor}}varphi(d) =sum_{d|{lfloor sqrt[3]n floor}}varphi(d)ig(n/d-(lfloor sqrt[3]n floor-1)/d ig)\ {lfloor sqrt[3]n floor}因子是不多的，只有sqrt[6]n个，所以复杂度也是这个 ]

[再看后面一个合式，我们可以将其转化后用反演来快速求和\ sum_{j=1}^{sqrt[3]n-1}sum_{i=j^3}^{(j+1)^3-1}(i,j) =sum_{j=1}^{sqrt[3]n-1}igg(sum_{i=1}^{(j+1)^3-1}(i,j)-sum_{i=1}^{j^3-1}(i,j) igg)\ 现在来求sum_{j=1}^{sqrt[3]n-1}sum_{i=1}^{f(j)}(i,j),其实可以把(i,j)展开做，但是考虑到欧拉替换\ 原式=sum_{j=1}^{sqrt[3]n-1}sum_{i=1}^{f(j)}sum_{d|i,d|j}varphi(d) =sum_{j=1}^{sqrt[3]n-1}sum_{d|j}varphi(d)frac{f(j)}{d} \ 注意因为i的上界是和j有关的，所以切换枚举次序时只能把d切换到i前，不能切换到j前，否则无法正确求出varphi(d)出现次数\ 原式=sum_{j=1}^{sqrt[3]n-1}sum_{d|j}varphi(d)igg(lfloor frac{(j+1)^3-1}{d} floor-lfloor frac{j^3-1}{d} floor igg) 显然有j=kd，所以不妨外层枚举d\ 原式=sum_{d=1}^{N=sqrt[3]n-1}varphi(d)sum_{k=1}^{frac{N}{d}}(3k^2d+3k+1) =sum_{d=1}^{N=sqrt[3]n-1}varphi(d)(3dsum_{k=1}^{frac{N}{d}}k^2+3sum_{k=1}^{frac{N}{d}}k+frac{N}{d}) \ 到此为止，我们可以预处理[1,1e7]范围的k^2,k的前缀和，varphi(i)的值,就可以O(1)求每个d的值 ]

#include <bits/stdc++.h>
const int N=1e7,XN=N+11,P=998244353;
int Add(int x,int const &y) {return (x+=y)>=P?x-P:x;}
int Minus(int x,int const &y) {return (x-=y)<0?x+P:x;}
int Mul(long long x,int const &y) {return x*y%P;}
int prime[XN],phi[XN],pcnt;
void Prep() {
	static bool notPrime[XN];
	phi[1]=1;
	for(int i=2;i<=N;++i) {
		if(!notPrime[i]) {
			prime[++pcnt]=i;
			phi[i]=i-1;
		}
		for(int j=1;j<=pcnt && i*prime[j]<=N;++j) {
			notPrime[i*prime[j]]=1;
			if(i%prime[j]==0) {
				phi[i*prime[j]]=phi[i]*prime[j];
				break;
			} else
				phi[i*prime[j]]=phi[i]*phi[prime[j]];
		}
	}
}

int Sum1(long long x) {return x*(x+1)/2%P;}
int Sum2(long long x) {return (x*(x+1))%(6ll*P)*(2*x+1)/6%P;}
int Calc(int a,__int128 L,__int128 R) {
	int res=0;
	for(int d=1;1LL*d*d<=a;++d)
		if(a%d==0) {
			res=Add(res,Mul((R/d-L/d)%P,phi[d]));
			if(a!=d*d)
				res=Add(res,Mul((R/(a/d)-L/(a/d))%P,phi[a/d]));
		}
	return res;
}

int main() {
	Prep();int T;fin>>T;
	while(T--) {
		__int128 n;fin>>n;
		if(n<=7) {fout<<n<<'
';} 
        else {
			int r;for(r=1;(__int128)(r+2)*(r+2)*(r+2)-1<=n;++r);
			int Ans=Calc(r+1,(__int128)(r+1)*(r+1)*(r+1)-1,n);
			for(int T=1;T<=r;++T)
				Ans=Add(Ans,Mul(phi[T],Add(Add(Mul(3*T,Sum2(r/T)),Mul(3,Sum1(r/T))),r/T)));
			fout<<Ans<<'
';
		}
	}
}