[莫比乌斯反演] BZOJ4816: [Sdoi2017]数字表格

来源:互联网 发布:手机底部导航栏源码 编辑:程序博客网 时间:2024/06/07 17:52

挺简单的。。。

i=1nj=1mfib[gcd(i,j)]=dfib[d]nimj[gcd(i,j)==d]

上面这个指数是经典的反演的问题,然后就变成:
dfib[d]kμ(k)ndkmdk

套路变形一下
T=1(d|Tfib[d]μ(Td))nTmT

中间这个f(T)=d|Tfib[d]μ(Td),直接暴力 nlogn 算。
T=1f(T)nTmT

分块求就好了。

这题我跑的很慢,于是网上看了一下发现一个小trick:(其实是我太菜之前不知道)
求一个数列的每个数的逆元和前缀积的逆元,可以先算 S1n,然后O(n) 推。
S1i=S1i+1ai+1  a1i=S1iSi1

#include<cstdio>#include<algorithm>using namespace std;typedef long long LL;const int maxn=1000005, N=1000000, MOD=1e9+7;int p[maxn],mu[maxn];bool vis[maxn];void get_mu(){    mu[1]=1;    for(int i=2;i<=N;i++){        if(!vis[i]) p[++p[0]]=i, mu[i]=-1;        for(int j=1;j<=p[0]&&(LL)i*p[j]<=N;j++){            vis[i*p[j]]=true;            if(i%p[j]==0){ mu[i*p[j]]=0; break; }            mu[i*p[j]]=-mu[i];        }    }}LL Pow(LL a,LL b){    if(b==-1) return Pow(a,MOD-2);    LL res=1;    for(;b;a=a*a%MOD,b>>=1) if(b&1) res=(res*a)%MOD;    return res; }int _test,n,m,fib[maxn],fib_n[maxn],f[maxn],f_n[maxn];LL ans;int f_mul(int L,int R){ return (LL)f[R]*f_n[L-1]%MOD; }int main(){    freopen("bzoj4816.in","r",stdin);    freopen("bzoj4816.out","w",stdout);    get_mu();    fib[0]=0; fib[1]=fib_n[1]=1; for(int i=2;i<=N;i++) fib[i]=(fib[i-1]+fib[i-2])%MOD, fib_n[i]=Pow(fib[i],MOD-2);    for(int i=1;i<=N;i++) f[i]=1;    for(int i=1;i<=N;i++)     for(int j=1;(LL)i*j<=N;j++){        int t=(mu[j]==0?1:(mu[j]==1?fib[i]:fib_n[i]));        f[i*j]=((LL)f[i*j]*t)%MOD;     }    f[0]=f_n[0]=1; for(int i=1;i<=N;i++) f[i]=((LL)f[i]*f[i-1])%MOD, f_n[i]=Pow(f[i],MOD-2);    scanf("%d",&_test);    while(_test--){        scanf("%d%d",&n,&m); int _min=min(n,m);        ans=1;        for(int i=1,nxt;i<=_min;i=nxt+1){            nxt=min(n/(n/i),m/(m/i));            ans=(ans*Pow(f_mul(i,nxt),(LL)(n/i)*(m/i)))%MOD;        }        printf("%lld\n",ans);    }    return 0;}
原创粉丝点击