BZOJ4816: [Sdoi2017]数字表格

题目描述

传送门

题目分析

发现就是要求:

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

其中(f(n))表示斐波那契数列第(n)项。

根据套路，枚举所有的(gcd)

[Ans=prod_{d=1}^{min(n,m)}prod_{i=1}^{n}prod_{j=1}^{m}[gcd(i,j)=d]f(d)+[gcd(i,j)!=d] ]

然后直接把(f(d))提出来

[Ans=prod_{d=1}^{n}f[d]^{sum_{i=1}^{n/d}sum_{j=1}^{m/d}[gcd(i,j)=1]} ]

然后着手开始化上面的式子：

[egin{align} & sum_{i=1}^{lfloorfrac{n}{d} floor}sum_{j=1}^{lfloorfrac{n}{d} floor}[gcd(i,j)=1]\ &=sum_{i=1}^{lfloorfrac{n}{d} floor}sum_{j=1}^{lfloorfrac{n}{d} floor}sum_{xmid gcd(i,j)}mu(x)\ &=sum_{x=1}^{min(lfloorfrac{n}{d} floor,lfloorfrac{m}{d} floor)}mu(x)sum_{i=1}^{lfloorfrac{n}{d} floor}sum_{j=1}^{lfloorfrac{n}{d} floor}[xmid i&&xmid j]\ &=sum_{x=1}^{min(lfloorfrac{n}{d} floor,lfloorfrac{m}{d} floor)}mu(x)sum_{i=1}^{lfloorfrac{n}{dx} floor}sum_{j=1}^{lfloorfrac{n}{dx} floor}1\ &=sum_{x=1}^{min(lfloorfrac{n}{d} floor,lfloorfrac{m}{d} floor)}mu(x)lfloorfrac{n}{dx} floorlfloorfrac{m}{dx} floor end{align} ]

然后令(T=dx)然后直接从整个式子中提出(T)

[Ans=prod_{T=1}^{n}prod_{dmid{T}}f(d)^{[n/T][m/T]mu(T/d)} ]

[ =prod_{T=1}^{n}(prod_{dmid{T}}f(d)^{mu(T/d)})^{[n/T][m/T]} ]

因为(n)足够小，所以这个里面式子可以直接暴力求出来。。。
然后外层整数分块，求解即可。

是代码呢

#include<iostream>
#include<cstdio>
#include<cstring>
#include<algorithm>
using namespace std;
const int MAXN=1e6+7;
#define ll long long
const int mo=1e9+7;
ll sumf[MAXN],fr[MAXN],f[MAXN],fi[MAXN][3],sum[MAXN],ans;
int mu[MAXN],prime[MAXN],n,m,mul[MAXN];
bool vis[MAXN];
inline ll pow(ll x,ll k)
{
    ll cnt=1;
    while(k){
        if(k&1) cnt*=x;
        cnt%=mo;x=x*x%mo;k>>=1;
    }
    return cnt;
}
inline void get_mu(int N)
{
    mu[1]=1;f[1]=1;
    for(int i=2;i<=N;i++){
        f[i]=(f[i-1]+f[i-2])%mo;
        if(!vis[i]){
            mu[i]=-1;
            prime[++prime[0]]=i;
        }
        for(int j=1;j<=prime[0];j++){
            if(i*prime[j]>N) break;
            vis[i*prime[j]]=1;
            if(i%prime[j]==0) break;
            mu[i*prime[j]]=-mu[i];
        }
    }
    for(int i=1;i<=N;i++) fi[i][0]=pow(f[i],mo-2),fi[i][1]=1,fi[i][2]=f[i];
    for(int i=1;i<=N;i++) mul[i]=1;
    for(int i=1;i<=N;i++)
        for(int j=i;j<=N;j+=i) mul[j]=1ll*mul[j]*fi[i][mu[j/i]+1]%mo;
    sumf[0]=1;
    for(int i=1;i<=N;i++) sumf[i]=1ll*sumf[i-1]*mul[i]%mo;
    fr[N]=pow(sumf[N],mo-2);
    for(int i=N-1;~i;i--) fr[i]=1ll*fr[i+1]*mul[i+1]%mo;
}
inline int read()
{
    int x=0,c=1;
    char ch=' ';
    while((ch<'0'||ch>'9')&&ch!='-')ch=getchar();
    while(ch=='-')c*=-1,ch=getchar();
    while(ch>='0'&&ch<='9')x=x*10+ch-'0',ch=getchar();
    return x*c;
}
int main()
{
    int T=read();
    get_mu(1000000);	
    while(T--){
        n=read(),m=read();
        ans=1;
        int mx=min(n,m);
        for(int l=1,r;l<=mx;l=r+1){
            r=min(n/(n/l),m/(m/l));
            ans=(1ll*ans*pow(1ll*sumf[r]*fr[l-1]%mo,1ll*(n/l)*(m/l)%(mo-1)))%mo;
        }
        printf("%lld
",ans);
    }
}