傳送門php
發現就是要求:
\[ Ans=\prod_{i=1}^{n}\prod_{j=1}^{m}f(gcd(i,j)) \]
其中\(f(n)\)表示斐波那契數列第\(n\)項。ios
根據套路,枚舉全部的\(gcd\)
\[ Ans=\prod_{d=1}^{min(n,m)}\prod_{i=1}^{n}\prod_{j=1}^{m}[gcd(i,j)=d]f(d)+[gcd(i,j)!=d] \]
而後直接把\(f(d)\)提出來
\[ Ans=\prod_{d=1}^{n}f[d]^{\sum_{i=1}^{n/d}\sum_{j=1}^{m/d}[gcd(i,j)=1]} \]
而後着手開始化上面的式子:spa
\[ \begin{align} &\ \ \ \ \sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{n}{d}\rfloor}[gcd(i,j)=1]\\ &=\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{x\mid gcd(i,j)}\mu(x)\\ &=\sum_{x=1}^{min(\lfloor\frac{n}{d}\rfloor,\lfloor\frac{m}{d}\rfloor)}\mu(x)\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{n}{d}\rfloor}[x\mid i\&\&x\mid j]\\ &=\sum_{x=1}^{min(\lfloor\frac{n}{d}\rfloor,\lfloor\frac{m}{d}\rfloor)}\mu(x)\sum_{i=1}^{\lfloor\frac{n}{dx}\rfloor}\sum_{j=1}^{\lfloor\frac{n}{dx}\rfloor}1\\ &=\sum_{x=1}^{min(\lfloor\frac{n}{d}\rfloor,\lfloor\frac{m}{d}\rfloor)}\mu(x)\lfloor\frac{n}{dx}\rfloor\lfloor\frac{m}{dx}\rfloor \end{align} \]code
而後令\(T=dx\)而後直接從整個式子中提出\(T\)
\[ Ans=\prod_{T=1}^{n}\prod_{d\mid{T}}f(d)^{[n/T][m/T]\mu(T/d)} \]
\[ \ \ \ \ \ \ \ \ \ \ =\prod_{T=1}^{n}(\prod_{d\mid{T}}f(d)^{\mu(T/d)})^{[n/T][m/T]} \]
由於\(n\)足夠小,因此這個裏面式子能夠直接暴力求出來。。。
而後外層整數分塊,求解便可。get
#include<iostream> #include<cstdio> #include<cstring> #include<algorithm> using namespace std; const int MAXN=1e6+7; #define ll long long const int mo=1e9+7; ll sumf[MAXN],fr[MAXN],f[MAXN],fi[MAXN][3],sum[MAXN],ans; int mu[MAXN],prime[MAXN],n,m,mul[MAXN]; bool vis[MAXN]; inline ll pow(ll x,ll k) { ll cnt=1; while(k){ if(k&1) cnt*=x; cnt%=mo;x=x*x%mo;k>>=1; } return cnt; } inline void get_mu(int N) { mu[1]=1;f[1]=1; for(int i=2;i<=N;i++){ f[i]=(f[i-1]+f[i-2])%mo; if(!vis[i]){ mu[i]=-1; prime[++prime[0]]=i; } for(int j=1;j<=prime[0];j++){ if(i*prime[j]>N) break; vis[i*prime[j]]=1; if(i%prime[j]==0) break; mu[i*prime[j]]=-mu[i]; } } for(int i=1;i<=N;i++) fi[i][0]=pow(f[i],mo-2),fi[i][1]=1,fi[i][2]=f[i]; for(int i=1;i<=N;i++) mul[i]=1; for(int i=1;i<=N;i++) for(int j=i;j<=N;j+=i) mul[j]=1ll*mul[j]*fi[i][mu[j/i]+1]%mo; sumf[0]=1; for(int i=1;i<=N;i++) sumf[i]=1ll*sumf[i-1]*mul[i]%mo; fr[N]=pow(sumf[N],mo-2); for(int i=N-1;~i;i--) fr[i]=1ll*fr[i+1]*mul[i+1]%mo; } inline int read() { int x=0,c=1; char ch=' '; while((ch<'0'||ch>'9')&&ch!='-')ch=getchar(); while(ch=='-')c*=-1,ch=getchar(); while(ch>='0'&&ch<='9')x=x*10+ch-'0',ch=getchar(); return x*c; } int main() { int T=read(); get_mu(1000000); while(T--){ n=read(),m=read(); ans=1; int mx=min(n,m); for(int l=1,r;l<=mx;l=r+1){ r=min(n/(n/l),m/(m/l)); ans=(1ll*ans*pow(1ll*sumf[r]*fr[l-1]%mo,1ll*(n/l)*(m/l)%(mo-1)))%mo; } printf("%lld\n",ans); } }