[SDOI2017]数字表格

Description:

求 $ \prod_{i=1}^n \prod_{j=1}^m f[gcd(i,j)]$

其中 \(f​\) 表示斐波那契数列

首先显然:

$ Ans = \prod_{d=1}^{min(n,m)} \prod_{i=1}^{ \lfloor \frac{n}{d} \rfloor} \prod _{j=1}^{\lfloor \frac{m}{d} \rfloor } f[d][gcd(i,j)==1] $

$ Ans = \prod_{d=1}{min(n,m)}f[d]^{} \mu(k) * \lfloor \frac{n}{kd} \rfloor \lfloor \frac{m}{kd} \rfloor} ​$

换元得:

\(Ans = \prod_{T=1}^n (\sum_{d|T}f[d]^{\mu(\frac{T}{d})})^{\lfloor \frac{n}{T} \rfloor\lfloor \frac{m}{T} \rfloor}​\)

\(Ans=\prod_{T=1}^n g(T)^{\lfloor \frac{n}{T} \rfloor\lfloor \frac{m}{T} \rfloor}\)

由于\(g(T)\)不能线性筛,故暴力枚举倍数

预处理出\(g(T)\)后,整除分块加快速幂即可

#include<bits/stdc++.h>
using namespace std;
const int mxn=1e6+5,mod=1e9+7;
int T,n,m,tot,mu[mxn],f[mxn],g[mxn],p[mxn],b[mxn],s[mxn],sum[mxn],vis[mxn],inv[mxn];

int qpow(int a,int b) {
    int ans=1,base=a;
    while(b) {
        if(b&1) ans=(1ll*ans*base)%mod;
        base=(1ll*base*base)%mod;
        b>>=1;
    }
    return ans;
}

void sieve(int lim) 
{
    mu[1]=1,f[1]=1; g[1]=1;
    for(int i=2;i<=lim;++i) {
        f[i]=(f[i-1]+f[i-2])%mod;
        g[i]=qpow(f[i],mod-2); s[i]=1;
        if(!vis[i]) p[++tot]=i,mu[i]=-1; 
        for(int j=1;j<=tot&&i*p[j]<=lim;++j) {
            vis[p[j]*i]=1;
            if(i%p[j]==0) break;
            else mu[p[j]*i]=-mu[i];
        }
    }
    for(int i=1;i<=lim;++i) 
        for(int j=i;j<=lim;j+=i)
            if(mu[i]==1) 
                s[j]=1ll*s[j]*f[j/i]%mod;
            else if(mu[i]==-1) 
                s[j]=1ll*s[j]*g[j/i]%mod;		
    inv[1]=g[1]; sum[1]=f[1];	sum[0]=inv[0]=1;		
    for(int i=2;i<=lim;++i) 
        b[i]=qpow(s[i],mod-2),sum[i]=1ll*sum[i-1]*s[i]%mod,inv[i]=1ll*inv[i-1]*b[i]%mod;
}

int main()
{
    scanf("%d",&T); sieve(1000000);
    while(T--) {
        scanf("%d%d",&n,&m); 
        if(n>m) swap(n,m); int ans=1;
        for(int l=1,r;l<=n;l=r+1) {	
            r=min(n/(n/l),m/(m/l));
            ans=(1ll*ans*qpow(1ll*sum[r]*inv[l-1]%mod,1ll*(n/l)*(m/l)%(mod-1))%mod);
        }
        printf("%d\n",ans);
    }
    return 0;
}
posted @ 2019-02-16 19:12  cloud_9  阅读(103)  评论(0编辑  收藏  举报