P6156 简单题/P6222 「P6156 简单题」加强版 题解

P6156 简单题

P6222 「P6156 简单题」加强版

写这份代码出了一堆bug，被 \(\color{\grey} \_\texttt{nobody}\_\) D没了/kk

原式

\[\begin{aligned}&=\sum_{i=1}^{n} \sum_{j=1}^{n} (i+j)^k \mu(\gcd(i,j))^2 \gcd(i,j)\\&=\sum_{d=1}^{n} \mu(d)^2 d\sum_{i=1}^{n} \sum_{j=1}^{n} (i+j)^k[gcd(i,j)==d]\\&=\sum_{d=1}^{n} \mu(d)^2 d\sum_{i=1}^{\frac{n}{d}} \sum_{j=1}^{\frac{n}{d}} (id+jd)^k[gcd(i,j)==1]\\&=\sum_{d=1}^{n} \mu(d)^2 d^{k+1}\sum_{i=1}^{\frac{n}{d}} \sum_{j=1}^{\frac{n}{d}} (i+j)^k [gcd(i,j)==1]\\&=\sum_{d=1}^{n} \mu(d)^2 d^{k+1} \sum_{i=1}^{\frac{n}{d}} \sum_{j=1}^{\frac{n}{d}} (i+j)^k \sum_{x|gcd(i,j)}\mu(x)\\&=\sum_{d=1}^{n} \mu(d)^2 d^{k+1} \sum_{x=1}^{\frac{n}{d}}\mu(x) \sum_{i=1}^{\frac{n}{dx}} \sum_{j=1}^{\frac{n}{dx}} (ix+jx)^k\\&=\sum_{d=1}^{n} \mu(d)^2 d^{k+1} \sum_{x=1}^{\frac{n}{d}}\mu(x)x^k \sum_{i=1}^{\frac{n}{dx}} \sum_{j=1}^{\frac{n}{dx}} (i+j)^k\\&=\sum_{d=1}^{n} \mu(d)^2 d^{k+1} \sum_{T=1,d|T}^{n}\mu(\dfrac{T}{d})(\dfrac{T}{d})^k \sum_{i=1}^{\frac{n}{T}} \sum_{j=1}^{\frac{n}{T}} (i+j)^k\\&=\sum_{T=1}^{n} \sum_{i=1}^{\frac{n}{T}} \sum_{j=1}^{\frac{n}{T}} (i+j)^k \sum_{d=1,d|T}^{n} \mu(\dfrac{T}{d}) (\dfrac{T}{d})^k \mu(d) d^k d\mu(d)\\&=\sum_{T=1}^{n} T^k \sum_{i=1}^{\frac{n}{T}} \sum_{j=1}^{\frac{n}{T}} (i+j)^k \sum_{d=1,d|T}^{n} d\mu^2(d) \mu(\dfrac{T}{d})\\\end{aligned} \]

\(T^k\) 非常显然的完全积性函数可以线性筛。

考虑如何快速求出 \(f(T)=\sum_{d=1,d|T}^{n} d\mu^2(d) \mu(\dfrac{T}{d})\)

\(\color{\grey} \_\texttt{nobody}\_\) 告诉我可以发掘 \(\mu\) 的性质。

那个式子是一堆积性函数狄利克雷卷积起来，还是积性函数。

线性筛当 \(p|T\) 时考虑进行以下分类讨论进行转移

\[f(p^k)=\begin{cases}1\ (k=0)\\1\times \mu^2(1)\times \mu(p)+p\times \mu(p)^2\times\mu(1)=p-1 \ (k=1)\\1\times \mu^2(1)\times \mu(p^2)+p^2\times \mu(p^2)^2\times\mu(1)+p\times \mu^2(p)\times \mu(p)=-p\ (k=2)\\0\ (k>2)\end{cases} \]

即

\[f(p^k)=\begin{cases}1\ (k=0)\\p-1 \ (k=1)\\-p\ (k=2)\\0\ (k>2)\end{cases} \]

于是可以线性筛了。

现在式子只剩

\[\sum_{i=1}^{n} \sum_{j=1}^{n} (i+j)^k \]

考虑记一个前缀和 \(s_n=\sum\limits_{i=1}^{n}i^k\)

那么原式 \(=\sum\limits_{i=n+1}^{2*n}s_i-\sum\limits_{i=1}^{n}s_i\) ，统计一下 \(s_i\) 的前缀和就可以 \(O(1)\) 了。

至此，每一个部分都已经预处理至 \(O(1)\) ，外层套个整除分块就可以单次 \(O(\sqrt{n})\) 回答询问，于是加强版也能过了。

简单版强样例：

Input
5000000 1000000000000000000
Output
305836999

简单版：

