zoukankan      html  css  js  c++  java
  • Jetpack[CSACADEMY]

    XXXIII.Jetpack[CSACADEMY]

    我们考虑先通过一些科技求出“一段长度为 \(2i\) 且相邻两位置差的绝对值为 \(1\) 且首尾都是 \(0\) 的序列的数量”,记其为 \(f_i\)

    大约的确可以列出奇奇怪怪的式子表示 \(f_i\) 然后使用奇奇怪怪的可以优化求值的过程,但假如大家对卡特兰数的模型足够熟习的话就应该能看出 \(f_i\) 即为卡特兰数第 \(i\) 项。

    显然,\(f_i\) 不能被直接使用。我们应使用的是“中间再无 \(0\) 的序列”,以避免重复计算。设其为 \(g_i\)

    一种方法是枚举 \(f_i\) 中第一次出现 \(0\) 的位置,然后列出由 \(f\)\(g\) 的卷积转移出 \(f\) 的式子再反推出 \(g\) 来,这样即可使用分治FFT解决。但我们发现,任何一组 \(g_i\),在删去左右两端的 \(1\) 后再全体减一,便可唯一对应到一个 \(f_{i-1}\)。故我们发现 \(g_i\) 即为卡特兰数第 \(i-1\) 项。

    然后我们重新令 \(f_i\) 表示长度为 \(i\) 的合法序列数量。于是有

    \[f_i=f_{i-1}+\sum\limits_{j=1}^{\min(m,i/2)}f_{i-2j}g_j \]

    本式可以直接上分治FFT解决。但是,因为模数是毒瘤的\(10^9+7\),而我的垃圾MTT写的极其丑陋,无法卡过。

    所以我们考虑直接套上XXVI.WD与积木中介绍的分治FFT转多项式求逆的方法,即可将复杂度优化到 \(O(n\log n)\)

    代码:

    #include<bits/stdc++.h>
    using namespace std;
    typedef long long ll;
    const int N=1<<20;
    int n,m,f[N],g[N],fac[N],inv[N];
    int ksm(int x,int y,int mod){
    	int rt=1;
    	for(;y;x=(1ll*x*x)%mod,y>>=1)if(y&1)rt=(1ll*rt*x)%mod;
    	return rt;
    }
    int INV(int x,int y){return ksm(x,y-2,y);}
    int rev[N],lim;
    struct Poly{
    	int mod,invlim;
    	const int G=3;
    	int ksm(int x,int y){
    		int rt=1;
    		for(;y;x=(1ll*x*x)%mod,y>>=1)if(y&1)rt=(1ll*rt*x)%mod;
    		return rt;
    	}
    	void NTT(int *a,int tp){
    		for(int i=0;i<lim;i++)if(i<rev[i])swap(a[i],a[rev[i]]);
    		for(int md=1;md<lim;md<<=1){
    			int rt=ksm(G,(mod-1)/(md<<1));
    			if(tp==-1)rt=ksm(rt,mod-2);
    			for(int stp=md<<1,pos=0;pos<lim;pos+=stp){
    				int w=1;
    				for(int i=0;i<md;i++,w=(1ll*w*rt)%mod){
    					int x=a[pos+i],y=(1ll*w*a[pos+md+i])%mod;
    					a[pos+i]=(x+y)%mod;
    					a[pos+md+i]=(x-y+mod)%mod;
    				}
    			}
    		}
    		if(tp==-1)for(int i=0;i<lim;i++)a[i]=(1ll*a[i]*invlim)%mod;
    	}
    	int A[N],B[N];
    	void mul(int *a,int *b){//using: Array A and B
    		invlim=INV(lim,mod);
    		for(int i=0;i<lim;i++)A[i]=B[i]=0;
    		for(int i=0;i<(lim>>1);i++)A[i]=a[i],B[i]=b[i];
    		NTT(A,1),NTT(B,1);
    		for(int i=0;i<lim;i++)A[i]=1ll*A[i]*B[i]%mod;
    		NTT(A,-1);
    	}
    }poly[3];
    const int mod0=998244353,mod1=1004535809,mod2=469762049;
    const ll mod01=1ll*mod0*mod1;
    const int inv0=INV(mod0,mod1),inv01=INV(mod01%mod2,mod2);
    const int p=1e9+7;
    void MTT(int *a,int *b,int *c,int LG){
    	lim=1<<LG;
    	for(int i=0;i<lim;i++)rev[i]=(rev[i>>1]>>1)|((i&1)<<(LG-1));
    	for(int i=0;i<3;i++)poly[i].mul(a,b);
    	for(int i=0;i<lim;i++){
    		ll x=1ll*(poly[1].A[i]-poly[0].A[i]+mod1)%mod1*inv0%mod1*mod0+poly[0].A[i];
    		c[i]=(1ll*(poly[2].A[i]-x%mod2+mod2)%mod2*inv01%mod2*(mod01%p)%p+x)%p;
    	}
    }
    int C[N];
    void INV(int *a,int *b,int LG){//using: Array C
    	b[0]=INV(a[0],p);
    	for(int k=1;k<=LG+1;k++){
    		MTT(b,a,C,k);
    		for(int i=0;i<(1<<k);i++)C[i]=(p-C[i])%p;
    		(C[0]+=2)%=p;
    		MTT(C,b,b,k);
    	}
    }
    int binom(int x,int y){return 1ll*fac[x]*inv[y]%p*inv[x-y]%p;}
    int main(){
    	scanf("%d%d",&n,&m),m=min(n>>1,m),poly[0].mod=mod0,poly[1].mod=mod1,poly[2].mod=mod2;
    	fac[0]=1;for(int i=1;i<=n;i++)fac[i]=1ll*fac[i-1]*i%p;
    	inv[n]=INV(fac[n],p);for(int i=n-1;i>=0;i--)inv[i]=1ll*inv[i+1]*(i+1)%p;
    	for(int i=0;i<m;i++)g[i+1]=1ll*binom(i<<1,i)*INV(i+1,p)%p;
    	m<<=1;for(int i=m;i;i--)if(i&1)g[i]=0;else g[i]=g[i>>1];
    	for(int i=1;i<=n;i++)(g[i]+=g[i-1])%=p;
    	(g[0]+=p-1)%=p;
    	for(int i=n;i>=0;i--)(g[i+1]+=g[i])%=p,g[i]=(p-g[i])%p;
    	int all=0;
    	while((1<<all)<=n)all++;
    	INV(g,f,all);
    //	for(int i=0;i<=n;i++)printf("%d ",f[i]);puts("");
    //	MTT(f,g,C,all+1);
    //	for(int i=0;i<=n;i++)printf("%d ",C[i]);puts("");
    	printf("%d\n",f[n]);
    	return 0;
    }
    

  • 相关阅读:
    「NOI2018」 你的名字
    「刷题笔记」杂题
    关于~
    「刷题笔记」网络流
    「考试」联赛模拟40-45,晚间小测4-9
    「考试」联赛模拟36-39,noip晚间小测2-3
    「刷题笔记」莫队
    「考试」CSP-S 2020
    「考试」noip模拟9,11,13
    「刷题笔记」概率与期望
  • 原文地址:https://www.cnblogs.com/Troverld/p/14610705.html
Copyright © 2011-2022 走看看