题目
题目描述
Doris刚刚学习了fibonacci数列。用$f[i]$表示数列的第$i$项，那么
$f[0]=0,f[1]=1,$
$f[n]=f[n-1]+f[n-2],n\geq 2$
Doris用老师的超级计算机生成了一个$n×m$的表格，
第$i$行第$j$列的格子中的数是$f[\gcd(i,j)]$，其中$\gcd(i,j)$表示$i,j$的最大公约数。
Doris的表格中共有$n×m$个数，她想知道这些数的乘积是多少。
答案对$10^9+7$取模。
输入输出格式
输入格式：
有多组测试数据。
第一个一个数$T$，表示数据组数。
接下来$T$行，每行两个数$n,m$
输出格式：
输出$T$行，第$i$行的数是第$i$组数据的结果
输入输出样例
输入样例#1：
3
2 3
4 5
6 7
输出样例#1：
1
6
960
说明
对$10\%$的数据，1\leq n,m\leq 100
对$30\%$的数据，$1\leq n,m\leq 1000$
另外存在$30\%$的数据，$T\leq 3$
对$100\%$的数据，$T\leq 1000,1\leq n,m\leq 10^6$
时间限制：5s
内存限制：128MB
题解
题目要求：$\prod_{i=1}^{n}\prod_{j=1}^{m}f[\gcd(i,j)]$
可以化简成：$\prod_{x=1}^{min(n,m)}f[x]^{\sum_{i=1}^{\left\lfloor\frac{n}{x}\right\rfloor}\sum_{j=1}^{\left\lfloor\frac{m}{x}\right\rfloor}[\gcd(i,j)]==1]}$
发现指数可以莫比乌斯反演一下：$\prod_{x=1}^{min(n,m)}f[x]^{\sum_{i=1}^{min(\left\lfloor\frac{n}{x}\right\rfloor,\left\lfloor\frac{m}{x}\right\rfloor)}\mu(i)\left\lfloor\frac{n}{ix}\right\rfloor\left\lfloor\frac{m}{ix}\right\rfloor}$
把$ix$提出来，即有：$\prod_{T=1}^{min(n,m)}\prod_{x\vert T}f[x]^{\mu(\frac{T}{x})\left\lfloor\frac{n}{T}\right\rfloor\left\lfloor\frac{m}{T}\right\rfloor}$
发现可以先处理出来$\prod_{x\vert T}f[x]^{\mu(\frac{T}{x})}$，搞个前缀积，然后整除分块带进去就行了，注意除的时候要求逆元
代码
# include<iostream>
# include<cstring>
# include<cstdio>
# include<algorithm>
# define LL long long
using namespace std;
const int MAX=1e6+1;
const LL mod=1e9+7;
int t,n,m,maxn,tot;
int u[MAX],pm[MAX];
LL f[2][MAX],A[MAX];
bool use[MAX];
LL pow(LL a,LL b)
{
LL ans=1;
for(;b;b>>=1,a=a*a%mod)
if(b&1) ans=ans*a%mod;
return ans;
}
int main()
{
u[1]=1,f[0][1]=f[1][1]=A[0]=A[1]=1;
for(int i=2;i<MAX;++i)
{
A[i]=1,f[0][i]=(f[0][i-1]+f[0][i-2])%mod,f[1][i]=pow(f[0][i],mod-2);
if(!use[i]) pm[++tot]=i,u[i]=-1;
for(int j=1;i*pm[j]<MAX;++j)
{
use[i*pm[j]]=1;
if(i%pm[j]==0)
{
u[i*pm[j]]=0;
break;
}
u[i*pm[j]]=-u[i];
}
}
for(int i=1;i<MAX;++i)
for(int j=1;j*i<MAX;++j)
if(u[i]==1) A[i*j]=A[i*j]*f[0][j]%mod;
else if(u[i]==-1) A[i*j]=A[i*j]*f[1][j]%mod;
for(int i=1;i<MAX;++i)
A[i]=A[i]*A[i-1]%mod;
scanf("%d",&t);
while(t--)
{
scanf("%d%d",&n,&m);
if(n>m) swap(n,m);
LL ans=1;
for(int l=1,r;l<=n;l=r+1)
r=min(n/(n/l),m/(m/l)),ans=ans*pow((A[r]*pow(A[l-1],mod-2))%mod,1ll*(n/l)*(m/l))%mod;
printf("%lld\n",ans);
}
return 0;
}
