分治NTT

首先是NTT的板子。

int cnt[N];
void NTT(int a[],int lim,bool type){
	for(int i=0;i<(1<<lim);i++)cnt[i]=(cnt[i>>1]>>1)|((i&1)<<(lim-1));
	for(int i=0;i<(1<<lim);i++)if(i<cnt[i])swap(a[i],a[cnt[i]]);
	for(int i=1;i<(1<<lim);i<<=1){
		int dd=qp(type?3:inv,(mod-1)/i>>1);
		for(int j=0;j<(1<<lim);j+=i<<1){
			int t=1;
			for(int k=0;k<i;k++){
				int x=t*a[i+j+k]%mod;
				a[i+j+k]=(a[j+k]-x+mod)%mod;
				a[j+k]=(a[j+k]+x)%mod;t=t*dd%mod;
			}
		}
	}
	if(type)return;int nowInv=qp(1<<lim,mod-2);
	for(int i=0;i<(1<<lim);i++)a[i]=a[i]*nowInv%mod;
}

然后是分治 NTT。分治 NTT 其实就是一个NTT的板子,然后通过分治的方法(类似cdq分治)来快速计算,复杂度是 O(Nlog2N),也很好理解。

#include<bits/stdc++.h>
//#define feyn
#define int long long
using namespace std;
const int N=1<<17;
const int mod=998244353;
const int inv=(mod+1)/3;
inline void read(int &wh){
    wh=0;int f=1;char w=getchar();
    while(w<'0'||w>'9'){if(w=='-')f=-1;w=getchar();}
    while(w<='9'&&w>='0'){wh=wh*10+w-'0';w=getchar();}
    wh*=f;return;
}
inline int qp(int s1,int s2){
	if(s2==1)return s1;int an=qp(s1,s2>>1);
	if(s2&1)return an*an%mod*s1%mod;else return an*an%mod;
}

int cnt[N];
void NTT(int a[],int lim,bool type){
	for(int i=0;i<(1<<lim);i++)cnt[i]=(cnt[i>>1]>>1)|((i&1)<<(lim-1));
	for(int i=0;i<(1<<lim);i++)if(i<cnt[i])swap(a[i],a[cnt[i]]);
	for(int i=1;i<(1<<lim);i<<=1){
		int dd=qp(type?3:inv,(mod-1)/i>>1);
		for(int j=0;j<(1<<lim);j+=i<<1){
			int t=1;
			for(int k=0;k<i;k++){
				int x=t*a[i+j+k]%mod;
				a[i+j+k]=(a[j+k]-x+mod)%mod;
				a[j+k]=(a[j+k]+x)%mod;t=t*dd%mod;
			}
		}
	}
	if(type)return;int nowInv=qp(1<<lim,mod-2);
	for(int i=0;i<(1<<lim);i++)a[i]=a[i]*nowInv%mod;
}

int m,f[N],g[N],a[N],b[N];
void solve(int l,int r,int lim){
	if(lim<=0)return;int mid=(l+r)>>1;
	solve(l,mid,lim-1);
	for(int i=0;i<(1<<lim);i++)a[i]=(l+i<=mid)?f[l+i]:0;
	for(int i=0;i<(1<<lim);i++)b[i]=g[i];
	NTT(a,lim,true);NTT(b,lim,true);
	for(int i=0;i<(1<<lim);i++)a[i]=a[i]*b[i]%mod;
	NTT(a,lim,false);
	for(int i=mid+1;i<=r;i++)f[i]=(f[i]+a[i-l])%mod;
	solve(mid+1,r,lim-1);
}

signed main(){
	
	#ifdef feyn
	freopen("in.txt","r",stdin);
	#endif
	
	read(m);int n=1;while((1<<n)<m)++n;
	for(int i=1;i<m;i++)read(g[i]);f[0]=1;
	solve(0,(1<<n)-1,n);
	for(int i=0;i<m;i++)printf("%lld ",f[i]);
	
	return 0;
}

