Bzoj4816 [Sdoi2017]数字表格

Time Limit: 50 Sec  Memory Limit: 128 MB
Submit: 646  Solved: 296

Description

Doris刚刚学习了fibonacci数列。用f[i]表示数列的第i项，那么 f[0]=0 f[1]=1 f
=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

HINT

 

Source

鸣谢infinityedge上传

 

数学问题 莫比乌斯反演 分块

推公式+分块

推出来的那个子函数不是积性函数，需要$O(nlogn)$ 筛出来

#include<iostream>
#include<algorithm>
#include<cstring>
#include<cstdio>
#include<cmath>
#define LL long long
using namespace std;
const int mxn=1000005;
const int mod=1e9+7;
int read(){
int x=0,f=1;char ch=getchar();
while(ch<'0' || ch>'9'){if(ch=='-')f=-1;ch=getchar();}
while(ch>='0' && ch<='9'){x=x*10+ch-'0';ch=getchar();}
return x*f;
}
int pri[mxn],mu[mxn],f[mxn],cnt=0;
LL g[mxn];
bool vis[mxn];
LL ksm(LL a,LL k){
LL res=1;
while(k){
if(k&1)res=res*a%mod;
a=a*a%mod;
k>>=1;
}
return res;
}
void init(){
f[1]=1;g[0]=1;mu[1]=1;
for(int i=2;i<mxn;i++)f[i]=((LL)f[i-1]+f[i-2])%mod;
for(int i=2;i<mxn;i++){
if(!vis[i]){pri[++cnt]=i;mu[i]=-1;}
for(int j=1;j<=cnt && (LL)pri[j]*i<mxn;j++){
int tmp=pri[j]*i;
vis[tmp]=1;
if(i%pri[j]==0){mu[tmp]=0;break;}
mu[tmp]=-mu[i];
}
}
for(int i=1;i<mxn;i++)g[i]=1;
for(int i=1;i<mxn;i++){
LL invi=ksm(f[i],mod-2);
for(int j=i,c=1;j<mxn;j+=i,c++){
if(!mu[c])continue;
if(mu[c]==1)g[j]=(LL)g[j]*f[i]%mod;
else g[j]=(LL)g[j]*invi%mod;
}
}
for(int i=2;i<mxn;i++)g[i]=(LL)g[i]*g[i-1]%mod;
return;
}
void calc(int n,int m){
if(n>m)swap(n,m);
LL ans=1;
for(int i=1,pos;i<=n;i=pos+1){
LL x=n/i,y=m/i;
pos=min(n/x,m/y);
LL inv=ksm(g[i-1],mod-2);
ans=ans*ksm(g[pos]*inv%mod,(LL)x*y%(mod-1))%mod;
}
printf("%lld\n",ans);
return;
}
int n,m;
int main(){
//  freopen("product.in","r",stdin);
//  freopen("product.out","w",stdout);
int i,j;
int T=read();
init();
while(T--){
n=read();m=read();
calc(n,m);
}
return 0;
}