Description
Doris刚刚学习了fibonacci数列。用f[i]表示数列的第i项，那么
f[0]=0
f[1]=1
f[n]=f[n-1]+f[n-2],n>=2
Doris用老师的超级计算机生成了一个n×m的表格，第i行第j列的格子中的数是f[gcd(i,j)]，其中gcd(i,j)表示i,
j的最大公约数。Doris的表格中共有n×m个数，她想知道这些数的乘积是多少。答案对10^9+7取模。
解题报告:
这题和之前的差不多.设(g[x])为(gcd==x)的数对个数，然后就是基本的反演了
(sum_{d=1}^nf[d]^{g[d]})
(g[d]=sum_{i|d}lfloor m/d
floorlfloor n/d
floormu(d/i))
这个明显可以数论分块，所以我们处理一下(g[d]^mu)的前缀积,对于除法再求个逆即可
#include <algorithm>
#include <iostream>
#include <cstdlib>
#include <cstring>
#include <cstdio>
#include <cmath>
#define RG register
#define il inline
#define iter iterator
#define Max(a,b) ((a)>(b)?(a):(b))
#define Min(a,b) ((a)<(b)?(a):(b))
using namespace std;
typedef long long ll;
const int N=1e6+5,mod=1e9+7;
int q[1005][2],lim,prime[N],num=0,pi[N];bool vis[N];
void solve(){
pi[1]=1;int to;
for(int i=2;i<=lim;i++){
if(!vis[i]){
prime[++num]=i;
pi[i]=-1;
}
for(int j=1;j<=num && i*prime[j]<=lim;j++){
to=i*prime[j];vis[to]=true;
if(i%prime[j])pi[to]=-pi[i];
else{
pi[to]=0;
break;
}
}
}
}
ll f[N],sum[N];
ll qm(ll x,ll k){
ll sum=1;
while(k){
if(k&1)sum*=x,sum%=mod;
x*=x;x%=mod;k>>=1;
}
return sum;
}
ll query(int n,int m){
int j;ll ret=1;
for(int i=1;i<=n;i=j+1){
j=min(n/(n/i),m/(m/i));
ret=ret*qm(sum[j]*qm(sum[i-1],mod-2)%mod,(ll)(n/i)*(m/i))%mod;
}
return ret;
}
ll ni[N],g[N];
void work()
{
int Q;cin>>Q;
for(int i=1;i<=Q;i++){
scanf("%d%d",&q[i][0],&q[i][1]);
if(q[i][0]>q[i][1])swap(q[i][0],q[i][1]);
lim=Max(lim,q[i][0]);
}
solve();
f[0]=0;f[1]=1;
for(int i=2;i<=lim;i++){
f[i]=f[i-1]+f[i-2];
if(f[i]>=mod)f[i]-=mod;
}
for(int i=0;i<=lim;i++)ni[i]=qm(f[i],mod-2);
for(int i=1;i<=lim;i++)g[i]=1;
for(int i=1;i<=lim;i++)
for(int j=i;j<=lim;j+=i)
if(pi[j/i]==1)g[j]=g[j]*f[i]%mod;
else if(pi[j/i]==-1)g[j]=g[j]*ni[i]%mod;
sum[0]=1;
for(int i=1;i<=lim;i++){
sum[i]=sum[i-1]*g[i]%mod;
}
for(int i=1;i<=Q;i++){
printf("%lld
",query(q[i][0],q[i][1]));
}
}
int main()
{
work();
return 0;
}