P4449 于神之怒加强版 题解

博客园同步

原题链接

简要题意：

给定一个固定的 $k$，$T$ 组询问求：

$$\sum _{i=1}^n\sum _{j=1}^m\gcd (i,j{)}^k$$

其中 $n,m,k\le 5\times 10^6$，$T\le 2\times 10^3$.

没什么好说的，没有部分分。根据 莫比乌斯反演性质 开始推式子！

$$\sum _{i=1}^n\sum _{j=1}^m\gcd (i,j{)}^k$$

$$=\sum _{i=1}^{\min (n,m)}{d}^k\sum _{i=1}^{\lfloor \frac{n}{d}\rfloor }\sum _{j=1}^{\lfloor \frac{n}{d}\rfloor }[\gcd (i,j)==1]$$

$$=\sum _{i=1}^{\min (n,m)}{d}^k\sum _{i=1}^{\lfloor \frac{n}{d}\rfloor }\sum _{j=1}^{\lfloor \frac{n}{d}\rfloor }\sum _{k|\gcd (i,j)}{\mu }_k$$

$$=\sum _{i=1}^{\min (n,m)}{d}^k\sum _{k=1}^{\lfloor \frac{\min (n,m)}{d}\rfloor }\lfloor \frac{n}{kd}\rfloor \lfloor \frac{m}{kd}\rfloor$$

$$=\sum _{i=1}^T\lfloor \frac{n}{T}\rfloor \lfloor \frac{m}{T}\rfloor \sum _{d|T}{d}^k\times {\mu }_{\frac{T}{d}}$$

显然 $O(T\sqrt{n})$ 的时间足以通过，前面的部分整除分块足够。但是对于后面这一块较难解决，则设：

$$g_x=\sum _{d|x}{d}^k\times {\mu }_{\frac{x}{d}}$$

通过观察可知这是个 积性函数。那么显然，对于：

$$x=\prod _{i=1}^q{p}_i^{a_i}$$

存在：

$$g_x=\prod _{i=1}^{q}g(p_i^{a_i})$$

$$=\prod _{i=1}^q(p_i^{k\times (a_i-1)}\times {\mu }_{p_i}+p_i^{k\times a_i}\times {\mu }_1)$$

$$=\prod _{i=1}^q{p}_i^{k\times (a_i-1)}\times (p_i^k-1)$$

首先，上面对于 $g$ 的几步需要解释。那个关于 $x$ 的式子实际上是 $x$ 质因数分解以后的结果，$p_i\in \text{prime}$.

然后第一步是线性筛积性函数的套路。

第二步是一个变形，考虑展开 $g(p_i^{a_i})$ 的过程.

第三步就是个乘法分配律，算出 $\mu$ 而已。

好了，看看这个式子，每个部分都可以筛。结束。

时间复杂度：$O(n+T\sqrt{n}+p\log k)$. $p$ 的值之后解释。

实际得分：$100pts$.

#pragma GCC optimize(2)
#include<bits/stdc++.h>
using namespace std;

typedef long long ll;
const int N=5e6+1;
const ll MOD=1e9+7;

inline int read(){char ch=getchar(); int f=1; while(ch<'0' || ch>'9') {if(ch=='-') f=-f; ch=getchar();}
	int x=0; while(ch>='0' && ch<='9') x=(x<<3)+(x<<1)+ch-'0',ch=getchar(); return x*f;}

inline void write(ll x) {
	if(x<0) {putchar('-');write(-x);return;}
	if(x<10) {putchar(char(x%10+'0'));return;}
	write(x/10);putchar(char(x%10+'0'));
}

int prime[N],mu[N];
ll g[N],f[N]; bool h[N];
int T,k,cnt=0; ll ans=0;

inline ll pw(ll x,ll y) {
	ll ans=1; while(y) {
		if(y&1) ans=ans*x%MOD;
		x=x*x%MOD; y>>=1;
	} return ans;
}

inline void Euler(int n) {
	f[1]=1; for(int i=2;i<=n;i++) {
		if(!h[i]) prime[++cnt]=i,g[cnt]=pw(i,k),f[i]=(g[cnt]-1+MOD)%MOD;
			//f 表示 分解质因数后后面的一堆东西
		for(int j=1;j<=cnt && i*prime[j]<=n;j++) {
			h[i*prime[j]]=1;
			if(i%prime[j]==0) {f[i*prime[j]]=1ll*f[i]*g[j]%MOD;break;}
			f[i*prime[j]]=1ll*f[i]*f[prime[j]]%MOD; //线性筛模板
		}
	}
	for(int i=1;i<=n;i++) f[i]=(f[i]+f[i-1])%MOD; //前缀和
}

int main() {
	T=read(); k=read(); Euler(N-1);
	while(T--) {
		int n=read(),m=read();
		ll ans=0; int r;
		for(int i=1;i<=min(n,m);i=r+1) { //整除分块
			r=min(n/(n/i),m/(m/i));
			ans=(ans+1ll*(f[r]-f[i-1]+MOD)*(n/i)%MOD*(m/i)%MOD)%MOD;
		} write(ans); putchar('\n');
	}
	return 0;
}

显然，你发现对于素数我们需要暴力计算 $p^k$. 那么，首先 $5\times 10^6$ 以内的素数有多少个呢？代码测试：

Link 代码

得到结果约 $3.4\times 10^5$，也即 $p\le 3.4\times 10^5$. 而 $\log k=\log (5\times 10^6)$，大概是 $22$，不会出问题，显然 $3.4\times 10^5\times 22$ 是 $10^7$ 级的，这也是程序较慢的原因（不过不影响通过，照样 $2.55s$）