#include<bits/stdc++.h>
using namespace std;
typedef long long LL;
typedef double db;
#define pb(x) push_back(x)
#define mkp(x,y) make_pair(x,y)
//#define getchar() (p1==p2&&(p2=(p1=buf)+fread(buf,1,1<<21,stdin),p1==p2)?EOF:*p1++)
//char buf[1<<21],*p1=buf,*p2=buf;
inline int read() {
	int x=0,f=1;char ch=getchar();
	while(!isdigit(ch)) {if(ch=='-')f=-1;ch=getchar();}
	while(isdigit(ch))x=x*10+(ch^48),ch=getchar();
	return x*f;
}
const int N=10000005;
const int mod = 998244353;
int n;LL k;
int pw[N],S[N],f[N], pri[N], cnt;
bool vis[N];
int qpow(int n,int k,int res=1){
	for(;k;k>>=1,n=1ll*n*n%mod)if(k&1)res=1ll*res*n%mod;
	return res;
}
void init(const int&N){
	pw[1]=1,f[1]=1;
	for(int i=2;i<=N;++i){
		if (!vis[i]) pri[++cnt] = i, pw[i] = qpow(i, k), f[i] = i - 1;
		for (int j = 1; j <= cnt && i * pri[j] <= N; ++ j) {
			vis[i * pri[j]] = 1, pw[i * pri[j]] = 1ll * pw[i] * pw[pri[j]] % mod;
			if (i % pri[j] == 0) {
				if(i / pri[j] % pri[j]) f[i * pri[j]] = -1ll * f[i / pri[j]] * pri[j] % mod;
				break;
			}
			pw[i * pri[j]] =1ll * pw[i] * pw[pri[j]] % mod, f [i * pri[j]] = 1ll * f[i] * f[pri[j]] % mod;
		}
	}
	for(int i=1;i<=N;++i)f[i]=(1ll*pw[i]*f[i]%mod+f[i-1])%mod,pw[i]=(pw[i]+pw[i-1])%mod;
	for(int i=1;i<=N;++i)S[i]=(S[i-1]+pw[i])%mod;
}
int query(int n){
	int res=0;
	for(int l=1,r;l<=n;l=r+1)
		r=n/(n/l),res=(res+1ll*(S[(n/l)*2]-S[n/l]*2)*(f[r]-f[l-1])%mod)%mod;
	return (res+mod)%mod;
}
signed main(){
	n=read(), scanf("%lld", &k), k %= mod - 1,init(n<<1);
	printf("%d\n",query(n));
}

加强版：

#include<bits/stdc++.h>
using namespace std;
typedef long long LL;
typedef double db;
typedef unsigned int unt;
#define pb(x) push_back(x)
#define mkp(x,y) make_pair(x,y)
//#define getchar() (p1==p2&&(p2=(p1=buf)+fread(buf,1,1<<21,stdin),p1==p2)?EOF:*p1++)
//char buf[1<<21],*p1=buf,*p2=buf;
inline unt read() {
	unt x=0;char ch=getchar();
	while(!isdigit(ch))ch=getchar();
	while(isdigit(ch))x=x*10+(ch^48),ch=getchar();
	return x;
}
const int N=20000005;
unt T,n,k;
unt pw[N],f[N], pri[1270610], cnt;
bool vis[N];
unt qpow(unt n,unt k,unt res=1){
	for(;k;k>>=1,n=n*n)if(k&1)res=res*n;
	return res;
}
void init(const unt&N){
	pw[1]=1,f[1]=1;
	for(unt i=2;i<=N;++i){
		if (!vis[i]) pri[++cnt] = i, pw[i] = qpow(i, k), f[i] = i - 1;
		for (unt j = 1; j <= cnt && i * pri[j] <= N; ++ j) {
			vis[i * pri[j]] = 1, pw[i * pri[j]] = pw[i] * pw[pri[j]];
			if (i % pri[j] == 0) {
				if(i / pri[j] % pri[j]) f[i * pri[j]] = - f[i / pri[j]] * pri[j];
				break;
			}
			pw[i * pri[j]] = pw[i] * pw[pri[j]], f [i * pri[j]] = f[i] * f[pri[j]];
		}
	}
	for(unt i=1;i<=N;++i)f[i]=pw[i]*f[i]+f[i-1],pw[i]=pw[i]+pw[i-1];
	for(unt i=1;i<=N;++i)pw[i]=pw[i-1]+pw[i];
}
unt query(unt n){
	unt res=0;
	for(unt l=1,r;l<=n;l=r+1)
		r=n/(n/l),res=(res+(pw[(n/l)<<1]-(pw[n/l]<<1))*(f[r]-f[l-1]));
	return res;
}
signed main(){
	T = read(), n=read(), k = read(),init(n<<1);
	while(T--)printf("%u\n",query(read()));
}

混合码风.jpg