BZOJ 4816[SDOI2017]数字表格(莫比乌斯反演)

题目链接

(Description)
(f_i)表示(fibonacci)数列第(i)项,求(prod_{i=1}^{n}prod_{j=1}^{m}f[gcd(i,j)])

(T<=10^3,n,m≤10^6)


(Solution)
再来推式子(默认(n<m))

[prod_{i=1}^{n}prod_{j=1}^mf[gcd(i,j)] ]

按照套路枚举(gcd)

[prod_{d=1}^nprod_{i=1}^{n/d}prod_{j=1}^{m/d}f[d][gcd(i,j)==1] ]

考虑每个(f[d])被乘了多少次

[prod_{d=1}^nf[d]^{sum_{i=1}^{n/d}sum_{j=1}^{m/d}[gcd(i,j)==1]} ]

指数非常熟悉

[sum_{i=1}^{n/d}sum_{j=1}^{m/d}[gcd(i,j)==1]=sum_{i=1}^{n/d}mu(i)lfloorfrac{n}{id} floorleftlfloorfrac{m}{id} ight floor ]

现在指数可以数论分块,对(f)做一个前缀积之后外层也可以数论分块,这就可以做到每次询问(O(nlogn))了。
但是这个式子还可以继续推,我们可以枚举(id)
(T=id)

[prod_{T=1}^nprod_{d|T}f[d]^{mu(T/d)lfloorfrac{n}{T} floorlfloorfrac{m}{T} floor} ]

根据数学运算法则就可以这样算

[prod_{T=1}^n(prod_{d|T}f[d]^{mu(T/d)})^{lfloorfrac{n}{T} floorlfloorfrac{m}{T} floor} ]

括号里面的式子可以预处理,乘方可以数论分块,这样就可以单词询问(O(sqrt n logn))了。

#include<complex>
#include<cstdio>
using namespace std;
const int mod=1e9+7;
const int N=1e6+7;
int Q,n,m,tot;
int prime[N],mu[N],f[N],g[N],pro[N];
bool check[N];
int qread()
{
	int x=0;
	char ch=getchar();
	while(ch<'0' || ch>'9')ch=getchar();
	while(ch>='0' && ch<='9'){x=x*10+ch-'0';ch=getchar();}
	return x;
}
int Fpow(long long b,int p)
{
	long long res=1;
	for(;p;p>>=1,b=b*b%mod)
		if(p&1)res=res*b%mod;
	return res;
}
void Init()
{
	check[1]=1;
	mu[1]=f[1]=pro[0]=pro[1]=g[1]=1;
	for(int i=2;i<N;i++)
	{
		if(!check[i])prime[++tot]=i,mu[i]=-1;
		for(int j=1;j<=tot && i*prime[j]<N;j++)
		{
			check[i*prime[j]]=1;
			if(i%prime[j])mu[i*prime[j]]=-mu[i];
			else break;
		}
		f[i]=(f[i-1]+f[i-2])%mod;pro[i]=1;
		g[i]=Fpow(f[i],mod-2);
	}
	for(int i=1;i<N;i++)
	{
		if(!mu[i])continue;
		for(int j=i;j<N;j+=i)
			pro[j]=1ll*pro[j]*(mu[i]==1?f[j/i]:g[j/i])%mod;
	}
	for(int i=2;i<N;i++)
		pro[i]=1ll*pro[i]*pro[i-1]%mod;
}
int main()
{
	Init();
	scanf("%d",&Q);
	while(Q--)
	{
		n=qread();m=qread();
		if(n>m)swap(n,m);
		int ans=1;
		for(int l=1,r;l<=n;l=r+1)
		{
			r=min(n/(n/l),m/(m/l));
			ans=1ll*ans*Fpow(1ll*pro[r]*Fpow(pro[l-1],mod-2)%mod,1ll*(n/l)*(m/l)%(mod-1))%mod;
		}
		printf("%d
",ans);
	}
	return 0;
}
原文地址:https://www.cnblogs.com/LeTri/p/10334467.html