[SDOI2017]数字表格

viii.[SDOI2017]数字表格

题意：求出

$\prod\limits_{i=1}^n\prod\limits_{j=1}^mf_{\gcd(i,j)}$，其中$f$是斐波那契数列。

就算是积，我们也一样能反演，只是反演到了指数头上。

$\begin{aligned}\prod_{i=1}^n\prod_{j=1}^mf_{\gcd(i,j)}& =\prod_{d=1}^{\min(n,m)}\prod_{i=1}^n\prod_{j=1}^mf_d[\gcd(i,j)=d]\\& =\prod_{d=1}^{\min(n,m)}f_d^{\sum_{i=1}^n\sum_{j=1}^m[\gcd(i,j)=d]}\end{aligned}$

观察它的指数$\sum_{i=1}^n\sum_{j=1}^m[\gcd(i,j)=d]$，正是我们这几题来最常见的反演形式。它等于$\sum_{i=1}^{n/d}\sum_{j=1}^{m/d}\sum_{x|i,x|j}\mu(x)=\sum_{x=1}^{\min(n/d,m/d)}\mu(x)\left\lfloor\dfrac{n}{dx}\right\rfloor\left\lfloor\dfrac{m}{dx}\right\rfloor$。

也就是说，

$\begin{aligned}&\prod_{i=1}^n\prod_{j=1}^mf_{\gcd(i,j)}\\=&\prod_{d=1}^{\min(n,m)}f_d^{\sum_{x=1}^{\min(n/d,m/d)}\mu(x)\left\lfloor\frac{n}{dx}\right\rfloor\left\lfloor\frac{m}{dx}\right\rfloor}\\=&\prod_{T=1}^{\min(n,m)}\prod_{d|T}(f_d)^{\mu(T/d)\left\lfloor n/T\right\rfloor\left\lfloor m/T\right\rfloor}\\=&\prod_{T=1}^{\min(n,m)}(\prod_{d|T}(f_d)^{\mu(T/d)})^{\left\lfloor n/T\right\rfloor\left\lfloor m/T\right\rfloor}\end{aligned}$

设$\prod_{d|T}(f_d)^{\mu(T/d)}=F(T)$，

则有$\prod\limits_{i=1}^n\prod\limits_{j=1}^mf_{\gcd(i,j)}=\prod\limits_{T=1}^{\min(n,m)}F(T)^{\left\lfloor n/T\right\rfloor\left\lfloor m/T\right\rfloor}$。

$F(T)$可以在$O(\sum\limits_{i=1}^n\dfrac{n}{i})\approx O(n\ln n)$时间内暴力算出。我们再对它求一个前缀积。这样就可以直接整除分块，在$O(\sqrt{n}\log n)$的时间内回答单次询问（$\log n$是快速幂复杂度）。

代码：

#include<bits/stdc++.h>
using namespace std;
const int N=1e6;
const int mod=1e9+7;
int n,m,T,mu[N+5],pri[N+5],fib[N+5],invfib[N+5],F[N+5],invF[N+5],res;
int ksm(int x,int y){
	int res=1;
	for(;y;x=(1ll*x*x)%mod,y>>=1)if(y&1)res=(1ll*res*x)%mod;
	return res;
}
void init(){
	mu[1]=1;
	for(int i=2;i<=N;i++){
		if(!pri[i])pri[++pri[0]]=i,mu[i]=-1;
		for(int j=1;j<=pri[0]&&i*pri[j]<=N;j++){
			pri[i*pri[j]]=true;
			if(!(i%pri[j]))break;
			mu[i*pri[j]]=-mu[i];
		}
	}
	fib[1]=invfib[1]=1;
	for(int i=2;i<=N;i++)fib[i]=(fib[i-1]+fib[i-2])%mod,invfib[i]=ksm(fib[i],mod-2);
	for(int i=0;i<=N;i++)F[i]=invF[i]=1;
	for(int i=1;i<=N;i++){
		if(!mu[i])continue;
		for(int j=i;j<=N;j+=i)F[j]=(1ll*F[j]*(mu[i]==1?fib[j/i]:invfib[j/i]))%mod;
	}
	for(int i=1;i<=N;i++)F[i]=(1ll*F[i]*F[i-1])%mod,invF[i]=ksm(F[i],mod-2);
}
int main(){
	scanf("%d",&T),init();
	while(T--){
		scanf("%d%d",&n,&m),res=1;
		for(int l=1,r;l<=min(n,m);l=r+1)r=min(n/(n/l),m/(m/l)),res=(1ll*res*ksm(1ll*F[r]*invF[l-1]%mod,1ll*(n/l)*(m/l)%(mod-1)))%mod;
		printf("%d\n",res);
	}
	return 0;
}