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[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

HINT

Source

鸣谢infinityedge上传

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

推公式+分块

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

 1 #include<iostream>
 2 #include<algorithm>
 3 #include<cstring>
 4 #include<cstdio>
 5 #include<cmath>
 6 #define LL long long
 7 using namespace std;
 8 const int mxn=1000005;
 9 const int mod=1e9+7;
10 int read(){
11     int x=0,f=1;char ch=getchar();
12     while(ch<'0' || ch>'9'){if(ch=='-')f=-1;ch=getchar();}
13     while(ch>='0' && ch<='9'){x=x*10+ch-'0';ch=getchar();}
14     return x*f;
15 }
16 int pri[mxn],mu[mxn],f[mxn],cnt=0;
17 LL g[mxn];
18 bool vis[mxn];
19 LL ksm(LL a,LL k){
20     LL res=1;
21     while(k){
22         if(k&1)res=res*a%mod;
23         a=a*a%mod;
24         k>>=1;
25     }
26     return res;
27 }
28 void init(){
29     f[1]=1;g[0]=1;mu[1]=1;
30     for(int i=2;i<mxn;i++)f[i]=((LL)f[i-1]+f[i-2])%mod;
31     for(int i=2;i<mxn;i++){
32         if(!vis[i]){pri[++cnt]=i;mu[i]=-1;}
33         for(int j=1;j<=cnt && (LL)pri[j]*i<mxn;j++){
34             int tmp=pri[j]*i;
35             vis[tmp]=1;
36             if(i%pri[j]==0){mu[tmp]=0;break;}
37             mu[tmp]=-mu[i];
38         }
39     }
40     for(int i=1;i<mxn;i++)g[i]=1;
41     for(int i=1;i<mxn;i++){
42         LL invi=ksm(f[i],mod-2);
43         for(int j=i,c=1;j<mxn;j+=i,c++){
44             if(!mu[c])continue;
45             if(mu[c]==1)g[j]=(LL)g[j]*f[i]%mod;
46             else g[j]=(LL)g[j]*invi%mod;
47         }
48     }
49     for(int i=2;i<mxn;i++)g[i]=(LL)g[i]*g[i-1]%mod;
50     return;
51 }
52 void calc(int n,int m){
53     if(n>m)swap(n,m);
54     LL ans=1;
55     for(int i=1,pos;i<=n;i=pos+1){
56         LL x=n/i,y=m/i;
57         pos=min(n/x,m/y);
58         LL inv=ksm(g[i-1],mod-2);
59         ans=ans*ksm(g[pos]*inv%mod,(LL)x*y%(mod-1))%mod;
60     }
61     printf("%lld
",ans);
62     return;
63 }
64 int n,m;
65 int main(){
66 //  freopen("product.in","r",stdin);
67 //  freopen("product.out","w",stdout);
68     int i,j;
69     int T=read();
70     init();
71     while(T--){
72         n=read();m=read();
73         calc(n,m);
74     }
75     return 0;
76 }