【BZOJ 4816】 4816: [Sdoi2017]数字表格 （莫比乌斯）

4816: [Sdoi2017]数字表格

Time Limit: 50 Sec  Memory Limit: 128 MB
Submit: 666  Solved: 312

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上传

 

 

【分析】

　　额。。。加法变乘法有时候还是忍不住求和。。

　　推式子。。

　　$$Ans=\Pi \Pi f(gcd(i,j))$$

　　$$=\Pi f(d)^{\sum 1 [gcd(i,j)==d]}$$

　　$$=\Pi f(d)^{\sum_{i=1}^{n/d}\sum_{j=1}^{m/d} 1 [gcd(i,j)==1]}$$

　　$$=\Pi f(d)^{\sum \mu(d')*\lfloor\dfrac{n}{d*d'}\rfloor \lfloor\dfrac{m}{d*d'}\rfloor}$$

　　【啊好辛苦

　　设D=d*d'

　　$$Ans=\Pi_{D}\Pi_{d|D} f(d)^{\lfloor\dfrac{n}{D}\rfloor \lfloor\dfrac{m}{D}\rfloor*\mu(\dfrac{D}{d})}$$

　　$$=\Pi \Pi (f(d)^{\mu(D/d)})^{\lfloor\dfrac{n}{D}\rfloor \lfloor\dfrac{m}{D}\rfloor}$$

　　设$g(D)=\Pi_{d|D}f(d)^{\mu(D/d)}$

　　这个$O(nlogn)$预处理，后面的D在询问时$\sqrt n $分块即可。

　　因为有快速幂，所以是$O(nlogn)+T*\sqrt n*logn$

 

 1 #include<cstdio>
 2 #include<cstdlib>
 3 #include<cstring>
 4 #include<iostream>
 5 #include<algorithm>
 6 using namespace std;
 7 #define Maxn 1000010
 8 #define Mod 1000000007
 9 
10 int g[Maxn],pri[Maxn],pl,mu[Maxn],f[Maxn];
11 int inv[Maxn],sm[Maxn];
12 bool vis[Maxn];
13 
14 int qpow(int x,int b)
15 {
16     int ans=1;
17     while(b)
18     {
19         if(b&1) ans=1LL*ans*x%Mod;
20         x=1LL*x*x%Mod;
21         b>>=1;
22     }
23     return ans;
24 }
25 
26 void init()
27 {
28     f[0]=0;f[1]=1;for(int i=2;i<=Maxn-10;i++) f[i]=(f[i-1]+f[i-2])%Mod;
29     for(int i=1;i<=Maxn-10;i++) inv[i]=qpow(f[i],Mod-2);
30     memset(vis,0,sizeof(vis));
31     mu[1]=1;
32     for(int i=2;i<=Maxn-10;i++)
33     {
34         if(!vis[i]) pri[++pl]=i,mu[i]=-1;
35         for(int j=1;j<=pl;j++)
36         {
37             if(i*pri[j]>Maxn-10) break;
38             vis[i*pri[j]]=1;
39             if(i%pri[j]==0) {mu[i*pri[j]]=0;break;}
40             mu[i*pri[j]]=-mu[i];
41         }
42     }
43     for(int i=0;i<=Maxn-10;i++) g[i]=1;
44     sm[0]=1;
45     int i;
46     for(i=1;i<=Maxn-10;i++)
47     {
48         for(int j=i;j<=Maxn-10;j+=i) if(mu[j/i]!=0)
49         {
50             if(mu[j/i]==-1) g[j]=1LL*g[j]*inv[i]%Mod;
51             else g[j]=1LL*g[j]*f[i]%Mod;
52         }
53         sm[i]=1LL*sm[i-1]*g[i]%Mod;
54     }
55 }
56 
57 
58 int main()
59 {
60     init();
61     int T;
62     scanf("%d",&T);
63     while(T--)
64     {
65         int n,m;
66         scanf("%d%d",&n,&m);
67         if(n>m) swap(n,m);
68         int ans=1;
69         for(int i=1;i<=n;)
70         {
71             int x=n/i,y=m/i,r1=n/x,r2=m/y;
72             r1=min(r1,r2);
73             ans=1LL*ans*qpow(1LL*sm[r1]*qpow(sm[i-1],Mod-2)%Mod,1LL*x*y%(Mod-1))%Mod;
74             i=r1+1;
75         }
76         printf("%d\n",ans);
77     }
78     return 0;
79 }
View Code

2017-04-26 20:21:03

posted @ 2017-04-26 20:16  konjak魔芋  阅读(266)  评论(0)  编辑  收藏  举报