bzoj4407:于神之怒加强版

传送门

这个题真的也是有点难度啊（应该是因为我太菜了）

\[ans=\sum_{i=1}^{n}\sum_{j=1}^{m}gcd(i,j)^k\\ \]

可以设

\[f(d)=\sum_{i=1}^{n}\sum_{j=1}^{m}[gcd(i,j)==d]d^k\\ ans=\sum_{d=1}^{min(n,m)}f(d)=\sum_{d=1}^{min(n,m)}\sum_{i=1}^{n}\sum_{j=1}^{m}[gcd(i,j)==d]d^k\\ ans=\sum_{d=1}^{min(n,m)}d^k\sum_{i=1}^{n}\sum_{j=1}^{m}[gcd(i,j)==d]\\ \]

设

\[F(d)=\sum_{i=1}^{n}\sum_{j=1}^{m}[gcd(i,j)==d]\\ g(x)=\sum_{x|d}F(d)=\sum_{i=1}^{n}\sum_{j=1}^{m}[x|gcd(i,j)]=\lfloor \frac{n}{x}\rfloor\lfloor \frac{m}{x}\rfloor\\ F(d)=\sum_{d|x}g(x)\mu(\frac{x}{d})=\sum_{d|x}\lfloor \frac{n}{x}\rfloor\lfloor \frac{m}{x}\rfloor\mu(\frac{x}{d})\\ \]

设\(T=\frac{x}{d}\)

\[F(d)=\sum_{T=1}^{min(n,m)/d}\lfloor \frac{n}{Td}\rfloor\lfloor \frac{m}{Td}\rfloor\mu(T)\\ \]

那么答案就是

\[ans=\sum_{d=1}^{min(n,m)}d^{k}F(d)\\ ans=\sum_{d=1}^{min(n,m)}d^{k}\sum_{T=1}^{min(n,m)/d}\lfloor \frac{n}{Td}\rfloor\lfloor \frac{m}{Td}\rfloor\mu(T)\\ \]

然后这个式子已经可以预处理\(d^k\)，然后数论分块来做了

不过会\(TLE\)，此时时间复杂度为\(O(nlogn+Tn)\)

所以我们还需要优化

考虑设\(k=Td\)

\[ans=\sum_{k=1}^{min(n,m)}\lfloor \frac{n}{k}\rfloor\lfloor \frac{m}{k}\rfloor\sum_{d|k}\mu(\frac{k}{d})d^k \]

考虑设\(s(d)=\sum_{d|k}\mu(\frac{k}{d})d^k\)

可以发现\(s\)是积性函数，因为\(\mu\)和\(d^k\)都是积性函数，积性函数和积性函数的狄利克雷卷积也是积性函数

所以\(s\)可以通过线筛求出

考虑如何线筛，对于\(a\)与\(b\)不互质的情况下，直接用\(s(ab)=s(a)*s(b)\)

当\(x\)为质数时，除了\(\mu(1)\)以及\(\mu(x)\)不为\(0\)，其余都为\(0\)，那么\(s(x)=\mu(1)*x^k+\mu(x)*1^k\)

考虑\(a\)有\(x\)个质因子\(prime_i\)，那么对于\(prime_i\)的贡献\(s(prime_{i}^{x})=\mu(1)*prime_i^x+\mu(prime_i)*prime_i^{x-1}\)

那么对于\(s(a*prime_i)\)，将相当于整体乘上一个\(prime_i^k\)

然后就可以筛了，我讲的有点不清楚，可以看代码理解一下

#include<cstdio>
#include<iostream>
#include<algorithm>
using namespace std;
void read(int &x) {
	char ch; bool ok;
	for(ok=0,ch=getchar(); !isdigit(ch); ch=getchar()) if(ch=='-') ok=1;
	for(x=0; isdigit(ch); x=x*10+ch-'0',ch=getchar()); if(ok) x=-x;
}
#define rg register
const int maxn=5e6+10,mod=1e9+7;bool vis[maxn];
int T,k,n,m,pri[maxn],tot,ans,f[maxn],g[maxn];
int mi(int a,int b)
{
	int ans=1;
	while(b)
	{
		if(b&1)ans=1ll*ans*a%mod;
		b>>=1,a=1ll*a*a%mod;
	}
	return ans;
}
void prepare()
{
	f[1]=1;
	for(rg int i=2;i<=5e6;i++)
	{
		if(!vis[i])pri[++tot]=i,g[i]=mi(i,k),f[i]=(g[i]-1+mod)%mod;
		for(rg int j=1;pri[j]*i<=5e6&&j<=tot;j++)
		{
			vis[pri[j]*i]=1;
			if(!(i%pri[j])){f[i*pri[j]]=1ll*f[i]*g[pri[j]]%mod;break;}
			else f[i*pri[j]]=1ll*f[i]*f[pri[j]]%mod;
		}
	}
	for(rg int i=1;i<=5e6;i++)f[i]=(f[i]+f[i-1])%mod;
}
int main()
{
	read(T),read(k),prepare();
	while(T--)
	{
		read(n),read(m);if(n>m)swap(n,m);ans=0;
		for(rg int i=1,j;i<=n;i=j+1)
		{
			j=min(n/(n/i),m/(m/i));
			ans=(ans+1ll*(f[j]-f[i-1]+mod)%mod*(n/i)%mod*(m/i)%mod)%mod;
		}
		printf("%d\n",ans);
	}
}