mobius反演。。。
∏ni=1∏mj=1fi[gcd(i,j)]
∏nk=1fi[k]∑ni=1∑mj=1[gcd(i,j)=k]
设f(d)=∑ni=1∑mj=1[gcd(i,j)=k] ,表示最大公约数为k的数对数
F(d)=?nd???md? 表示公约数为k的数对数
根据莫比乌斯反演的公式f(d)=∑d|nμ(nd)?F(n)
所以式子可以变成
但是实际上式子还能进一步的化简,设T=i?k
设h(T)=∏d|nfi[d]μ(Td),h(T)可以预处理,那么回答询问的时间复杂度就是O(√n+√m)。
#include<iostream> #include<cstdio> #include<cstdlib> #include<cstring> #include<algorithm> using namespace std; #define lon long long #define mod 1000000007 #define N 1000100 lon sum[N],f[N],inv[N],g[N]; int mu[N],pri[N],tot=0,n,m,vis[N]; lon ksm(lon x,int y){ lon res=1; while(y){ if(y&1)res=(res*x)%mod; x=(x*x)%mod;y>>=1; } return res; } void ycl(){ f[0]=0;f[1]=1; for(int i=2;i<N;i++)f[i]=(f[i-1]+f[i-2])%mod; for(int i=1;i<N;i++)inv[i]=ksm(f[i],mod-2); mu[1]=1; for(int i=2;i<N;i++){ if(!vis[i]){ vis[i]=1;pri[++tot]=i; mu[i]=-1; } for(int j=1;((j<=tot)&&(i*pri[j]<N));j++){ vis[i*pri[j]]=1; if(i%pri[j]==0){mu[i*pri[j]]=0;break;} mu[i*pri[j]]=-mu[i]; } } for(int i=1;i<N;i++)g[i]=1; for(int i=1;i<N;i++)//g[i] for(int j=i;j<N;j+=i){ if(mu[j/i]==1)g[j]=(g[j]*f[i])%mod; else if(mu[j/i]==-1)g[j]=(g[j]*inv[i])%mod; } sum[0]=1; for(int i=1;i<N;i++) sum[i]=(sum[i-1]*g[i])%mod; } int main(){ lon ans,res,in; ycl(); int t;scanf("%d",&t); while(t--){ scanf("%d%d",&n,&m); if(n>m)swap(n,m);//n==min(n,m) ans=1;int last; for(int i=1;i<=n;i=last+1){ last=min(n/(n/i),m/(m/i)); in=ksm(sum[i-1],mod-2); res=(sum[last]*in)%mod; ans=(ans*ksm(ksm(res,n/i),m/i))%mod; } cout<<ans<<endl; } return 0; }
原文:http://www.cnblogs.com/harden/p/6701708.html