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取模。
Input
有多组测试数据。
第一个一个数T，表示数据组数。
接下来T行，每行两个数n,m
T<=1000,1<=n,m<=10^6
Output
输出T行，第i行的数是第i组数据的结果
Sample Input
3
2 3
4 5
6 7
Sample Output
1
6
960
Solution
莫比乌斯反演以及作死地推式子，开始：
[ans=prod_{i=1}^nprod_{j=1}^mf[gcd(i,j)]
]
[displaystyle =prod_{d=1}^nf[d]^{sum_{i=1}^{lfloor frac{n}{d}
floor}mu(i)lfloor frac{n}{id}
floor lfloor frac{m}{id}
floor}
]
套路设(T=id)
[ans=prod_{d=1}^nf[d]^{sum_{d|T}^nmu(frac{T}{d})lfloor frac{n}{T}
floorlfloor frac{m}{T}
floor}
]
[=prod_{d=1}^nprod_{d|T}^nf[d]^{mu(frac{T}{d})lfloor frac{n}{T}
floorlfloor frac{m}{T}
floor}
]
[=prod_{T=1}^nprod_{d|T}f[d]^{mu(frac{T}{d})lfloor frac{n}{T}
floorlfloor frac{m}{T}
floor}
]
[=prod_{T=1}^n(prod_{d|T}f[d]^{mu(frac{T}{d})})^{lfloor frac{n}{T}
floorlfloor frac{m}{T}
floor}
]
从第二步到第三步，枚举方式之所以能够改变是因为两者的实际意义是一样的。
想一想是为什么。
(prod_{d=1}^nprod_{d|T}^n)含义是枚举数对((d,T))，保证(T)是(d)的倍数且(T)的值不超过(n)
(prod_{T=1}^nprod_{d|T})难道不是一样的吗？
回到题目，只要括号里面的东西能求前缀和（前缀积？。。）括号外面的东西显然可以整除分块
那么就考虑括号里面的东西
要考虑吗。。
线性筛筛完后暴力求。。。复杂度(O(nlnn))可以接受
然后？
就做完了
#include<bits/stdc++.h>
#define ll long long
const int MAXN=1000000+10,Mod=1e9+7;
int cnt,prime[MAXN],vis[MAXN],mu[MAXN];
ll g[MAXN],f[MAXN],fn[MAXN],F[MAXN];
template<typename T> inline void read(T &x)
{
T data=0,w=1;
char ch=0;
while(ch!='-'&&(ch<'0'||ch>'9'))ch=getchar();
if(ch=='-')w=-1,ch=getchar();
while(ch>='0'&&ch<='9')data=((T)data<<3)+((T)data<<1)+(ch^'0'),ch=getchar();
x=data*w;
}
template<typename T> inline void write(T x,char c=' ')
{
if(x<0)putchar('-'),x=-x;
if(x>9)write(x/10);
putchar(x%10+'0');
if(c!=' ')putchar(c);
}
template<typename T> inline void chkmin(T &x,T y){x=(y<x?y:x);}
template<typename T> inline void chkmax(T &x,T y){x=(y>x?y:x);}
template<typename T> inline T min(T x,T y){return x<y?x:y;}
template<typename T> inline T max(T x,T y){return x>y?x:y;}
inline ll qexp(ll a,ll b)
{
ll res=1;
while(b)
{
if(b&1)res=res*a%Mod;
a=a*a%Mod;
b>>=1;
}
return res;
}
inline void init()
{
memset(vis,1,sizeof(vis));
vis[0]=vis[1]=0;
mu[1]=1;
g[0]=g[1]=1;
for(register int i=2;i<MAXN;++i)
{
g[i]=1;
if(vis[i])
{
prime[++cnt]=i;
mu[i]=-1;
}
for(register int j=1;j<=cnt&&i*prime[j]<MAXN;++j)
{
vis[i*prime[j]]=0;
if(i%prime[j])mu[i*prime[j]]=-mu[i];
else break;
}
}
f[1]=fn[1]=1;
for(register int i=2;i<MAXN;++i)
{
f[i]=(f[i-1]+f[i-2])%Mod;
fn[i]=qexp(f[i],Mod-2);
}
for(register int i=1;i<MAXN;++i)
if(mu[i]!=0)
for(register int j=i;j<MAXN;j+=i)(g[j]*=(mu[i]==1?f[j/i]:fn[j/i]))%=Mod;
for(register int i=1;i<MAXN;++i)(g[i]*=g[i-1])%=Mod;
}
inline ll solve(ll n,ll m)
{
ll res=1;
if(n>m)std::swap(n,m);
for(register int i=1;;)
{
if(i>n)break;
int j=min(n/(n/i),m/(m/i));
(res*=qexp(g[j]*qexp(g[i-1],Mod-2)%Mod,(n/i)*(m/i)%(Mod-1)))%=Mod;
i=j+1;
}
return (res+Mod)%Mod;
}
int main()
{
init();
int T;
read(T);
while(T--)
{
ll n,m;
read(n);read(m);
write(solve(n,m),'
');
}
return 0;
}