BZOJ:4816: [Sdoi2017]数字表格
4816: [Sdoi2017]数字表格
Time Limit: 50 Sec Memory Limit: 128 MBSubmit: 501 Solved: 222
[Submit][Status][Discuss]
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
2 3
4 5
6 7
Sample Output
1
6
960
6
960
算是基础数论吧……
然后我就由于过于自信而gg,又WA又T
推一下式子很容易得到
$$ans=\prod_{d=1}^{n}\prod_{i|d}F_i^{\mu(\frac{d}{i})*\left\lfloor\frac{n}{d}\right\rfloor*\left\lfloor\frac{m}{d}\right\rfloor}$$
一开始我就拿这个直接分段上$O(n^{\frac{3}{4}})$然后顺利地T掉了……
把$\prod_{i|d}F_i^{\mu(\frac{d}{i})}$预处理出来就好了。
#include<cstdio> #include<algorithm> #define MN 1000001 using namespace std; int read_p,read_ca; inline int read(){ read_p=0;read_ca=getchar(); while(read_ca<'0'||read_ca>'9') read_ca=getchar(); while(read_ca>='0'&&read_ca<='9') read_p=read_p*10+read_ca-48,read_ca=getchar(); return read_p; } const int MOD=1e9+7; int T,p[MN],num=0,mu[MN],F[MN],I[MN],n,m,mmh,S[MN]; bool bo[MN]; inline void M(int &x){while(x>=MOD)x-=MOD;while(x<0)x+=MOD;} inline void _M(int &x){while(x>=MOD-1)x-=MOD-1;while(x<0)x+=MOD-1;} inline int min(int x,int y){return x<y?x:y;} inline int mi(int x,int y){ int mmh=1; if (y<0) y+=MOD-1; for (;y;x=1LL*x*x%MOD,y>>=1) if (y&1) mmh=1LL*mmh*x%MOD; return mmh; } int main(){ register int i,j,k; mu[1]=1; for (i=2;i<MN;i++){ if (!bo[i]) p[++num]=i,mu[i]=-1; for (j=1;j<=num&&i*p[j]<MN;j++) if (bo[i*p[j]]=1,i%p[j])mu[i*p[j]]=-mu[i];else{mu[i*p[j]]=0;break;} } F[0]=0;F[1]=1;I[0]=I[1]=1;S[1]=1; for (i=2;i<MN;i++) M(F[i]=F[i-1]+F[i-2]),I[i]=mi(F[i],MOD-2),S[i]=1; for (i=1;i<MN;i++) if (mu[i]) for (j=i,k=1;j<MN;j+=i,k++) S[j]=1LL*S[j]*(mu[i]==1?F[k]:I[k])%MOD;S[0]=1; for (i=1;i<MN;i++) S[i]=1LL*S[i]*S[i-1]%MOD; T=read(); while(T--){ n=read();m=read(); if (n>m) swap(n,m); mmh=1; for (i=1;i<=n;i=j+1){ j=min(n/(n/i),m/(m/i)); mmh=1LL*mmh*mi(1LL*S[j]*mi(S[i-1],MOD-2)%MOD,1LL*(n/i)*(m/i)%(MOD-1))%MOD; } printf("%d\n",mmh); } }