[bzoj4816] [Sdoi2017]数字表格

来自FallDream的博客，未经允许，请勿转载，谢谢。

---

题意：T组数据，求$\prod_{i=1}^{n}\prod_{j=1}^{m}f[gcd(i,j)]$，其中$f[i]$是斐波那契数列的第i项。T<=1000 n,m<=10^6

题解：让n<=m,考虑枚举斐波那契数列每一项被计算了多少次，列出式子

$$Ans=\prod_{k=1}^{n}f[k]^{\sum_{i=1}^{\lfloor\frac{n}{k}\rfloor}\sum_{m=1}^{\lfloor\frac{n}{k}\rfloor}\sum_{d|i,d|j}\mu(d)}$$

老套路，把miu拆出来

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

令$T=kd$，得到

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

发现$\prod_{d|T}{f[d]^{\mu(\frac{T}{d})}}$可以预处理然后前缀积，而$\lfloor\frac{n}{T}\rfloor$和$\lfloor\frac{m}{T}\rfloor$只有根号种取值，所以这道题就做完了

复杂度$O(T\sqrt{n}logn)$

#include<cstdio>
#include<iostream>
#define ll long long
#define MN 1000000
#define mod 1000000007
using namespace std;
inline int read()
{
    int x=0;char ch=getchar();
    while(ch<'0'||ch>'9')ch=getchar();
    while(ch>='0'&&ch<='9'){x=x*10+ch-'0'; ch=getchar();}
    return x;
}

int mu[MN+5],f[MN+5],g[MN+5],s[MN],num=0,n,m,ans,inv[MN+5];
bool b[MN+5];

int pow(int x,int k)
{
    if(k==-1) return inv[x];
    int sum=1;
    for(;k;k>>=1,x=1LL*x*x%mod)
        if(k&1) sum=1LL*sum*x%mod;
    return sum; 
}

int main()
{    
    mu[1]=f[1]=g[1]=g[0]=1;
    for(int i=2;i<=MN;i++)f[i]=(f[i-1]+f[i-2])%mod,g[i]=1;
    for(int i=2;i<=MN;i++)
    {
        if(!b[i]) s[++num]=i,mu[i]=-1;
        for(int j=1;s[j]*i<=MN;j++)
        {
            b[s[j]*i]=1;
            if(i%s[j]==0) break;
            else mu[s[j]*i]=-mu[i];
        }
    }

    for(int i=1;i<=MN;i++) inv[i]=pow(f[i],mod-2);    
    for(int i=1;i<=MN;i++)
        for(int j=i;j<=MN;j+=i)
            g[j]=1LL*g[j]*(mu[j/i]==-1?inv[i]:pow(f[i],mu[j/i]))%mod;
    for(int i=2;i<=MN;i++) g[i]=1LL*g[i]*g[i-1]%mod;
    for(int T=read();T;T--)
    {
        n=read();m=read();ans=1;
        if(n>m)swap(n,m);
        for(int i=1,last;i<=n;i=last+1)
        {
            last=min(n/(n/i),m/(m/i));
            ans=1LL*ans*pow(1LL*g[last]*pow(g[i-1],mod-2)%mod,1LL*(n/i)*(m/i)%(mod-1))%mod;
        } 
        printf("%d\n",ans);
    }
    return 0;
}