然后是两个例题。

P4841

fx 代表大小为 x 的图的答案,gx 代表不要求联通的方案。于是显然有:

fx=i=1x1(x1i)figxi=i=1x1x!fi(i1)!gxi(xi)!

然后就是分治 NTT 的板子了。

#include<bits/stdc++.h>
//#define feyn
#define int long long 
using namespace std;
const int N=1<<18;
const int mod=1004535809;
const int inv=(mod+1)/3;
inline void read(int &wh){
    wh=0;int f=1;char w=getchar();
    while(w<'0'||w>'9'){if(w=='-')f=-1;w=getchar();}
    while(w<='9'&&w>='0'){wh=wh*10+w-'0';w=getchar();}
    wh*=f;return;
}
inline int qp(int s1,int s2){
	if(s2==0)return 1;if(s2==1)return s1;int an=qp(s1,s2>>1);
	if(s2&1)return an*an%mod*s1%mod;else return an*an%mod;
}

int cnt[N];
void NTT(int a[],int lim,bool type){
	for(int i=0;i<(1<<lim);i++)cnt[i]=(cnt[i>>1]>>1)|((i&1)<<(lim-1));
	for(int i=0;i<(1<<lim);i++)if(i<cnt[i])swap(a[i],a[cnt[i]]);
	for(int i=1;i<(1<<lim);i<<=1){
		int dd=qp(type?3:inv,(mod-1)/i>>1);
		for(int j=0;j<(1<<lim);j+=i<<1){
			int t=1;
			for(int k=0;k<i;k++){
				int x=a[i+j+k]*t%mod;
				a[i+j+k]=(a[j+k]-x+mod)%mod;
				a[j+k]=(a[j+k]+x)%mod;t=t*dd%mod;
			}
		}
	}
	if(type)return;int nowInv=qp((1<<lim),mod-2);
	for(int i=0;i<(1<<lim);i++)a[i]=a[i]*nowInv%mod;
}

int m,p[N]={1},q[N]={1},f[N]={1},g[N],a[N],b[N];

void solve(int l,int r,int lim){
	if(lim==0)return f[l]=((qp(2,l*(l-1)/2)-f[l])%mod+mod)%mod,void();
	int mid=(l+r)>>1;solve(l,mid,lim-1);
	
	for(int i=l;i<=r;i++)a[i-l]=i<=mid?f[i]*q[i-1]%mod:0;
	for(int i=0;i<(1<<lim);i++)b[i]=g[i]*q[i]%mod;
	NTT(a,lim,true);NTT(b,lim,true);
	for(int i=0;i<(1<<lim);i++)a[i]=a[i]*b[i]%mod;
	NTT(a,lim,false);
	for(int i=mid+1;i<=r;i++)f[i]=(f[i]+a[i-l]*p[i-1])%mod;
	solve(mid+1,r,lim-1);
}

signed main(){
	
	#ifdef feyn
	freopen("in.txt","r",stdin);
	#endif
	
	read(m);int lim=0;while((1<<lim)<=m)++lim;
	for(int i=1;i<(1<<lim);i++)p[i]=p[i-1]*i%mod,q[i]=qp(p[i],mod-2);
	for(int i=1;i<(1<<lim);i++)g[i]=qp(2,i*(i-1)/2);
	solve(0,(1<<lim)-1,lim);printf("%lld\n",f[m]);
	
	return 0;
}
posted @   Feynn  阅读(373)  评论(0编辑  收藏  举报
相关博文:
阅读排行:
· 10年+ .NET Coder 心语 ── 封装的思维:从隐藏、稳定开始理解其本质意义
· 地球OL攻略 —— 某应届生求职总结
· 提示词工程——AI应用必不可少的技术
· Open-Sora 2.0 重磅开源!
· 周边上新:园子的第一款马克杯温暖上架
点击右上角即可分享
微信分享提示