反演题,推式子么=。=
$\prod\limits_{d=1}^{min(n,m)}\prod\limits_{i=1}^n\prod\limits_{j=1}^m[gcd(i,j)==d]fib[d]$
把$fib[d]$前提,前面的连乘就跑到指数上去了
$\prod\limits_{d=1}^{min(n,m)}fib[d]^{\sum\limits_{i=1}^n\sum\limits_{j=1}^m[gcd(i,j)==d]}$
开始反演那坨指数,等等这玩意不是做过么=。=
$\sum\limits_{i=1}^n\sum\limits_{j=1}^m[gcd(i,j)==d]$
$\sum\limits_{i=1}^{\left\lfloor\frac{n}{d}\right\rfloor}\sum\limits_{j=1}^{\left\lfloor\frac{m}{d}\right\rfloor}[gcd(i,j)==1]$
$\sum\limits_{i=1}^{min(\left\lfloor\frac{n}{d}\right\rfloor,\left\lfloor\frac{m}{d}\right\rfloor)}μ(i)\left\lfloor\frac{n}{id}\right\rfloor\left\lfloor\frac{m}{id}\right\rfloor$
于是把$id$捉出来,在原来的整个式子里枚举$id$(不是那个$id$,都懂)
$\prod\limits_{k=1}^{min(n,m)}(\prod\limits_{d|k}fib[d]^{μ(\frac{k}{d})})^{\left\lfloor\frac{n}{k}\right\rfloor\left\lfloor\frac{m}{k}\right\rfloor}$
停,可以做了
对于$\prod\limits_{d|k}fib[d]^{μ(\frac{k}{d})}$,预处理,大力把每个数乘到倍数上去,复杂度$O(n\log n)$
对于$\left\lfloor\frac{n}{k}\right\rfloor\left\lfloor\frac{m}{k}\right\rfloor$这个指数,可以数论分块,这样再加个快速幂每次回答复杂度就是$O(\sqrt n\log mod)$了,可能有点卡常?我倒是一次过了
1 #include<cstdio> 2 #include<cstring> 3 #include<algorithm> 4 using namespace std; 5 const int N=1000005,M=1000000,mod=1e9+7; 6 int fib[N],ifb[N],pfb[N],ipf[N]; 7 int pri[N],npr[N],mul[N]; 8 int T,n,m,x,y,mn,cnt,ans; 9 int qpow(int x,int k) 10 { 11 if(k==1) return x; 12 int tmp=qpow(x,k/2); 13 return k%2?1ll*tmp*tmp%mod*x%mod:1ll*tmp*tmp%mod; 14 } 15 void exGCD(int a,int b,int &x,int &y) 16 { 17 if(!b) {x=1,y=0; return ;} 18 exGCD(b,a%b,y,x); y-=a/b*x; 19 } 20 int Inv(int b) 21 { 22 exGCD(b,mod,x,y); 23 return (x+mod)%mod; 24 } 25 void prework() 26 { 27 register int i,j; 28 npr[1]=true,mul[1]=1,fib[1]=1,ifb[1]=1,pfb[1]=1; 29 for(i=2;i<=M;i++) 30 { 31 fib[i]=(fib[i-1]+fib[i-2])%mod; 32 ifb[i]=Inv(fib[i]),pfb[i]=1; 33 if(!npr[i]) 34 pri[++cnt]=i,mul[i]=-1; 35 for(j=1;j<=cnt&&i*pri[j]<=M;j++) 36 { 37 npr[i*pri[j]]=true; 38 if(i%pri[j]==0) break; 39 else mul[i*pri[j]]=-mul[i]; 40 } 41 } 42 for(i=1;i<=M;i++) 43 for(j=i;j<=M;j+=i) 44 if(mul[j/i]) pfb[j]=1ll*pfb[j]*((~mul[j/i])?fib[i]:ifb[i])%mod; 45 pfb[0]=ipf[0]=1; 46 for(i=1;i<=M;i++) 47 { 48 ipf[i]=Inv(pfb[i]); 49 pfb[i]=1ll*pfb[i-1]*pfb[i]%mod; 50 ipf[i]=1ll*ipf[i-1]*ipf[i]%mod; 51 } 52 } 53 int main() 54 { 55 register int i,j; 56 scanf("%d",&T),prework(); 57 while(T--) 58 { 59 scanf("%d%d",&n,&m); 60 mn=min(n,m),ans=1; 61 for(i=1;i<=mn;i=j+1) 62 { 63 j=min(n/(n/i),m/(m/i)); 64 ans=1ll*ans*qpow(1ll*pfb[j]*ipf[i-1]%mod,1ll*(n/i)*(m/i)%(mod-1))%mod; 65 } 66 printf("%d\n",ans); 67 } 68 return 0; 69 }