「SDOI2017」数字表格

「SDOI2017」数字表格

链接：https://www.luogu.com.cn/problem/P3704

---

显然答案是：

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

化简考虑枚举 \(f\) 的下标：

\[\begin{aligned} &=\prod_d\prod_{i=1}^{\lfloor\frac{n}{d}\rfloor}\prod_{j=1}^{\lfloor\frac{m}{d}\rfloor}f_d[gcd(i,j)=1] \\ &=\prod_df_d^{\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{d}\rfloor}\sum_{k|gcd(i,j)}\mu(k)} \\ &=\prod_df_d^{\sum_{k=1}^{\min(\lfloor\frac{n}{d}\rfloor\lfloor\frac{m}{d}\rfloor)}\mu(k)\times\lfloor\frac{n}{dk}\rfloor\lfloor\frac{m}{dk}\rfloor} \end{aligned} \]

这个式子我们可以外层数论分块内层数论分块在 \(O(n)\) 的时间内算出，但是显然不够。继续优化式子。

设 \(T=k\times d\)。将 \(T\) 提出，枚举 \(T\)。

那么原式为：

\[\prod_{T=1}(\prod_{d|T}f_d^{\mu(\frac{T}{d})})^{\lfloor\frac{n}{T}\rfloor\lfloor\frac{m}{T}\rfloor} \]

而对于 \(\prod_{d|T}f_d^{\mu(\frac{T}{d})}\) 我们可以枚举倍数在 \(O(n\ln n)\) 的时间内求出。

外层再用数论分块优化，算上求逆元，单次查询是 \(\sqrt{n}\log MOD\) 的。

总复杂度是 \(O(T\sqrt{n}\log MOD+n\ln n\log MOD)\)。

代码如下：

#include<bits/stdc++.h>
#define ll long long
using namespace std;
const int MOD = 1e9+7;
const int MAXN = 1e6+5;
const int MX = 1e6;
int n,m,prime[MAXN],p[MAXN],pc;
ll f[MAXN],Sum[MAXN],mu[MAXN],g[MAXN],inv[MAXN];
ll qpw(ll x,ll b)
{
	ll now=1,tmp=x,ans=1;
	while(now<=b)
	{
		if(now&b) ans=ans*tmp%MOD;
		tmp=tmp*tmp%MOD;
		now<<=1;
	}
	return ans;
}
int main()
{
	int T;scanf("%d",&T);
	f[0]=0;f[1]=1;
	for(int i=2;i<=MX;++i) f[i]=(f[i-1]+f[i-2])%MOD;
	for(int i=0;i<=MX;++i) inv[i]=qpw(f[i],MOD-2);
	mu[1]=1;
	for(int i=2;i<=MX;++i)
	{
		if(!p[i]) prime[++pc]=i,mu[i]=-1;
		for(int j=1;j<=pc&&1ll*prime[j]*i<=MX;++j)
		{
			p[i*prime[j]]=1;
			if(i%prime[j]==0)
			{
				mu[i*prime[j]]=0;
				break;
			}
			mu[i*prime[j]]=-mu[i];
		}
	}
	for(int i=0;i<=MX;++i) g[i]=1;
	for(int i=1;i<=MX;++i)
		for(int j=i;j<=MX;j+=i)
			if(mu[j/i]!=0) g[j]=g[j]*(mu[j/i]==1?f[i]:inv[i])%MOD;
	for(int i=2;i<=MX;++i) g[i]=g[i]*g[i-1]%MOD;
	while(T--)
	{
		scanf("%d %d",&n,&m);
		ll Ans=1;
		for(ll l=1,r;l<=min(n,m);l=r+1)
		{
			r=min(n/(n/l),m/(m/l));
			Ans=Ans*qpw(g[r]*qpw(g[l-1],MOD-2)%MOD,(n/l)*(m/l))%MOD;
		}
		printf("%lld\n",Ans);
	}
	return 0;
}